The analyses of interior penalty discontinuous Galerkin methods of any order k for solving elliptic and parabolic problems with Dirac line sources are presented. For the steady state case, we prove convergence of the method by deriving a priori error estimates in the L 2 norm and in weighted energy norms. In addition, we prove almost optimal local error estimates in the energy norm for any approximation order. Further, almost optimal local error estimates in the L 2 norm are obtained for the case of piecewise linear approximations whereas suboptimal error bounds in the L 2 norm are shown for any polynomial degree. For the time-dependent case, convergence of semi-discrete and of backward Euler fully discrete scheme is established by proving error estimates in L 2 in time and in space. Numerical results for the elliptic problem are added to support the theoretical results.
more »
« less
A discontinuous Galerkin method for sequences of earthquakes and aseismic slip on multiple faults using unstructured curvilinear grids
SUMMARY Physics-based simulations provide a path to overcome the lack of observational data hampering a holistic understanding of earthquake faulting and crustal deformation across the vastly varying space–time scales governing the seismic cycle. However, simulations of sequences of earthquakes and aseismic slip (SEAS) including the complex geometries and heterogeneities of the subsurface are challenging. We present a symmetric interior penalty discontinuous Galerkin (SIPG) method to perform SEAS simulations accounting for the aforementioned challenges. Due to the discontinuous nature of the approximation, the spatial discretization natively provides a means to impose boundary and interface conditions. The method accommodates 2-D and 3-D domains, is of arbitrary order, handles subelement variations in material properties and supports isoparametric elements, that is, high-order representations of the exterior boundaries, interior material interfaces and embedded faults. We provide an open-source reference implementation, Tandem, that utilizes highly efficient kernels for evaluating the SIPG linear and bilinear forms, is inherently parallel and well suited to perform high-resolution simulations on large-scale distributed memory architectures. Additional flexibility and efficiency is provided by optionally defining the displacement evaluation via a discrete Green’s function approach, exploiting advantages of both the boundary integral and volumetric methods. The optional discrete Green’s functions are evaluated once in a pre-computation stage using algorithmically optimal and scalable sparse parallel solvers and pre-conditioners. We illustrate the characteristics of the SIPG formulation via an extensive suite of verification problems (analytic, manufactured and code comparison) for elastostatic and quasi-dynamic problems. Our verification suite demonstrates that high-order convergence of the discrete solution can be achieved in space and time and highlights the benefits of using a high-order representation of the displacement, material properties and geometries. We apply Tandem to realistic demonstration models consisting of a 2-D SEAS multifault scenario on a shallowly dipping normal fault with four curved splay faults, and a 3-D intersecting multifault scenario of elastostatic instantaneous displacement of the 2019 Ridgecrest, CA, earthquake sequence. We exploit the curvilinear geometry representation in both application examples and elucidate the importance of accurate stress (or displacement gradient) representation on-fault. This study entails several methodological novelties. We derive a sharp bound on the smallest value of the SIPG penalty ensuring stability for isotropic, elastic materials; define a new flux to incorporate embedded faults in a standard SIPG scheme; employ a hybrid multilevel pre-conditioner for the discrete elasticity problem; and demonstrate that curvilinear elements are specifically beneficial for volumetric SEAS simulations. We show that our method can be applied for solving interesting geophysical problems using massively parallel computing. Finally, this is the first time a discontinuous Galerkin method is published for the numerical simulations of SEAS, opening new avenues to pursue extreme scale 3-D SEAS simulations in the future.
more »
« less
- NSF-PAR ID:
- 10447649
- Date Published:
- Journal Name:
- Geophysical Journal International
- Volume:
- 233
- Issue:
- 1
- ISSN:
- 0956-540X
- Page Range / eLocation ID:
- 586 to 626
- Format(s):
- Medium: X
- Sponsoring Org:
- National Science Foundation
More Like this
-
-
more » « less
SUMMARY We introduce a new finite-element (FE) based computational framework to solve forward and inverse elastic deformation problems for earthquake faulting via the adjoint method. Based on two advanced computational libraries, FEniCS and hIPPYlib for the forward and inverse problems, respectively, this framework is flexible, transparent and easily extensible. We represent a fault discontinuity through a mixed FE elasticity formulation, which approximates the stress with higher order accuracy and exposes the prescribed slip explicitly in the variational form without using conventional split node and decomposition discrete approaches. This also allows the first order optimality condition, that is the vanishing of the gradient, to be expressed in continuous form, which leads to consistent discretizations of all field variables, including the slip. We show comparisons with the standard, pure displacement formulation and a model containing an in-plane mode II crack, whose slip is prescribed via the split node technique. We demonstrate the potential of this new computational framework by performing a linear coseismic slip inversion through adjoint-based optimization methods, without requiring computation of elastic Green’s functions. Specifically, we consider a penalized least squares formulation, which in a Bayesian setting—under the assumption of Gaussian noise and prior—reflects the negative log of the posterior distribution. The comparison of the inversion results with a standard, linear inverse theory approach based on Okada’s solutions shows analogous results. Preliminary uncertainties are estimated via eigenvalue analysis of the Hessian of the penalized least squares objective function. Our implementation is fully open-source and Jupyter notebooks to reproduce our results are provided. The extension to a fully Bayesian framework for detailed uncertainty quantification and non-linear inversions, including for heterogeneous media earthquake problems, will be analysed in a forthcoming paper.
-
more » « less
Abstract We generated dense, high‐resolution 3‐D ground displacement maps for the 2016 MW7.8 Kaikōura, New Zealand earthquake—the most geometrically and kinematically complex rupture yet recorded—from stereo WorldView optical satellite imagery using a new methodology that combines subpixel image correlation with a ray‐tracing approach. Our analysis reveals fundamental new details of near‐field displacement patterns, which cannot easily be obtained through other methods. From our detailed correlation maps, we measured fault slip in 3‐D along 19 faults at 500‐m spacing. Minimum resolvable horizontal slip is ~0.1 m, and vertical is ~0.5 m. Net slip measurements range from <1 to ~12 m. System‐level kinematic analysis shows that slip on faults north of the Hope fault was oriented primarily subparallel to the Pacific‐Australian plate motion direction. In contrast, slip on faults to the south was primarily at high angle to the plate motion and secondarily parallel to plate motion. Fault kinematics are in some locations consistent with long‐term uplift patterns, but inconsistent in others. Deformation within the Seaward Kaikōura Range may indicate an attempt by the plate boundary fault system to geometrically simplify. Comparison of published field measurements along the Kekerengu fault with our correlation‐derived measurements reveals that ~36% of surface displacement was accommodated as distributed off‐fault deformation when considering only field measurements of discrete slip. Comparatively, field measurements that project previously linear features (e.g., fence lines) into the fault over apertures >5–100 m capture nearly all (~90%) of the surface deformation.
-
An H(div)-conforming finite element method for the Biot’s consolidation mo- del is developed, with displacements and fluid velocity approximated by elements from BDM_k space. The use of H(div)-conforming elements for flow variables ensures the local mass conservation. In the H(div)-conforming approximation of displacement, the tan- gential components are discretised in the interior penalty discontinuous Galerkin frame- work,and the normal components across the element interfaces are continuous. Having introduced a spatial discretisation, we develop a semi-discrete scheme and a fully dis- crete scheme,prove their unique solvability and establish optimal error estimates for each variable.more » « less
-
Abstract An interior penalty discontinuous Galerkin method is devised to approximate minimizers of a linear folding model by discontinuous isoparametric finite element functions that account for an approximation of a folding arc. The numerical analysis of the discrete model includes an a priori error estimate in case of an accurate representation of the folding curve by the isoparametric mesh. Additional estimates show that geometric consistency errors may be controlled separately if the folding arc is approximated by piecewise polynomial curves. Various numerical experiments are carried out to validate the a priori error estimate for the folding model.more » « less
- Free Publicly Accessible Full Text
-
Accepted Manuscript1.0
- Journal Article:
- https://doi.org/10.1093/gji/ggac467
