Matrices coming from elliptic Partial Differential Equations have been shown to have a low-rank property that can be efficiently exploited in multifrontal solvers to provide a substantial reduction of their complexity. Among the possible low-rank formats, the Block Low-Rank format (BLR) is easy to use in a general purpose multifrontal solver and its potential compared to standard (full-rank) solvers has been demonstrated. Recently, new variants have been introduced and it was proved that they can further reduce the complexity but their performance has never been analyzed. In this article, we present a multithreaded BLR factorization and analyze its efficiency and scalability in shared-memory multicore environments. We identify the challenges posed by the use of BLR approximations in multifrontal solvers and put forward several algorithmic variants of the BLR factorization that overcome these challenges by improving its efficiency and scalability. We illustrate the performance analysis of the BLR multifrontal factorization with numerical experiments on a large set of problems coming from a variety of real-life applications.

/pdf/performance-and-scalability-of-the-block-low-rank-3orlgawx2u.pdf

Performance and Scalability of the Block Low-Rank Multifrontal Factorization on Multicore Architectures

Geophysical inversion attempts to estimate the distribution of physical properties in the Earth's interior from observations collected at or above the surface. Inverse problems are commonly posed as least-squares optimization problems in high-dimensional parameter spaces. Existing approaches are largely based on deterministic gradient-based methods, which are limited by nonlinearity and nonuniqueness of the inverse problem. Probabilistic inversion methods, despite their great potential in uncertainty quantification, still remain a formidable computational task. In this paper, I explore the potential of deep learning methods for electromagnetic inversion. This approach does not require calculation of the gradient and provides results instantaneously. Deep neural networks based on fully convolutional architecture are trained on large synthetic datasets obtained by full 3-D simulations. The performance of the method is demonstrated on models of strong practical relevance representing an onshore controlled source electromagnetic CO2 monitoring scenario. The pre-trained networks can reliably estimate the position and lateral dimensions of the anomalies, as well as their resistivity properties. Several fully convolutional network architectures are compared in terms of their accuracy, generalization, and cost of training. Examples with different survey geometry and noise levels confirm the feasibility of the deep learning inversion, opening the possibility to estimate the subsurface resistivity distribution in real time.

Deep learning electromagnetic inversion with convolutional neural networks

S U M M A R Y We present a numerical algorithm for 3-D electromagnetic (EM) simulations in conducting media with general electric anisotropy. The algorithm is based on the finite-difference discretization of frequency-domain Maxwell’s equations on a Lebedev grid, in which all components of the electric field are collocated but half a spatial step staggered with respect to the magnetic field components, which also are collocated. This leads to a system of linear equations that is solved using a stabilized biconjugate gradient method with a multigrid preconditioner. We validate the accuracy of the numerical results for layered and 3-D tilted transverse isotropic (TTI) earth models representing typical scenarios used in the marine controlled-source EM method. It is then demonstrated that not taking into account the full anisotropy of the conductivity tensor can lead to misleading inversion results. For synthetic data corresponding to a 3-D model with a TTI anticlinal structure, a standard vertical transverse isotropic (VTI) inversion is not able to image a resistor, while for a 3-D model with a TTI synclinal structure it produces a false resistive anomaly. However, if the VTI forward solver used in the inversion is replaced by the proposed TTI solver with perfect knowledge of the strike and dip of the dipping structures, the resulting resistivity images become consistent with the true models.

Fully anisotropic 3-D EM modelling on a Lebedev grid with a multigrid pre-conditioner

A meaningful solution to an inversion problem should be composed of the preferred inversion model and its uncertainty and resolution estimates. The model uncertainty estimate describes an equivalent model domain in which each model generates responses which fit the observed data to within a threshold value. The model resolution matrix measures to what extent the unknown true solution maps into the preferred solution. However, most current geophysical electromagnetic (also gravity, magnetic and seismic) inversion studies only offer the preferred inversion model and ignore model uncertainty and resolution estimates, which makes the reliability of the preferred inversion model questionable. This may be caused by the fact that the computation and analysis of an inversion model depend on multiple factors, such as the misfit or objective function, the accuracy of the forward solvers, data coverage and noise, values of trade-off parameters, the initial model, the reference model and the model constraints. Depending on the particular method selected, large computational costs ensue. In this review, we first try to cover linearised model analysis tools such as the sensitivity matrix, the model resolution matrix and the model covariance matrix also providing a partially nonlinear description of the equivalent model domain based on pseudo-hyperellipsoids. Linearised model analysis tools can offer quantitative measures. In particular, the model resolution and covariance matrices measure how far the preferred inversion model is from the true model and how uncertainty in the measurements maps into model uncertainty. We also cover nonlinear model analysis tools including changes to the preferred inversion model (nonlinear sensitivity tests), modifications of the data set (using bootstrap re-sampling and generalised cross-validation), modifications of data uncertainty, variations of model constraints (including changes to the trade-off parameter, reference model and matrix regularisation operator), the edgehog method, most-squares inversion and global searching algorithms. These nonlinear model analysis tools try to explore larger parts of the model domain than linearised model analysis and, hence, may assemble a more comprehensive equivalent model domain. Then, to overcome the bottleneck of computational cost in model analysis, we present several practical algorithms to accelerate the computation. Here, we emphasise linearised model analysis, as efficient computation of nonlinear model uncertainty and resolution estimates is mainly determined by fast forward and inversion solvers. In the last part of our review, we present applications of model analysis to models computed from individual and joint inversions of electromagnetic data; we also describe optimal survey design and inversion grid design as important applications of model analysis. The currently available model uncertainty and resolution analyses are mainly for 1D and 2D problems due to the limitations in computational cost. With significant enhancements of computing power, 3D model analyses are expected to be increasingly used and to help analyse and establish confidence in 3D inversion models.

/pdf/uncertainty-and-resolution-analysis-of-2d-and-3d-inversion-48y3cwprd5.pdf

Uncertainty and Resolution Analysis of 2D and 3D Inversion Models Computed from Geophysical Electromagnetic Data

Nous nous interessons a l'utilisation d'approximations de rang faible pour 
reduire le cout des solveurs creux directs multifrontaux. Parmi les differents formats
matriciels qui ont ete proposes pour exploiter la propriete de rang faible dans les
solveurs multifrontaux, nous nous concentrons sur le format Block Low-Rank (BLR)
dont la simplicite et la flexibilite permettent de l'utiliser facilement dans un solveur
multifrontal algebrique et generaliste. Nous presentons differentes variantes de la
factorisation BLR, selon comment les mises a jour de rang faible sont effectuees, et
comment le pivotage numerique est gere.
D'abord, nous etudions la complexite theorique du format BLR qui, contrairement a d'autres formats comme les formats hierarchiques, etait inconnue jusqu'a
present. Nous prouvons que la complexite theorique de la factorisation 
multifrontale BLR est asymptotiquement inferieure a celle du solveur de rang plein. Nous
montrons ensuite comment les variantes BLR peuvent encore reduire cette complexite. 
Nous etayons nos bornes de complexite par une etude experimentale.
Apres avoir montre que les solveurs multifrontaux BLR peuvent atteindre une
faible complexite, nous nous interessons au probleme de la convertir en gains de
performance reels sur les architectures modernes. Nous presentons d'abord une
factorisation BLR multithreadee, et analysons sa performance dans des environnements
 multicœurs a memoire partagee. Nous montrons que les variantes BLR sont
cruciales pour exploiter efficacement les machines multicœurs en ameliorant l'intensite
 arithmetique et la scalabilite de la factorisation. Nous considerons ensuite
a la factorisation BLR sur des architectures a memoire distribuee.
Les algorithmes presentes dans cette these ont ete implementes dans le solveur
MUMPS. Nous illustrons l'utilisation de notre approche dans trois applications 
industrielles provenant des geosciences et de la mecanique des structures. Nous 
comparons egalement notre solveur avec STRUMPACK, base sur des approximations
Hierarchically Semi-Separable. Nous concluons cette these en rapportant un resultat
 sur un probleme de tres grande taille (130 millions d'inconnues) qui illustre les
futurs defis poses par le passage a l'echelle des solveurs multifrontaux BLR.

Solveurs multifrontaux exploitant des blocs de rang faible : complexité, performance et parallélisme

We have developed an efficient numerical scheme for fast multimodel 3D electromagnetic simulations by applying a Schur complement approach to a frequency-domain finite-difference method. The scheme is based on direct solvers and developed with constrained inversion algorithms in view. Such algorithms normally need many forward modeling jobs with different resistivities for the target zone and/or background formation. We geometricallydivide the computational domain intotwo subdomains: an anomalous subdomain, the resistivities of which were permitted to change, and a background subdomain, having fixed resistivities.The system matrix ispartiallyfactorizedby precomputing a Schur complement to eliminate unknowns associated with the background subdomain. The Schur complement system is then solved to compute fields inside the anomalous subdomain. Finally, the background subdomain fields are computed using inexpensive local substitutions. For each successive simulation, only the relatively small Schur complement system has to be solved, which results in significant computational savings. We applied this approach to two moderately sized 3D problems in marine controlled-source electromagnetic modeling: (1) a deepwater model in which the resistivities of the seawater and the air layer were kept fixed and (2) a model in which focused inversion was performed in a scenario in which the resistivities of the background formation, the air layer, and the seawater were known. We found a significant reduction of the modeling time in inversion that depended on the relative sizes of the constrained and unconstrained volumes: the smaller the unconstrained volume, the larger the savings. Specifically, for a focused inversion of the Troll oil field in the North Sea, the gain amounted up to 80% of the total modeling time.

Fast multimodel finite-difference controlled-source electromagnetic simulations based on a Schur complement approach

We put forward the idea of using a Block Low-Rank (BLR) multifrontal direct solver to efficiently solve the linear systems of equations arising from a finite-difference discretization of the frequency-domain Maxwell equations for 3-D electromagnetic (EM) problems. The solver uses a low-rank representation for the off-diagonal blocks of the intermediate dense matrices arising in the multifrontal method to reduce the computational load. A numerical threshold, the so-called BLR threshold, controlling the accuracy of low-rank representations was optimized by balancing errors in the computed EM fields against savings in floating point operations (flops). Simulations were carried out over large-scale 3-D resistivity models representing typical scenarios for marine controlled-source EM surveys, and in particular the SEG SEAM model which contains an irregular salt body. The flop count, size of factor matrices and elapsed run time for matrix factorization are reduced dramatically by using BLR representations and can go down to, respectively, 10, 30 and 40 per cent of their full-rank values for our largest system with N = 20.6 million unknowns. The reductions are almost independent of the number of MPI tasks and threads at least up to 90 × 10 = 900 cores. The BLR savings increase for larger systems, which reduces the factorization flop complexity from O(N2) for the full-rank solver to O(Nm) with m = 1.4–1.6. The BLR savings are significantly larger for deep-water environments that exclude the highly resistive air layer from the computational domain. A study in a scenario where simulations are required at multiple source locations shows that the BLR solver can become competitive in comparison to iterative solvers as an engine for 3-D controlled-source electromagnetic Gauss–Newton inversion that requires forward modelling for a few thousand right-hand sides.

https://hal.inria.fr/hal-01672952/document

Large-scale 3D EM modeling with a Block Low-Rank multifrontal direct solver

S U M M A R Y We have developed an efficient numerical scheme for 3-D electromagnetic (EM) simulations using an exponential finite-difference (FD) method with non-uniform grids. The method uses the set of exponential basis functions {1, exp[±(νx x + νy y + νz z)]}, where the exponents νx , νy and νz must be chosen carefully depending on the simulation frequency and local node conductivity. The method achieves an approximation of the oscillatory and exponentially decaying EM fields that is better than that obtained via the low-degree polynomial fitting from standard FDs—and hence also leads to more accurate results. An important property of the exponential FD method is that it tends to the standard FD method when the exponents νx , νy and νz tend to zero. We applied the standard and exponential FD methods to three marine controlled-source EM modelling scenarios: deep-water, shallow-water and intermediate water depth. For the deep-water scenario, we found that the proposed exponential FD method gave two to three times more accurate results as compared to the standard FD method on the same grid. For the shallow-water and intermediate water depth scenarios, the exponential FD method improved the accuracy of the upgoing fields; it gave 2–2.5 times more accurate results for the upgoing fields than the standard FD method on the same grid. Consequently, the method can achieve the same accuracy with a coarser grid and hence is faster than the standard FD method, as demonstrated using a frequency-domain iterative solver.

Piyoosh Jaysaval

Papers

Fully anisotropic 3-D EM modelling on a Lebedev grid with a multigrid pre-conditioner

Fast multimodel finite-difference controlled-source electromagnetic simulations based on a Schur complement approach

Large-scale 3D EM modeling with a Block Low-Rank multifrontal direct solver

Efficient 3-D controlled-source electromagnetic modelling using an exponential finite-difference method

Geophysical Methods for Stratigraphic Identification