Structural Similarity in Joint Inverse Problems

A novel 3D simultaneous joint inversion scheme for gravity and seismic travel time data is developed to solve for near-surface complex velocity distributions. The method incorporates industry-standard gravity and travel time inversion techniques while the joint inversion problem is solved by the introduction of various regularization functions about the model such as a-priori parameter distribution information, solution-space bounds, structural similarity via cross-gradient constraints and rock physics relations. The effectiveness of our joint inversion is demonstrated against a synthetic model representing a complex pattern of near-surface anomalies incorporating low and high velocity and density bodies. Results demonstrate the superiority of our approach where the shallow anomalies are better reconstructed by the joint inversion rather than that obtained by the single-domain inversions. The developed algorithm is tested with real data from Saudi Arabia acquired over a wadi structure. The results show a significant uplift of the time stack using the seismic-gravity joint inversion velocity model. The developed methodology is part of a multi-geophysics platform for near-surface velocity model building in complex geology scenarios.

SummaryIn geophysical studies the problem of joint inversion of multiple experimental data sets obtained by different methods is conventionally considered as a scalar one. Namely, a solution is found by minimization of linear combination of functions describing the fit of the values predicted from the model to each set of data. In the present paper we demonstrate that this standard approach is not always justified and propose to consider a joint inversion problem as a multiobjective optimization problem (MOP), for which the misfit function is a vector. The method is based on analysis of two types of solutions to MOP considered in the space of misfit functions (objective space). The first one is a set of complete optimal solutions that minimize all the components of a vector misfit function simultaneously. The second one is a set of Pareto optimal solutions, or trade-off solutions, for which it is not possible to decrease any component of the vector misfit function without increasing at least one other. We investigate connection between the standard formulation of a joint inversion problem and the multiobjective formulation and demonstrate that the standard formulation is a particular case of scalarization of a multiobjective problem using a weighted sum of component misfit functions (objectives). We illustrate the multiobjective approach with a non-linear problem of the joint inversion of shear wave splitting parameters and longitudinal wave residuals. Using synthetic data and real data from three passive seismic experiments, we demonstrate that random noise in the data and inexact model parametrization destroy the complete optimal solution, which degenerates into a fairly large Pareto set. As a result, non-uniqueness of the problem of joint inversion increases. If the random noise in the data is the only source of uncertainty, the Pareto set expands around the true solution in the objective space. In this case the ‘ideal point’ method of scalarization of multiobjective problems can be used. If the uncertainty is due to inexact model parametrization, the Pareto set in the objective space deviates strongly from the true solution. In this case all scalarization methods fail to find the solution close to the true one and a change of model parametrization is necessary.

Inverts permittivity and conductivity with structural constraint in GPR FWI based on truncated Newton method

Harmonic pumping tests consist in stimulating an aquifer by the means of hydraulic stimulations at some discrete frequencies. The inverse problem consisting in retrieving the hydraulic properties is inherently ill posed and is usually underdetermined when considering the number of well head data available in field conditions. To better constrain this inverse problem, we add self‐potential data recorded at the ground surface to the head data. The self‐potential method is a passive geophysical method. Its signals are generated by the groundwater flow through an electrokinetic coupling. We showed using a 3‐D saturated unconfined synthetic aquifer that the self‐potential method significantly improves the results of the harmonic hydraulic tomography. The hydroelectric forward problem is obtained by solving first the Richards equation, describing the groundwater flow, and then using the result in an electrical Poisson equation describing the self‐potential problem. The joint inversion problem is solved using a reduction model based on the principal component geostatistical approach. In this method, the large prior covariance matrix is truncated and replaced by its low‐rank approximation, allowing thus for notable computational time and storage savings. Three test cases are studied, to assess the validity of our approach. In the first test, we show that when the number of harmonic stimulations is low, combining the harmonic hydraulic and self‐potential data does not improve the inversion results. In the second test where enough harmonic stimulations are performed, a significant improvement of the hydraulic parameters is observed. In the last synthetic test, we show that the electrical conductivity field required to invert the self‐potential data can be determined with enough accuracy using an electrical resistivity tomography survey using the same electrodes configuration as used for the self‐potential investigation.

Minimum-structure, or Occam’s, style of inversion deals with the fundamental non-uniqueness of the inverse problem by finding the simplest Earth model that reproduces the observations. As an additional consequence of this approach, minimum-structure inversion is also reliable and robust. Because of these reasons, it has been extensively utilized in mineral and petroleum exploration problems, and lithospheric studies. The method has been adapted and extended in many ways to obtain more reliable and realistic models of the Earth’s subsurface. Joint inversion of geophysical data-sets is one of the most important extensions of minimum-structure inversion. This method can reduce the non-uniqueness of the inverse problem by combing two, or more, different geophysical data-sets in a single inverse problem. Different geophysical methods have different sensitivity to different physical properties, hence, it is hoped that the null space for one type of data can be spanned by the other.Joint inversion algorithms can be divided into two main categories, structural-based and petrophysical-based joint inversion methods, depending on the coupling measure used between the physical property models. We have adopted the fuzzy c-mean (FCM) clustering technique which is a petrophysical-based method to jointly invert MT and gravity data-sets. The optimization of this method is not as challenging as for structural-based approaches. We have also performed constrained FCM clustering for independent MT and gravity inversions to compare the constructed models of this method with the joint inversion, and independent MT and gravity inversions. The FCM clustering method makes effective use of statistical petrophysical data which may exist in complex geological structures, or can be anticipated, to encourage the inverted physical property values to move towards the a priori petrophysical data as target clusters.The capabilities of the joint and constrained FCM clustering inversion are evaluated on synthetic and real examples. The constructed density and conductivity models from the joint inversion have a more plausible representation of the true model’s geometry and have a reasonable range of the recovered physical property values compared to the independent constrained FCM clustering technique and independent unconstrained MT and gravity inversions.Keywords: Fuzzy c-mean clustering, Gravity, Joint inverse problem, MT, Unstructured tetrahedral mesh

Stochastic optimization is key to efficient inversion in PDE-constrained optimization. Using ‘simultaneous shots’, or random superposition of source terms, works very well in simple acquisition geometries where all sources see all receivers, but this rarely occurs in practice.We develop an approach that interpolates data to an ideal acquisition geometry while solving the inverse problem using simultaneous shots. The approach is formulated as a joint inverse problem, combining ideas from low-rank interpolation with full-waveform inversion. Results using synthetic experiments illustrate the flexibility, efficiency and potential of the approach, especially in comparison to naive workflows where data are first interpolated and then inverted by simultaneous shots. The new approach effectively uses the wave equation as a prior on the data, and can handle high subsampling (85% missing data), noisy data, and both random and periodic sampling. These results suggest that the joint is promising for realistic acquisition scenarios.

History matching of 4D seismic and well production data has been developed recently for reservoir characterization. This is a data driven optimization process to derive reservoir model parameters and constitutes a joint inverse problem. Considering the inherent non‐uniqueness in the inverse problem and the unique feature of Bayesian inference in data integration and uncertainty analysis, this joint inverse problem is formulated in a Bayesian framework and solved stochastically by reconstructing the posterior probability density (PPD) surface using a new multi‐scale Markov Chain Monte Carlo (MCMC) algorithm. In this new MCMC method, a technique of multi‐scaling is used to take advantage of the benefits from both the fine scale model and the coarse scale model. Although the coarse scale does not provide reliable information about the model, it helps speed up the convergence of the fine scale model to a good estimation, and, by exchanging information between the fine and coarse scales, works like a regularization operator to smooth the fine scale model and make it more realistic. The resulting PPD samples can also be used to quantify corresponding uncertainties in order to facilitate risk assessment associated with reservoir decision making and management. We use a numerical example to demonstrate how we derive reservoir's static and dynamic properties as well as quantify uncertainties in a Bayesian framework using the multi‐scale MCMC algorithm.

Rayleigh wave dispersion curves can be inverted to retrieve subsurface seismic velocity profiles. The inverse problem is ill‐posed, nonlinear and poorly conditioned, necessitating the application of global optimization methods. We present the application of the multi‐objective grey wolf optimization algorithm to perform joint inversion of the phase velocity dispersion curves corresponding to the fundamental and higher order modes of Rayleigh waves to obtain shear (S‐) and primary (P‐) wave velocity profiles. Multi‐objective grey wolf optimization is an extension of the grey wolf optimization algorithm for application to multi‐objective optimization and can be adapted to solve joint inversion problems. We compare the joint inversion results obtained from the multi‐objective grey wolf optimizer with those obtained from Markov chain Monte Carlo and fundamental mode inversion using the grey wolf optimizer on synthetic examples. The errors associated with phase velocity measurements are simulated by adding frequency‐dependent noise, with a higher level of noise added to the phase velocities corresponding to lower frequencies as compared to the higher frequencies. In the multimode joint inversion problem, the multi‐objective grey wolf optimizer gives a suite of solutions corresponding to each model parameter. The suite of S‐wave and P‐wave velocity profiles estimated from the multi‐objective grey wolf optimizer matches closely with the true model for the synthetic case studies even in the presence of noise. However, the suite of solutions has a greater spread for the last few layers, qualitatively indicating a higher degree of uncertainty in the predicted model parameter. The uncertainty in the solution for the deeper layers is a consequence of the uncertainty in the phase velocity at lower frequencies. We demonstrate the efficacy of the algorithm on recorded data from a shallow seismic survey conducted at the Indian Institute of Technology Bombay. The results from the multi‐objective grey wolf optimizer are in close agreement with those from Markov chain Monte Carlo, and the depth of investigation is found to be greater in comparison to results from refraction traveltime inversion.

In carbon capture and storage projects, and in unconventional plays, microseismic monitoring and optical fiber are critical components of the measure, monitor, and verify value chain. A velocity model is required to estimate source location (hypocenter), source parameters, and source mechanism of a detected microseismic event. Incorrect event locations are often the result of an inaccurate knowledge of the velocity model. We propose a new method to simultaneously invert for the hypocenter and the velocity model to provide a robust long-term microseismic monitoring workflow. Such problem has been studied in several areas of seismology over the last few decades. However, those studies focusing on large-scale earthquakes have remained of limited interest to reservoir-scale applications such as short-term and long-term microseismic monitoring and induced seismic monitoring. In such domains, the integration of sonic log-derived information into the joint inversion problem is critical as the scale is fundamentally different. Our algorithm respects the resolution of sonic measurement, while it calibrates the wavelength where microseismic data have sensitivity: the number of unknown parameters in the velocity model is decoupled from the number of layers included in the model. Therefore, we can solve for the velocity and event location inversion problem in a stable manner while respecting the resolution of the initial velocity model. In the present article, we introduce the science and technique behind the simultaneous inversion for the hypocenter and the velocity model and share case study applications based on a synthetic dataset and a real monitoring campaign.

Marchenko focusing and imaging are novel methods for correctly handling multiple scattered energy while processing seismic data. However, strict requirements in the acquisition geometry, specifically the colocation of sources and receivers as well as dense and regular sampling, currently constrain their practical applicability. We have reformulated the Marchenko equations to handle the case in which there are gaps in the source geometry while the receiver sampling remains regular (or the opposite, by means of reciprocity). Using synthetic data based on a velocity model that produces strong interbed multiples, we test different solvers for the newly formulated inversion problem, and we compare these results to results obtained by applying standard Marchenko inversion to a previously reconstructed data set. When using the unreconstructed data set, the ability of the Marchenko equations to retrace multiple reflected energy deteriorates. Sparsity-promoting Marchenko inversion, while improving the appearance of focusing functions, barely decreases multiple leakage in gathers and does not visibly improve the final image when compared to standard least-squares inversion. On the other hand, reconstructing the wavefield in advance restores the proper functioning of the Marchenko methods. Further, we test a joint inversion technique designed for time-lapse data with nonrepeated geometries and originally intended to be solved using sparsity-promoting inversion. Motivated by our previous results, we compare images produced by this method to images produced by solving the same joint inversion problem without sparsity constraint. We find that the joint inversion alone hardly improves the resulting images but, when combined with the sparsity constraint, it leads to better noise and multiple suppression and produces a clean time-lapse image. Overall, none of the results from sparsity-promoting inversion techniques match the results obtained when reconstructing the wavefield in advance. We show that this can be explained by the comparatively slow convergence rate of the sparsity-promoting Marchenko inversion.

Muography is a relatively new geophysical imaging method that uses muons to provide estimates of average densities along particular lines of sight. Muography can only see above the horizontal elevation of the detector and it is therefore attractive to attempt a joint inversion of muography data with gravity data, which is also responsive to density but generally requires combination with another geophysical data set to overcome issues related to non-uniqueness and poor depth resolution. Some previous work has investigated this joint inverse problem and demonstrated the potential improvements to be gained by jointly inverting muography and gravity data. However, there has yet to be a thorough investigation of different numerical approaches for formulating the joint inverse problem. Particularly important is how to account for the fact that even though the two data types are sensitive to the same physical quantity, density, they respond through different response functions. Moreover, the two measurements are affected by different systematic uncertainties that are difficult to model. In this work, we considered an approximation where the two density quantities, inferred from the two data types, can be related by an unknown scalar offset. We considered various existing and new joint inversion methods that might solve this problem and we applied them to a synthetic volcano imaging scenario based on the Puy de Dome volcano in the Central Massif region of France. We used unstructured meshes in our modelling to adequately honour the significant topography in that scenario. Our experiments indicated that the most successful joint inversion method for this type of geological scenario was one in which the data misfit function is reformulated to automatically determine the best-fitting offset following a least-squares minimization argument. However, other approaches showed merit and we suggest several of the investigated methods be applied and compared for any specific joint inversion scenario.

Particle swarm optimization (PSO) is a relatively new global optimization approach inspired by the social behavior of birds and fishes. Although this approach has proven to provide excellent convergence rates in different optimization problems, it has seldom been used in geophysical inversion. Here, we propose a PSO-based inversion strategy to jointly invert GPR and P-wave seismic traveltimes from co-located crosshole experiments. Using a synthetic data example, we demonstrate the potential of our approach. Comparing our results to the input models as well as to velocity models found by separately inverting the data using a standard linearized inversion approach, illustrates the benefits of using an efficient global optimization approach for such a joint inversion problem. These include a straightforward appraisal of uncertainty, non-uniqueness, and resolution issues as well as the possibility of an improved and more objective interpretation.

Wave attenuation is an important physical property of hydrocarbon-bearing sediments that is rarely taken into account in site characterization with seismic data. We present a 1D viscoelastic waveform inversion scheme for determining the quality factor [Formula: see text] from the normal-incidence surface seismic and zero-offset vertical seismic profile (VSP) data simultaneously. The joint inversion problem is solved by the damped least-squares method, and the inversion result is successful using synthetic data. The effects of initial model thickness, [Formula: see text] value, and the existence of noise were studied through a synthetic example. By extracting all of the information contained in the waveforms, the waveform inversion of the seismic reflection and transmission data becomes a powerful tool for estimating [Formula: see text]. For a more comprehensive image of [Formula: see text], the tomographic inversion of [Formula: see text] is applied to the walkaway VSP and prestack surface seismic data, using the waveform inversion result as the initial model. Results from applying the method to a real seismic line and zero-offset VSP data from the Nanyang oilfield, central China, indicate that [Formula: see text] from the tomographic inversion of reflection and transmission data contains useful information on medium properties, which can aid in reservoir appraisal.

Two related formulations are proposed for target-oriented joint least-squares migration/inversion of time-lapse seismic data sets. Time-lapse seismic images can be degraded when reservoir overburden is complex or when acquisition geometries significantly differ, because the migration operator does not compensate for the resulting amplitude and phase distortions. Under these circumstances, time-lapse amplitudes are poor indicators of production-related changes in reservoir properties. To correct for such image degradation, time-lapse imaging is posed as joint inverse problems that utilize concatenations of target-oriented approximations to the linear least-squares imaging Hessian. In both formulations, spatial and temporal constraints ensure inversion stability and geologically consistent time-lapse images. Using two numerical time-lapse data sets, we confirmed that these formulations can attenuate illumination artifacts caused by complex overburden or geometry differences, and that they yield better-quality images than obtainable with migration.

Geophysical monitoring of geologic carbon sequestration is critical for risk assessment during and after carbon dioxide (CO2) injection. Integration of multiple geophysical measurements is a promising approach to achieve high‐resolution reservoir monitoring. However, joint inversion of large geophysical data is challenging due to high computational costs and difficulties in effectively incorporating measurements from different sources and with different resolutions. This study develops a differentiable physics model for large‐scale joint inverse problems with reparameterization of model variables by neural networks and implementation of a differentiable programming approach of the forward model. The proposed physics‐informed neural network model is completely differentiable and thus enables end‐to‐end training with automatic differentiation for multi‐objective optimization by multiphysics data assimilation. The application to the Sleipner benchmark model demonstrates that the proposed method is effective in estimation of reservoir properties from seismic and resistivity data and shows promising results for CO2 storage monitoring. Moreover, the global parameters that are assumed to be uncertain in the rock‐physics model are accurately quantified by integration of a Bayesian neural network.