Direct and Inverse Problems for Models with Uncertain Data

Modelling real-world physical systems always includes uncertainties (at least implicitly) regarding the measurement and determination of the data of the model. This for instance involves material properties, loading and also domains. We are concerned with infinite dimensional space dependent systems which are described by PDE. Typical applications are the simulation of mechanical systems and of groundwater flow in porous media. Since the solutions directly depend on the data, the measurement of this dependence in probabilistic terms is required in order to asses the reliability of the calculated solution. By this, it becomes feasible to deduce and quantify the solution uncertainty depending on data uncertainties.


Illustration of forward propagation of uncertainties and inverse parameter identification.

The numerical treatment of these problems usually requires a computational effort which is significantly larger than for deterministic problems and which may easily become intractable. Thus, the development and analysis of efficient model reduction techniques for the numerical evaluation is important to enable the treatment of relevant applications. Depending on the data model, we employ modern sampling methods (multilevel Monte-Carlo) as well as projection methods based on the polynomial chaos (stochastic collocation and stochastic FEM).

Depending on the chosen problem and data model, different methods are employed which either perform a discretization of the function space or determine statistics of functionals based on samples. In particular, these are
  • projection methods based on polynomial chaos (stochastic Galerkin FEM and, in wider sense, also collocation und regression)
  • modern sampling methods (multilevel and Quasi Monte Carlo)
  • methods based on stochastic differential equations (SDEs)

Adaptive Functional Representations and Modern Compression Techniques

The evaluation of forward or inverse stochastic problems often can be accellerated significantly (in terms of convergence rates) by employing functional representations. Moreover, such a representation allows for a numerical analysis, which resembles and extends concepts known from deterministic PDEs. With discretizations in generalized polynomial chaos as e.g. used for stochastic Galerkin methods, (reliable) a posteriori error estimates can be computed. These iteratively lead to problem-adapted solution spaces with optimal convergence.

Since the discrete algebraic systems exhibit a high dimensionality, model reduction techniques are often inevitable. Apart from the adaptation of the discrete space, we employ modern hierarchical tensor methods for the compression of the operators and the parametric solutions. This is tightly related to a possible representation as low-rank manifold. The tensor train (TT) format has proven advantageous for this.

Applications for the developed methods can e.g. be found with sampling-free Bayesian inversion, topology optimization under uncertainties, and the determination of effective material models in case of media with multiscale properties.

Mesh adaptivity with oscillating coefficient field.

SDE- and Sampling-based Methods for Random PDEs

As an alternative to Monte Carlo methods, we investigate methods which exploit the equivalence of PDEs with random data and stochastic differential equations (Feynman-Kac). Such highly parallelizable approaches typically require the reconstruction of a global solution representation. However, they allow for a complete separate and local control of all discretization parameters. To be more specific, pointwise solution realizations in the physical domain are determined by some appropriate numerical method (Euler-Maruyama). Subsequently these are either used for a global or local polynomial regression or an interpolation on a mesh.

Regression approach for the pointwise solution of an SDE equivalent to the random PDE.

Publications

  Articles in Refereed Journals

  • M. Eigel, Ch. Merdon, J. Neumann, An adaptive multilevel Monte--Carlo method with stochastic bounds for quantities of interest in groundwater flow with uncertain data, SIAM/ASA Journal on Uncertainty Quantification, 4 (2016) pp. 1219--1245.
    Abstract
    The focus of this work is the introduction of some computable a posteriori error control to the popular multilevel Monte Carlo sampling for PDE with stochastic data. We are especially interested in applications in the geosciences such as groundwater flow with rather rough stochastic fields for the conductive permeability. With a spatial discretisation based on finite elements, a goal functional is defined which encodes the quantity of interest. The devised goal-oriented error estimator enables to determine guaranteed a posteriori error bounds for this quantity. In particular, it allows for the adaptive refinement of the mesh hierarchy used in the multilevel Monte Carlo simulation. In addition to controlling the deterministic error, we also suggest how to treat the stochastic error in probability. Numerical experiments illustrate the performance of the presented adaptive algorithm for a posteriori error control in multilevel Monte Carlo methods. These include a localised goal with problem-adapted meshes and a slit domain example. The latter demonstrates the refinement of regions with low solution regularity based on an inexpensive explicit error estimator in the multilevel algorithm.

  • M. Eigel, Ch. Merdon, Local equilibration error estimators for guaranteed error control in adaptive stochastic higher-order Galerkin finite element methods, SIAM/ASA Journal on Uncertainty Quantification, 4 (2016) pp. 1372--1397.
    Abstract
    Equilibration error estimators have been shown to commonly lead to very accurate guaranteed error bounds in the a posteriori error control of finite element methods for second order elliptic equations. Here, we extend previous results by the design of equilibrated fluxes for higher-order finite element methods with nonconstant coefficients and illustrate the favourable performance of different variants of the error estimator within two deterministic benchmark settings. After the introduction of the respective parametric problem with stochastic coefficients and the stochastic Galerkin FEM discretisation, a novel a posteriori error estimator for the stochastic error in the energy norm is devised. The error estimation is based on the stochastic residual and its decomposition into approximation residuals and a truncation error of the stochastic discretisation. Importantly, by using the derived deterministic equilibration techniques for the approximation residuals, the computable error bound is guaranteed for the considered class of problems. An adaptive algorithm allows the simultaneous refinement of the deterministic mesh and the stochastic discretisation in anisotropic Legendre polynomial chaos. Several stochastic benchmark problems illustrate the efficiency of the adaptive process.

  • F. Lanzara, V. Maz'ya, G. Schmidt, Approximation of solutions to multidimensional parabolic equations by approximate approximations, Applied and Computational Harmonic Analysis. Time-Frequency and Time-Scale Analysis, Wavelets, Numerical Algorithms, and Applications, 41 (2016) pp. 749--767.

  • M. Eigel, C.J. Gittelson, Ch. Schwab, E. Zander, A convergent adaptive stochastic Galerkin finite element method with quasi-optimal spatial meshes, ESAIM: Mathematical Modelling and Numerical Analysis, 49 (2015) pp. 1367--1398.
    Abstract
    We analyze a-posteriori error estimation and adaptive refinement algorithms for stochastic Galerkin Finite Element methods for countably-parametric elliptic boundary value problems. A residual error estimator which separates the effects of gpc-Galerkin discretization in parameter space and of the Finite Element discretization in physical space in energy norm is established. It is proved that the adaptive algorithm converges, and to this end we establish a contraction property satisfied by its iterates. It is shown that the sequences of triangulations which are produced by the algorithm in the FE discretization of the active gpc coefficients are asymptotically optimal. Numerical experiments illustrate the theoretical results.

  • F. Lanzara, G. Schmidt, On the computation of high-dimensional potentials of advection-diffusion operators, Mathematika. A Journal of Pure and Applied Mathematics, 61 (2015) pp. 309--327.

  • W. Giese, M. Eigel, S. Westerheide, Ch. Engwer, E. Klipp, Influence of cell shape, inhomogeneities and diffusion barriers in cell polarization models, Physical Biology, 12 (2015) pp. 066014/1--066014/18.
    Abstract
    In silico experiments bear the potential to further the understanding of biological transport processes by allowing a systematic modification of any spatial property and providing immediate simulation results for the chosen models. We consider cell polarization and spatial reorganization of membrane proteins which are fundamental for cell division, chemotaxis and morphogenesis. Our computational study is motivated by mating and budding processes of S. cerevisiae. In these processes a key player during the initial phase of polarization is the GTPase Cdc42 which occurs in an active membrane-bound form and an inactive cytosolic form. We use partial differential equations to describe the membrane-cytosol shuttling of Cdc42 during budding as well as mating of yeast. The membrane is modeled as a thin layer that only allows lateral diffusion and the cytosol is modeled as a volume. We investigate how cell shape and diffusion barriers like septin structures or bud scars influence Cdc42 cluster formation and subsequent polarization of the yeast cell. Since the details of the binding kinetics of cytosolic proteins to the membrane are still controversial, we employ two conceptual models which assume different binding kinetics. An extensive set of in silico experiments with different modeling hypotheses illustrate the qualitative dependence of cell polarization on local membrane curvature, cell size and inhomogeneities on the membrane and in the cytosol. We examine that spatial inhomogenities essentially determine the location of Cdc42 cluster formation and spatial properties are crucial for the realistic description of the polarization process in cells. In particular, our computer simulations suggest that diffusion barriers are essential for the yeast cell to grow a protrusion.

  • TH. Arnold, A. Rathsfeld, Reflection of plane waves by rough surfaces in the sense of Born approximation, Mathematical Methods in the Applied Sciences, 37 (2014) pp. 2091--2111.
    Abstract
    The topic of the present paper is the reflection of electromagnetic plane waves by rough surfaces, i.e., by smooth and bounded perturbations of planar faces. Moreover, the contrast between the cover material and the substrate beneath the rough surface is supposed to be low. In this case, a modification of Stearns' formula based on Born approximation and Fourier techniques is derived for a special class of surfaces. This class contains the graphs of functions if the interface function is a radially modulated almost periodic function. For the Born formula to converge, a sufficient and almost necessary condition is given. A further technical condition is defined, which guarantees the existence of the corresponding far field of the Born approximation. This far field contains plane waves, far-field terms like those for bounded scatterers, and, additionally, a new type of terms. The derived formulas can be used for the fast numerical computations of far fields and for the statistics of random rough surfaces.

  • M. Eigel, C. Gittelson, Ch. Schwab, E. Zander, Adaptive stochastic Galerkin FEM, Computer Methods in Applied Mechanics and Engineering, 270 (2014) pp. 247--269.

  • F. Lanzara, V. Maz'ya, G. Schmidt, Fast cubature of volume potentials over rectangular domains by approximate approximations, Applied and Computational Harmonic Analysis. Time-Frequency and Time-Scale Analysis, Wavelets, Numerical Algorithms, and Applications, 36 (2014) pp. 167--182.
    Abstract
    In the present paper we study high-order cubature formulas for the computation of advection-diffusion potentials over boxes. By using the basis functions introduced in the theory of approximate approximations, the cubature of a potential is reduced to the quadrature of one dimensional integrals. For densities with separated approximation, we derive a tensor product representation of the integral operator which admits efficient cubature procedures in very high dimensions. Numerical tests show that these formulas are accurate and provide approximation of order O(h6) up to dimension 108.

  Preprints, Reports, Technical Reports

  • M. Eigel, M. Marschall, R. Schneider, Bayesian inversion with a hierarchical tensor representation, Preprint no. 2363, WIAS, Berlin, 2016, DOI 10.20347/WIAS.PREPRINT.2363 .
    Abstract, PDF (744 kByte)
    The statistical Bayesian approach is a natural setting to resolve the ill-posedness of inverse problems by assigning probability densities to the considered calibration parameters. Based on a parametric deterministic representation of the forward model, a sampling-free approach to Bayesian inversion with an explicit representation of the parameter densities is developed. The approximation of the involved randomness inevitably leads to several high dimensional expressions, which are often tackled with classical sampling methods such as MCMC. To speed up these methods, the use of a surrogate model is beneficial since it allows for faster evaluation with respect to calibration parameters. However, the inherently slow convergence can not be remedied by this. As an alternative, a complete functional treatment of the inverse problem is feasible as demonstrated in this work, with functional representations of the parametric forward solution as well as the probability densities of the calibration parameters, determined by Bayesian inversion. The proposed sampling-free approach is discussed in the context of hierarchical tensor representations, which are employed for the adaptive evaluation of a random PDE (the forward problem) in generalized chaos polynomials and the subsequent high-dimensional quadrature of the log-likelihood. This modern compression technique alleviates the curse of dimensionality by hierarchical subspace approximations of the involved low rank (solution) manifolds. All required computations can be carried out efficiently in the low-rank format. A priori convergence is examined, considering all approximations that occur in the method. Numerical experiments demonstrate the performance and verify the theoretical results.

  • M. Eigel, J. Neumann, R. Schneider, S. Wolf, Stochastic topology optimisation with hierarchical tensor reconstruction, Preprint no. 2362, WIAS, Berlin, 2016, DOI 10.20347/WIAS.PREPRINT.2362 .
    Abstract, PDF (8552 kByte)
    A novel approach for risk-averse structural topology optimization under uncertainties is presented which takes into account random material properties and random forces. For the distribution of material, a phase field approach is employed which allows for arbitrary topological changes during optimization. The state equation is assumed to be a high-dimensional PDE parametrized in a (finite) set of random variables. For the examined case, linearized elasticity with a parametric elasticity tensor is used. Instead of an optimization with respect to the expectation of the involved random fields, for practical purposes it is important to design structures which are also robust in case of events that are not the most frequent. As a common risk-aware measure, the Conditional Value at Risk (CVaR) is used in the cost functional during the minimization procedure. Since the treatment of such high-dimensional problems is a numerically challenging task, a representation in the modern hierarchical tensor train format is proposed. In order to obtain this highly efficient representation of the solution of the random state equation, a tensor completion algorithm is employed which only required the pointwise evaluation of solution realizations. The new method is illustrated with numerical examples and compared with a classical Monte Carlo sampling approach.

  • F. Lanzara, V. Maz'ya, G. Schmidt, A fast solution method for time dependent multidimensional Schrödinger equations, Preprint no. 2279, WIAS, Berlin, 2016.
    Abstract, PDF (3794 kByte)
    In this paper we propose fast solution methods for the Cauchy problem for the multidimensional Schrödinger equation. Our approach is based on the approximation of the data by the basis functions introduced in the theory of approximate approximations. We obtain high order approximations also in higher dimensions up to a small saturation error, which is negligible in computations, and we prove error estimates in mixed Lebesgue spaces for the inhomogeneous equation. The proposed method is very efficient in high dimensions if the densities allow separated representations. We illustrate the efficiency of the procedure on different examples, up to approximation order 6 and space dimension 200.

  • F. Anker, Ch. Bayer, M. Eigel, J. Neumann, J.G.M. Schoenmakers, Adaptive SDE based interpolation for random PDEs, Preprint no. 2200, WIAS, Berlin, 2015.
    Abstract, PDF (1343 kByte)
    A numerical method for the fully adaptive sampling and interpolation of PDE with random data is presented. It is based on the idea that the solution of the PDE with stochastic data can be represented as conditional expectation of a functional of a corresponding stochastic differential equation (SDE). The physical domain is decomposed subject to a non-uniform grid and a classical Euler scheme is employed to approximately solve the SDE at grid vertices. Interpolation with a conforming finite element basis is employed to reconstruct a global solution of the problem. An a posteriori error estimator is introduced which provides a measure of the different error contributions. This facilitates the formulation of an adaptive algorithm to control the overall error by either reducing the stochastic error by locally evaluating more samples, or the approximation error by locally refining the underlying mesh. Numerical examples illustrate the performance of the presented novel method.

  • F. Anker, Ch. Bayer, M. Eigel, M. Ladkau, J. Neumann, J.G.M. Schoenmakers, SDE based regression for random PDEs, Preprint no. 2192, WIAS, Berlin, 2015.
    Abstract, PDF (1752 kByte)
    A simulation based method for the numerical solution of PDE with random coefficients is presented. By the Feynman-Kac formula, the solution can be represented as conditional expectation of a functional of a corresponding stochastic differential equation driven by independent noise. A time discretization of the SDE for a set of points in the domain and a subsequent Monte Carlo regression lead to an approximation of the global solution of the random PDE. We provide an initial error and complexity analysis of the proposed method along with numerical examples illustrating its behaviour.

  • M. Eigel, M. Pfeffer, R. Schneider, Adaptive stochastic Galerkin FEM with hierarchical tensor representations, Preprint no. 2153, WIAS, Berlin, 2015.
    Abstract, PDF (587 kByte)
    The solution of PDE with stochastic data commonly leads to very high-dimensional algebraic problems, e.g. when multiplicative noise is present. The Stochastic Galerkin FEM considered in this paper then suffers from the curse of dimensionality. This is directly related to the number of random variables required for an adequate representation of the random fields included in the PDE. With the presented new approach, we circumvent this major complexity obstacle by combining two highly efficient model reduction strategies, namely a modern low-rank tensor representation in the tensor train format of the problem and a refinement algorithm on the basis of a posteriori error estimates to adaptively adjust the different employed discretizations. The adaptive adjustment includes the refinement of the FE mesh based on a residual estimator, the problem-adapted stochastic discretization in anisotropic Legendre Wiener chaos and the successive increase of the tensor rank. Computable a posteriori error estimators are derived for all error terms emanating from the discretizations and the iterative solution with a preconditioned ALS scheme of the problem. Strikingly, it is possible to exploit the tensor structure of the problem to evaluate all error terms very efficiently. A set of benchmark problems illustrates the performance of the adaptive algorithm with higher-order FE. Moreover, the influence of the tensor rank on the approximation quality is investigated.

  Talks, Poster

  • M. Eigel, Adaptive stochastic Galerkin FEM with hierarchical tensor representations, Advances in Uncertainty Quantification Methods, Algorithms and Applications (UQAW 2016), January 5 - 10, 2016, King Abdullah University of Science and Technology (KAUST), Thuwal, Saudi Arabia, January 8, 2016.

  • M. Eigel, Adaptive stochastic Galerkin FEM with hierarchical tensor representations, 15th Conference on the Mathematics of Finite Elements and Applications (Brunel MAFELAP 2016), Minisymposium ``Uncertainty Quantification Using Stochastic PDEs and Finite Elements'', June 14 - 17, 2016, Brunel University London, Uxbridge, UK, June 14, 2016.

  • M. Eigel, Adaptive stochastic Galerkin FEM with hierarchical tensor representations, Joint Annual Meeting of DMV and GAMM, Section 18 ``Numerical Methods of Differential Equations'', March 7 - 11, 2016, Technische Universität Braunschweig, March 10, 2016.

  • M. Eigel, Bayesian inversion using hierarchical tensor approximations, SIAM Conference on Uncertainty Quantification, Minisymposium 67 ``Bayesian Inversion and Low-rank Approximation (Part II)'', April 5 - 8, 2016, Lausanne, Switzerland, April 6, 2016.

  • M. Eigel, Some aspects of adaptive random PDEs, Oberseminar, Rheinisch-Westfälische Technische Hochschule Aachen, Institut für Geometrie und Praktische Mathematik, July 21, 2016.

  • J. Neumann, Adaptive SDE based sampling for random PDE, SIAM Conference on Uncertainty Quantification, Minisymposium 142 ``Error Estimation and Adaptive Methods for Uncertainty Quantification in Computational Sciences -- Part II'', April 5 - 8, 2016, Lausanne, Switzerland, April 8, 2016.

  • J. Neumann, The phase field approach for topology optimization under uncertainties, ZIB Computational Medicine and Numerical Mathematics Seminar, Konrad-Zuse-Zentrum für Informationstechnik Berlin, August 25, 2016.

  • J. Pellerin, RINGMesh: A programming library for geological model meshes, The 17th annual conference of the International Association for Mathematical Geosciences, September 5 - 13, 2015, Freiberg, September 8, 2015.

  • M. Eigel, Adaptive stochastic Galerkin FEM with hierarchical tensor representations, 2nd GAMM AGUQ Workshop on Uncertainty Quantification, September 10 - 11, 2015, Chemnitz, September 10, 2015.

  • M. Eigel, Fully adaptive higher-order stochastic Galerkin FEM in low-rank tensor representation, International Conference on Scientific Computation And Differential Equations (SciCADE 2015), September 14 - 18, 2015, Universität Potsdam, September 15, 2015.

  • M. Eigel, Guaranteed error bounds for adaptive stochastic Galerkin FEM, Technische Universität Braunschweig, Institut für Wissenschaftliches Rechnen, April 1, 2015.

  • M. Eigel, Stochastic adaptive FEM, Forschungsseminar Numerische Mathematik, Humboldt-Universität zu Berlin, Institut für Mathematik, January 28, 2015.

  • CH. Bayer, SDE based regression for random PDEs, Direct and Inverse Problems for PDEs with Random Coefficients, WIAS Berlin, November 13, 2015.

  • M. Eigel, A posteriori error control in stochastic FEM and MLMC, 27th Chemnitz FEM Symposium 2014, September 22 - 24, 2014, September 24, 2014.

  • M. Eigel, Adaptive spectral methods for stochastic optimisation problems, Technische Universität Berlin, Institut für Mathematik, May 22, 2014.

  • M. Eigel, Adaptive stochastic FEM, Universität Heidelberg, Interdisziplinäres Zentrum für Wissenschaftliches Rechnen (IWR), June 5, 2014.

  • M. Eigel, Guaranteed a posteriori error control with adaptive stochastic Galerkin FEM, SIAM Conference on Uncertainty Quantification (UQ14), March 31 - April 3, 2014, Savannah, USA, April 1, 2014.

  • M. Ladkau, Brownian motion approach for spatial PDEs with stochastic data, International Workshop ``Advances in Optimization and Statistics'', May 15 - 16, 2014, Russian Academy of Sciences, Institute of Information Transmission Problems (Kharkevich Institute), Moscow, May 16, 2014.

  • J. Neumann, A posteriori error estimators for problems with uncertain data, Norddeutsches Kolloquium über Angewandte Analysis und Numerische Mathematik (NoKo), Christian-Albrechts-Universität zu Kiel, May 10, 2014.

  • J. Neumann, Stochastic bounds for quantities of interest in groundwater flow with uncertain data, Université Paris-Sud, Laboratoire d'Analyse Numérique, Orsay, France, October 9, 2014.