We seek to quantify bulk viscoelastic flow, afterslip, and locking, within a rheological framework that is consistent over the entire earthquake cycle. We address this using an ensemble smoother. We construct a 2D finite element seismic cycle model with a power-law rheology in the asthenosphere. A priori information, such as a realistic temperature field and a coseismic slip distribution, is integrated into the model. Model pre-stresses are initialized during repeated earthquake cycles wherein the accumulated slip deficit is released entirely. We tailor the last earthquake to match the observed co-seismic slip of the 2011 Tohoku earthquake. The heterogeneous rheology structure is derived from the temperature field and experimental flow laws. Additionally, we simulate afterslip using a thin viscoelastic shear zone. We focus on constraining power-law flow parameters for the asthenosphere and the shear zone. We assimilate 3D GEONET GNSS displacement time series acquired before and after the 2011 Tohoku earthquake. Power-law viscosity parameters are successfully retrieved for all domains. The data require separate viscoelastic domains in the mantle wedge above and below ~50 km depth. The sub-slab asthenosphere has viscoelastic properties that are distinctly different from the mantle wedge. The trade-off between the power-law activation energy and water fugacity hinders their individual estimation. The wedge viscosity is >10^19 Pa·s during the interseismic phase. Postseismic afterslip and bulk viscoelastic relaxation can be individually resolved from the surface deformation data. Afterslip is substantial between 40-50 km depth and extends to 80 km depth. Bulk viscoelastic relaxation in the wedge concentrates above 150 km depth with viscosities <10^18 Pa·s. Landward motion of the near-trench region occurs during the early postseismic period without the need for a separate low-viscosity channel below the slab.@en

Bayesian-based data assimilation methods integrate observational data into geophysical forward models to obtain the temporal evolution of an improved state vector, including its uncertainties. We explore the potential of a variant, a particle method, to estimate mechanical parameters of the overriding plate during the interseismic period. Here we assimilate vertical surface displacements into an elementary flexural model to estimate the elastic thickness of the overriding plate, and the locations and magnitudes of line loads acting on the overriding plate to produce flexure. Assimilation of synthetic observations sampled from a different forward model than is used in the particle method, reveal that synthetic seafloor data within 150 km from the trench are required to properly constrain parameters for long wavelength solutions of the upper plate (i.e. wavelength ∼500 km). Assimilation of synthetic observations sampled from the same flexural model used in the particle method shows remarkable convergence towards the true parameters with synthetic on-land data only for short to intermediate wavelength solutions (i.e. wavelengths between ∼100 and 300 km). In real-data assimilation experiments we assign representation errors due to discrepancies between our incorrect or incomplete physical model and the data. When assimilating continental data prior to the 2011 Mw Tohoku-Oki earthquake (1997-2000), an unrealistically low effective elastic plate thickness for Tohoku of ∼5-7 km is estimated. Our synthetic experiments suggest that improvements to the physical forward model, such as the inclusion of a slab, a megathrust interface and viscoelasticity of the mantle, including accurate seafloor data, and additional geodetic observations, may refine our estimates of the effective elastic plate thickness. Overall, we demonstrate the potential of using the particle method to constrain geodynamic parameters by providing constraints on parameters and corresponding uncertainty values. Using the particle method, we provide insights into the data network sensitivity and identify parameter trade-offs.

@en

Our ability to forecast earthquakes and slow slip events is hampered by limited information on the current state of stress on faults. Ensemble data assimilation methods permit estimating the state by combining physics-based models and observations, while considering their uncertainties. We use an ensemble Kalman filter (EnKF) to estimate shear stresses, slip rates and the state θ acting on a fault point governed by rate-and-state friction embedded in a 1-D elastic medium. We test the effectiveness of data assimilation by conducting perfect model experiments. We assimilate noised shear-stress and velocity synthetic values acquired at a small distance to the fault. The assimilation of uncertain shear stress observations improves in particular the estimates of shear stress on fault segments hosting slow slip events, while assimilating observations of velocity improves their slip-rate estimation. Both types of observations help equally well to better estimate the state θ. For earthquakes, the shear stress observations improve the estimation of shear stress, slip rates and the state θ, whereas the velocity observations improve in particular the slip-rate estimation. Data assimilation significantly improves the estimates of the temporal occurrence of slow slip events and to a large extent also of earthquakes. Rapid and abrupt changes in velocity and shear stress during earthquakes lead to non-Gaussian priors for subsequent assimilation steps, which breaks the assumption of Gaussian priors of the EnKF. In spite of this, the EnKF still provides estimates that are unexpectedly close to the true evolution. In fact, the forecastability for earthquakes for the same alarm duration is very similar to slow slip events, having a very low miss rate with an alarm duration of just 10 per cent of the recurrence interval of the events. These results confirm that data assimilation is a promising approach for the combination of uncertain physics and indirect, noisy observations for the forecasting of both slow slip events and earthquakes.

@en

The feasibility of physics-based forecasting of earthquakes depends on how well models can be calibrated to represent earthquake scenarios given uncertainties in both models and data. We investigate whether data assimilation can estimate current and future fault states, i.e., slip rate and shear stress, in the presence of a bias in the friction parameter. We perform state estimation as well as combined state-parameter estimation using a sequential-importance resampling particle filter in a zero-dimensional (0D) generalization of the Burridge–Knopoff spring–block model with rate-and-state friction. Minor changes in the friction parameter ϵ can lead to different state trajectories and earthquake characteristics. The performance of data assimilation with respect to estimating the fault state in the presence of a parameter bias in ϵ depends on the magnitude of the bias. A small parameter bias in ϵ (+3 %) can be compensated for very well using state estimation (R2 = 0.99), whereas an intermediate bias (−14 %) can only be partly compensated for using state estimation (R2 = 0.47). When increasing particle spread by accounting for model error and an additional resampling step, R2 increases to 0.61. However, when there is a large bias (−43 %) in ϵ, only state-parameter estimation can fully account for the parameter bias (R2 = 0.97). Thus, simultaneous state and parameter estimation effectively separates the error contributions from friction and shear stress to correctly estimate the current and future shear stress and slip rate. This illustrates the potential of data assimilation for the estimation of earthquake sequences and provides insight into its application in other nonlinear processes with uncertain parameters.@en

Geodetic observations of vertical land motion following a megathrust earthquake are key to a better understanding of processes and parameters controlling the dynamics at subduction margins. The relative contributions of dominant drivers during the postseismic phase, such as viscoelastic relaxation, afterslip and relocking, remain difficult to estimate individually and are often derived at the end of an observation period. Data assimilation can provide the means to estimate model parameters by combining physical models with observations whilst taking their uncertainties into account. We use Bayesian inference in the form of an ensemble smoother to estimate geodynamic parameters during the postseismic phase of the megathrust earthquake cycle. The ensemble smoother uses a Monte Carlo approach to represent the probability density distribution (pdf) of model states with a finite number of realizations. Prior estimates of the imperfect physical model are combined with the likelihood of noisy observations to estimate the posterior pdf of model parameters. We discuss a synthetic data experiment where observations are sampled from a 3D finite element model with noise added to represent errors in the data. With a smoother, observations at all time steps are assimilated in one go, to ensure consistency between estimated parameters and model outputs. We assimilate vertical and horizontal surface displacements into a 2D finite element viscoelastic earthquake cycle model with a power-law rheology. We incorporate heterogeneity into the viscoelastic structure and estimate the viscosity field based on a heterogeneous temperature field and experimental flow laws. Preliminary results show that model parameters, such as the extent of the cold nose, maximum depth of afterslip and flow law parameters can be recovered remarkably well by assimilating synthetic on- and offshore surface observations. We will apply the smoother to postseismic surface displacements of the Tohoku 2011 earthquake to estimate model parameters and driving processes across spatio-temporal domains.@en

Gas production in the Groningen field in the Netherlands reactivates faults and induces earthquakes, triggering severe societal unrest. Laboratory experiments indicate the relevant lithologies are velocity-strengthening [1], which favors stable creeping rather than instabilities and thus in theory prohibits earthquake nucleation. This makes the observed seismicity hard to explain. However, slide-hold-slide experiments on the same rocks reveal healing: stress drop and energy release of the subsequent slip after reactivation increase as healing time increases [2]. Given that the Groningen region has been tectonically inactive for millions of years, we will show that the healing history is a key factor for the nucleation of induced earthquakes.

We simulate earthquake sequences on a normal fault governed by rate-and-state friction, crosscutting a depleting gas reservoir. We reproduced the healing behavior observed in the laboratory and extended healing time up to millions of years. Our simulations suggest induced seismicity occur after fault reactivation, no matter whether the frictional behavior is velocity-weakening (VW) or -strengthening (VS). Unlike under VW friction where larger earthquakes with longer recurrence intervals occur as production continues, no subsequent earthquakes are expected under VS friction. No theoretically predictable nucleation size is measured: nucleation zone keeps expanding simultaneously when shear stress is dropping. Rupture propagation is mostly restricted to within the reservoir. Its reduced possibility of propagating into the under- and overburden reduces energy release and thus moment magnitude. The more VS the fault is, the lower the maximum slip rate will be, if under the same healing time. In other words, for a larger (a-b), a longer healing time is required before a large earthquake can happen. Therefore, if the healing history is well constrained, we can have an estimate of the seismic hazard of the coming event before it occurs. However, in order to have an accurate hazard estimate, dynamic weakening mechanisms activated at slip rate of cm/s may not be ignored. Adding flash heating shortens the required healing time by at least eight orders of magnitude in one of our representative models. Parameter space is explored systematically to build analytical expressions of earthquake magnitude as a function of healing time and frictional parameters. This study will help to understand in which tectonic histories, induced seismicity can still occur under VS friction. This is important for forecasting induced seismicity and minimizing the hazard in similar lithologies across northwestern Europe and in subsurface production and storage sites worldwide.

@en