Ellipsoids (v1.0): 3D Magnetic modelling of ellipsoidal bodies

. A considerable amount of literature has been published on the magnetic modelling of uniformly magnetized ellipsoids since the second half of the nineteenth century. Ellipsoids have ﬂexibility to represent a wide range of geometrical forms, are the only known bodies which can be uniformly magnetized in the presence of a uniform inducing ﬁeld and are the only


Introduction
Based on the mathematical theory of the magnetic induction developed by Poisson (1824), Maxwell (1873) affirmed that, if U is the gravitational potential produced by any body with uniform density ρ and arbitrary shape at a point (x, y, z), then − ∂U ∂x is the magnetic scalar potential produced at the same point by the same body if it has a uniform magnetization oriented along x with intensity ρ.Maxwell (1873) generalized this idea as a way of determining the magnetic scalar potential produced by any uniformly magnetized body in a given direction.By presuming that this uniform magnetization is due to induction and that it is proportional to the resulting magnetic field (intensity) inside the body, he postulated that the resulting field must also be uniform and parallel to the magnetization.This uniformity is due to the fact that the resulting field is defined as the negative gradient of the magnetic scalar potential.As a consequence of this uniformity, the gravitational potential U at points within the body must be a quadratic function of the spatial coordinates.Apparently, Maxwell (1873) was the first one to postulate that ellipsoids are the only finite bodies having a gravitational potential which satisfies this property and hence can be uniformly magnetized in the presence of a uniform inducing magnetic field.This property can be extended to other bodies defined as limiting cases of an ellipsoid (e.g., spheres, elliptic cylinders), however all the remaining non-ellipsoidal bodies cannot be uniformly magnetized in the presence of a uniform inducing field.
Another particularity of ellipsoids is that they are the only bodies which enable an analytical computation of its selfdemagnetization.The self-demagnetization contributes to decrease the magnitude of the magnetization along the shortest axes of a body.It is a function of the body shape and gives rise to shape anisotropy (Uyeda et al., 1963;Thompson and Oldfield, 1986;Dunlop and Özdemir, 1997;Clark and Emerson, 1999;Tauxe, 2003).It is well-established in the literature that the selfdemagnetization can be neglected if the body has a susceptibility lower than 0.1 SI (Emerson et al., 1985;Clark et al., 1986;Eskola and Tervo, 1980;Guo et al., 1998Guo et al., , 2001;;Purss and Cull, 2005;Hillan and Foss, 2013;Austin et al., 2014;Clark, 2014).
On the other hand, neglecting the self-demagnetization in geological bodies with high susceptibilities (> 0.1 SI) may strongly mislead the interpretation obtained from magnetic methods.This limiting value, however, seems to be determined empirically and, so far, there has been little discussion about how it was determined.Farrar (1979) demonstrated the importance of the ellipsoidal model in taking into account the self-demagnetization and determining reliable drilling directions on the Tennant Creek field, Australia.Posteriorly, Hoschke (1991) also showed how the ellipsoidal model proved to be highly successful in locating and defining ironstone bodies in the Tennant Creek field.Clark (2000) provides a good discussion about the influence of the self-demagnetization in magnetic interpretation of the Osborne copper-gold deposit, Australia.This deposit is hosted by ironstone bodies that have very high susceptibility.According to Clark (2000), neglecting the effects of self-demagnetisation led errors of ≈ 55 • in the interpreted dip.Recently, Austin et al. (2014) used magnetic modelling and rock property measurements to show that, contrary to previous interpretations, the magnetization of the Candelaria iron oxide copper-gold deposit, Chile, is not dominated by the induced component.Rather, the deposit has a relatively weak remanent magnetization and is strongly affected by self-demagnetization.These examples show the importance of the self-demagnetization and the ellipsoidal model in producing trustworthy geological models of high-susceptibility orebodies, which may save significant cost associated with drilling.
A vast literature about the magnetic modelling of ellipsoidal bodies was developed in which are to be found the names of many researchers.Nevertheless, interest in this subject has not yet died out, as is evidenced by a list of modern papers in this field.Besides, the geoscientific community lacks an : a ::: free : easy-to-use tool to simulate the magnetic field produced by uniformly magnetized ellipsoids.Such a tool could prove to be useful either for teaching and researching geophysics.
In this work, we present a review of the vast literature about the magnetic modelling of ellipsoidal bodies and a theoretical discussion about the determination of the isotropic susceptibility value above which the self-demagnetization must be taken into consideration.We propose an alternative way of determining this value based on the body shape and the maximum relative error allowed in the resultant magnetization.This alternative approach is validated by the results obtained with numerical simulations.We also provide a set of routines to model the magnetic field produced by ellipsoids.The routines are written in Python language as part of the Fatiando a Terra (Uieda et al., 2013), which is an open-source library for modelling and inversion in geophysics.We attempt to use the best practices of continuous integration, documentation, unit-testing, and version-control for the purpose of providing a reliable and easy-to-use code.
The magnetic modelling of an ellipsoidal body is commonly performed in a particular Cartesian coordinate system that is aligned with the body semi-axes and has the origin coincident with the body centre (Fig. 1b).For convenience, we denominate this particular coordinate system as local coordinate system.The relationship between the Cartesian coordinates (x, ỹ, z) of a point in a local coordinate system and the Cartesian coordinates (x, y, z) of the same point in the main system is given by: where r = [ x ỹ z ] , r and r c are defined in Eq. 1 and the matrix V (Eq.2) is defined according to the ellipsoid type.
Subsequently, quantities referred to the local coordinate system (Fig. 1b) are represented with the simbol "∼".

Theoretical background
Consider a magnetized ellipsoid immersed in a uniform inducing magnetic field H 0 (in Am −1 ) given by where • denotes the Euclidean norm (or 2-norm) and D and I are respectively, the declination and inclination of the localgeomagnetic field in the main coordinate system (Fig. 1a).This field represents the main component of the Earth's magnetic field, which is usually assumed to be generated by the Earth's liquid core.In the absence of conduction currents, the total magnetic field H(r) at the position r (Eq. 1) of a point referred to the main coordinate system is defined as follows (Sharma, 1966;Eskola and Tervo, 1980;Reitz et al., 1992;Stratton, 2007): where the second term is the negative gradient of the magnetic scalar potential V (r) given by: In this equation, r = [ x y z ] is the position vector of a point located within the volume ϑ, the integral is conducted over the variables x , y and, z and M(r ) is the magnetization vector (in Am −1 ).Eq. 11 is valid anywhere, independently if the position vector r represents a point located inside or outside the magnetized body (DuBois, 1896;Reitz et al., 1992;Stratton, 2007).
By using the magnetization M defined by Eq. 12, the total magnetic field H(r) (Eq.10) can be rewritten as follows: where N(r) is a symmetrical matrix whose ij-element n ij (r) is given by r 1 = x, r 2 = y, r 3 = z are the elements of the position vector r (Eq.1), and Notice that the scalar function f (r) (Eq.17) is proportional to the gravitational potential that would be produced by the ellipsoidal body with volume ϑ if it had a uniform density equal to the inverse of the gravitational constant.It can be shown that the elements n ij (r) are finite whether r is a point within or without the volume ϑ (Peirce, 1902;Webster, 1904).The matrix N(r) (Eq.15) is called depolarization tensor (Solivérez, 1981(Solivérez, , 2008)).
The following part of this paper moves on to describe the magnetic field H(r) (Eq.15) at points located both within and without the volume ϑ of the ellipsoidal body.However, the mathematical developments are conveniently performed in the local coordinate system (Fig. 1b) related to the respective ellipsoidal body.

Coordinate transformation
To continue our description of the magnetic modelling of ellipsoidal bodies, it is convenient to perform two important coordinate transformations.The first one transforms the scalar function f (r) (Eq.17) from the main coordinate system (Fig. 1a) into a new scalar function f (r) referred to the local coordinate system (Fig. 1b).The function f (r) was first presented by Dirichlet ( 1839) to describe the gravitational potential produced by homogeneous ellipsoids.Posteriorly, several authors also deduced and used this function for describing the magnetic and gravitational fields produced by triaxial, prolate, and oblate ellipsoids (Maxwell, 1873;Thomson and Tait, 1879;DuBois, 1896;Peirce, 1902;Webster, 1904;Kellogg, 1929;Stoner, 1945;Osborn, 1945;Peake and Davy, 1953;Macmillan, 1958;Chang, 1961;Lowes, 1974;Clark et al., 1986;Tejedor et al., 1995;Stratton, 2007).
It is convenient to use f † (r) and f ‡ (r) to define the function f (r) evaluated, respectively, at points r inside and outside the volume ϑ of the ellipsoidal body.The scalar function f † (r) is given by where This function represents the gravitational potential that would be produced by the ellipsoidal body at points located within its volume ϑ if it had a uniform density equal to the inverse of the gravitational constant.Notice that, in this case, the gravitational potential is a quadratic function of the spatial coordinates x, ỹ, and z, which supported the Maxwell's (1873) postulate about uniformly magnetized ellipsoids.In a similar way, the function f ‡ (r) is given by where R(u) is defined by Eq. 19 and the parameter λ is defined according to the ellipsoid type as a function of the spatial coordinates x, ỹ, and z (see Appendix B).For readers interested in additional information about the parameter λ, we recommend Webster (1904, p. 234), Kellogg (1929, p. 184) and Clark et al. (1986).
The second important coordinate transformation is defined with respect to Eq. 15.By properly using the orthogonality of matrix V (Eq.2), the magnetic field H(r) (Eq.15) can be transformed from the main coordinate system (Fig. 1a) to the local coordinate system (Fig. 1b) as follows: V H(r) where the superscript "∼" denotes quantities referred to the respective local coordinate system.
In Eq.21, the transformed depolarization tensor Ñ(r) is calculated as a function of the original depolarization tensor N(r) (Eq.15).In this case, the elements of Ñ(r) are calculated as a function of the second derivatives of the function f (r) (Eq.17), which is defined in the main coordinate system (Fig. 1a).It can be shown (Appendix A), however, that the elements ñij (r) of Ñ(r) can also be calculated as follows: where r1 = x, r2 = ỹ, and r3 = z are the elements of the transformed vector r (Eq.8) and f (r) is given by Eq. 18 or 20, depending if r represents a point located within or without the volume ϑ of the ellipsoidal body.

Depolarization tensor Ñ †
Let Ñ † be the transformed depolarization tensor calculated for the case in which r (Eq.8) represents a point located outside ::::: inside the ellipsoidal body.In this case, the elements of Ñ † are calculated according to Eq. 22, with f (r) given by f † (r) (Eq.18).As we have already pointed out, the f † (r) (Eq.18) is a quadratic function of the spatial coordinates x, ỹ and z.
Consequently, the elements ñ † ij , i = 1, 2, 3, j = 1, 2, 3, of Ñ † do not depend on the elements of the transformed position vector r (Eq.8).Besides, the off-diagonal elements are zero and the diagonal elements are given by (Stoner, 1945): where R(u) is defined by Eq. 19 and e 1 = a, e 2 = b, and e 3 = c.These elements are commonly known as demagnetizing factors and are defined according to the ellipsoid type.Here, we calculate the demagnetizing factors in the SI system.Consequently, they satisfies the condition ñ † 11 + ñ † 22 + ñ † 33 = 1, independently of the ellipsoid type.It is worth stressing that, according to Eq. 23, the demagnetizing factors ñ † ii are constants defined by the ellipsoid semi-axes a, b, and c.
Notice that, according to Eqs.21 and A7, where Ñ † is a diagonal matrix and V (Eq.2) is an orthogonal matrix.This equation shows that, for the particular case in which r and consequently r represent a point inside the volume ϑ of the ellipsoid, the elements ñ † ii of Ñ † represent the eigenvalues while the columns of V represent the eigenvectors of the original depolarization tensor N(r).

Triaxial ellipsoids
For triaxial ellipsoids (e.g., a > b > c), the demagnetizing factors obtained by solving Eq. 23 are given by: where and

Prolate ellipsoids
For prolate (e.g., a > b = c) ellipsoids, the functions g i (Eq.37) are given by: and where g 3 = g 2 .These formulas can be obtained by properly manipulating those presented by (Emerson et al., 1985).

Oblate ellipsoids
For oblate (e.g., a < b = c) ellipsoids, the functions g i (Eq.37) are given by: and where g 3 = g 2 .Similarly to the case of prolate ellipsoid shown previously, these formulas can be obtained by properly manipulating those presented by (Emerson et al., 1985).

Internal magnetic field and magnetization
By considering r as a point within the volume ϑ of the ellipsoid and using the Maxwell's postulate about the uniformity of the magnetic field H(r) inside ellipsoidal bodies, we can use Eq.21 for defining the resultant uniform magnetic field H † inside the ellipsoidal body as follows: where I is the identity matrix and Ñ † is defined in the previous section.
Let us pre-multiply the uniform internal field H † (Eq.45) by the transformed susceptibility tensor K (Eq.21) to obtain where M represents the transformed magnetization, as can be easily verified by using Eqs. 12 and 21.The matrix identity used for obtaining the second line of Eq. 46 is given by Searle (1982, p. 151).
Equation 46 can be easily generalized for the case in which the ellipsoid has also a uniform remanent magnetization MR .
Let us first consider that the uniform remanent magnetization satisfies the condition HA = K−1 MR , where HA represents a hypothetical uniform ancient field.Then, if we assume that H0 , in Eqs.45 and 46, is in fact the sum of the inducing magnetic field H0 and the hypothetical ancient field HA , we obtain the following generalized equation where Despite of the coordinate system transformation represented by the matrix V (Eq.2), Eq. 47 is consistent with that given by Clark et al. (1986, Eq. 38).It shows the combined effect of the anisotropy of magnetic susceptibility (AMS) and the shape anisotropy.The AMS is represented by the susceptibility tensor K (Eq.13) and reflects the preferred orientation of the magnetic minerals forming the body.The susceptibility tensor appears in Eq. 47, defined in the main coordinate system (Fig. 1a), and Eq.48, defined in the local coordinate system (Fig. 1b).The shape anisotropy is represented, in Eq. 47, by the depolarization tensor Ñ † and reflects the self-demagnetization associated to the body shape.Notice that the resultant magnetization M (Eq.47) does not necessarily have the same direction as the inducing field H 0 (Eq.9).The angular difference between the resultant magnetization and the inducing field depends on the combined effect of the anisotropy of magnetic susceptibility and the shape anisotropy.
For the particular case in which the susceptibility is isotropic, the susceptibility tensor is defined according to Eq. 14.In this case, the magnetization M (Eqs.12 and 47), referred to the main coordinate system (Fig. 1a), and the matrix Λ (Eq.48) can be rewritten as follows: and Despite the coordinate transformation represented by matrix V (Eq.2), this equation is in perfect agreement with those presented by Guo et al. (2001, Eqs. 13-15).The first term, depending on the inducing field H 0 (Eq.9), represents the induced magnetization whereas the term depending on M R is the remanent magnetization.Equation 49 reveals that, as pointed out by many authors (e.g., Maxwell, 1873;DuBois, 1896;Stoner, 1945;Clark et al., 1986;Stratton, 2007), the induced magnetization opposes the inducing field if it is parallel to an ellipsoid axis, independently of the ellipsoid type.Otherwise, the magnetization is not necessarily parallel to the inducing field.If we additionally consider that χ << 1, the matrix Λ (Eq.50) approaches to the identity and the magnetization M (Eq.49) can be approximated by: which is the classical equation describing the resultant magnetization in applied geophysics (Blakely, 1996, p. 89).Notice that, in this particular case, the induced magnetization is parallel to the inducing field H 0 (Eq.9), whether it is parallel to an ellipsoid axis or not.Usually, Eq. 51 is considered a good approximation for χ ≤ 0.1 SI.Although this value has been largely used in the literature, there have been few empirical and/or theoretical investigations about it.

Relationship between χ and the relative error in M
In the case of isotropic susceptibility, the resultant magnetization M (Eq.49) may be determined by solving the following linear system: where, according to Eq. 50, As we have already pointed out, the approximated magnetization M (Eq.51) represents the particular case in which the matrix Λ (Eq.50), and consequently the matrix Λ −1 (Eq.53), are close to the identity.
Consider a perturbed matrix δΛ −1 given by and, similarly, a perturbed magnetization vector δM given by By using these two equations, we may rewrite that of the approximated magnetization M (Eq.51) as follows: Now, by subtracting the true magnetization M (Eq.52) from this linear system (Eq.56) and rearranging the terms, we obtain the following linear system for the perturbed magnetization δM (Eq.55): By using the concept of vector norm and its corresponding operator norm (Demmel, 1997;Golub and Loan, 2013), we may use Eq.57 to write the following inequality: where δM and M denote Euclidean norms (or 2-norms) and the term δΛ −1 denotes the matrix 2-norm of δΛ −1 .By using Eqs.53 and 54 and the orthogonal invariance of the matrix 2-norm (Demmel, 1997;Golub and Loan, 2013), we define δΛ −1 as follows: where ñ † max is the demagnetization factor associated to the shortest ellipsoid semi-axis.For a triaxial ellipsoid, ñ † max ≡ ñ †

33
(Eq. 27), for a prolate ellipsoid, ñ † max ≡ ñ † 22 (Eq.31), and, for an oblate ellipsoid, ñ † max ≡ ñ † 11 (Eq.32).It is worth stressing that, independently of the ellipsoid type, ñ † max is a scalar function of the ellipsoid semi-axes.In Eq. 58, the ratio δM M −1 represents the relative error in the approximated magnetization M (Eq.51) with respect to the true magnetization M (Eqs.49 and 52).Given a target relative error and an ellipsoid with given semi-axes, we may use the inequality represented by Eq. 59 to define which represents the maximum isotropic susceptibility that the ellipsoidal body can assume in order o guarantee a relative error lower than or equal to .For isotropic susceptibilities greater than χ max , there is no guarantee that the relative error in the approximated magnetization M (Eq.51) with respect to the true magnetization M (Eqs.49 and 52) is lower than or equal to .The geoscientific community has been using χ max = 0.1 SI as a limit value for neglecting the self-demagnetization and, consequently, use magnetization M (Eq.51) as a good approximation of the true magnetization M (Eqs.49 and 52).Equation 60, on the other hand, defines χ max as a function of the ellipsoid semi-axes, according to a user-specified relative error .

External magnetic field and total-field anomaly
The magnetic field ∆H(r) produced by an ellipsoid at external points is calculated from Eqs. 21 and 47 as the difference between the resultant field H(r) and the inducing field H 0 : where Ñ ‡ (r) is the transformed depolarization tensor whose elements ñ ‡ ii (r) and ñ ‡ ij (r) are defined, respectively, by Eqs.34 and 35.∆H(r) represents the magnetic field produced by a uniformly magnetized body located in the crust.Equation 65gives the magnetic field (in A m −1 ) produced by an ellipsoid.However, in geophysics, the most widely used field is the magnetic induction (in nT).Fortunately, this conversion can be easily done by multiplying Eq. 65 by k m = 10 9 µ 0 , where µ 0 represents the magnetic constant (in H m −1 ).For geophysical applications, it is preferable to calculate the total-field anomaly produced by the magnetic sources.This scalar quantity is given by (Blakely, 1996): where B 0 = k m H 0 and ∆B(r) = k m ∆H(r), with H 0 and ∆H(r) defined, respectively, by Eqs. 9 and 65.In practical situations, however, B 0 >> ∆B(r) and, consequently, the following approximation is valid (Blakely, 1996): 3 Computational implementation and reproducibility The code is implemented in the Python language, by using the NumPy and SciPy libraries (van der Walt et al., 2011), as part of the open-source source library Fatiando a Terra (Uieda et al., 2013).Our code is very modular and has a test suite formed by a considerable number of assertions, unit tests, doc tests, and integration tests.We refer the readers interested in best practices for scientific computing to Wilson et al. (2014).
The numerical simulations presented here were generated with the Jupyter Notebook (http://jupyter.org),which is a web application that allows creating and sharing documents that contain live code, equations, visualizations and explanatory text.
Besides using Fatiando a Terra (Uieda et al., 2013), the numerical simulations use the NumPy library (van der Walt et al., 2011) to perform numerical computations and the Matplotlib library (Hunter, 2007) to plot the results and generate figures.
The Jupyter Notebooks used to produce all the results presented here are available in a repository on GitHub (https://github.com/pinga-lab/magnetic-ellipsoid).The larger the variable u, the larger the resulting semi-axes a, b, and c, but the smaller the relative difference between them.
Table 2 shows the parameters defining a synthetic orebody which is based on that presented by Farrar (1979) to represent the Warrego orebody.Figure 5 shows the total-field anomaly ∆T (r) (Eq.67) produced by the synthetic body on a regular grid of 100 × 100 points at a constant vertical coordinate z = 0 m.The total-field anomaly varies from ≈ −176 :::::: ≈ −71 : nT to ≈ 237 ::::: ≈ 482 : nT, resulting in a peak-to-peak amplitude of ≈ 413 ::::: ≈ 553 nT, and was calculated by using the true magnetization M defined in Eqs.49 and 52.
We have calculated the difference between the total-field anomaly ∆T (r) (Eq.67) calculated with the true magnetization M (Eqs.49 and 52) and that calculated with the approximated magnetization M (Eq.51).The differences were calculated by using the synthetic body defined in Tab. 2, but with three different isotropic susceptibilities.Figures 6a, 6b and 6c show the differences calculated by using, respectively, isotropic susceptibilities χ = 1.69 SI (Tab.2), χ 1 = 0.1 SI and χ 2 = 0.116 SI.
As expected, the differences calculated by using the higher isotropic susceptibility (Fig. 6a) are very large.The peak-to-peak amplitude is ≈ 141 :::: ≈ 40 : nT and represents ≈ 34% ::::: ≈ 8% of the peak-to-peak amplitude of the total-field anomaly shown in Fig. 5.This result reinforces that the use of the approximated magnetization M (Eq.51) may negatively impact the magnetic modelling of bodies with high isotropic susceptibility.
Finally, Fig. 6c shows the differences calculated by using χ 2 = 0.116 SI.This value was calculated by using Eq.60 with a target relative error = 8% and the ñ † max defined by Eq. 27.By using this isotropic susceptibility, it is expected that the calculated relative error δM M −1 (Eq.58) in the magnetization be lower than or equal to the target relative error = 8%.
By considering the functions f (r) (Eq.17) and f (r) evaluated at the same point, but on different coordinate systems, we have: , j = 1, 2, 3 , which, from Eq. A1, can be given by ∂f (r) Now, by deriving ∂f (r) ∂rj (Eq.A2) with respect to the ith element r i of the position vector r (Eq.1), we obtain: ) where F(r) is a 3 × 3 matrix whose ij-th element is ∂ 2 f (r) ∂ ri ∂ rj .From Eq. A3, we obtain  Table 2. Parameters defining a synthetic orebody.This model is based on the that presented by Farrar (1979)  .Difference between the total-field anomaly calculated with the approximated magnetization M (Eq.51) and with the true magnetization M (Eqs.49 and 52).The total-field anomalies are : in : nT :: and :::: were : calculated with Eq. 67, on a regular grid of 100 × 100 points, at the constant vertical coordinate z = 0 m.The differences are produced by the synthetic orebody defined in Tab. 2, but with different isotropic susceptibilities: (a) the isotropic susceptibility defined in Tab. 2, (b) an isotropic susceptibility χ = 0.1 SI, and (c) an isotropic susceptibility χ = 0.116 SI.This last value was calculated with Eq. 60, by using = 8%.

Figure 1 .Figure 5 .
Figure 1.Schematic representation of the coordinate systems used to represent an ellipsoidal body.(a) Main coordinate system with axes x, y, and z pointing to North, East, ans down, respectively.The dark gray plane contains the centre (xc, yc, zc) (white circle) and two unit vectors, u and w, defining two semi-axes of the ellipsoidal body.For triaxial and prolate ellipsoids, u and w define, respectively, the semiaxes a and b.For oblate ellipsoids, u and w define the semi-axes b and c, respectively.The strike direction is defined by the intersection of the dark gray plane and the horizontal plane (represented in light gray), which contains the x and y axes.The angle ε between "minus :: the x" :::: -axis and the strike direction is called strike.The angle ζ between the horizontal plane and the dark gray plane is called dip.The line containing the unit vector u defines the plunge direction.The angle η between the strike direction and the plunge direction ::: line :::::::: containing :: the ::: unit ::::: vector :: u : is called rake.The projection of the plunge direction :: this :::: line on the horizontal plane ::: (not :::::: shown) is called plunge azimuth direction ::: dip ::::::: direction :::::::::::::::::::::::::::::::::::::::: (Pollard and Fletcher, 2005; Allmendinger et al., 2012).(b) Local coordinate system with origin at the ellipsoid centre (xc, yc, zc) (black dot) and axes defined by unit vectors v1, v2, and v3.These unit vectors define the semi-axes a, b, and c of triaxial, prolate, and oblate ellipsoids in the same way.For triaxial and prolate ellipsoids, the unit vectors u and w shown in (a) coincide with v1 and v2, respectively.For oblate ellipsoids, the unit vectors u and w shown in (a) coincide with v2 and v3, respectively.
Figure6.Difference between the total-field anomaly calculated with the approximated magnetization M (Eq.51) and with the true magneti-