Numerical modelling of post-seismic rupture propagation after the Sumatra 26.12.2004 earthquake constrained by GRACE gravity data

In the last decades, the development of the surface and satellite geodetic and geophysical observations brought a new insights into the seismic cycle, documenting new features of inter-, co-, and post-seismic processes. In particular since 2002 satellite mission GRACE provides monthly models of the global gravity field with unprecedented accuracy showing temporal variations of the Earth's gravity field, including those caused by mass redistribution associated with earthquake processes. When combined with GPS measurements, these new data have allowed to assess the relative importance of afterslip and viscoelastic relaxation after the Sumatra 26.12.2004 earthquake. Indeed the observed post-seismic crustal displacements were fitted well by a viscoelastic relaxation model assuming Burgers body rheology for the asthenosphere (60–220 km deep) with a transient viscosity as low as 4 × 10^17 Pas and constant ~ 10^19 Pas steady state viscosity in the 60–660-km depth range. However, even the low-viscosity asthenosphere provides the amplitude of strain which gravity effect does not exceed 50 per cent of the GRACE gravity variations, thus additional localized slip of about 1 m was suggested at downdip extension of the coseismic rupture. Post-seismic slip at coseismic rupture or its downdip extension has been suggested by several authors but the mechanism of the post-seismic fault propagation has never been investigated numerically. Depth and size of localized slip area as well as rate and time decay during the post-seismic stage were either assigned a priory or estimated by fitting real geodesy or gravity data. In this paper we investigate post-seismic rupture propagation by modelling two consequent stages. First, we run a long-term, geodynamic simulation to self-consistently produce the initial stress and temperature distribution. At the second stage, we simulate a seismic cycle using results of the first step as initial conditions. The second short-term simulation involves three substeps, including additional stress accumulation after part of the subduction channel was locked; spontaneous coseismic slip; formation and development of damage zones producing afterslip. During the last substep post-seismic stress leads to gradual ~1 m slip localized at three faults around ~100-km downdip extension of the coseismic rupture. We used the displacement field caused by the slip to calculate pressure and density variations and to simulate gravity field variations. Wavelength of calculated gravity anomaly fits well to that of the real data and its amplitude provides about 60 per cent of the observed GRACE anomaly. Importantly, the surface displacements caused by the estimated afterslip are much smaller than those registered by GPS networks. As a result cumulative effect of Burgers rheology viscoelastic relaxation (which explains measured GPS displacements and about a half of gravity variations) plus post-seismic slip predicted by damage rheology model (which causes much smaller surface displacements but provides another half of the GRACE gravity variations) fits well to both sets of the real data. Hence, the presented numerical modelling based on damage rheology supports the process of post-seismic downdip rupture propagation previously hypothesized from the GRACE gravity data.