Stewart Greenhalgh

A New Vector Waveform Inversion Algorithm for Simultaneous Updating of Conductivity and Permittivity Parameters From Combination Crosshole/Borehole-to-Surface GPR Data

We have developed a new full-waveform groundpenetrating radar (GPR) multicomponent inversion scheme for imaging the shallow subsurface using arbitrary recording configurations. It yields significantly higher resolution images than conventional tomographic techniques based on first-arrival times and pulse amplitudes. The inversion is formulated as a nonlinear least squares problem in which the misfit between observed and modeled data is minimized. The full-waveform modeling is implemented by means of a finite-difference time-domain solution of Maxwells equations. We derive here an iterative gradient method in which the steepest descent direction, used to update iteratively the permittivity and conductivity distributions in an optimal way, is found by cross-correlating the forward vector wavefield and the backward-propagated vectorial residual wavefield. The formulation of the solution is given in a very general, albeit compact and elegant, fashion. Each iteration step of our inversion scheme requires several calculations of propagating wavefields. Novel features of the scheme compared to previous full-waveform GPR inversions are as follows: 1) The permittivity and conductivity distributions are updated simultaneously (rather than consecutively) at each iterative step using improved gradient and step length formulations; 2) the scheme is able to exploit the full vector wavefield; and 3) various data sets/survey types (e.g., crosshole and borehole-to-surface) can be individually or jointly inverted. Several synthetic examples involving both homogeneous and layered stochastic background models with embedded anomalous inclusions demonstrate the superiority of the new scheme over previous approaches.

Zonation for 3D aquifer characterization based on joint inversions of multimethod crosshole geophysical data

Predictive groundwater modeling requires accurate information about aquifer characteristics. Geophysical imaging is a powerful tool for delineating aquifer properties at an appropriate scale and resolution, but it suffers from problems of ambiguity. One way to overcome such limitations is to adopt a simultaneous multitechnique inversion strategy. We have developed a methodology for aquifer characterization based on structural joint inversion of multiple geophysical data sets followed by clustering to form zones and subsequent inversion for zonal parameters. Joint inversions based on cross-gradient structural constraints require less restrictive assumptions than, say, applying predefined petrophysical relationships and generally yield superior results. This approach has, for the first time, been applied to three geophysical data types in three dimensions. A classification scheme using maximum likelihood estimation is used to determine theparameters of a Gaussian mixture model that defines zonal geometries ...

3D crosshole ERT for aquifer characterization and monitoring of infiltrating river water

The hydrogeological properties and responses of a productive aquifer in northeastern Switzerland are investigated. For this purpose, 3D crosshole electrical resistivity tomography (ERT) is used to define the main lithological structures within the aquifer (through static inversion) and to monitor the water infiltration from an adjacent river. During precipitation events and subsequent river flooding, the river water resistivity increases. As a consequence, the electrical characteristics of the infiltrating water can be used as a natural tracer to delineate preferential flow paths and flow velocities. The focus is primarily on the experiment installation, data collection strategy, and the structural characterization of the site and a brief overview of the ERT monitoring results. The monitoring system comprises 18 boreholes each equipped with 10 electrodes straddling the entire thickness of the gravel aquifer. A multichannel resistivity system programmed to cycle through various four-point electrode configurations of the 180 electrodes in a rolling sequence allows for the measurement of approximately 15,500 apparent resistivity values every 7 h on a continuous basis. The 3D static ERT inversion of data acquired under stable hydrological conditions provides a base model for future time-lapse inversion studies and the means to investigate the resolving capability of our acquisition scheme. In particular, it enables definition of the main lithological structures within the aquifer. The final ERT static model delineates a relatively high-resistivity, low-porosity, intermediate-depth layer throughout the investigated aquifer volume that is consistent with results from well logging and seismic and radar tomography models. The next step will be to define and implement an appropriate time-lapse ERT inversion scheme using the river water as a natural tracer. The main challenge will be to separate the superposed time-varying effects of water table height, temperature, and salinity variations associated with the infiltrating water.

Solutions, algorithms and inter-relations for local minimization search geophysical inversion

This paper presents in a general form the most popular local minimization search solutions for geophysical inverse problems—the Tikhonov regularization solutions, the smoothest model solutions and the subspace solutions, from which the inter-relationships between these solutions are revealed. For the Tikhonov regularization solution, a variety of forms exist—the general iterative formula, the iterative linearized scheme, the Levenberg–Marquardt version, the conjugate gradient solver (CGS) and local-search quadratic approximation CGS. It is shown here that the first three solutions are equivalent and are just a specified form of the gradient solution. The local-search quadratic approximation CGS is shown to be a more general form from which these three solutions emerge. The smoothest model solution (Occams inversion) is a subset of the Tikhonov regularization solutions, for which a practical algorithm generally comprises two parts: (1) determine the subset of the Tikhonov regularization solutions which satisfy the desired tolerance of data fit, (2) choose the solution which best fits the data in the subset. The solution procedure is actually an adaptive Tikhonov regularization solution. The subspace solution is mathematically no more than a dimension-decreased transform of model parameterization in the Tikhonov regularization solution or the smoothest model solution. The crucial step is to choose the basis vectors which span the whole model space in the parameterization transform. A simple form of the transform is developed which yields greater flexibility for geophysical inverse applications.

'Shortest path' ray tracing for most general 2D/3D anisotropic media

This paper presents a simple method for seismic ray tracing in a general anisotropic medium, which may include complex structures and compound materials, such as water, isotropic and anisotropic rocks, fine-layers and parallel small cracked blocks. The anisotropy may be defined by up to 21 density-normalized elastic moduli which vary with spatial position. The method presented is a direct extension of the irregular network shortest path method for an isotropic solid. For this extension, we first apply analytic solutions of the wave velocities (phase velocity and group velocity) for a general anisotropic medium as a transform or mapping operator to convert the elastic-moduli-described medium into the direction-dependent group-velocity models for the three independent wave modes (qP, qS1, qS2). We then utilize the shortest path method to trace raypaths through such group-velocity models for the three modes. We also give an alternative derivation of Fermats principle of stationary time in anisotropic media. With this method, the travel times and ray paths of the first arrivals emanating from a source to multiple receivers can be simultaneously obtained for the three modes. Some 2D/3D numerical experiments are performed to show the accuracy and applicability of the method. From these results, one can see that the method may be applied to kinematic modelling and inversion in 2D and 3D seismic or seismological applications.

Nonlinear traveltime inversion scheme for crosshole seismic tomography in tilted transversely isotropic media

Crosshole seismic tomography often is applied to image the velocity structure of an interwell medium. If the rocks are anisotropic, the tomographic technique must be adapted to the complex situation; otherwise, it leads to a false interpretation. We propose a nonlinear kinematic inversion method for crosshole seismic tomography in composite transversely isotropic media with known dipping symmetry axes. This method is based on a new version of the first-order traveltime perturbation equation. It directly uses the derivative of the phase velocity rather than the eigenvectors of the body-wave modes to overcome the singularity problem for application to the two quasi-shear waves. We applied an iterative nonlinear solver incorporating our kinematic ray-tracing scheme and directly compute the Jacobian matrix in an arbitrary reference medium. This reconstructs the five elastic moduli or Thomsen parameters from the first-arrival traveltimes of the three seismic body waves (qP, qSV, qSH) in strongly and weakly anisotropic media. We conducted three synthetic experiments that involve determining anisotropic parameters for a homogeneous rock, reconstructing a fault embedded in a strongly anisotropic background, and imaging a complicated four-layer model containing a small channel and a buried dipping interface. We compared results of our nonlinear inversion method with isotropic tomography and the traditional linear anisotropic inversion scheme, which showed the capability and superiority of the new scheme for crosshole tomographic imaging.

On the computation of elastic wave group velocities for a general anisotropic medium

This paper deals exclusively with the computation of the group velocities for the three wave modes (qP, qS1, qS2) in a general anisotropic medium, which may involve up to 21 density-normalized elastic moduli. We tackled the shear-wave singularity problem through two independent approaches: (1) an eigenvalue method, and (2) an eigenvector method. In the former, we derive analytic formulae and introduce an approximation for the directional derivative of the phase velocity at the singularity points. In the latter, we develop two simple schemes to find the eigenvectors of the quasi-shear waves at the singularity points. Computational experiments have been conducted to show the merits and validity of both approaches. Furthermore, the numerical results demonstrate that both methods produce consistent and satisfactory results for any degree of anisotropic media, notwithstanding the possible discrepancy between the specific solution for a general TI medium and the general solution for arbitrary anisotropy. The former is more suitable than the general solutions for such media, because it achieves complete polarization discrimination of the qSV and qSH modes. For more complex forms on anisotropy, e.g. orthorhombic, the general solutions yield mixed versions of the two quasi-shear waves qS1 and qS2, which have many singularity points in the phase velocity space, but overall recover the true modes.

Comparison of source-independent methods of elastic waveform inversion

In this paper, we investigate several source-independent methods of nonlinear full-waveform inversion of multicomponent elastic-wave data. This includes iterative estimation of source signature (IES), standard trace normalization (STN), and average trace normalization (ATN) inversion methods. All are based on the finite-element method in the frequency domain. One synthetic elastic crosshole model is used to compare the recovered images with all these methods as well as the known source signature (KSS) inversion method. The numerical experiments show that the IES method is superior to both STN and ATN methods in two-component, elastic-wave inversion in the frequency domain when the source signature is unknown. The STN and ATN methods have limitations associated with near-zero amplitudes (or polarity reversals) in traces from one of the components, which destroy the energy balance in the normalized traces and cause a loss of frequency information. But the ATN method is somewhat superior to the STN method in suppressing random noise and improving stability, as the developed formulas and the numerical experiments show. We suggest the IES method as a practical procedure for multicomponent seismic inversion.

Precision analysis of first-break times in grid models

Rectangular grid velocity models and their derivatives are widely used in geophysical inversion techniques. Specifically, seismic tomographic reconstruction techniques, whether they be based on raypath methods (Bregman et al., 1989; Moser, 1991; Schneider et al., 1992; Cao and Greenhalgh, 1993; Zhou, 1993) or full wave equation methods (Vidale, 1990; Qin and Schuster, 1993; Cao and Greenhalgh, 1994) for calculating synthetic arrival times, involve propagation through a grid model. Likewise, migration of seismic reflection data, using asymptotic ray theory or finite difference/pseudospectral methods (Stolt and Benson, 1986; Zhe and Greenhalgh, 1997) involve assigning traveltimes to upward and downward propagating waves at every grid point in the model. The traveltimes in both cases depend on the grid specification. However, the precision level of such numerical models and their dependence on the model parameters is often unknown. In this paper, we describe a two‐dimensional velocity model and derive an err...

Multichannel, full waveform and flexible electrode combination resistivity-imaging system

The development of a new multichannel and multielectrode, full waveform resistivity acquisition system is presented. This system, which can be built in any electronics workshop without major difficulties, allows for large data collection capacity and flexibility but without the sacrifice of channels, as in most other systems that use a preselected array. A superior interpretation quality is enabled by the increased data and the monitoring of noise through full waveform display. The proposed hardware is robust because no electronics need to be installed in the cable. Therefore, it is suitable for 2D or 3D surface or crosshole tomography. A major feature of the data collection strategy is to capture all of the electrode combinations subject to having one current electrode in each cable. This is efficient because there is no extra cost in doing this. It is an automatic process; one does not have to make any prior decisions as to what electrode combinations to use. The data sets collected with this strategy are suitable for an inversion process involving a large amount of data because the number of voltage values measured using a couple of current electrodes is maximized. Synthetic modeling shows that the proposed approach allows for better imaging of the subsurface than with conventional arrays and the results are comparable to results obtained with an optimized experimental design strategy. A field experiment illustrates the effectiveness of the presented strategy. The inverted image of the subsurface is well-correlated with the available borehole information.