Table of contents

Volume 16

Number 5, October 2000

Previous issue Next issue

SPECIAL SECTION: ELECTROMAGNETIC IMAGING AND INVERSION OF THE EARTH'S SUBSURFACE (INVITED PAPERS)

001

and

In the past couple of decades, there have been considerable advances in the efficiency, robustness and efficacy of inversion methods in almost every field. This is apparent in the many journal articles published on this subject, the ever-growing number of specialized sessions organized in major conferences, and in the several books and proceedings that have appeared in the past few years. Inverse scattering continues to be a primary area of research to a large number of scientists and practitioners.

Many developments have been driven by several new applications and some old ones, such as mathematical physics, atmospheric sciences, geophysical prospecting, quantum mechanics, remote sensing, underwater acoustics, nondestructive testing and evaluation, medical imaging, to mention only a few. These developments have been applied to a variety of data that include direct current (DC), monochromatic, multi-frequency and transients.

The majority of the effort has gone into improving on the robustness and speed of the various inversion techniques. A considerable effort has also gone into addressing the nonuniqueness and ill-conditioning that are inherent in almost all inverse scattering problems which arise mainly as a consequence of the inadequacy and redundancy in measurements. Many researchers have attempted to quantify these uncertainties by resorting to probabilistic methods.

This special section is dedicated to the imaging and inversion of the Earth's subsurface by electromagnetic means with special emphasis, though nonexclusive, toward oil and mineral exploration, evaluation and monitoring. The section includes 14 invited papers that cover a wide range of topics related to this theme. This involves 30 authors from 21 academic, governmental and industrial research groups worldwide. The subjects covered in this section encompass theoretical and computational developments as well as experimental results. In addition to the original and specialized contributions presented in these papers, many authors have included an informative review of the most recent advances in the particular subject covered in their articles.

The papers in this special section are listed not in any particular order but, nevertheless, we attempted to group them according to the broad similar subjects they cover. Next, we give highlights of the various contributions.

de Hoop sheds some insight on three closely related topics: the first deals with how to accurately model transient electromagnetic fields in strongly heterogeneous and dispersive media by means of domain-integrated field equations; the second topic deals with how to preserve causality in both the direct and the parameter estimation problems; and the third is on the nonlinear inversion problem cast in terms of the contrast-source wavefield formulation.

Maurer, Boerner and Curtis start with the premise that the potential usefulness of the presently developed multidimensional electromagnetic inversion methods in geophysics is hindered by the deficiencies of critical information in the data. They review a number of linear and nonlinear design strategies for geophysical surveys that attempt to provide maximum information content within a given budget.

Dorn, Miller and Rappaport present a two-step shape reconstruction method for cross-well EM imaging. The first step recovers equivalent sources from which an initial permittivity distribution is inferred. The second step combines the 'level set technique' for representing shapes of the reconstructed objects and an 'adjoint field technique' for solving the nonlinear inverse problem.

Yagle develops a discrete layer stripping algorithm for the solution of a two-dimensional inverse conductivity problem. The fully nonlinear solution involves the transformation of the elliptic problem into an hyperbolic one similar to a one-dimensional Schrödinger equation but with a time-varying potential.

Bloemenkamp and van den Berg reconstruct two-dimensional dielectric objects buried in a dielectric half space from transient electric fields. The approach is based on the 'contrast source inversion method', originally developed for time-harmonic waves, after reformulating it directly in the time domain.

Saillard, Vincent and Micolau address the problem of reconstructing the boundary and inverting the permittivity of a bounded homogeneous cylindrical object surrounded by a randomly distributed set of point-like scatterers in a cross-well tomography setup.

Hjelt and Pirttijärvi have further studied the thin conducting plate model which has proven to be versatile in describing various relevant geological features in crystalline bedrocks. Using measurements in both the frequency and time domains they show how to assess by simple methods the interpretability of the various model parameters that define the plate.

Bourgeois, Suignard and Perrusson focus on mineral exploration where a highly conductive bounded orebody embedded in a much less conductive host Earth is to be appraised. They study a number of simple inverse models for the rapid and robust interpretation of diffusive electromagnetic response resulting from the interaction of the orebody with a primary source and observed along nearby boreholes.

Haber, Ascher and Oldenburg consider nonlinear optimization techniques where the forward problems are computed by using the discretized differential forms of Maxwell's equations. By employing such a sparse matrix formulation, the authors demonstrate how to carry out the full Newton iteration with only a modest additional cost as compared to the quasi-Newton method.

Zhang and Liu aim at the reconstruction of an axisymmetric conductivity distribution using single-well electromagnetic induction measurements by means of a two-step linear inversion scheme. The approach is based on a fast Fourier and Hankel transform which is preconditioned by the extended Born approximation.

Zhdanov and Hursan use the quasi-analytical approximation to significantly speed-up the three-dimensional inversion of low-frequency electromagnetic data in connection with geophysical exploration and prospecting. The authors show a number of examples, using multi-frequency, surface-to-surface configurations.

Alumbaugh applies linearized and nonlinear inversion techniques in a carefully set Bayesian framework to both synthetic and oil-field real data in two-dimensional cross-well scenarios. Using a Monte Carlo estimation technique and quadratic programming, the author is able to obtain estimates of parameter uncertainty and cross-correlations in the reconstructed images, which allows him to compare the performances of linearized and nonlinearized techniques.

Malinverno and Torres-Verdín introduce a Bayesian inference approach to map the boundary (saturation front) of the region swept by water around a water injection well in an oil reservoir. The solution methodology is applied to direct current data synthesized by an array of electrodes. The authors investigate the link between uncertainties in the location of the front and lack of background information such as thicknesses and resistivities of the reservoir layers as well as the presence of noise in the data.

Newman and Hoversten examine several strategies for solving two- and three-dimensional geophysical inverse problems which they illustrate using a number of experimental data acquired in cross-well configurations in connection with environmental and enhanced oil recovery applications. The analysis is pursued within a partial differential equation framework of the electromagnetic fields.

We hope that the exceptional list of original contributions included in this special section will give the interested reader a fairly good sample of the type of approaches and techniques that are successfully applied to the electromagnetic prospecting of the Earth's subsurface.

Our sincere thanks go to the Honorary Editors, Professors F Natterer and F A Grünbaum for giving us the unique opportunity to act as Guest Editors to this special section of Inverse Problems. The strong commitment and patience of the authors at each stage of the process are very much appreciated. Finally, an expression of gratitude goes to Ms Elaine Longden-Chapman the Publisher of Inverse Problems for her meticulous help in coordinating the many efforts that went into preparing this special section.

1083

In this contribution some `areas for exploration in electromagnetic modelling and inversion' are discussed, the focus being on applications in the evaluation and monitoring of subsurface fossil energy reservoirs in the Earth. Three topics receive attention: (1) the modelling of electromagnetic wavefields in strongly heterogeneous media (such as rock), for which the domain-integrated field equations method is presented as a tool, (2) the feasibility of carrying out inversion in the time Laplace-transformed domain at a set of real, positive, values of the transform parameter to which Lerch's uniqueness theorem of the one-sided, causal Laplace transformation applies, and (3) a general operator formalism that is based on the iterative decrease in the norm of the mismatch in data equation and object equation as they occur in the contrast source formulation of the inversion problem.

1097

, and

Acquiring information on the Earth's electric and magnetic properties is a critical task in many geophysical applications. Since electromagnetic (EM) geophysical methods are based on nonlinear relationships between observed data and subsurface parameters, designing experiments that provide the maximum information content within a given budget can be quite difficult. Using examples from direct-current electrical and frequency-domain EM applications, we review four approaches to quantitative experimental design. Repeated forward modelling is effective in feasibility studies, but may be cumbersome and time-consuming for studying complete data and model spaces. Examining Fréchet derivatives provides more insights into sensitivity to perturbations of model parameters, but only in the linear space around the trial model and without easily accounting for combinations of model parameters. A related sensitivity measure, the data importance function, expresses the influence each data point has on determining the final inversion model. It considers simultaneously all model parameters, but provides no information on the relative position of the individual points in the data space. Furthermore, it tends to be biased towards well resolved parts of the model space. Some of the restrictions of these three methods are overcome by the fourth approach, statistical experimental design. This robust survey planning method, which is based on global optimization algorithms, can be customized for individual needs. It can be used to optimize the survey layout for a particular subsurface structure and is an appropriate procedure for nonlinear experimental design in which ranges of subsurface models are considered simultaneously.

1119

, and

A two-step shape reconstruction method for electromagnetic (EM) tomography is presented which uses adjoint fields and level sets. The inhomogeneous background permittivity distribution and the values of the permittivities in some penetrable obstacles are assumed to be known, and the number, sizes, shapes, and locations of these obstacles have to be reconstructed given noisy limited-view EM data. The main application we address in the paper is the imaging and monitoring of pollutant plumes in environmental cleanup sites based on cross-borehole EM data. The first step of the reconstruction scheme makes use of an inverse scattering solver which recovers equivalent scattering sources for a number of experiments, and then calculates from these an approximation for the permittivity distribution in the medium. The second step uses this result as an initial guess for solving the shape reconstruction problem. A key point in this second step is the fusion of the `level set technique' for representing the shapes of the reconstructed obstacles, and an `adjoint field technique' for solving the nonlinear inverse problem. In each step, a forward and an adjoint Helmholtz problem are solved based on the permittivity distribution which corresponds to the latest best guess for the representing level set function. A correction for this level set function is then calculated directly by combining the results of these two runs. Numerical experiments are presented which show that the derived method is able to recover one or more objects with nontrivial shapes given noisy cross-borehole EM data.

1157

We develop a discrete layer-stripping algorithm for the 2D inverse conductivity problem. Unlike previous algorithms, this algorithm transforms the problem into a time-varying 1D Schrödinger equation inverse scattering problem, discretizes this problem and then solves the discrete problem exactly. This approach has three advantages: (i) the poor conditioning inherent in the problem is concentrated in the solution of a linear integral transform at the beginning of the problem, to which standard regularization techniques may be applied and (ii) feasibility conditions on the transformed data are obtained, satisfaction of which ensures that (iii) the solution of the discrete nonlinear inverse scattering problem is exact and stable. Other contributions include solution of discrete Schrödinger equation inverse potential problems with time-varying potentials by both layer-stripping algorithms and solution of nested systems of equations which amount to a time-varying discrete version of the Gel'fand-Levitan equation. An analytic and numerical example is supplied to demonstrate the operation of the algorithm.

1173

and

This paper describes a time-domain method for reconstructing the permittivity of a two-dimensional object buried in a half-space from knowledge of how the object reflects known incident radiation. Inspired by the contrast source inversion method for monochromatic waves and multi-frequency waves, we have explored the possibilities to formulate the contrast source inversion method directly in the time domain. Therefore, first the formulation of the forward problem is sought in terms of a space-time domain integral equation, and then, for a known permittivity, the space-time field solution is obtained iteratively, using a conjugate gradient space-time fast fourier technique. This conjugate gradient type of iteration is the key to arrive at an iterative algorithm for the inverse problem, where in each iteration, both the unknown contrast source and the unknown permittivity are updated from a minimization of a cost functional of two terms, one containing a normalized norm of the discrepancies between the modelled scattered field and the measured scattered field data, and the other containing the normalized norm of the errors in the satisfaction of the space-time domain integral equation in the object domain. Some preliminary numerical results show the effectiveness of the reconstruction algorithms.

1195

, and

This paper deals with the reconstruction of the boundary and the permittivity of homogeneous objects from synthetic data. The method is based on a rigorous integral method and a conjugate gradient algorithm requiring the solution of only two direct problems at each step. To model some noise resulting from the inhomogeneity of the background, this object is surrounded by a set of point-like scatterers, arbitrarily distributed between the two parallel lines of measurement. Since convergence of the algorithm strongly depends on the initial guess, a method to detect a target in a random medium has been adapted. Statistical numerical results for reconstruction in a crosswell configuration are also presented.

1209

and

A thin conducting plate has proven to be a versatile model in describing various relevant geologic features in crystalline bedrock. During the last few years we have developed several inversion codes for some commonly used electromagnetic (EM) measurements that are based on this model. Both when developing and when applying inversion codes, it is necessary to have a good understanding of how the properties of the model parameters affect the model response. We have studied the parameters of the plate model by using the misfit mapping technique and singular value decomposition. These studies include both time and frequency domain EM systems. We highlight the major differences between a plate model in free space and in conductive surroundings.

1225

, and

We have developed a number of simple inverse models for rapid geometric interpretation of the low-frequency (diffusive) electromagnetic field scattered by a bounded 3D heterogeneity in a conductive background. These models, based on the nature of the interaction of the source primary field with the anomalous body in its environment, consist mainly of magnetic and electric dipoles and their combinations, and of their respective finite-size equivalents, i.e. closed and open current circuits (filaments). These `equivalent-source' models have been validated on synthetic data calculated with `exact' 3D modelling codes for different anomalous targets of canonical shapes (parallelepipeds and spheres). In this paper we present the practical application of our dipolar models in mineral exploration and discuss the reliability of each kind of information derived from them. In particular it is seen that the 3D attitude of the target is always well retrieved, and that the 3D location of its centre is accurate to within a few per cent of its distance to the measuring profile. Though this specific inverse problem is generally ill posed due to limited spatial coverage of the observation grid - a common characteristic in geophysics - it is shown that in the near field of the target the cases of non-uniqueness are scarce and controllable, as opposed to what is expected in the far field.

1263

, and

This paper considers optimization techniques for the solution of nonlinear inverse problems where the forward problems, like those encountered in electromagnetics, are modelled by differential equations. Such problems are often solved by utilizing a Gauss-Newton method in which the forward model constraints are implicitly incorporated. Variants of Newton's method which use second-derivative information are rarely employed because their perceived disadvantage in computational cost per step offsets their potential benefits of faster convergence. In this paper we show that, by formulating the inversion as a constrained or unconstrained optimization problem, and by employing sparse matrix techniques, we can carry out variants of sequential quadratic programming and the full Newton iteration with only a modest additional cost. By working with the differential equation explicitly we are able to relate the constrained and the unconstrained formulations and discuss the advantages of each. To make the comparisons meaningful we adopt the same global optimization strategy for all inversions. As an illustration, we focus upon a 1D electromagnetic (EM) example simulating a magnetotelluric survey. This problem is sufficiently rich that it illuminates most of the computational complexities that are prevalent in multi-source inverse problems and we therefore describe its solution process in detail. The numerical results illustrate that variants of Newton's method which utilize second-derivative information can produce a solution in fewer iterations and, in some cases where the data contain significant noise, requiring fewer floating point operations than Gauss-Newton techniques. Although further research is required, we believe that the variants proposed here will have a significant impact on developing practical solutions to large-scale 3D EM inverse problems.

1281

and

We invert for the axisymmetric conductivity distribution from borehole electromagnetic induction measurements using a two-step linear inversion method based on a fast Fourier and Hankel transform (FFHT) enhanced extended Born approximation (EBA). In this method, the inverse problem is first cast as an under-determined linear least-norm problem for the induced electric current density; from the solution of this induced current density, the unknown conductivity distribution is then obtained by solving an over-determined linear problem using the newly developed, FFHT-enhanced EBA. Numerical results show that this inverse method is applicable to a very high conductivity contrast. It is a natural extension of the original two-step linear inversion method of Torres-Verdin and Habashy to axisymmetric media. In the first step, the CPU time costs O(N2). In the second step, the CPU time costs O(Nlog 2N) where N is the number of unknowns. Because of the FFHT algorithm, this inverse method is actually more efficient than the conventional, brute-force first-order Born approximation.

1297

and

In this paper we address one of the most challenging problems of electromagnetic (EM) geophysical methods: three-dimensional (3D) inversion of EM data over inhomogeneous geological formations. The difficulties in the solution of this problem are two-fold. On the one hand, 3D EM forward modelling is an extremely complicated and time-consuming mathematical problem itself. On the other hand, the inversion is an unstable and ambiguous problem. To overcome these difficulties we suggest using, for forward modelling, the new quasi-analytical (QA) approximation developed recently by Zhdanov et al (Zhdanov M S, Dmitriev V I, Fang S and Hursan G 1999 Geophysics at press). It is based on ideas similar to those developed by Habashy et al (Habashy T M, Groom R W and Spies B R 1993 J. Geophys. Res.98 1759-75) for a localized nonlinear approximation, and by Zhdanov and Fang (Zhdanov M S and Fang S 1996a Geophysics61 646-65) for a quasi-linear approximation. We assume that the anomalous electrical field within an inhomogeneous domain is linearly proportional to the background (normal) field through a scalar electrical reflectivity coefficient, which is a function of the background geoelectrical cross-section and the background EM field only. This approach leads to construction of the QA expressions for an anomalous EM field and for the Frechet derivative operator of a forward problem, which simplifies dramatically the forward modelling and inversion. To obtain a stable solution of a 3D inverse problem we apply the regularization method based on using a focusing stabilizing functional introduced by Portniaguine and Zhdanov (Portniaguine O and Zhdanov M S 1999 Geophysics64 874-87). This stabilizer helps generate a sharp and focused image of anomalous conductivity distribution. The inversion is based on the re-weighted regularized conjugate gradient method.

1323

Linearized and nonlinear techniques are presented for determining estimates of parameter uncertainty within a two-dimensional iterative Born scheme. The scheme employs low frequency (<100 kHz) magnetic dipole sources in one well, and uses measurements of the vertical magnetic field in a second well to invert for the electrical conductivity distribution between the two boreholes. For computational efficiency a localized nonlinear approximation is employed to compute the sensitivity matrix. Parameter variance estimates are determined using an iterative Monte Carlo technique that assumes the data contain measurement noise, and that constraint assumptions imposed on the model are in error. The a posteriori model covariance matrix is determined statistically for the linearized technique by rerunning the last iteration of the nonlinear inversion N times, each time adding random errors to the data and constraints. The nonlinear approach involves rerunning the full inversion N times. Two oil field examples from California indicate that the linearized approach produces the same general pattern in the uncertainty estimates as the nonlinear estimation process. However, the linearized estimates are smaller in magnitude and show less spatial variation compared to the full nonlinear estimates, and the deviation between the two techniques increases as the contrast between the maximum and minimum conductivities within the inversion domain becomes greater.

1343

and

This paper describes a method to infer the geometry of the saturation front (the boundary of the region swept by water) around a water injection well in an oil reservoir. Information on how injected water displaces oil in different layers is essential to best design a sweep of the reservoir. We use a Bayesian inference approach, which combines DC electrical measurements from an electrode array in the injection well with background information on the thicknesses and resistivities of each reservoir layer. The solution is the posterior probability density function of the saturation front location. We estimate its characteristics by a combination of nonlinear optimization and Monte Carlo sampling. Using synthetic data, we show how to jointly invert measurements acquired at different times after the beginning of injection and how the uncertainty in the inferred location of the saturation front depends on uncertainties in the background information available and on the amount of noise in the data.

1357

and

We analyse a variety of solution strategies for nonlinear two-dimensional (2D) and three-dimensional (3D) electromagnetic (EM) inverse problems. To impose a realistic parameterization on the problem, the finite-difference equations arising from Maxwell equations are employed in the forward problem. Krylov subspace methods are then used to solve the resulting linear systems. Because of the efficiencies of the Krylov methods, they are used as the computational kernel for solving 2D and 3D inverse problems, where multiple solutions of the forward problem are required. We derive relations for computing the full Hessian matrix and functional gradient that are needed for computing the model update, via the Newton iteration. Different strategies are then discussed for economically carrying out the Newton iteration for 2D and 3D problems, including the incorporation of constraints necessary to stabilize the inversion process. Two case histories utilizing EM inversion are presented. These include inversion of cross-well data for monitoring electrical conductivity changes arising from an enhanced oil recovery project and the usefulness of cross-well EM methods to characterize the transport pathways for contaminants in the subsurface.

PAPERS

1377

We give a constructive method for extending time domain data for the inverse scattering problem for the one-dimensional wave equation. We show that a reflection operator on L2(-T,T) with T finite is essentially a Hankel operator and then modify the Nehari extension of the kernel of the reflection operator to obtain a reflection operator on L2(Bbb R) that is consistent with conservation of energy. This extension result allows frequency domain techniques to be used when the time domain data are only available for finite time, and we demonstrate this by using the frequency domain characterization of reflection coefficients to obtain a new proof of the characterization of reflection operators on L2(-T,T). We also give an a priori estimate for the operator norm of the reflection operator on L2(-T,T) and use the theory of Toeplitz operators to show how the singular values of these reflection operators are related to the reflection coefficient.

1405

Iterative algorithms for image reconstruction often involve minimizing some cost function h(x) that measures the degree of agreement between the measured data and a theoretical parametrized model. In addition, one may wish to have x satisfy certain constraints. It is usually the case that the cost function is the sum of simpler functions:

h(x) = ∑i = 1Ihi(x).

Partitioning the set {i = 1,...,I} as the union of the disjoint sets Bn,n = 1,...,N, we let

hn(x) = ∑iBnhi(x).

The method presented here is block iterative, in the sense that at each step only the gradient of a single hn(x) is employed. Convergence can be significantly accelerated, compared to that of the single-block (N = 1) method, through the use of appropriately chosen scaling factors.

The algorithm is an interior point method, in the sense that the images xk + 1 obtained at each step of the iteration satisfy the desired constraints. Here the constraints are imposed by having the next iterate xk + 1 satisfy the gradient equation

F(xk + 1) = ∇F(xk)-tnhn(xk),

for appropriate scalars tn, where the convex function F is defined and differentiable only on vectors satisfying the constraints.

Special cases of the algorithm that apply to tomographic image reconstruction, and permit inclusion of upper and lower bounds on individual pixels, are presented. The focus here is on the development of the underlying convergence theory of the algorithm. Behaviour of special cases has been considered elsewhere.

1421

, , and

A simple numerical algorithm for the reconstruction of micropore size distribution (MPSD) dealing with equilibrium adsorption isotherm data is presented. The mathematical model is based on the theory of volume filling of micropores (TVFMP) and the regularization method was used in the numerical algorithm for solving the TVFMP inverse problem. The two different approaches to inversion of the initial integral equation, namely, the MPSD reconstruction in terms of both micropore volume and micropore surface area, are analysed. The various exponent values in the local Dubinin-Astakhov isotherm and two types of smoothing functionals are analysed with the aim of studying their influence on the results of the MPSD reconstruction. On the basis of numerical experiments with the various model distributions, a simple scheme for choosing the regularization's optimal parameter is proposed. This scheme allows one to use observed data in calculations, obtained during only one experimental run, that essentially decrease the time of a cycle; for example an adsorption experiment - a MPSD reconstruction. The N2 adsorption isotherm for the properly characterized molecular sieve, carbon HGS-638, was used for testing the proposed technique. The results of MPSD reconstruction calculated in terms of the surface area model are in good agreement with the known data for this sample. The present technique can be applied to analysis of a wide class of micropore materials, such as zeolites, molecular sieves and carbon sorbents.

1441

, and

Under certain conditions, starting from the Riemann inversion formula, which gives an explicit representation of the inverse Laplace transform in the complex form, we derive an integral equation, of convolution type, whose solution is the inverse Laplace transform function.

This formula can be used if the Laplace transform has a finite number of singularities, located everywhere in the complex plane, and provided that their corresponding residues are known. It only requires the knowledge of the Laplace transform function on the real negative axis. Preliminary numerical experiments illustrating the reliability of the inversion algorithm are described.

1457

In this paper, the iteratively regularized Gauss-Newton method is applied to compute the stable solutions of nonlinear inverse problems. For the numerical realization, the discretization of this method is considered and the iterative solution is used to approximate the exact solution. An a posteriori rule is suggested to choose the stopping index of iteration and the convergence and rates of convergence are also derived under certain conditions. An example from the parameter identification of a differential equation problem is given to illustrate the required conditions.

1477

, and

It has been observed that sugar crystals growing in solution exhibit growth rate dispersion, that is, variation in growth rate from one crystal to the next. We consider the problem of estimating the distribution of growth rates in batch-grown crystals, given only samples of their sizes at a number of fixed points in time. The problem can be expressed as a tomographic image reconstruction problem, in which we try to reconstruct the joint density of initial size and growth from a set of marginal densities obtained by integrating the joint density in a number of different directions.

1487

, , and

This paper discusses the electrical impedance tomography (EIT) problem: electric currents are injected into a body with unknown electromagnetic properties through a set of contact electrodes. The corresponding voltages that are needed to maintain these currents are measured. The objective is to estimate the unknown resistivity, or more generally the impedivity distribution of the body based on this information. The most commonly used method to tackle this problem in practice is to use gradient-based local linearizations. We give a proof for the differentiability of the electrode boundary data with respect to the resistivity distribution and the contact impedances. Due to the ill-posedness of the problem, regularization has to be employed. In this paper, we consider the EIT problem in the framework of Bayesian statistics, where the inverse problem is recast into a form of statistical inference. The problem is to estimate the posterior distribution of the unknown parameters conditioned on measurement data. From the posterior density, various estimates for the resistivity distribution can be calculated as well as a posteriori uncertainties. The search of the maximum a posteriori estimate is typically an optimization problem, while the conditional expectation is computed by integrating the variable with respect to the posterior probability distribution. In practice, especially when the dimension of the parameter space is large, this integration must be done by Monte Carlo methods such as the Markov chain Monte Carlo (MCMC) integration. These methods can also be used for calculation of a posteriori uncertainties for the estimators. In this paper, we concentrate on MCMC integration methods. In particular, we demonstrate by numerical examples the statistical approach when the prior densities are non-differentiable, such as the prior penalizing the total variation or the L1 norm of the resistivity.

1523

In this paper we study the regularization of linear and nonlinear ill-posed operator equations by projection onto finite-dimensional spaces with a posteriori chosen space dimension. We show that this results in a regularization method, i.e. a stable solution method for the ill-posed problem, that converges to an exact solution as the data noise level goes to zero, with optimal rates under additional regularity conditions.

1541

and

We consider the inverse problem to determine the shape of a local perturbation of a perfectly conducting plate from a knowledge of the far-field pattern of the scattering of TM polarized time-harmonic electromagnetic waves by reformulating it as an inverse scattering problem for a planar domain with corners. For its approximate solution we propose a regularized Newton iteration scheme. For a foundation of Newton type methods we establish the Fréchet differentiability of the solution to the scattering problem with respect to the boundary and investigate the injectivity of the linearized mapping. Some numerical examples of the feasibility of the method are presented. For the sake of completeness, the first part of the paper outlines the solution of the direct scattering problem via an integral equation of the first kind including the numerical solution.

1561

In the inversion of data using the relaxation-spectrum model of electrical and rheological phenomena and processes, an important goal is to obtain a comprehensive and accurate approximation of the distribution of relaxation times (DRT) associated with the phenomenon or process under investigation. The estimation of such a DRT poses theoretical as well as experimental challenges. In the latter, it is often necessary to measure many decades of frequency-response data before sufficient information has been collected to allow an adequate estimate of the structure of the required spectrum to be calculated. At the theoretical level, among other things, one must solve the inverse problem of the relaxation-spectrum estimation, along with the frequent need to work with many decades of data. Some consequences of these two aspects of the problem are examined in this paper.

In order to address this problem, accurate (or noisy), wide-range frequency-response data sets derived from the Kohlrausch-Williams-Watts (KWW) continuous DRT function (with a value of its β0 parameter of 0.5) are inverted numerically using a complex nonlinear least-squares method with variable τ free parameters and variable quadrature weighting, a method superior to Tikhonov regularization for data of the present type when inversion accuracy need not be sacrificed for increased resolution. The relative errors in the resulting DRT point estimates are investigated for various frequency ranges for both frequency-response data derived from the usual KWW DRT and from such a DRT which is abruptly cut off at its low τ end. Although the inversions are ill-posed, DRT errors were found to be very small over appreciable τ ranges, but the effects of cutoff of the range of the DRT and of a limited frequency range led to rapid relative-error increases at the ends of the DRT τ range, and allowed separation and quantification of the effects of these two limitations.

Important resolution differences are illustrated between the results of inversions of accurate and of noisy data by the present method and by Tikhonov regularization. Finally, the temporal response was calculated very simply from frequency-response inversion DRT estimates and compared with exact stretched-exponential response function values associated with the KWW dispersion model. Extremely small relative errors were found, but they nevertheless showed the effects of the above limitations at short and long times. Transformation to the time domain of the kind illustrated here for wide-range data is far simpler, more accurate, and more convenient than Fourier transformation.

1585

, and

Electrical resistance tomography, the nonlinear inverse problem of retrieving the resistivity of an unknown body from steady-state boundary measurements, is addressed by a quadratic model that approximates the operator mapping the unknown parameters into boundary data. The forward problem, i.e. the computation of the measurements for a given resistivity, is solved by averaging the solutions, computed by a finite element method, of the two complementary formulations. Numerical results show that the complementary solutions, for the cases under investigation, are affected by approximately opposite numerical errors. Therefore, their average gives a good estimate of the solution even for coarse meshes. The overhead resulting from the need to compute two solutions is balanced by the selection of a coarse mesh. In particular, upper and lower bounds to the discretization errors are computed by the complementary solutions and used to optimize the finite element mesh. The retrieval of the resistivity is stated as the search for the minimizer of an error functional related to the quadratic model. The algebraic properties of quadratic models are exploited to estimate the number of linearly independent equations and to obtain an upper bound to the maximum number of unknown parameters, which guarantee an error functional without local minima.