Project A3 • Adaptive implicit space-time discretization for wave equations

Principal investigators

  Prof. Dr. Willy Dörfler (7/2015 - )
  Prof. Dr. Christian Wieners (7/2015 - )

Project summary

For acoustic, elastic, and electromagnetic wave equations, we developed fully adaptive discretizations in space and time with integrated error control and efficient implicit solution methods with, at least numerically validated (almost) optimal complexity. The system can have still millions of degrees of freedom. Therefore it must be parallelized on several computational cores and solved using a space-time multigrid preconditioner.

One motivation for developing space-time methods is the design of modern computer facilities with an enormous number of processor cores, where the parallel realization of conventional methods becomes inefficient. Since these machines allow a fully implicit space-time approach, new parallel solution techniques are required to solve the huge linear systems.

Linear wave equations

  • Acoustic waves: find pressure \(p\) and velocity \(\boldsymbol{v}\) with \[\partial_t p - \kappa \operatorname{div} \boldsymbol{v} = 0 \,, \qquad \rho\partial_t \boldsymbol{v} - \nabla p = \boldsymbol{f}\] for permeability \(\kappa\) and density \(\rho\).

  • Elastic waves: find stress \(\boldsymbol{\sigma}\) and velocity \(\boldsymbol{v}\) with \[\partial_t \boldsymbol{\sigma} - \boldsymbol{C} \boldsymbol{\varepsilon}(\boldsymbol{v}) = \boldsymbol{0} \,, \qquad \rho\partial_t \boldsymbol{v} - \boldsymbol{\operatorname{div}} \boldsymbol{\sigma} = \boldsymbol{f}\] with strain rate \(\boldsymbol{\varepsilon}(\boldsymbol{v}) = \operatorname{sym} (\mathrm{D} \boldsymbol{v})\) and the constitutive relation \[\boldsymbol{C}\boldsymbol{\varepsilon} = 2\mu\boldsymbol{\varepsilon} + \lambda\operatorname{trace}(\boldsymbol{\varepsilon})\textbf{1}\,.\]

We also realized space-time methods for electromagnetic waves [[Fin16]] and for viscoacoustic and viscoelastic waves.

Discretization

We consider three different discretizations.

  • (symmetric) Space-time discontinuous Petrov-Galerkin discretization (ST-DPG)
    The Discontinuous Petrov–Galerkin Method is a methodology to construct robust discretization schemes for linear PDEs. Due to its guaranteed stability independently of the mesh, it is well suited for adaptive applications. Moreover, as a Least-Squares type method, the resulting linear system is symmetric and positive definite which can be exploited by the iterative solver.

  • Space-time discontinuous Petrov-Galerkin discretization (ST-PGDG)
    The discretization uses the trial space $V_h \subset \mathrm{H}^1(0,1;\mathrm{L}_2(\Omega;\mathbb{R}^m))$, which is discontinuous in space but continuous in time. The test space $W_h \subset \mathrm{L}_2((0,T) \times \Omega ;\mathbb{R}^m))$ is discontinuous in space and time. On a tensor product space-time mesh with a fixed mesh in space and a time series $0 = t_0 < \dots < t_N=T$ with a local selection of polynomial degrees in space and time $p_R$ and $q_R$ in every cell, we set for the local test space $W_{h,R}=\mathbb{P}_{q_R-1}((t_{n-1},t_n); \mathbb{R}^m) \otimes \mathbb{P}_{p_R}(K;\mathbb{R}^m)$. Then, the local ansatz spaces $V_{h,R} = V_h|_R$ take the form \begin{align*} V_{h,R} = \{\mathbf{v}_{h,R} \in \mathrm{L}_2(R,\mathbb{R}^m) \colon& \mathbf{v}_{h,R}(t,\mathbf{x}) = \frac{t_n-t}{t_n-t_{n-1}} \mathbf{v}_{h,R}(t_{n-1},\mathbf{x})+ \frac{t_n-t}{t_n-t_{n-1}}\mathbf{w}_{h,R}(t,\mathbf{x}), \\ &\mathbf{v}_{h,R} \in V_h|_{[0,t_{n-1}]}, \mathbf{w}_{h,R} \in W_{h,R}, (t,\mathbf{x}) \in R = I_n \times K\} \end{align*} The global test space is discontinuous in space and time and defined by \begin{align*} W_h = \{\mathbf{w}_h \in \mathrm{L}_2((0,T);\mathrm{H}) \colon \mathbf{w}_{h,R} = \mathbf{w}_{h}|_R \} \end{align*} The global ansatz space $V_h$ is only discontinuous in space but continuous in time. More details are given in [DFWZ19].

  • Space-time discontinuous Galerkin discretization (ST-DGDG)
    Similar to above, this discretization also uses a tensor product space-time mesh. We have each cell $R = (t_{n-1},t_n) \times K$ the polynomial degrees in space and time $(q_R, p_q)$ for its local ansatz and test space. Hence, we have the following local polynomial space $V_{h,R}=\mathbb{P}_{q_R}((t_{n-1},t_n); \mathbb{R}^m) \otimes \mathbb{P}_{p_R}(K;\mathbb R^m) \subset \mathrm{L}_2((0,T) \times \Omega ;\mathbb{R}^m))$. Globally, this yields a discontinuous space in space and time.

Numerical analysis

Using the variational setting, we can show inf-sup stability and a priori bounds for the first two discretizations.

The bilinear form $b_h(\cdot,\cdot)=(L_h \cdot,\cdot)_{0,Q}$ is inf-sup stable in $V_h \times W_h$ with $\beta=C_1T^{-1}$, i.e., \begin{align*} \operatorname{sup}_{\mathbf{w_h} \in W_h\setminus \{0\} } b_h(\mathbf{v_h},\mathbf{w_h}) \geq \beta \|\mathbf{v_h}\|_{V_h},\quad \mathbf{v_h} \in V_h. \end{align*} A direct result is that for given $\mathbf{f}\in L_2(Q;\mathbb{R}^J)$ a unique solution $\mathbf{u_h}\in V_h$ exists solving \begin{align*} (L_h\mathbf{u_h},\mathbf{w_h})_{0,Q}=(f,\mathbf{w_h})_{0,Q}, \quad \mathbf{w_h} \in W_h \end{align*} and satisfying the a priori bound $\|\mathbf{w_h}\|_{V_h} \leq C_2 T\|M_h^{-1}\Pi_h \mathbf{f}\|_W$.

  • Numerical analysis for ST-DPG

    An extensive numerical study in [Ern17] demonstrates that the theoretically predicted convergence rates are achieved in practice for various benchmark problems in one and two spatial dimensions. While the analysis of DPG in [Ern17] and [EW18b] is restricted to homogeneous materials in space, [EW18a] extends the arguments to the case of inhomogeneous material distributions.
  • Numerical analysis for ST-PGDG

    The discrete operator $L_h$ is obtained by Galerkin approximation and upwind flux. For the proof we refer to [Fin16], [DFW16] and [DFWZ19].
  • Numerical analysis for ST-DGDG

    Using a mesh-dependent DG-norm $\|\cdot\|_{h,DG}$, we obtain the following convergence result for a sufficiently smooth solution $\mathbf{u} \in H^s(Q,\mathbb{R})$ with $s \geq 1$ and $s \leq \min_{n,K}\{p_{n,K}, q_{n,K}\}+1$, we have \begin{align*} \|\mathbf{u} - \mathbf{u}_h\|_{h,\mathrm{DG}} \leq C h^{s-1/2}\|D^s\mathbf{u}\|_Q + CTh^{-1/2}\|M_h^{-1/2}(M_h-M) \partial_t \mathbf{u}\|_Q \end{align*} For a proof and the construction of a residual error estimator for this setting as well as numerical experiments we refer to [CDW23].

Numerical experiments

  • Numerical experiment with ST-DPG

    Here, we also demonstrate that the space-time DPG method can be applied successfully to an application driven benchmark problem that has been inspired by the Marmousi benchmark. For this benchmark, we use the DPG method with tensor-product polynomials of order \(k\) in each space-time cell and tensor-product elements of order \(k+1\) on the space-time interfaces.

    The Marmousi model is a two space dimensions seismic model devised by the Institut Français du Petrole to test seismic data inversion.

    In contrast to the \(pq\)-adaptive algorithm used with ST-PGDG, we decided to reduce the computational effort by truncating the space-time mesh.

    The next animation shows, that we don't lose any informations by truncating the mesh. We also verify the results by comparing them with ST-PGDG and a classical time stepping scheme.

  • Numerical experiment with ST-PGDG

    The second example illustrates seismic tunnel exploration: An artificially generated surface wave (at \(x_\text{mid}\)) in the tunnel propagates into the solid and the reflected waves are measured in a certain region (marked red in the figure). We solve the acoustic wave equation in two space dimensions using a space-time discontinuous Galerkin discretization. In space we consider discontinuous polynomials over quadrilaterals of degree \(p\) and continuous polynomials of degree \(q\) in time (for details see [DFWZ19]). On the right we see a sketch of the computational domain \(\Omega\). The region of interest is marked in red. Starting with a finite volume discretization \((p=0, q=1)\) we perform four adaptive refinement steps using a dual weighted goal-oriented error estimator.
    The animation shows the space-time solution being sliced in time. On the bottom is the polynomial degree with 0 (blue) to 4 (red). The adaptive algorithm is controlled by a goal functional with support in the region of interest.
    We compared the uniform and adaptive refinement in the acoustic case observe that the adaptive algorithm saves 70% of the degrees of freedom while achieving the same accuracy compared with uniform refinement.
     
    ref-step #DoF (effort) GMRES steps with MG-PC \(\vartriangle\,E_{\text{ex}}\) in %
      uniform refinement
    r=1 534 528 (00%) 7 1.85
    r=2 2 405 376 (00%) 13 0.22
    r=3 6 414 336 (00%) 19 0.49
    r=4 13 363 200 (00%) 27 0.25
      adaptive refinement
    r=0 133 632 (00%) 5 91.01
    r=1 291 411 (55%) 7 1.27
    r=2 819 279 (34%) 13 0.58
    r=3 1 875 753 (29%) 20 0.56
    r=4 3 594 969 (27%) 28 0.38

    Acoustic wave: Uniform vs. adaptive refinement on \(44\,544 =928 \times 48\) space-time cells distributed on 64 processor cores. The error \(\vartriangle E_\text{ex}(\mathbf{u_h}) = |E(\mathbf{u_h}) -E_\text{ex}|\) of the goal functional is approximately estimated with respect to a linear extrapolation of the uniform results.

    The parallel scaling behavior of the parallel multilevel preconditioner is tested for different numbers of processes. On mesh level 4 we have 2 850 816 space-time cells, a linear discretization in space and time results in 34 209 792 DoFs for the acoustic case. The computing time for solving this huge linear system system with the parallel multigrid method scales nearly optimal.

    Tunnel Experiment: Strong scaling for \(\sim34\) Mio. DoFs
    computed on the bwUniCluster.

    For further explanations can be found in [DFWZ19] and a introduction into the used software in [Wie16] and Summer school SpaceTime.

  • Numerical experiment with ST-DPG

    We consider a model problem in order to verify our proofs given in [CDW23]. With a smooth initial value and a discontinuous material distribution we obtain the expected convergence rate in the DG-Norm.

    Further, a Riemann problem was considered with a discontinuous initial value and constant material parameters. We observe convergence with a reduced rate, since the above smoothness assumption is not satisfied.

Publications

  1. , , and . Space-time discontinuous Galerkin methods for weak solutions of hyperbolic linear symmetric Friedrichs systems. J. Sci. Comput., 94:27, December . URL https://doi.org/10.1007/s10915-022-02076-3. Online first. [preprint] [files] [bibtex]

  2. . The heterogeneous multiscale method for dispersive Maxwell systems. Multiscale Model. Simul., 20(2):769–797, June . URL https://doi.org/10.1137/21M1443960. [preprint] [bibtex]

  3. , , , , , and . Visco-acoustic full waveform inversion: From a DG forward solver to a Newton-CG inverse solver. Comput. Math. Appl., 100:126–140, October . URL https://doi.org/10.1016/j.camwa.2021.09.001. [preprint] [bibtex]

  4. , , and . Parallel space-time solutions for the linear visco-acoustic and visco-elastic wave equation. In W. E. Nagel, D. H. Kröner, and M. M. Resch, editors, High Performance Computing in Science and Engineering '19, pages 589–599, May . Springer, Cham. [preprint] [bibtex]

  5. and . The parallel finite element system M++ with integrated multilevel preconditioning and multilevel Monte Carlo methods. Comput. Math. Appl., 81:391–406, January . URL https://doi.org/10.1016/j.camwa.2020.03.004. [preprint] [bibtex]

  6. , , and . Space-time discontinuous Galerkin methods for linear hyperbolic systems and the application to the forward problem in seismic imaging. In R. Klöfkorn, E. Keilegavlen, F. A. Radu, and J. Fuhrmann, editors, Finite Volumes for Complex Applications IX - Methods, Theoretical Aspects, Examples. FVCA 2020, volume 323 of Springer Proceedings in Mathematics & Statistics, pages 447–485, June . Springer, Cham. [bibtex]

  7. , , , and . Parallel adaptive discontinuous Galerkin discretizations in space and time for linear elastic and acoustic waves. In U. Langer and O. Steinbach, editors, Space-Time Methods: Applications to Partial Differential Equations, volume 25 of Radon Series on Computational and Applied Mathematics, chapter 2, pages 61–88. De Gruyter, Berlin/Boston, September . [preprint] [bibtex]

  8. and . A space-time discontinuous Petrov–Galerkin method for acoustic waves. In U. Langer and O. Steinbach, editors, Space-Time Methods: Applications to Partial Differential Equations, volume 25 of Radon Series on Computational and Applied Mathematics, chapter 3, pages 89–116. De Gruyter, Berlin/Boston, September . [preprint] [bibtex]

  9. and . Space-time discontinuous Petrov–Galerkin methods for linear wave equations in heterogeneous media. Comput. Methods Appl. Math., 19(3):465–481, July . URL https://doi.org/10.1515/cmam-2018-0190. [preprint] [bibtex]

  10. and . Parallel inelastic heterogeneous multi-scale simulations. In S. Diebels and S. Rjasanow, editors, Multi-scale simulation of composite materials — Results from the project MuSiKo, pages 57–96. Springer, Berlin, February . [bibtex]

  11. , , and . Space-time discontinuous Galerkin discretizations for linear first-order hyperbolic evolution systems. Comput. Methods Appl. Math., 16(3):409–428, July . URL https://doi.org/10.1515/cmam-2016-0015. [preprint] [bibtex]

Preprints

  1. , , and . A fully parallelized and budgeted multi-level Monte Carlo method and the application to acoustic waves. CRC 1173 Preprint 2023/18, Karlsruhe Institute of Technology, July . [bibtex]

  2. . A space-time Petrov–Galerkin method for linear wave equations (Lecture notes for Zurich Summer School 2016: Numerical Methods for Wave Propagation). CRC 1173 Preprint 2018/19, Karlsruhe Institute of Technology, August . [bibtex]

Theses

  1. . A fully parallelized and budgeted multi-level Monte Carlo framework for partial differential equations. PhD thesis, Karlsruhe Institute of Technology (KIT), February . [bibtex]

  2. . Numerical homogenization of time-dependent Maxwell's equations with dispersion effects. PhD thesis, Karlsruhe Institute of Technology (KIT), November . [bibtex]

  3. . A parallel and adaptive space-time discontinuous Galerkin method for visco-elastic and visco-acoustic waves. PhD thesis, Karlsruhe Institute of Technology (KIT), November . [bibtex]

  4. . Space-time methods for acoustic waves with applications to full waveform inversion. PhD thesis, Karlsruhe Institute of Technology (KIT), December . [bibtex]

  5. . A parallel and adaptive space-time method for Maxwell's equations. PhD thesis, Karlsruhe Institute of Technology (KIT), July . [bibtex]

Former staff
Name Title Function
Dr. Postdoctoral and doctoral researcher
Dr. Postdoctoral and doctoral researcher
Dr. Postdoctoral and doctoral researcher
Dr. Postdoctoral and doctoral researcher