Forward Solver Research Articles

The nth-Order Comprehensive Adjoint Sensitivity Analysis Methodology for Response-Coupled Forward/Adjoint Linear Systems (nth-CASAM-L): I. Mathematical Framework

This work presents the mathematical framework of the nth-Order Comprehensive Adjoint Sensitivity Analysis Methodology for Response-Coupled Forward/Adjoint Linear Systems (abbreviated as “nth-CASAM-L”), which is conceived for obtaining the exact expressions of arbitrarily-high-order (nth-order) sensitivities of a generic system response with respect to all of the parameters (including boundary and initial conditions) underlying the respective forward/adjoint systems. Since many of the most important responses for linear systems involve the solutions of both the forward and the adjoint linear models that correspond to the respective physical system, the sensitivity analysis of such responses makes it necessary to treat linear systems in their own right, rather than treating them as particular cases of nonlinear systems. This is in contradistinction to responses for nonlinear systems, which can depend only on the forward functions, since nonlinear operators do not admit bona-fide adjoint operators (only a linearized form of a nonlinear operator admits an adjoint operator). The nth-CASAM-L determines the exact expression of arbitrarily-high order sensitivities of responses to the parameters underlying both the forward and adjoint models of a nonlinear system, thus enable the most efficient and accurate computation of such sensitivities. The mathematical framework underlying the nth-CASAM is developed in linearly increasing higher-dimensional Hilbert spaces, as opposed to the exponentially increasing “parameter-dimensional” spaces in which response sensitivities are computed by other methods, thus providing the basis for overcoming the “curse of dimensionality” in sensitivity analysis and all other fields (uncertainty quantification, predictive modeling, etc.) which need such sensitivities. In particular, for a scalar-valued valued response associated with a nonlinear model comprising TP parameters, the 1st-−CASAM-L requires 1 additional large-scale adjoint computation (as opposed to TP large-scale computations, as required by other methods) for computing exactly all of the 1st-−order response sensitivities. All of the (mixed) 2nd-order sensitivities are computed exactly by the 2nd-CASAM-L in at most TP computations, as opposed to TP(TP + 1)/2 computations required by all other methods, and so on. For every lower-order sensitivity of interest, the nth-CASAM-L computes the “TP next-higher-order” sensitivities in one adjoint computation performed in a linearly increasing higher-dimensional Hilbert space. Very importantly, the nth-CASAM-L computes the higher-level adjoint functions using the same forward and adjoint solvers (i.e., computer codes) as used for solving the original forward and adjoint systems, thus requiring relatively minor additional software development for computing the various-order sensitivities.

Factorization method with one plane wave: from model-driven and data-driven perspectives

The factorization method provides a necessary and sufficient condition for characterizing the shape and position of an unknown scatterer by using far-field patterns of infinitely many time-harmonic plane waves at a fixed frequency (which are also called the multistatic data response matrix). This paper is concerned with the factorization method with a single far-field pattern to recover an arbitrary convex polygonal scatterer/source. Its one-wave version relies on the absence of analytical continuation of the scattered/radiated wave-fields in corner domains. It can be regarded as a domain-defined sampling method and does not require forward solvers. In this paper we provide a rigorous mathematical justification of the one-wave factorization method and present some preliminary numerical examples. In particular, the proposed method can be interpreted as a model-driven and data-driven imaging scheme, and it shows how to incorporate a priori knowledge about the unknown target into the test scatterers for the purpose of detecting obstacles/sources with specific features.

Numerical calculation of extremal Steklov eigenvalues in 3D and 4D

We develop a numerical method for solving shape optimization of functionals involving Steklov eigenvalues and apply it to the problem of maximization of the k-th Steklov eigenvalue, under volume constraint. A similar study in the planar case was already addressed in the literature using the boundary integral equation method. Here we extend that study to the 3D and 4D cases, using the Method of Fundamental Solutions as the forward solver.

Multiscale modeling and simulation methods for electromagnetic and multiphysics problems

Multiscale modeling and simulation methods for electromagnetic and multiphysics problems

Anatomical atlas of the upper part of the human head for electroencephalography and bioimpedance applications

Objective. The objective of this work is to develop a 4D (3D+T) statistical anatomical atlas of the electrical properties of the upper part of the human head for cerebral electrophysiology and bioimpedance applications. Approach. The atlas was constructed based on 3D magnetic resonance images (MRI) of 107 human individuals and comprises the electrical properties of the main internal structures and can be adjusted for specific electrical frequencies. T1w+T2w MRI images were used to segment the main structures of the head while angiography MRI was used to segment the main arteries. The proposed atlas also comprises a time-varying model of arterial brain circulation, based on the solution of the Navier–Stokes equation in the main arteries and their vascular territories. Main results. High-resolution, multi-frequency and time-varying anatomical atlases of resistivity, conductivity and relative permittivity were created and evaluated using a forward problem solver for EIT. The atlas was successfully used to simulate electrical impedance tomography measurements indicating the necessity of signal-to-noise between 100 and 125 dB to identify vascular changes due to the cardiac cycle, corroborating previous studies. The source code of the atlas and solver are freely available to download. Significance. Volume conductor problems in cerebral electrophysiology and bioimpedance do not have analytical solutions for nontrivial geometries and require a 3D model of the head and its electrical properties for solving the associated PDEs numerically. Ideally, the model should be made with patient-specific information. In clinical practice, this is not always the case and an average head model is often used. Also, the electrical properties of the tissues might not be completely known due to natural variability. Anatomical atlases are important tools for in silico studies on cerebral circulation and electrophysiology that require statistically consistent data, e.g. machine learning, sensitivity analyses, and as a benchmark to test inverse problem solvers.

A Reconstruction Scheme Using an Average Conductivity Element for Pitting Corrosion Defects in Steam Generator Tubes From Eddy Current Testing Signals

Pitting corrosion is a typical defect in steam generator tubes of nuclear power plants (NPPs) which may cause leak of coolant water in the first circulation loop of an NPP. Different from a crack, the pitting defect has small hole-like geometrical feature, which makes the conventional inversion scheme for crack profile reconstruction using eddy current testing (ECT) signals unsuitable. In this article, a numerical scheme to reconstruct the pitting corrosion defect from ECT signals is proposed and validated, including a new efficient ECT signal forward solver and an iterative inversion algorithm. For efficient calculation of the ECT signals due to a pitting hole, defect model based on a finite element of average conductivity is proposed. The model is validated by comparing its results with those of the conventional finite element method and boundary element method (FEM-BEM) hybrid code of $A$ – $\phi $ formulation using a detailed mesh for a flat-bottom hole. To reconstruct the defect profile such as the depth and diameter, the particle swarm optimization algorithm was adopted to solve the inversion problem iteratively. Finally, ECT experiments were conducted for test piece tubes with artificial pitting defects of different sizes and the defect profile parameters were reconstructed from the measured ECT signals with the proposed numerical schemes. Both the depths and the diameters of the artificial pitting defects, i.e., flat-bottom holes, were properly reconstructed referring to the true defect sizes.

Verification method of Monte Carlo codes for transport processes with arbitrary accuracy

In this work, we present a robust and powerful method for the verification, with arbitrary accuracy, of Monte Carlo codes for simulating random walks in complex media. Such random walks are typical of photon propagation in turbid media, scattering of particles, i.e., neutrons in a nuclear reactor or animal/humans’ migration. Among the numerous applications, Monte Carlo method is also considered a gold standard for numerically “solving” the scalar radiative transport equation even in complex geometries and distributions of the optical properties. In this work, we apply the verification method to a Monte Carlo code which is a forward problem solver extensively used for typical applications in the field of tissue optics. The method is based on the well-known law of average path length invariance when the entrance of the entities/particles in a medium obeys to a simple cosine law, i.e., Lambertian entrance, and annihilation of particles inside the medium is absent. By using this law we achieve two important points: (1) the invariance of the average path length guarantees that the expected value is known regardless of the complexity of the medium; (2) the accuracy of a Monte Carlo code can be assessed by simple statistical tests. We will show that we can reach an arbitrary accuracy of the estimated average pathlength as the number of simulated trajectories increases. The method can be applied in complete generality versus the scattering and geometrical properties of the medium, as well as in presence of refractive index mismatches in the optical case. In particular, this verification method is reliable to detect inaccuracies in the treatment of boundaries of finite media. The results presented in this paper, obtained by a standard computer machine, show a verification of our Monte Carlo code up to the sixth decimal digit. We discuss how this method can provide a fundamental tool for the verification of Monte Carlo codes in the geometry of interest, without resorting to simpler geometries and uniform distribution of the scattering properties.

Using deep generative neural networks to account for model errors in Markov chain Monte Carlo inversion

SUMMARY Most geophysical inverse problems are non-linear and rely upon numerical forward solvers involving discretization and simplified representations of the underlying physics. As a result, forward modelling errors are inevitable. In practice, such model errors tend to be either completely ignored, which leads to biased and over-confident inversion results, or only partly taken into account using restrictive Gaussian assumptions. Here, we rely on deep generative neural networks to learn problem-specific low-dimensional probabilistic representations of the discrepancy between high-fidelity and low-fidelity forward solvers. These representations are then used to probabilistically invert for the model error jointly with the target geophysical property field, using the computationally cheap, low-fidelity forward solver. To this end, we combine a Markov chain Monte Carlo (MCMC) inversion algorithm with a trained convolutional neural network of the spatial generative adversarial network (SGAN) type, whereby at each MCMC step, the simulated low-fidelity forward response is corrected using a proposed model-error realization. Considering the crosshole ground-penetrating radar traveltime tomography inverse problem, we train SGAN networks on traveltime discrepancy images between: (1) curved-ray (high fidelity) and straight-ray (low fidelity) forward solvers; and (2) finite-difference-time-domain (high fidelity) and straight-ray (low fidelity) forward solvers. We demonstrate that the SGAN is able to learn the spatial statistics of the model error and that suitable representations of both the subsurface model and model error can be recovered by MCMC. In comparison with inversion results obtained when model errors are either ignored or approximated by a Gaussian distribution, we find that our method has lower posterior parameter bias and better explains the observed traveltime data. Our method is most advantageous when high-fidelity forward solvers involve heavy computational costs and the Gaussian assumption of model errors is inappropriate. Unstable MCMC convergence due to non-linearities introduced by our method remain a challenge to be addressed in future work.

Pseudo-Likelihood Estimation for Parameters of Stochastic Time-Fractional Diffusion Equations

Although stochastic fractional partial differential equations have received increasing attention in the last decade, the parameter estimation of these equations has been seldom reported in literature. In this paper, we propose a pseudo-likelihood approach to estimating the parameters of stochastic time-fractional diffusion equations, whose forward solver has been investigated very recently by Gunzburger, Li, and Wang (2019). Our approach can accurately recover the fractional order, diffusion coefficient, as well as noise magnitude given the discrete observation data corresponding to only one realization of driving noise. When only partial data is available, our approach can also attain acceptable results for intermediate sparsity of observation.

Visco-acoustic full waveform inversion: From a DG forward solver to a Newton-CG inverse solver

Full waveform inversion (FWI) entails the ill-posed reconstruction of material parameters (such as sound speed and attenuation) from measurements of complete wave fields (full seismograms). In this paper we present a novel framework for FWI in the visco-acoustic regime. The new framework is based on a new elegant derivation of the system of state and adjoint PDEs which are approximated by the discontinuous Galerkin (DG) method. The inverse problem is then solved by the well established regularization scheme CG-REGINN which has not yet been applied in the context of FWI. For the DG discretization we provide a preconditioner for the efficient computation of the time steps by GMRES which yields optimal convergence estimates in space and time and which is confirmed by numerical tests. The inverse solver expresses the required Fréchet derivative and its adjoint in the DG setting. Successful reconstructions in a simplified cross-well setting serve as a proof of concept for our framework and demonstrate the applicability of our new combination of DG method and inverse solver. Some of the inversion experiments use seismograms generated by an independent finite difference time domain forward solver to avoid inverse crime.

Boundary element fast multipole method for modeling electrical brain stimulation with voltage and current electrodes

Objective. To formulate, validate, and apply an alternative to the finite element method (FEM) high-resolution modeling technique for electrical brain stimulation—the boundary element fast multipole method (BEM-FMM). To include practical electrode models for both surface and embedded electrodes. Approach. Integral equations of the boundary element method in terms of surface charge density are combined with a general-purpose fast multipole method and are expanded for voltage, shunt, current, and floating electrodes. The solution of coupled and properly weighted/preconditioned integral equations is accompanied by enforcing global conservation laws: charge conservation law and Kirchhoff’s current law. Main results. A sub-percent accuracy is reported as compared to the analytical solutions and simple validation geometries. Comparison to FEM considering realistic head models resulted in relative differences of the electric field magnitude in the range of 3%–6% or less. Quantities that contain higher order spatial derivatives, such as the activating function, are determined with a higher accuracy and a faster speed as compared to the FEM. The method can be easily combined with existing head modeling pipelines such as headreco or mri2mesh. Significance. The BEM-FMM does not rely on a volumetric mesh and is therefore particularly suitable for modeling some mesoscale problems with submillimeter (and possibly finer) resolution with high accuracy at moderate computational cost. Utilizing Helmholtz reciprocity principle makes it possible to expand the method to a solution of EEG forward problems with a very large number of cortical dipoles.

Infinite rotational motion generation and analysis of a spherical parallel manipulator with coaxial input axes

This paper presents a novel framework for the analysis of a 3-RRR spherical parallel manipulator with coaxial input axes (coaxial SPM) with the focus on its infinite rotational motion capabilities and its effects on the manipulator’s characteristics. The framework consists of three phases. At first, an approach for obtaining unique solutions of forward and inverse kinematics problems is introduced for setting up univocal relation between coaxial SPM’s input joint positions and orientation of its end-effector. At the second phase, a method for generating infinite rotational motions of an end-effector is formulated. The third phase outlines numerical computation procedures of the coaxial SPM’s workspaces in the joint and Cartesian spaces, excluding singularity configurations and mechanical collisions of SPM’s links during infinite rotational motion. A 3D design model and an experimental prototype of the coaxial SPM is presented and utilized for numerical analysis and experimental verification of the presented framework supplemented by an accompanying video demonstration.

Solution of large-scale 3D controlled-source electromagnetic modeling problem using efficient iterative solvers

With geophysical surveys evolving from traditional 2D to 3D models, the large volume of data adds challenges to inversion, especially when aiming to resolve complex 3D structures. An iterative forward solver for the controlled-source electromagnetic (CSEM) method requires less memory than that for a direct solver; however, it is not easy to iteratively solve an ill-conditioned linear system of equations arising from finite-element discretization of Maxwell’s equations. To solve this problem, we have developed efficient and robust iterative solvers for frequency- and time-domain CSEM modeling problems. For the frequency-domain problem, we first transform the linear system into its equivalent real-number format, and then introduce an optimal block-diagonal preconditioner. Because the condition number of the preconditioned linear equation system has an upper bound of [Formula: see text], we can achieve fast solution convergence when applying a flexible generalized minimum residual solver. Applying the block preconditioner further results in solving two smaller linear systems with the same coefficient matrix. For the time-domain problem, we first discretize the partial differential equation for the electric field in time using an unconditionally stable backward Euler scheme. We then solve the resulting linear equation system iteratively at each time step. After the spatial discretization in the frequency domain, or space-time discretization in the time domain, we exploit the conjugate-gradient solver with auxiliary-space preconditioners derived from the Hiptmair-Xu decomposition to solve these real linear systems. Finally, we check the efficiency and effectiveness of our iterative methods by simulating complex CSEM models. The most significant advantage of our approach is that the iterative solvers we adopt have almost the same accuracy and robustness as direct solvers but require much less memory, rendering them more suitable for large-scale 3D CSEM forward modeling and inversion.

Two-dimensional magnetotelluric finite element modeling by a hybrid Helmholtz-curl formulae system

The magnetotelluric (MT) method is an important geophysical electromagnetic induction exploration tool, with proper inversion routines, we can image the electrical conductivity structures of the Earth. A reliable inversion routine needs to iteratively call forward modeling solvers which solve the electromagnetic induction equations in the Earth. In this study, a novel hybrid Helmholtz-curl formulae system is presented for two-dimensional (2D) MT finite element modeling. In the air space with small values of conductivities, the well-behaved Helmholtz equation is used to represent the propagation of the electromagnetic field. In the Earth domain with positive values of conductivities, the well-behaved curl-curl electromagnetic equation is used to represent the inductive and propagating behaviors of the electromagnetic field. On the air-Earth interface, with arbitrary surface topographies, continuity conditions of the electromagnetic fields are enforced to guarantee the solution's uniqueness. Several synthetic models with both simple and complex geo-electrical structures incorporating topographies are used to verify and test the performance of the newly developed hybrid Helmholtz-curl formulae system. Numerical tests show that the hybrid Helmholtz-curl formulae system is numerically stable when using low conductivity values in the air space and using low signal frequencies. The hybrid Helmholtz-curl formulae system is an alternative algorithm to develop 2D MT forward modeling solvers.

Discrete Dipole Approximation-Based Microwave Tomography for Fast Breast Cancer Imaging.

This paper describes a fast microwave tomography reconstruction algorithm based on the two-dimensional discrete dipole approximation. Synthetic data from a finite-element based solver and experimental data from a microwave imaging system are used to reconstruct images and to validate the algorithm. The microwave measurement system consists of 16 monopole antennas immersed in a tank filled with lossy coupling liquid and a vector network analyzer. The low-profile antennas and lossy nature of system make the discrete dipole approximation an ideal forward solver in the image reconstructions. The results show that the algorithm can readily reconstruct a 2D plane of a cylindrical phantom. The proposed forward solver combined with the nodal adjoint method for computing the Jacobian matrix enables the algorithm to reconstruct an image within 6 seconds. This implementation provides a significant time savings and reduced memory requirements and is a dramatic improvement over previous implementations.

Source Localization of EEG Brainwaves Activities via Mother Wavelets Families for SWT Decomposition.

A Brain-Computer Interface (BCI) is a system used to communicate with an external world through the brain activity. The brain activity is measured by electroencephalography (EEG) signal and then processed by a BCI system. EEG source reconstruction could be a way to improve the accuracy of EEG classification in EEG based brain-computer interface (BCI). The source localization of the human brain activities can be an important resource for the recognition of the cognitive state, medical disorders, and a better understanding of the brain in general. In this study, we have compared 51 mother wavelets taken from 7 different wavelet families, which are applied to a Stationary Wavelet Transform (SWT) decomposition of an EEG signal. This process includes Haar, Symlets, Daubechies, Coiflets, Discrete Meyer, Biorthogonal, and reverse Biorthogonal wavelet families in extracting five different brainwave subbands for source localization. For this process, we used the Independent Component Analysis (ICA) for feature extraction followed by the Boundary Element Model (BEM) and the Equivalent Current Dipole (ECD) for the forward and inverse problem solutions. The evaluation results in investigating the optimal mother wavelet for source localization eventually identified the sym20 mother wavelet as the best choice followed by bior6.8 and coif5.

Forward modelling of gravity data on unstructured grids using an adaptive mimetic finite-difference method

One of the techniques to optimize the gradient-based inversion of geophysical gravity data is to deal with the large and dense Jacobian matrices implicitly. This is made possible by solving the normal systems iteratively and to use pseudo-forward problems to calculate the effect of the Jacobian matrices. For this, numerical forward solvers should be used instead of analytical solvers. In this study, we propose mimetic finite-difference schemes for the solution of the forward gravity problem on unstructured grids. While the mimetic finite-difference method shares many traits with the finite-element and finite-volume methods, it has the advantage of naturally accommodating grids with arbitrary polyhedral elements. We use cell-based and vertex-based mimetic finite-difference schemes to solve a diffusion problem in mixed form, and Poisson's problem, respectively. The cell-based scheme is solved for the gravitational potential and attraction at the centroids and the facets of the elements, respectively, and the vertex-based scheme is solved for the potential at the vertices of the elements. We compare these mimetic finite-difference schemes with cell-based and vertex-based finite-volume schemes and a vertex-based finite-element scheme, in terms of accuracy of the potential. We then utilize the versatility of the mimetic finite-difference method, in accommodating arbitrary elements, and develop an adaptive mesh refinement tool. The mesh adaptivity is performed by an iterative h-refinement where goal-oriented error estimates are used to mark the elements for refinement. The marked elements are decomposed into new tetrahedra by regular subdivision. Since arbitrary polyhedra are naturally permitted in the mimetic finite-difference method, the added nodes are not regarded as hanging nodes, and therefore, any modification of the scheme or extra refinement is avoided. This characteristic allows preserving the quality of the initial grid and simplifying the refinement procedure. We demonstrate the practicality of the presented mesh adaptivity for the forward modelling of gravity and gravity gradient data using realistic geological models. We show that the cell-based mimetic scheme is more efficient than the vertex-based scheme in terms of computational resources. We also demonstrate that in the presence of irregular interfaces in the input domain, the regular refinement approach, presented in this study, is more reliable than a standard mesh regeneration technique.

Approaching geoscientific inverse problems with vector-to-image domain transfer networks

We present vec2pix, a deep neural network designed to predict categorical or continuous 2D subsurface property fields from one-dimensional measurement data (e.g., time series), thereby offering an alternative approach to solve inverse problems. The performance of the method is investigated through two types of synthetic inverse problems: (a) a crosshole ground penetrating radar (GPR) tomography experiment with GPR travel times being used to infer a 2D velocity field, and (2) a multi-well pumping experiment within an unconfined aquifer with time series of transient hydraulic heads being used to retrieve a 2D hydraulic conductivity field. For each type of problem, both a multi-Gaussian and a binary channelized subsurface domain with long-range connectivity are considered. Using a training set of 20,000 examples (implying as many forward model evaluations), the method is found to recover a 2D model that is in much closer agreement with the true model than the closest training model in the forward-simulated data space. Further testing with smaller training sample sizes shows only a moderate reduction in performance when using 5000 training examples only. Although these 5000 to 20,000 forward runs can be performed offline in parallel, the associated computational expense may still be prohibitive for very demanding forward solvers. In addition, re-training is required for each new measurement configuration. Even if the recovered models are visually close to the true ones, the data misfits associated with their forward responses are generally larger than the noise level used to contaminate the true data. Uncertainty of the inverse solution is partially assessed using deep ensembles, in which the network is trained repeatedly with random initialization. Overall, this study advances understanding of how to use deep learning to infer subsurface models from indirect measurement data. More work is needed to evaluate the suitability of our proposed approach for real data, for which model errors and uncertainty in the geologic scenario will bring additional complexities.

Semi-Analytical Solutions for Barrier and American Options Written on a Time-Dependent Ornstein–Uhlenbeck Process

In this article, we develop semi-analytical solutions for the barrier (perhaps, time-dependent) and American options written on the underlying stock that follows a time-dependent Ornstein–Uhlenbeck process with a lognormal drift. Semi-analytical means that given the time-dependent interest rate, continuous dividend and volatility functions, one need to solve a linear (for the barrier option) or nonlinear (for the American option) Volterra equation of the second kind (or a Fredholm equation of the first kind). After that, the option prices in all cases are presented as one-dimensional integrals of combination of the preceding solutions and Jacobi theta functions. We also demonstrate that computationally our method is more efficient than the backward finite difference method traditionally used for solving these problems, and can be as efficient as the forward finite difference solver while providing better accuracy and stability. TOPICS:Derivatives, options, statistical methods Key Findings ▪ For the first time the method of generalized integral transform, invented in physics for solving an initial-boundary value parabolic problem at [0, y(t)] with a moving boundary [y(t)], is applied to finance. ▪ Using this method, pricing of barrier and American options, where the underlying follows a time-dependent OU process (the Bachelier model with drift) are solved in a semi-analytical form. ▪ It is demonstrated that computationally this method is more efficient than the backward and even forward finite difference method traditionally used for solving these problems whereas providing better accuracy and stability.

Two-dimensional Bayesian inversion of magnetotelluric data using trans-dimensional Gaussian processes.

Bayesian inversion of electromagnetic data produces crucial uncertainty information on inferred subsurface resistivity. Due to their high computational cost, however, Bayesian inverse methods have largely been restricted to computationally expedient 1-D resistivity models. In this study, we successfully demonstrate, for the first time, a fully 2-D, trans-dimensional Bayesian inversion of magnetotelluric (MT) data. We render this problem tractable from a computational standpoint by using a stochastic interpolation algorithm known as a Gaussian process (GP) to achieve a parsimonious parametrization of the model vis-a-vis the dense parameter grids used in numerical forward modelling codes. The GP links a trans-dimensional, parallel tempered Markov chain Monte Carlo sampler, which explores the parsimonious model space, to MARE2DEM, an adaptive finite element forward solver. MARE2DEM computes the model response using a dense parameter mesh with resistivity assigned via the GP model. We demonstrate the new trans-dimensional GP sampler by inverting both synthetic and field MT data for 2-D models of electrical resistivity, with the field data example converging within 10d on 148 cores, a non-negligible but tractable computational cost. For a field data inversion, our algorithm achieves a parameter reduction of over 32× compared to the fixed parameter grid used for the MARE2DEM regularized inversion. Resistivity probability distributions computed from the ensemble of models produced by the inversion yield credible intervals and interquartile plots that quantitatively show the non-linear 2-D uncertainty in model structure. This uncertainty could then be propagated to other physical properties that impact resistivity including bulk composition, porosity and pore-fluid content.