BZPEERReading, trust and collaboration for research
Menu

Discover

GraphFollow connections around a paper, topic or saved view.

Workspaces

HomePick up where your research flow left off.InboxHandle requests, replies and notifications from one place.CollectionsCurate reading lists, evidence sets and shareable dossiers.ResearchSee your private provenance graph, trail and imports.CollaborationMove from discovery into focused discussion rooms.PublishDraft, preview and publish full-text research articles.

Network

ExpertsFind specialists worth following or asking for help.CommunitiesJoin research circles organized around shared work.PartnersSee external organizations active around the same topics.

Opportunities

ListingsBrowse roles, calls and collaborations with clearer fit signals.

Account

Log inReturn to your reading and collaboration workspace.Create accountStart saving work, following people and joining teams.
Log inCreate account

Researcher profile

Saeed Vatankhah

Saeed Vatankhah contributes to research discovery and scholarly infrastructure.

ResearcherAffiliation not importedOpen to collaborate
math.NANumerical Analysisphysics.geo-phComputational Engineering, Finance, and ScienceMethodologyphysics.comp-phphysics.data-an

Trust snapshot

Quick read

Trust 21 - EmergingVerification L1Unclaimed author
14works
0followers
7topics
4close collaborators

Actions

Decide how to stay connected

Follow researcher0

Identity and collaboration

How to connect with this researcher

Claiming links this public author record to a researcher profile and unlocks direct collaboration workflows.

Log in to claim

Direct collaboration

Open a focused conversation when the fit is right

Claim this author entity first to unlock direct invitations.

Research graph

See the researcher in context

Open full explorer

Inspect adjacent work, topics, institutions and collaborators without jumping out to a separate graph page.

Building this graph slice

BZPEER is loading the nearby papers, people, topics and institutions for this page.

Published work

14 published item(s)

preprint2021arXiv

Large-scale focusing joint inversion of gravity and magnetic data with Gramian constraint

A fast algorithm for the large-scale joint inversion of gravity and magnetic data is developed. It uses a nonlinear Gramian constraint to impose correlation between density and susceptibility of 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 nonlinear constraint is imposed using two separate Lagrange parameters, one for each model domain. This combined approach provides more similarity between the reconstructed models. 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. The sensitivity matrices exhibit a block Toeplitz Toeplitz block structure for each depth layer of the model domain. 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 nonlinear 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 nonlinear joint inversion algorithm, with $L_p$-norm stabilisers, for the reconstruction of large model domains on a standard laptop computer. It is demonstrated, that the p=1 choice provides sparse reconstructed solutions with sharp boundaries, and $p=2$ provides smooth and blurred models. Gravity and magnetic data obtained over an area in northwest of Mesoproterozoic St. Francois Terrane, southeast of Missouri, USA are inverted.

0 citations0 reviewsNo code linkedOpen access
preprint2020arXiv

A fast methodology for large-scale focusing inversion of gravity and magnetic data using the structured model matrix and the $2D$ fast Fourier transform

Focusing inversion of potential field data for the recovery of sparse subsurface structures from surface measurement data on a uniform grid is discussed. For the uniform grid the model sensitivity matrices exhibit block Toeplitz Toeplitz block structure, by blocks for each depth layer of the subsurface. Then, through embedding in circulant matrices, all forward operations with the sensitivity matrix, or its transpose, are realized using the fast two dimensional Fourier transform. Simulations demonstrate that this fast inversion algorithm can be implemented on standard desktop computers with sufficient memory for storage of volumes up to size $n \approx 1M$. The linear systems of equations arising in the focusing inversion algorithm are solved using either Golub Kahan bidiagonalization or randomized singular value decomposition algorithms in which all matrix operations with the sensitivity matrix are implemented using the fast Fourier transform. These two algorithms are contrasted for efficiency for large-scale problems with respect to the sizes of the projected subspaces adopted for the solutions of the linear systems. The presented results confirm earlier studies that the randomized algorithms are to be preferred for the inversion of gravity data, and that it is sufficient to use projected spaces of size approximately $m/8$, for data sets of size $m$. In contrast, the Golub Kahan bidiagonalization leads to more efficient implementations for the inversion of magnetic data sets, and it is again sufficient to use projected spaces of size approximately $m/8$. Moreover, it is sufficient to use projected spaces of size $m/20$ when $m$ is large, $m \approx 50000$, to reconstruct volumes with $n \approx 1M$. Simulations support the presented conclusions and are verified on the inversion of a practical magnetic data set that is obtained over the Wuskwatim Lake region in Manitoba, Canada.

0 citations0 reviewsNo code linkedOpen access
preprint2020arXiv

Generalized L$_p$-norm joint inversion of gravity and magnetic data using cross-gradient constraint

A generalized unifying approach for $L_{p}$-norm joint inversion of gravity and magnetic data using the cross-gradient constraint is presented. The presented framework incorporates stabilizers that use $L_{0}$, $L_{1}$, and $L_{2}$-norms of the model parameters, and/or the gradient of the model parameters. Furthermore, 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 that provide, for example, depth weighting and imposition of hard constraint data. The developed efficient algorithm can, therefore, be employed to provide physically-relevant smooth, sparse, or blocky target(s) which are relevant to the geophysical community. Here, the nonlinear objective function, that describes the inclusion of all stabilizing terms and the fit to data measurements, 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. The general framework is then validated for three-dimensional synthetic models for both sparse and smooth reconstructions, and the results are compared with those of individual gravity and magnetic inversions. It is demonstrated that the presented joint inversion algorithm is practical and significantly improves reconstructed models obtained by independent inversion.

0 citations0 reviewsNo code linkedOpen access
preprint2019arXiv

A Tutorial and Open Source Software for the Efficient Evaluation of Gravity and Magnetic Kernels

Fast computation of three-dimensional gravity and magnetic forward models is considered. Measurement data is assumed to be obtained on a uniform grid which is staggered with respect to the discretization of the parameter volume. Then, the resulting kernel sensitivity matrices exhibit block-Toeplitz Toeplitz-block (BTTB) structure. These matrices are symmetric for the gravity problem but non-symmetric for the magnetic problem. In each case, the structure facilitates fast forward computation using two-dimensional fast Fourier transforms. The construction of the kernel matrices and the application of the transform for fast forward multiplication, for each problem, is carefully described. But, for purposes of comparison with the transform approach, the generation of the unique entries that define a given kernel matrix is also explained. It is also demonstrated how the matrices, and hence transforms, are adjusted when padding around the volume domain is introduced. The transform algorithms for fast forward matrix multiplication with the sensitivity matrix and its transpose, without the direct construction of the relevant matrices, are presented. Numerical experiments demonstrate the significant reduction in computation time that is achieved using the transform implementation. Moreover, it becomes feasible, both in terms of reduced memory requirements and computational time, to implement the transform algorithms for large three-dimensional volumes. All presented algorithms, including with variable padding, are coded for optimal memory, storage and computation as an open source MATLAB code which can be adapted for any convolution kernel which generates a BTTB matrix. This work, therefore, provides a general tool for the efficient simulation of gravity and magnetic field data, as well as any formulation which admits a sensitivity matrix with the required structure.

0 citations0 reviewsNo code linkedOpen access
preprint2019arXiv

Convergence of Regularization Parameters for Solutions Using the Filtered Truncated Singular Value Decomposition

The truncated singular value decomposition may be used to find the solution of linear discrete ill-posed problems in conjunction with Tikhonov regularization and requires the estimation of a regularization parameter that balances between the sizes of the fit to data function and the regularization term. The unbiased predictive risk estimator is one suggested method for finding the regularization parameter when the noise in the measurements is normally distributed with known variance. In this paper we provide an algorithm using the unbiased predictive risk estimator that automatically finds both the regularization parameter and the number of terms to use from the singular value decomposition. Underlying the algorithm is a new result that proves that the regularization parameter converges with the number of terms from the singular value decomposition. For the analysis it is sufficient to assume that the discrete Picard condition is satisfied for exact data and that noise completely contaminates the measured data coefficients for a sufficiently large number of terms, dependent on both the noise level and the degree of ill-posedness of the system. A lower bound for the regularization parameter is provided leading to a computationally efficient algorithm. Supporting results are compared with those obtained using the method of generalized cross validation. Simulations for two-dimensional examples verify the theoretical analysis and the effectiveness of the algorithm for increasing noise levels, and demonstrate that the relative reconstruction errors obtained using the truncated singular value decomposition are less than those obtained using the singular value decomposition. This is a pre-print of an article published in BIT Numerical Mathematics. The final authenticated version is available online at: https://doi.org/10.1007%2Fs10543-019-00762-7.

0 citations0 reviewsNo code linkedOpen access
preprint2019arXiv

Improving the use of the randomized singular value decomposition for the inversion of gravity and magnetic data

The large-scale focusing inversion of gravity and magnetic potential field data using $L_1$-norm regularization is considered. The use of the randomized singular value decomposition methodology facilitates tackling the computational challenge that arises in the solution of these large-scale inverse problems. As such the powerful randomized singular value decomposition is used for the numerical solution of all linear systems required in the algorithm. A comprehensive comparison of the developed methodology for the inversion of magnetic and gravity data is presented. These results indicate that there is generally an important difference between the gravity and magnetic inversion problems. Specifically, the randomized singular value decomposition is dependent on the generation of a rank $q$ approximation to the underlying model matrix, and the results demonstrate that $q$ needs to be larger, for equivalent problem sizes, for the magnetic problem as compared to the gravity problem. Without a relatively large $q$ the dominant singular values of the magnetic model matrix are not well-approximated. The comparison also shows how the use of the power iteration embedded within the randomized algorithm is used to improve the quality of the resulting dominant subspace approximation, especially in magnetic inversion, yielding acceptable approximations for smaller choices of $q$. The price to pay is the trade-off between approximation accuracy and computational cost. The algorithm is applied for the inversion of magnetic data obtained over a portion of the Wuskwatim Lake region in Manitoba, Canada

0 citations0 reviewsNo code linkedOpen access
preprint2018arXiv

IGUG: A MATLAB package for $3$D inversion of gravity data using graph theory

The open source MATLAB package IGUG for 3D inversion of gravity data is presented. It is based on methodology that was introduced by Bijani et al (2015), in which a homogeneous subsurface body is modeled by an ensemble of simple point masses. The model parameters are the Cartesian coordinates of the point masses and their total mass. Associating the point masses to the vertices of a weighted full graph with weights computed by the Euclidean pairwise distances separating vertices, Kruskal's algorithm is used to find the minimum spanning tree of the graph. Stabilization is achieved using an equidistance function that restricts the spatial distribution of the masses and favors a homogeneous subsurface structure. A regularization parameter $λ$ is introduced to balance the two terms of the objective function, and reasonable physically-relevant bound constraints are imposed on the model parameters. Then the bound constrained objective function is solved using a genetic algorithm. A new diagnostic approach is presented for determining a suitable choice for $λ$, requiring a limited number of solutions for a small set of $λ$ without using the L-curve as suggested by Bijani et al. Simulations for synthetic examples demonstrate the efficiency and effectiveness of the implementation of the algorithm. It is verified that the constraints on the model parameters are not restrictive. Included in the package is the script GMD.m which is used for generating synthetic data and for putting measurement data in the format required for the inversion implemented by IGUG.m The script Diagnostic_Results.m is included for analyzing and visualizing the results. The software can be used to verify the simulations and the analysis of real data that is presented here, The real data set uses gravity data from the Mobrun ore body, north east of Noranda, Quebec, Canada.

0 citations0 reviewsNo code linkedOpen access
preprint2017arXiv

3-D Projected $ L_{1}$ inversion of gravity data

Sparse inversion of gravity data based on $L_1$-norm regularization is discussed. An iteratively reweighted least squares algorithm is used to solve the problem. At each iteration the solution of a linear system of equations and the determination of a suitable regularization parameter are considered. The LSQR iteration is used to project the system of equations onto a smaller subspace that inherits the ill-conditioning of the full space problem. We show that the gravity kernel is only mildly to moderately ill-conditioned. Thus, while the dominant spectrum of the projected problem accurately approximates the dominant spectrum of the full space problem, the entire spectrum of the projected problem inherits the ill-conditioning of the full problem. Consequently, determining the regularization parameter based on the entire spectrum of the projected problem necessarily over compensates for the non-dominant portion of the spectrum and leads to inaccurate approximations for the full-space solution. In contrast, finding the regularization parameter using a truncated singular space of the projected operator is efficient and effective. Simulations for synthetic examples with noise demonstrate the approach using the method of unbiased predictive risk estimation for the truncated projected spectrum. The method is used on gravity data from the Mobrun ore body, northeast of Noranda, Quebec, Canada. The $3$-D reconstructed model is in agreement with known drill-hole information.

0 citations0 reviewsNo code linkedOpen access
preprint2017arXiv

A fast algorithm for regularized focused 3-D inversion of gravity data using the randomized SVD

A fast algorithm for solving the under-determined 3-D linear gravity inverse problem based on the randomized singular value decomposition (RSVD) is developed. The algorithm combines an iteratively reweighted approach for $L_1$-norm regularization with the RSVD methodology in which the large scale linear system at each iteration is replaced with a much smaller linear system. Although the optimal choice for the low rank approximation of the system matrix with m rows is q=m, acceptable results are achievable with q<<m. In contrast to the use of the LSQR algorithm for the solution of the linear systems at each iteration, the singular values generated using the RSVD yield a good approximation of the dominant singular values of the large scale system matrix. The regularization parameter found for the small system at each iteration is thus dependent on the dominant singular values of the large scale system matrix and appropriately regularizes the dominant singular space of the large scale problem. The results achieved are comparable with those obtained using the LSQR algorithm for solving each linear system, but are obtained at reduced computational cost. The method has been tested on synthetic models along with the real gravity data from the Morro do Engenho complex from central Brazil.

0 citations0 reviewsNo code linkedOpen access
preprint2017arXiv

Total variation regularization of the $3$-D gravity inverse problem using a randomized generalized singular value decomposition

We present a fast algorithm for the total variation regularization of the $3$-D gravity inverse problem. Through imposition of the total variation regularization, subsurface structures presenting with sharp discontinuities are preserved better than when using a conventional minimum-structure inversion. The associated problem formulation for the regularization is non linear but can be solved using an iteratively reweighted least squares algorithm. For small scale problems the regularized least squares problem at each iteration can be solved using the generalized singular value decomposition. This is not feasible for large scale problems. Instead we introduce the use of a randomized generalized singular value decomposition in order to reduce the dimensions of the problem and provide an effective and efficient solution technique. For further efficiency an alternating direction algorithm is used to implement the total variation weighting operator within the iteratively reweighted least squares algorithm. Presented results for synthetic examples demonstrate that the novel randomized decomposition provides good accuracy for reduced computational and memory demands as compared to use of classical approaches.

0 citations0 reviewsNo code linkedOpen access
preprint2016arXiv

Hybrid and iteratively reweighted regularization by unbiased predictive risk and weighted GCV for projected systems

Tikhonov regularization for projected solutions of large-scale ill-posed problems is considered. The Golub-Kahan iterative bidiagonalization is used to project the problem onto a subspace and regularization then applied to find a subspace approximation to the full problem. Determination of the regularization parameter for the projected problem by unbiased predictive risk estimation, generalized cross validation and discrepancy principle techniques is investigated. It is shown that the regularized parameter obtained by the unbiased predictive risk estimator can provide a good estimate for that to be used for a full problem which is moderately to severely ill-posed. A similar analysis provides the weight parameter for the weighted generalized cross validation such that the approach is also useful in these cases, and also explains why the generalized cross validation without weighting is not always useful. All results are independent of whether systems are over or underdetermined. Numerical simulations for standard one dimensional test problems and two dimensional data, for both image restoration and tomographic image reconstruction, support the analysis and validate the techniques. The size of the projected problem is found using an extension of a noise revealing function for the projected problem Hn\u etynková, Ple\u singer, and Strako\u s, [\textit{BIT Numerical Mathematics} {\bf 49} (2009), 4 pp. 669-696.]. Furthermore, an iteratively reweighted regularization approach for edge preserving regularization is extended for projected systems, providing stabilization of the solutions of the projected systems and reducing dependence on the determination of the size of the projected subspace.

0 citations0 reviewsNo code linkedOpen access
preprint2014arXiv

Application of the $χ^2$ principle and unbiased predictive risk estimator for determining the regularization parameter in 3D focusing gravity inversion

The $χ^2$ principle and the unbiased predictive risk estimator are used to determine optimal regularization parameters in the context of 3D focusing gravity inversion with the minimum support stabilizer. At each iteration of the focusing inversion the minimum support stabilizer is determined and then the fidelity term is updated using the standard form transformation. Solution of the resulting Tikhonov functional is found efficiently using the singular value decomposition of the transformed model matrix, which also provides for efficient determination of the updated regularization parameter each step. Experimental 3D simulations using synthetic data of a dipping dike and a cube anomaly demonstrate that both parameter estimation techniques outperform the Morozov discrepancy principle for determining the regularization parameter. Smaller relative errors of the reconstructed models are obtained with fewer iterations. Data acquired over the Gotvand dam site in the south-west of Iran are used to validate use of the methods for inversion of practical data and provide good estimates of anomalous structures within the subsurface.

0 citations0 reviewsNo code linkedOpen access
preprint2014arXiv

Automatic estimation of the regularization parameter in 2-D focusing gravity inversion: an application to the Safo manganese mine in northwest of Iran

We investigate the use of Tikhonov regularization with the minimum support stabilizer for underdetermined 2-D inversion of gravity data. This stabilizer produces models with non-smooth properties which is useful for identifying geologic structures with sharp boundaries. A very important aspect of using Tikhonov regularization is the choice of the regularization parameter that controls the trade off between the data fidelity and the stabilizing functional. The L-curve and generalized cross validation techniques, which only require the relative sizes of the uncertainties in the observations are considered. Both criteria are applied in an iterative process for which at each iteration a value for regularization parameter is estimated. Suitable values for the regularization parameter are successfully determined in both cases for synthetic but practically relevant examples. Whenever the geologic situation permits, it is easier and more efficient to model the subsurface with a 2-D algorithm, rather than to apply a full 3-D approach. Then, because the problem is not large it is appropriate to use the generalized singular value decomposition for solving the problem efficiently. The method is applied on a profile of gravity data acquired over the Safo mining camp in Maku-Iran, which is well known for manganese ores. The presented results demonstrate success in reconstructing the geometry and density distribution of the subsurface source.

0 citations0 reviewsNo code linkedOpen access
preprint2014arXiv

Regularization Parameter Estimation for Underdetermined problems by the $χ^2$ principle with application to $2D$ focusing gravity inversion

The $χ^2$-principle generalizes the Morozov discrepancy principle (MDP) to the augmented residual of the Tikhonov regularized least squares problem. Weighting of the data fidelity by a known Gaussian noise distribution on the measured data, when the regularization term is weighted by unknown inverse covariance information on the model parameters, the minimum of the Tikhonov functional is a random variable following a $χ^2$-distribution with $m+p-n$ degrees of freedom, model matrix $G:$ $m \times n$ and regularizer $L:$ $p\times n$. It is proved that the result holds also for $m<n$ when $m+p\ge n$. A Newton root-finding algorithm is used to find the regularization parameter $α$ which yields the optimal inverse covariance weighting in the case of a white noise assumption on the mapped model data. It is implemented for small-scale problems using the generalized singular value decomposition. Numerical results verify the algorithm for the case of regularizers approximating zero to second order derivative approximations, contrasted with the methods of generalized cross validation and unbiased predictive risk estimation. The inversion of underdetermined $2D$ focusing gravity data produces models with non-smooth properties, for which typical implementations in this field use the iterative minimum support (MS) stabilizer and both regularizer and regularizing parameter are updated each iteration. For a simulated data set with noise, the regularization parameter estimation methods for underdetermined data sets are used in this iterative framework, also contrasted with the L-curve and MDP. Experiments demonstrate efficiency and robustness of the $χ^2$-principle, moreover the L-curve and MDP are generally outperformed. Furthermore, the MS is of general use for the $χ^2$-principle when implemented without the knowledge of a mean value of the model.

0 citations0 reviewsNo code linkedOpen access