Recovering the Magnetic Image of Mars from Satellite Observations

One of the possible approaches to reconstructing the map of the distribution of magnetization parameters in the crust of Mars from the data of the Mars MAVEN orbiter mission is considered. Possible ways of increasing the accuracy of reconstruction of the magnetic image of Mars are discussed.


Introduction
The study of the magnetic fields of the planets is one of the ways to obtain information about the internal structure of the planet and its evolution. Access to information about the magnetic fields of planets became possible due to the emergence and development of interplanetary missions. The first results of measuring the magnetic field of Mars were obtained in 1965 by the US mission Mariner 4, which discovered the absence of the global magnetosphere of Mars. In the 1970s, Soviet Mars 2, 3, 5 and later Phobos 2 detected a rather weak magnetic field of about 60 nT in the vicinity of the equator and 120 nT in the vicinity of the pole [1]. In 1996, the MGS (Mars Global Surveyor) mission with a magnetometer reflectometer onboard was sent to Mars [2]. Thanks to this device, data were obtained on the magnetic field of Mars at various altitudes above the surface of Mars. These data made it possible to solve inverse problems of restoration of such physical parameters as magnetization [3,4]. It was found that the crust of Mars is quite strongly magnetized in some areas. From this it was concluded that, although Mars now does not have a global magnetic field, it may have had an active magnetic dynamo earlier [5]. In this regard, modeling of the remanent magnetic field is an important task for studying the deep structure of Mars and for testing models of the magnetic dynamo of Mars in the past, as stated in the literature, see [6][7][8][9][10][11][12][13][14][15]. In the latter case, this allows the creators of magnetic dynamo models to verify their models and examine what happens after it disappears,whether the distribution of remanent magnetization predicted by these models in the crust of Mars coincides with the values observed today. The experimental data of the MGS mission made it possible to obtain a primary estimate of the global distribution of magnetic field sources in the planet's crust [6,16].
The most recent research is based on experimental observation data [17] of the MAVEN mission (NASA's Mars MAVEN orbiter [18] that began work in 2014). A variety of methods developed for studying the Earth's magnetic field can be applied to studying the magnetic field of Mars (see, for example, [3,4]). Previously widespread approaches to modeling local (the intensely magnetized southern highlands with strong magnetic anomalies reaching 1500 nT at 200 km altitude) [7][8][9] and global regions magnetic field distributions induced by remanent magnetism of the Martian crust were divided into the "spherical harmonic expansions" models [10,11] and the "equivalent source" models [12,13]. Current crustal magnetic field models [14,15] consider the MAVEN magnetometer dataset at an altitude of about 135 km [2]. This is due to the fact that the data obtained by the abovementioned MGS magnetometer are mainly distributed in an orbit with an altitude of 370-430 km, some data are available at an altitude of 90-170 km. Therefore, in these studies, only a part of the data from the MAVEN mission was used in order to compare the results obtained on the basis of the data of the MAVEN magnetometer with the results obtained earlier on the basis of the data of the MGS magnetometer.
The structure of this work is as follows. In Section 2, we describe a model based on restoring magnetization parameters from traditional magnetic data that are total magnetic intensity (TMI) [3,4]. Further, in Section 3 we use this model to reconstruct the equivalent distribution of magnetization parameters in the near-surface layer (crust) of Mars. This approach is similar to the approach proposed in [19]-the method of "sweeping" the sources of fields from a multidimensional region of space to its boundary. On the example of studying the gravitational field of Mars, it was shown [20,21] that the distribution of two-dimensional sources on a plane under the surface of Mars resembles the field itself in outline. In [20], a simulation of anomalous gravity masses distribution in the sourthwestern part of Elysium Planitia under the surface was performed on the base of the modified S-approximation method. The same technique was used to restore the equivalent magnetic sources in [22]. The results obtained in these works are based on the so-called linear integral representation method proposed in [23,24] and then evolved in [25][26][27][28][29][30]. Therefore, the approach we are considering to restore the equivalent distribution of magnetic field sources is encouraging. The proposed algorithm was used to process the magnetic data of the MAVEN mission at all altitudes. The obtained results demonstrate the effectiveness of the proposed approach. Section 4 discusses the possibilities of increasing the accuracy of the retrieval of the unknown magnetization parameters in the case of future interplanetary missions using instruments that allow measuring magnetic gradient tensor data [31][32][33].

General Problem Statement
The equation describing the magnetic field induction B f ield induced by an object with a magnetization distribution M(r) and a localization in a domain V has the form [34] where |r − r s | = (x − x s ) 2 + (y − y s ) 2 + (z − z s ) 2 is a distance between the point r s = (x s , y s , z s ), at which the triaxial sensor s measuring magnetic field induction B f ield is located, and the point r = (x, y, z) of the domain V, at which a magnetic source with a total magnetic moment per unit volume M(r) is placed, µ 0 is a permeability in vacuum. The inverse problem is to determine the magnetic moment density M(r), r ∈ V, by using measured magnetic field induction B f ield (r s ) at points r s , s = 1, S.

Remark 1.
Note that this problem statement is physically defined. In the present formulation, it is necessary to restore one vector function from the results of measurements of one vector function as well. Taking into account the fact that each component of a vector function is a scalar function, it requires reconstructing three scalar functions from the results of measurements of three scalar functions. This formulation leads to a system of three scalar equations with three unknown functions.

Reconstruction of the Equivalent Distribution of Magnetization Parameters in the Near-Surface Layer (Crust) of Mars
For convenience, we rewrite Equation (1) as where Here, for shorthand, we use the notation: In addition to the default Cartesian coordinate system, we also introduce a spherical coordinate system (see Figure 1): Here, R is the average radius of Mars. Based on the assumption that the region V is a near-surface spherical layer of Mars of depth h and taking into account (3), Equation (2) can be rewritten as:

Remark 2.
Here B ≡ B f ield and K ≡ K TMI . The introduction of additional notation is due to the fact that the matrix K (and, accordingly, the vector-function B) can be extended by taking into account additional physical information in the formulation of the problem. This will be done at the end of Section 4.

Let us introduce grids
Furthermore, for simplicity, we will proceed from the assumption that the spherical layer is thin enough. As a consequence, with respect to the variable ρ we introduce a grid with only one node ρ h = R − h 2 . As a result, approximating the integrals in (4) using the midpoint rule, we obtain We take into account that (1) measurements are made for all s = 1, S, (2) B and M re vector-functions (in particular, M = M x i + M y j + M z k). Thus, we obtain a system with 3 × S equations (which correspond to the measurement of three components of the vectorfunction B at S points) with 3 × N ϕ × N θ unknowns (which correspond to the grid values of the three components of the vector-function M on the introduced grid Φ N ϕ × Θ N θ ).
To these equations we add the following natural physical conditions.

1.
Matching condition along one of the meridians. This condition means that the magnetic image must be 2π-periodic in the variable ϕ: These conditions give additional 3 × N θ equations.

2.
Matching condition at the South Pole. This means that all the grid values of the components of vector-function M must match when θ = θ N θ : These conditions give additional 3 × (N ϕ − 1) equations.

3.
Matching condition at the North Pole.
Thus, taking into account additional physical conditions, we obtain a system of linear algebraic equations consisting of 3S + 3N θ + 6N ϕ − 6 equations, each of which contains 3N ϕ N θ unknown (the grid values of the components of the vector-function M). Note that the resulting system for a sufficiently large S is overdetermined only from the mathematical point of view. From a physical point of view, the resulting system is still definite.
This system of linear algebraic equations can be written in matrix form Here, the vector M contains the grid values of three components of the unknown vector function M, the first 3S components of the vector on the right side of B contain the results of experimental measurements of three components of the vector function B, and subsequent components of the vector B are zeros.

Using a Regularizing Algorithm
When processing experimental data, instead of exactly known vector B and matrix A their approximate values B δ and A h are usually known, such that B δ − B E δ, A − A h E→E h. This is due both to errors in the measurements of the magnetic field (they add an error in the values of the components of the vector B) and to errors in the accuracy of determining the position of the satellite relative to Mars (they add an error in the values of the components of the matrix A). Under the described conditions, the problem is ill-posed. To solve it, it is necessary to construct a regularizing algorithm. We will use an algorithm based on minimization of the Tikhonov functional [35] For any α > 0 there is a unique extremal of the Tikhonov functional M α η , η = {δ, h}, realizing the minimum of F α [M]. The algorithm of the generalized residual principle [35] can be used to select the regularization parameter. Then, choosing the parameter α = α(η) according to the generalized residual principle The method of conjugate gradients is used as a method for minimizing the Tikhonov functional.

Using the Conjugate Gradient Method
Let M (k) be a minimizing sequence, p (k) and q (k) be auxiliary vectors, and p (0) = 0, M (1) -any admissible point. Then, the formulas of the conjugate gradient method for finding the element M (3N ϕ N θ ) , realizing the minimum of the functional (6) can be written as

Remark 4.
It should be noted that during numerical experiments we define α = 0, M (1) = 0 and use the iteration number k as the regularization parameter. In this case, the criterion for terminating the iterative process is consistent with the error in specifying the input data by means of the condition [36] When numerically searching for the minimum of the functional (6) it is possible to effectively use multiprocessor systems, the specifics of working with these are detailed in [37][38][39].

Experimental Data Processing Results
The data [17] of observations of NASA's Mars MAVEN orbiter [18] were taken as experimental data. The observation area was a set of points at which the magnetic induction of Mars was measured on the 135th day of 2020. Every 32nd point was taken from the data file [17], starting with the first and ending with the 86400th (see Figure 2). Thus, S = 2700. For calculations, grids with the number of intervals N ϕ = 200 and N θ = 200 were used. The grid data made it possible to reconstruct the qualitative distribution of the magnetization parameters over the surface of Mars. The calculations used the values for the average radius of Mars R = 3,389,500 (m) and the thickness of the near-surface layer h = 1000 (m). The results are shown in Figure 3 in polar coordinates. Note that changes in the value of h will lead to a proportional change in the values of the components of the vector M, but the picture of its normalized value shown in Figure 3 will remain unchanged.

Discussion
In recent years, with the development of advanced technology, acquisition of the full tensor gradient magnetic data became available. Much research has shown the advantages of magnetic gradient tensor (MGT) surveys as compared to the conventional total magnetic intensity (TMI) surveys [40][41][42][43][44]. The main conclusion is that better inversion results can be obtained with full MGT data. Thus, we see this approach as promising in the future, since the use of solely magnetic data does not always give satisfactory results. Some work notes their insufficiency and the importance of using data on the magnetic field at low altitudes, which makes it possible to increase the accuracy of retrieving the desired magnetic parameters when using the TMI model.
We define full tensor magnetic gradient B tensor , which unlike the magnetic induction B f ield (that has only three components) has nine components and can be written in the following matrix form [31]: actually, we have only five different components of the tensor matrix.
The diagonal elements and non-diagonal elements of tensor matrix B tensor have the form (details of the derivation of these formulas are presented in [31]): Here, for each character from the set {i, j} one of the characters must be substituted among {x, y, z}.
By adding these equations to Equation (2), we obtain the following system of two 3D Fredholm integral equations of the 1st kind: where Thus, using the full magnetic gradient tensor in Equation (4) the structure of B and K will change. Now B ≡ [B f ield B tenzor ] T and K ≡ [K TMI K MGT ] T .
As a result, we obtain a physically overdetermined problem, which consists of determining three unknown scalar functions from the results of experimental measurements of the other eight scalar functions. The accuracy of recovering the unknown functions from this formulation of the inverse problem will be much higher [31].

Conclusions
The method developed by the authors for determining the magnetic properties of rocks from satellite data can be successfully applied in the process of complex interpretation of data (the so-called joint inversion) on the physical fields of Mars and its topography. The distribution of the magnetization parameters equivalent to the external field, found using the algorithm presented in the article, makes it possible to recreate a qualitative picture of the internal structure of the planet [45][46][47][48][49]. This is due to the fact that the magnetic permeability and magnetic susceptibility are closely related to such magnetic distributions. In the future, it is planned to test the methodology when solving problems in various statements, both simpler and more complicated. It is very important to emphasize that we are processing a dataset all over Mars, i.e., we solve the problem of interpretation in the global version. The method described in the article provides a high quality interpretation of the data on the magnetic field of Mars, which were obtained from satellites flying at different heights above the surface of the Red Planet. To take into account the nuances of mathematical formulations, a more detailed study of the local features of the magnetic field of the Red Planet is required. We hope that, as information comes from the surface of the planet, we will be able to improve the accuracy and reliability of our method. Building an analytical model of the magnetic field of Mars from satellite data is a difficult problem that researchers from different countries are trying to solve. Verification of certain theoretical studies can be carried out only on the basis of experimental measurement data. At the same time, the difficulties that theorists face when solving interpretation problems can serve as an incentive for technologists and developers of space equipment when planning an experiment.