This article was originally published on September 16, 2020 in the Geodynamics Blog of the European Geophysical Union by and . It is reproduced here with permission.

Advanced geodynamic models of giant earthquakes
Though giant earthquakes are disastrous, they provide essential information to investigate earthquake physics. In this week’s news and views, Thyagarajulu Gollapalli, a PhD student jointly from the Monash University and the Indian Institute of Technology, Bombay, discusses our present understanding …
Link to the original article at

Earthquakes with magnitude more than or equal to 8.5 (Mw >= 8.5) are known as giant earthquakes. A finer classification puts them into two further categories, great earthquakes (8.5 <= Mw < 9.0) and mega earthquakes (Mw >= 9.0). These earthquakes are so huge that they can be compared with the energy released from a few trillion kilograms of explosives. In this century, we have already witnessed two such giant earthquakes, e.g., the Sumatra 2004 Mw 9.2 and the Tohoku 2011 Mw 9.0. Jointly, they have crushed near about half a million human lives and triggered economic damage of more than $360 billion. Despite such large-scale socio-economic disaster, giant earthquakes are crucial for geoscientists as they enhance our limited understanding of earthquake physics. In this article, I will be discussing various concepts related to giant earthquakes: what are the present understanding of physics of such events, our existing knowledge gap and how we can fill those in near future.

Present understanding of giant earthquakes

Along the subduction zone, cold and denser plate sinks into the mantle. The contact surface between the subducting and overriding plate is called the “megathrust fault” (Figure 1). These are the only possible locations for giant earthquakes. They preferably occur at some shallow to an intermediate depth of the megathrust fault, where shallow level earthquakes along these faults can initiate tsunami.  Apart from the subduction zones, continental collision zones (e.g., Indo-Eurasia collision zone) have a strong potential to produce mega or great earthquakes [Dal Zilio et al., 2019].

Figure 1: Megathrust fault zone within a generalised subduction zone, highlighting the different segments of fault. Modified after [Bilek and Lay, 2018].

From the observational points we have understood that, a) the giant earthquakes occur in clusters (Figure 2), b) they have very long rupture lengths around 1000 km [Shearer and Bürgmann, 2010], c) they tend to occur within a short interval of the period in the same segment of megathrust (Figure 2).

Figure 2: Earthquakes Mw >= 8.0 on the subduction zone megathrust from the USGS seismic catalog.

The physical factors that control such large rupture length within a short time and intensify their magnitude have remained enigmatic. Previous studies suggested that parameters like convergence rate and lithosphere age [Ruff and Kanamori, 1980], sediment thickness at the trench [Ruff, 1989], state of the upper plate strain [Heuret et al., 2012] and the roughness of the seafloor [van Rijsingen et al., 2018]  play important role in controlling occurrence of giant earthquakes. Apart from these statistical analysis, most giant earthquakes have been modelled using the plate kinematics ignoring the complete dynamics of subduction systems. In the classical kinematic model of Ruff and Kanamori  [1980], they found a strong correlation between giant earthquake occurrence with young and fast-converging lithosphere. Although, such correlation is absent for recent Sumatra 2004 and Tohoku 2011 mega earthquakes. In a global geodynamic model, Conrad et al. [2004] concluded that slabs are likely to be detached from their subducting plate, which in turn reduce overall slab pull to produce a giant earthquake. However, their models did not incorporate detailed geometry of subducting slabs and thus failed to match the present-day slab pull force arising from slab depth variation along Sumatra trench. Slab pull force for Sumatra 2004 event estimated from intra-plate stress analysis was higher compared to surrounding region [Gordon and Houseman, 2015; Sandiford et al., 2005]. All these previous efforts to decipher the physics of giant earthquakes have remained isolated from each other and an integrated approach involving the subduction dynamics and plate kinematics is still lacking. In our advanced geodynamic models we incorporate the associated forces, rheologies and accurate geometry of the subduction system to mimic more realistic condition.

Figure 3: Relevant forces acting on the convergent margin that control tectonics and seismicity. Red and green arrows represent slab pull and mantle flow tractions. The blue curve represents interface strength. Source

Advanced geodynamic models

Three main forces are closely interconnected along the megathrust faults. The most prominent of them is the slab pull force, which is represented by the red arrow in figure 3. The magnitude of slab pull force depends on the density contrast between the slab and surrounding mantle, gravitational acceleration, slab depth inside the mantle, and slab thickness. While descending the slab is coupled with the mantle and exerts traction (shown in green in Figure 3) on both the subducting and the overriding plates. Another critical parameter that appears in the dynamic analysis is the strength of the interface or cohesion that depends on the properties of the plate.

One such important megathrust fault is located along the Indian Ocean trench. Here, the Indo-Australian oceanic plate subducts beneath the overriding Sunda plate. It has been recorded to release high seismic energy in the recent past. Observations revealed that there is a highly variable plate coupling along this megathrust (Figure 4a). The Sumatra segment of the trench hosted three closely spaced earthquakes of magnitude Mw >= 8.5 within three years rupturing ~1800 km of the plate boundary [McCaffrey, 2009; Subarya et al., 2006]. Spatial distribution of such events and their tendency to occur within a short interval are at odds with our current understanding of subduction zone seismicity.

Figure 4: Seismotectonic map of the Southeast Asian subduction zone. (a) Rupture patches of 2004, 2005, and 2007 are plotted along the Sunda trench with three distinct colours []. Along the convergent margins slab depth in the upper mantle is plotted [Hayes et al., 2018]. Black arrows represent the present-day velocities of Indian (IN), Capricorn (CP), and Australian (AU) plates with respect to the Sunda plate from NNR-MORVEL56 in no net rotation reference frame [Argus et al., 2011]. (b) Top view of the 3-D model with upper plate and base of the ocean lithosphere constructed from crustal ages. The longitudinal and latitudinal extent of the 3-D model is the same as the map in (a). (c) Three-dimensional view of the model with the geometry of the slab immersed in the upper mantle and lower mantle.

Presently, I am developing a novel 3-D spherical numerical models of the whole southeast Asian subduction zone to address the enigmas of giant earthquakes. In my models, I emphasize on the motion of converging lithosphere, slab geometry and forces arising from mantle flow.  I reconstructed the configuration of the Sunda margin fetching relevant datasets from slab2, crustal age and tomography in the lower mantle which are proxy for slab pull, ridges push, and large-scale mantle tractions (Figure 4b, c). This setup allows us to capture the complete dynamics of the subduction system. Recent numerical models of subduction zones have greatly advanced our understanding of tectonic forces that are generated at the convergent margins and how they drive plate convergence, deformation and episodic earthquakes [Sobolev and Muldashev, 2017; van Zelst et al., 2019]. However, most of these models did not consider the third dimension of the subduction zone. To critically analyse and understand the earthquake physics behind 1800 km rupture along the Indo-Australian plate boundary, I also consider property variations along the strike of the subducting plate.

My research aims to merge cutting edge numerical models to simulate the fundamental force balance of subduction and plate convergence to understand tectonic stress loading along megathrust. Preliminary results show that force balance at the convergent margin controls the coupling and stress evolution between subducting and overriding plates. These are the key physical factors that control the rupture length and periodicity of giant earthquakes. Though, the exact nature of interaction between the slab buoyancy, ridge push force, plate-coupling and deformation mechanism is yet to be quantified. There is also a plan to incorporate seismic cycles in my 3-D self-consistent dynamic models to constrain the time scale related to giant earthquakes: a compelling need to protect the human lives living near the megathrust regions.

Note: Advanced numerical geodynamics models are constructed using Underworld2 framework, an open-source tool which utilizes PIC-FEM approach. A big shout-out to the underworld developers and maintainers for their continual help and support.


Argus, D. F., R. G. Gordon, and C. DeMets (2011), Geologically current motion of 56 plates relative to the no‐net‐rotation reference frame, Geochemistry, Geophysics, Geosystems, 12(11).

Bilek, S. L., and T. Lay (2018), Subduction zone megathrust earthquakes, Geosphere, 14(4), 1468-1500.

Conrad, C. P., S. Bilek, and C. Lithgow-Bertelloni (2004), Great earthquakes and slab pull: interaction between seismic coupling and plate–slab coupling, Earth and Planetary Science Letters, 218(1-2), 109-122.

Dal Zilio, L., Y. van Dinther, T. Gerya, and J.-P. Avouac (2019), Bimodal seismicity in the Himalaya controlled by fault friction and geometry, Nature communications, 10(1), 48.

Gordon, R. G., and G. A. Houseman (2015), Deformation of Indian Ocean lithosphere: Evidence for a highly nonlinear rheological law, Journal of Geophysical Research: Solid Earth, 120(6), 4434-4449.

Hayes, G. P., G. L. Moore, D. E. Portner, M. Hearne, H. Flamme, M. Furtney, and G. M. Smoczyk (2018), Slab2, a comprehensive subduction zone geometry model, Science, eaat4723.

Heuret, A., C. Conrad, F. Funiciello, S. Lallemand, and L. Sandri (2012), Relation between subduction megathrust earthquakes, trench sediment thickness and upper plate strain, Geophysical Research Letters, 39(5).

McCaffrey, R. (2009), The tectonic framework of the Sumatran subduction zone, Annual Review of Earth and Planetary Sciences, 37, 345-366.

Ruff, L., and H. Kanamori (1980), Seismicity and the subduction process, Physics of the Earth and Planetary interiors, 23(3), 240-252.

Ruff, L. J. (1989), Do trench sediments affect great earthquake occurrence in subduction zones?, in Subduction Zones Part II, edited, pp. 263-282, Springer.

Sandiford, M., D. Coblentz, and W. P. Schellart (2005), Evaluating slab-plate coupling in the Indo-Australian plate, Geology, 33(2), 113-116.

Shearer, P., and R. Bürgmann (2010), Lessons learned from the 2004 Sumatra-Andaman megathrust rupture, Annual Review of Earth and Planetary Sciences, 38, 103-131.

Sobolev, S. V., and I. A. Muldashev (2017), Modeling seismic cycles of great megathrust earthquakes across the scales with focus at postseismic phase, Geochemistry, Geophysics, Geosystems, 18(12), 4387-4408.

van Rijsingen, E., S. Lallemand, M. Peyret, D. Arcay, A. Heuret, F. Funiciello, and F. Corbi (2018), How subduction interface roughness influences the occurrence of large interplate earthquakes, Geochemistry, Geophysics, Geosystems, 19(8), 2342-2370.

van Zelst, I., S. Wollherr, A.-A. Gabriel, E. H. Madden, and Y. van Dinther (2019), Modeling megathrust earthquakes across scales: One‐way coupling from geodynamics and seismic cycles to dynamic rupture, Journal of Geophysical Research: Solid Earth, 124(11), 11414-11446.