A fast algorithm for the large-scale joint inversion of gravity and magnetic data is developed. The algorithm uses a non-linear Gramian constraint to impose correlation between the density and susceptibility of the reconstructed models. The global objective function is formulated in the space of the weighted parameters, but the Gramian constraint is implemented in the original space, and the non-linear constraint is imposed using two separate Lagrange parameters, one for each model domain. It is significant that this combined approach, using the two spaces provides more similarity between the reconstructed models. Moreover, it is shown theoretically that the gradient for the use of the unweighted space is not a scalar multiple of that used for the weighted space, and hence cannot be accounted for by adjusting the Lagrange parameters. It is assumed that the measured data are obtained on a uniform grid and that a consistent regular discretization of the volume domain is imposed. Then, the sensitivity matrices exhibit a block-Toeplitz-Toeplitz-block structure for each depth layer of the model domain, and both forward and transpose operations with the matrices can be implemented efficiently using two dimensional fast Fourier transforms. This makes it feasible to solve for large scale problems with respect to both computational costs and memory demands, and to solve the non-linear problem by applying iterative methods that rely only on matrix–vector multiplications. As such, the use of the regularized reweighted conjugate gradient algorithm, in conjunction with the structure of the sensitivity matrices, leads to a fast methodology for large-scale joint inversion of geophysical data sets. Numerical simulations demonstrate that it is possible to apply a non-linear joint inversion algorithm, with Lp-norm stabilisers, for the reconstruction of large model domains on a standard laptop computer. It is demonstrated, that while the p = 1 choice provides sparse reconstructed solutions with sharp boundaries, it is also possible to use p = 2 in order to provide smooth and blurred models. The methodology is used for inverting gravity and magnetic data obtained over an area in northwest of Mesoproterozoic St Francois Terrane, southeast of Missouri, USA.
Within the iron metallogenic province of southeast Missouri, USA, there are several mines that contain not only economic iron resources, magnetite and/or hematite, but also contain rare earth elements, copper and gold. An area including three major deposits, Pea Ridge, Bourbon and Kratz Spring, was selected for detailed modelling for the upper crustal magnetic susceptibility and density structures. For the study area, ground gravity and high-resolution airborne magnetic and gravity gradiometry data sets are available. An efficient and novel joint inversion algorithm for the simultaneous inversion of these multiple data sets is presented. The Gramian coupling constraint is used to correlate the reconstructed density and magnetic susceptibility models. The implementation relies on the structures of the sensitivity matrices and an efficient minimization algorithm to achieve significant reductions in the memory requirements and computational costs. Consequently, it is feasible to use a laptop computer for the inversion of multiple data sets, each containing thousands of data points, for the recovery of models on the study area, each including approximately one million model parameters. This is the first time that these multiple data sets have been simultaneously inverted for this area. The L1-norm stabilizer is used to provide compact and focused images of the ore deposits. For contrast, independent inversions of each data set are also discussed. In general, our results provide new insights about the concealed ore deposits in the Mesoproterozoic basement rocks of southeast Missouri. Both short- and long-wavelength anomalies exist in the recovered models; these provide a high-resolution image of the subsurface. The geometry and physical properties of the known deposits are determined very well. Additionally, some unknown concealed deposits are revealed; these could be economically valuable and should be considered in future geophysical and geological investigations.
more » « less- NSF-PAR ID:
- 10435637
- Publisher / Repository:
- Oxford University Press
- Date Published:
- Journal Name:
- Geophysical Journal International
- Volume:
- 235
- Issue:
- 2
- ISSN:
- 0956-540X
- Page Range / eLocation ID:
- p. 1064-1085
- Format(s):
- Medium: X
- Sponsoring Org:
- National Science Foundation
More Like this
-
more » « less
SUMMARY -
An efficient algorithm for the Lp -norm joint inversion of gravity and magnetic data using the cross-gradient constraint is presented. The presented framework incorporates stabilizers that use Lp -norms ( 0≤p≤2 ) of the model parameters, and/or the gradient of the model parameters. The formulation is developed from standard approaches for independent inversion of single data sets, and, thus, also facilitates the inclusion of necessary model and data weighting matrices, for example, depth weighting and hard constraint matrices. Using the block Toeplitz Toeplitz block structure of the underlying sensitivity matrices for gravity and magnetic models, when data are obtained on a uniform grid, the blocks for each layer of the depth are embedded in block circulant circulant block matrices. Then, all operations with these matrices are implemented efficiently using 2-D fast Fourier transforms, with a significant reduction in storage requirements. The nonlinear global objective function is minimized iteratively by imposing stationarity on the linear equation that results from applying linearization of the objective function about a starting model. To numerically solve the resulting linear system, at each iteration, the conjugate gradient algorithm is used. This is improved for large scale problems by the introduction of an algorithm in which updates for the magnetic and gravity parameter models are alternated at each iteration, further reducing total computational cost and storage requirements. Numerical results using a complicated 3-D synthetic model and real data sets obtained over the Galinge iron-ore deposit in the Qinghai province, north-west (NW) of China, demonstrate the efficiency of the presented algorithm.more » « less
-
null (Ed.)SUMMARY We discuss the focusing inversion of potential field data for the recovery of sparse subsurface structures from surface measurement data on a uniform grid. For the uniform grid, the model sensitivity matrices have a block Toeplitz Toeplitz block structure for each block of columns related to a fixed depth layer of the subsurface. Then, all forward operations with the sensitivity matrix, or its transpose, are performed using the 2-D fast Fourier transform. Simulations are provided to show that the implementation of the focusing inversion algorithm using the fast Fourier transform is efficient, and that the algorithm can be realized on standard desktop computers with sufficient memory for storage of volumes up to size n ≈ 106. The linear systems of equations arising in the focusing inversion algorithm are solved using either Golub–Kahan bidiagonalization or randomized singular value decomposition algorithms. These two algorithms are contrasted for their efficiency when used to solve large-scale problems with respect to the sizes of the projected subspaces adopted for the solutions of the linear systems. The results confirm earlier studies that the randomized algorithms are to be preferred for the inversion of gravity data, and for data sets of size m it is sufficient to use projected spaces of size approximately m/8. For the inversion of magnetic data sets, we show that it is more efficient to use the Golub–Kahan bidiagonalization, and that it is again sufficient to use projected spaces of size approximately m/8. Simulations support the presented conclusions and are verified for the inversion of a magnetic data set obtained over the Wuskwatim Lake region in Manitoba, Canada.more » « less
-
more » « less
SUMMARY Gravity surveys constitute an important method for investigating the Earth's interior based on density contrasts related to Earth material differentials. Because lithology depends on the environment and the period of formation, there are generally clear boundaries between rocks with different lithologies. Inversions with convex functions for approximating the L0 norm are used to detect boundaries in reconstructed models. Optimizations can easily be found because of the convex transformations; however, the volume of the reconstructed model depends on the weighting parameter and the density constraint rather than the model sparsity. To determine and adapt the modelling size, a novel non-convex framework for gravity inversion is proposed. The proposed optimization aims to directly reduce the L0 norm of the density matrix. An improved iterative hard thresholding algorithm is developed to linearly reduce the L0 penalty during the inner iteration. Accordingly, it is possible to determine the modelling scale during the iteration and achieve an expected scale for the reconstructed model. Both simple and complex model experiments demonstrate that the proposed method efficiently reconstructs models. In addition, granites formed during the Yanshanian and Indosinian periods in the Nanling region, China, are reconstructed according to the modelling size evaluated in agreement with the magnetotelluric profile and density statistics of rock samples. The known ores occur at the contact zones between the sedimentary rocks and the reconstructed Yanshanian granites. The ore-forming bodies, periods, and processes are identified, providing guidance for further deep resource exploration in the study area.
-
Societal demand for critical metals used in the high-tech and green industries has led to an increased interest in REEs associated with ore deposits. Hydrothermal mineralization of monazite (CePO4) in various REE deposits can display significant variations in REE mineralogy and rock chemistry. Monazite displays textural and REE compositional variations, such as those observed in the giant Bayan Obo REE deposit in China and in the Pea Ridge iron-oxide-apatite (IOA) deposit in Missouri. The coupling of compositional variations of monazite with thermodynamic modeling of fluid-rock interaction processes may provide a useful vectoring tool in these ore deposits. However, interpreting these geochemical fingerprints requires building an internally consistent thermodynamic dataset for REE minerals and their relevant aqueous complexes. In this study, a series of hydrothermal solubility experiments were carried out using synthetic monazite crystals (i.e., LaPO4, PrPO4, NdPO4, and EuPO4) to assess the consistency of reported mineral calorimetric data and the thermodynamic data of the aqueous REE complexes. The solubility experiments were conducted in aqueous HClO4-H3PO4–bearing solutions at temperatures between 100° and 250°C and at saturated water vapor pressure. Equilibrium constants (Ks0) for the dissolution reaction of monazite end members were retrieved as a function of temperature and extrapolated to standard conditions of 25°C and 1 bar. Results indicate significant differences between the new solubility data and those reported in the literature. We demonstrate the impact of these new thermodynamic data in a series of fluid-rock interaction models using the GEMS code package (http://gems.web.psi.ch) and the MINES thermodynamic database (http://tdb.mines.edu). The simulated monazite stability can be correlated to field observations and allows for the prediction of the behavior of REE in hydrothermal fluids and their association to alteration observed in ore deposits.more » « less
