The “Graviquake” model, proposed in 2015 as an alternative to the elastic dislocation model, posits that normal faults are passive features dominated by coseismic gravitational collapse into a dilated crustal wedge, and that normal faulting is fundamentally distinct from strike‐slip and reverse faulting. Developed using finite‐element modeling before the 2016 central Apennines earthquake sequence, the model was revamped based on interpreted Differential Interferometric Synthetic Aperture Radar data from these events and used as evidence for a gravitational collapse episode. However, this interpretation relies on miscalculated elevation changes and is not corroborated by independent geophysical and seismological observations. Our analysis exposes fundamental flaws in the Graviquake model. By assuming that faults are passive players, it underrepresents the dynamic role of strain accumulation and release in rocks adjacent to faults. The hypothesized rapid expulsion of overpressurized fluids appears inconsistent with observed diffusion rates and lacks supporting seismological evidence. Part of the uplifted–subsided volume imbalance is likely an artifact arising from data processing, and in part is a transient effect due to the delayed response of the lower crust. Moment tensor analyses detect no isotropic components indicative of gravitational collapse, and observed ground motion and stress‐drop levels remain fully consistent with elastic dislocation theory. In addition, finite‐element modeling of normal faulting replicates observed surface deformation without invoking a collapsing wedge. The Graviquake model proposes a representation of normal‐faulting mechanics that differs significantly from established models and observations. Gravity does play a role in normal faulting, but the elastic dislocation theory remains the definitive framework of fault mechanics. Reinterpreting the 2016 earthquakes as a cascade of gravitational episodes, based on incorrect data processing and modeling, fails to substantiate the Graviquake hypothesis. Persistence in advocating this model could mislead seismic hazard assessment and undermine our understanding of normal faulting.

KEY POINTS

  • We critically evaluate the Graviquake model.

  • The Graviquake model is not physically viable.

  • Persistence in pushing this model undermines our understanding of normal faulting.

Gravity plays an integral role in deformation of the Earth’s crust, with forces such as ridge push and slab pull (Forsyth and Uyeda, 1975) leading to continental collisions, thrusting, subduction, as well as rifting and transform boundaries. Faulting accommodates much of this deformation in the brittle crust, often generating earthquakes.

Doglioni et al. (2015) have posited a specific role of gravity regarding extensional normal faulting, describing isolated normal faults as essentially gravity driven; they explicitly invoked the collapse of dilatated crust underlying normal faults as a necessary feature for their inception and progressive dislocation. More specifically, these investigators contend that “… in crustal extensional settings, gravity is the main energy source for hanging‐wall faults collapsing. Gravitational potential is about 100 times larger than the observed magnitude, far more than enough to explain the earthquake”. They also maintain that “… normal faults have a different mechanism of energy accumulation and dissipation (Graviquakes) with respect to other tectonic settings (strike slip and contractional), where elastic energy allows motion even against gravity”.

From the very title of their initial article, Doglioni et al. (2015) referred to these special earthquakes as “Graviquakes” (as they act “pro‐gravity”); conversely, standard strike‐slip faulting (“iso‐gravity”) and reverse‐/thrust‐faulting (“counter‐gravity”) earthquakes, which are caused by the accumulation of primarily elastic energy, were referred to as “Elastoquakes” (Fig. 1).

Based on these assumptions, Doglioni et al. (2015) moved on to state that the physics of normal‐faulting earthquakes is distinct from that of strike‐slip and compressional settings, and a companion article by Petricca et al. (2015) explored in detail the physics of the Graviquake model. Both articles highlighted important potential implications for seismic hazard assessment. For example, Doglioni et al. (2015) stated that “… the steeper the normal fault, the larger is the vertical displacement and the larger is the seismic energy released…”. This statement, however, is contradicted by the 28 December 1908 M 7.1 Messina Straits normal‐faulting earthquake (hereafter, M is used to indicate moment magnitudes); the largest and most catastrophic Italian earthquake of the instrumental era was caused by up to 5 m of slip over a 29°–42° dipping fault plane (e.g., Pino et al., 2009).

Ten years after its first publication, and having been invoked and discussed in a number of articles, the Graviquake model is still poorly constrained. In this study, we conduct basic analytical and numerical tests to assess its feasibility.

Although the following text focuses primarily on the observational consequences of the Graviquake model, it is relevant to report here the theoretical findings of Dahlen (1977), a renowned theoretical geophysicist who investigated the energy balance in earthquake faulting. Dahlen stated that “The need to consider a gravitational contribution to the energy release, especially if permanent elevation changes in the vicinity of the earthquake focus have occurred, has, on the other hand, often been noted, although seldom quantitatively explored. A not uncommon misconception appears to be that the gravitational energy contribution can be adequately accounted for by adding the first‐order work done against gravitational body forces to the classical expression ΔEclassical, for the energy change in the absence of gravity. We have shown this ‘ad hoc’ procedure to be completely unfounded, as it amounts essentially to having counted the relatively enormous change in gravitational potential energy twice. Both the total energy change ΔE and the seismic energy Es can be expressed in terms of quantities only on the fault surface of an earthquake, and in this form both ΔE and Es, and therefore also the seismic efficiency η, are explicitly independent of the rotation and the self‐gravitation of the Earth.”. Furthermore, Savage and Walsh (1978) reached similar conclusions; they showed that the principal energy exchange occurs between the gravitational field and the initial elastic energy stored in the medium; the resulting released energy is simply the small difference between the large gravitational and elastic energy changes. This does not imply that the gravitational field is not important because (1) it sets vertical and normal stresses, (2) influences fault stability, and (3) governs the hydrostatic pressure while balancing the effective stress and the shear stress controlling the accumulated elastic strain energy and its release at failure.

On 24 August 2016, a year after the publication of the article by Doglioni et al. (2015), an M 6.1 shock struck the town of Amatrice and other nearby villages, in the central Apennines, about 100 km northeast of Rome. It caused extensive damage and about 300 casualties, marking the onset of a dominantly normal‐faulting sequence that extended for nearly 80 km along strike in one of the most earthquake‐prone districts of the entire Italian peninsula (Fig. 2). For the following discussion, we note that all the moment magnitudes used for constructing our dataset of corner frequency versus moment couples are listed in the moment tensor (MT) catalog determined by Herrmann, Benz, and Ammon (2011 and updates). This MT catalog is similar to that published routinely by Istituto Nazionale di Geofisica e Vulcanologia (INGV; Scognamiglio et al., 2006).

On 26 October, an M 5.9 shock struck near Visso, about 35 km north of Amatrice, and on 30 October, an additional M 6.5 shock struck the historical town of Norcia and its surroundings, halfway between the first two largest shocks. The sequence continued into 2017 with a sequence of five shocks in the M range of 5.3–5.5 in the Campotosto Lake‐Capitignano area, about 10–15 km southeast of Amatrice. The sequence, and specifically the 30 October shock, gave Doglioni and coworkers a convenient and very well‐recorded case for testing the Graviquake model.

DInSAR observations supplied by the Japanese Advanced Land Observation Satellite (ALOS2) for the 30 October shock were initially analyzed by Valerio et al. (2018). Later on, Bignami et al. (2019) argued that “…the coseismically subsided hanging‐wall volume is about 0.12  km3, whereas the adjacent uplifted volumes amounted to only 0.016  km3. Thus, the subsided volume was ∼7.5 times larger than the uplifted one…”. They concluded that this imbalance between uplifted and subsided crustal volumes is systematic in normal‐faulting earthquakes and provides evidence for the existence of Graviquakes, events driven by gravity, not by the elastic dislocation mechanism.

Segall and Heimisson (2019) responded to Bignami et al. (2019) with a short note specifically devoted to refuting any violation of the elastic dislocation theory, and hence the need to invoke the hypothesis of the volume collapse. They first highlighted the work of Ward (1986), which suggested that the crustal volumes affected by a strong earthquake do not necessarily balance out between uplift and subsidence immediately after the event. Rather, this balance might only be attained over time, due to the viscoelastic behavior of the deepest portions of the Earth’s crust. This process is observed in the central Basin and Range province, in which long‐term geodetic analyses show persistent postseismic uplift that is fit by finite‐element modeling of the viscoelastic response (Thompson and Parsons, 2016). Therefore, it is normal to observe a temporary imbalance after an earthquake, initially as excess subsidence in extensional regimes and excess uplift in compressional regimes, simply due to the rheological properties of the Earth’s crust. Concerning these properties, we remark that the proposed viscosity of the deeper crust in the Apennines (Dalla Via et al., 2005) is similar to that estimated for seismic regions in the United States (e.g., Pollitz et al., 1998, 2000, 2001). In other words, the process equalizes postseismically, following viscoelastic uplift of the footwall.

Based on studies by Okada (1985) and Ward (1986), Segall and Heimisson (2019) derived a new equation to estimate the expected coseismic crustal volume imbalance. For a Poisson’s ratio of 0.25 and a dip angle of 45°, they found “… a net volume decrease of 0.08  km3, which is within 20% of the Bignami et al. (2019) estimate of 0.10  km3.”

They also criticized the methodology adopted by Bignami et al. (2019) to perform their volume calculations. The elastic dislocation theory dictates that the subsidence caused by a shallow normal‐faulting earthquake is concentrated in an elongated depression parallel to the seismogenic fault and slightly longer than the fault itself; ∼30 km in the case of the 30 October 2016 earthquake, which caused subsidence of up to one meter (Bignami et al., 2019; see their fig. 5). Nevertheless, Global Positioning System (GPS) data indicated that uplift affected a much larger area, extending as far as the Adriatic coast; a roughly 80 × 80 km region, in which absolute uplift values decrease gradually from a maximum of ∼12 cm to zero (see fig. 5 in Bignami et al., 2019). Segall and Heimisson (2019) argued that Bignami et al. (2019) miscalculated the volumes because they masked out areas where DInSAR data showed coseismic elevation changes between +3 cm and −3 cm (see fig. 5 in Bignami et al., 2019). This exclusion was somehow inevitable due to the 24 cm wavelength of the ALOS2 satellite sensor, which makes it virtually impossible to detect movements smaller than about 3 cm. For this reason, Bignami et al. (2019) stated that 3 cm represents “a realistic error range for the estimated co‐seismic displacement”: a honest admission of a strong bias that may have led to significant underestimations of the uplifted volume within the full dislocation field captured by DInSAR data.

According to Segall and Heimisson (2019), had the volumes been calculated correctly—that is, considering the large area where uplift could not be resolved due to the limitations of the technique—the DInSAR data from the 30 October earthquake would be fully compatible with the predictions of classical elastic dislocation theory after taking into account the caveats discussed by Ward in his 1986 article.

Segall and Heimisson concluded as follows: “We do not claim that the elastic dislocation model is unique. Occam’s razor, however, suggests that a simpler, well‐tested theory (elastic dislocation theory) should be preferred.”. The same conclusion, that is, that ruptures on both normal and reverse faults are consistent with elastic rebound theory, was reached through numerical modeling by Simpson (2024), who stated that: “Despite differences in energetics, model earthquakes on both reverse and normal faults are consistent with elastic rebound theory. Although normal‐fault earthquakes in models are accompanied by significant drops in gravitational potential energy, this occurs after and is caused by the release of elastic strain energy.”

Doglioni and his colleagues replied to Segall and Heimisson (2019) with an article authored by Bignami et al. (2020). Their response overlooked the two main arguments brought forward by Segall and Heimisson (2019)—the role of viscoelastic relaxation of the deeper crust and the lack of resolving power of DInSAR data—and contended that using the formulas proposed by Okada (1985) was appropriate also for a Graviquake, because “… coseismic slip affects an elastic medium.” They also maintained that “… the Segall and Heimisson (2019) formulation is therefore too simplified and much less precise than InSAR data…”; but the drastic reduction of their computational domain (i.e., the masked areas above) inevitably returns substantially underestimated uplift. The issues raised by Segall and Heimisson (2019) were also not addressed by Doglioni (2024), whose study strongly supports the Graviquake model.

It is important to acknowledge that estimates of uplift and subsidence from DInSAR measurements can be sensitive to the specific data processing methods employed. Factors such as the choice of reference areas, the correction for atmospheric phase delays, the spatial integration domain, and the temporal stacking strategy can all influence the resulting deformation fields (e.g., Pepe et al., 2005; Biggs et al., 2007). Alternative processing approaches, such as the application of atmospheric phase screens or longer time‐series analyses, may reduce uncertainties and recover weak uplift signals that are otherwise obscured by noise or decorrelation effects.

Although no DInSAR studies to date have specifically tested the volumetric collapse scenario proposed by the Graviquake model, the available Interferometric Synthetic Aperture Radar (InSAR) and Global Navigation Satellite Systems (GNSS) analyses of the 2016 central Apennines earthquake sequence (e.g., Cheloni et al., 2017) consistently support interpretations based on elastic dislocation models without invoking significant net volumetric collapse. Future studies applying refined InSAR time‐series analysis techniques may further improve the quantification of the uplift–subsidence balance, but the existing body of observations remains consistent with elastic fault slip as the dominant mechanism.

Following Segall and Heimisson (2019) and Simpson (2024) studies, we feel it is mandatory to reevaluate the physics of the Graviquake model from multiple perspectives that have not yet been explored; particularly from a purely seismological one, that so far has been overlooked by all the articles that have dealt with this new model.

The model

The idea behind the Graviquake model is that the differential stress on a normal fault is substantially lower than that on a reverse fault (Sibson, 1974; Mandl, 1988; Kohlstedt et al., 1995). This is suggested to be one of the reasons why normal‐faulting earthquakes systematically exhibit a smaller maximum magnitude (Neely and Stein, 2021) and smaller peak ground accelerations (McGarr, 1984) than reverse‐ and thrust‐faulting earthquakes. In what follows, we are going to take the Graviquake model literally, fully embracing the reasoning and the parameters adopted by its authors to see whether any contradictions or shortcomings arise.

Doglioni et al. (2015) describe the cornerstones of the Graviquake model as: “…during the interseismic period, while the lower crust shears steadily, the brittle upper crust is locked and a dilating wedge is inferred. The width of this triangle is here hypothetically imaged to affect an antithetic section to the locked fault, say 3 km thick.”. The extensional stress growing beyond a strength threshold “…determines the collapse of the normal‐fault hanging wall, dissipating mostly gravitational potential energy, being the elastic component a minor factor in the fall, if any” (the Graviquake model is summarized in the cartoon of our Fig. 1, modified from fig. 1 of Bignami et al., 2020).

A vertical cross section of the severely dilated crustal wedge envisioned by Doglioni et al. (2015) and by Petricca et al. (2015) is characterized by a triangular shape in cross section (Fig. 1, upper‐ and lower‐right frames). Before starting to slide on the fault plane and collapsing, however, this dilated wedge must support the weight of the entire hanging wall. For this reason, the cracks that permeate the dilated wedge need to be filled with fluids kept at lithostatic pressure at all depths (see fig. 1 in Doglioni et al., 2015). A lower fluid pressure (e.g., hydrostatic) would imply unsustainable stresses within the hanging wall. The authors of the Graviquake model did not specify what would be the source for the fluids filling the dilatant wedge, although a deep origin can be hypothesized without losing generality.

The abrupt collapse of the fluid‐saturated dilated wedge (described as “instantaneous” by Bignami et al., 2019) must cause a similarly quick migration of a very large volume of overpressurized fluids. Following Townend and Zoback (2000), we remark that the easiest way for such fluids to escape would be by leaking along the fault plane itself. In fact, Townend and Zoback demonstrated that it is the network of the critically stressed faults that represents the plumbing system of the Earth’s crust, and it is responsible for keeping the crust strong by limiting the pore‐fluid pressure near hydrostatic values. However, no massive fluid expulsion was observed along the normal faults that caused the most recent Apennines earthquakes (e.g., the 6 April 2009, Mw 6.1 L’Aquila earthquake, and the multiple mainshocks of the 2016 Amatrice–Visso–Norcia–Capitignano earthquake sequence). What scientists did observe was a large postseismic groundwater surplus between August 2016 and November 2017, estimated in the range between 4  108  and  5  108  m3 by Mastrorillo et al. (2020). According to these investigators, however, only part of this surplus can be explained by an increase of crustal permeability (e.g., Roeloffs, 1998; Brodsky et al., 2003; Manga and Brodsky, 2006; Manga et al., 2012); they attributed most of the remainder to the disruption in the upgradient portion of the aquifer caused by the earthquake fault rupture, promoting an increase of the hydraulic gradient which in turn determined a permanent increase of the global groundwater flow. In other words, most of the observed surplus would be the result of permanent modifications of the aquifers’ structure, and in particular, a result of the eastward shift of the piezometric divide caused by seismogenic fault slip.

As advocated by Doglioni et al. (2015), and more radically by Bignami et al. (2019), the Graviquake model prescribes a coseismic “free fall” of the hanging wall, which must cause the sudden expulsion of the fluids that permeate the fractures of the dilated wedge. Such a massive fluid displacement must occur as a diffusion process, which is expected to leave a clear signature in the recorded seismograms because its duration must be comparable with one of the time intervals characterizing the seismic source: the earthquake rise time, 0.5–3.0 s, or its rupture duration, ∼5–10 s, for an M 6.0–6.5 event (Heaton, 1990; Sommerville, 1995; Mai and Beroza, 2000). Finally, we note that the surplus of groundwater estimated by Mastrorillo et al. (2020) is between 4 and 5 times larger than the supposedly missing volume of 108  m3 estimated by Bignami et al. (2019).

As for the duration of the coseismic “free fall” of the crustal block hanging over the dilatant wedge, Doglioni et al. (2015) hypothesized a duration of 10–100 s. For the coseismic rupture of an M 6.0–6.5 event, however, a duration within such a range would produce only a “slow earthquake” (e.g., Ide et al., 2007). A much thinner dilated wedge than that specified by Doglioni et al. (2015) would allow for a much shorter diffusion time, but it would raise further issues because the crack density of the dilated wedge would become several times larger than that hypothesized by the same investigators to allow for the required fault slip.

For the sudden (“instantaneous”) collapse, Bignami et al. (2019) invoked a nonelastic, gravitational earthquake mechanism, hypothesizing that in the dilatant wedge “… once rocks are stressed above their yield stress during the interseismic stage, fractures occur, decreasing or eliminating the rock’s elasticity.”

Three alternative scenarios may be envisioned for the behavior of the hanging‐wall wedge:

  1. It collapses slowly, causing a “slow earthquake” with no high‐frequency energy radiation;

  2. It collapses quickly, raising issues about the permeability values needed for the diffusion of the fluids saturating the dilated wedge. Note, however, that a sudden collapse would radiate seismic waves, producing a substantial isotropic component in the MT solutions, as commonly observed in underground mine collapses (Mayeda and Walter, 1996; Ford et al., 2008): see the subsequent section on investigating the seismic MTs of the primary events in the sequence;

  3. The collapse is slow and decoupled from fast, elastically driven fault slip, causing no emission of seismic waves from the isotropic source.

Allthough scenarios 1 and 2 describe a Graviquake, in which the collapse and the sliding on the fault are coupled, scenario 3 describes a regular earthquake generated by an elastic shear dislocation (Elastoquake). Both Graviquake scenarios appear inconsistent with empirical observations, but the third scenario is also unacceptable, for it relies on the existence of a dilated subvertical wedge for which no empirical evidence has been produced.

Hydraulic behavior

According to Malagnini et al. (2012), following the 6 April 2009 M 6.0 L’Aquila earthquake, a postmainshock diffusion process dominated the seismicity in the Campotosto Lake area, 20–30 km north of the epicentral area, for approximately ten days. The process caused a series of three M ∼ 5 earthquakes, all having roughly the same centroid depth (7 km) and presumably being driven by northbound diffusion of pressurized fluids.

Interestingly, the fault system lying below Campotosto Lake ruptured again on 18 January 2017, in the latest phase of the 2016–2017 central Apennines sequence, causing a series of five M 5.0–5.5 events that marked a southbound diffusion episode. Malagnini et al. (2022) demonstrated the diffusive nature of triggering leading to these earthquakes. The occurrence of these fluid‐induced earthquake sequences confirms the hydraulic continuity of the whole fault system, favoring fluid migration.

The findings mentioned earlier demonstrate that, although episodes of fluid diffusions are very common in the central Apennines, their duration (from several hours to several days) is in direct contrast with Bignami et al. (2019) allegation that the closure of the dilated wedge during a Graviquake is instantaneous. Again, in what follows, we take Bignami et al.’s assertions literally and test if any contradictions arise.

First, we estimate the characteristic diffusion time for the fastest fluid expulsions (diffusion) that we have evidence for, and assume that such diffusion is the result of the coseismic “free fall” envisioned by Doglioni and colleagues. Then, based on the estimated duration, we infer the expected seismic radiation at high frequency to see if it matches the calculated corner frequencies and the observed ground motion. The diffusivity value used here for the central Apennines is based on results from Malagnini et al. (2022, see their fig. 8), who found that the fastest, massive diffusion episode was triggered by the Amatrice mainshock and was characterized by an extremely high diffusivity parameter:
From (1), we can retrieve an extremely high value of permeability:
We point out that the extreme diffusivity/permeability values listed in (1) and (2) are the highest ones among those estimated by Malagnini et al. (2022) during the multiple diffusion episodes that characterized the Amatrice–Visso–Norcia earthquake sequence. As such, they produce the most favorable scenario for the validation of a critical aspect of the application of the Graviquake model in the same region analyzed by Bignami et al. (2019). The described scenario applies only to the Amatrice mainshock, however, the other mainshocks provided evidence of much slower diffusion episodes. The value of diffusivity for the causative fault of the Amatrice earthquake (D=2·103  m2/s) can be compared, for example, with the two independent estimates of postseismic diffusivity along the Campotosto Lake fault: D=60±10  m2/s, provided by Malagnini et al. (2012) (corresponding to a diffusivity 4  1013κ6  1013m2), and the value of postseismic diffusivity D=50  m2/s, calculated by Malagnini et al. (2022). Brodsky and Saffer (2020) estimate the range of typical fault zone diffusivity values to be between 101  and  105  m2/s.
To obtain the characteristic time for the diffusion process, we use equation (3) of Townend and Zoback (2000):
in which τ is the characteristic time for a diffusive process (recall that Doglioni and coworkers report τ=10100  s), l is the characteristic width of the collapsing wedge (here we use l1·103  m), and D is hydraulic diffusivity.

The dilatant wedge was originally envisioned to be ∼3 km thick (Doglioni et al., 2015). In what follows, we use one‐third of such thickness (1 km) as an arbitrary choice rather than the original one of 3 km. This smaller thickness of the collapsing wedge would make our scenario more favorable to the Graviquake model.

Based on the diffusivity obtained for the Amatrice earthquake, and on our arbitrary choice of the thickness of the dilated wedge, we estimate the following diffusion time,
which is incompatible with the generation of the high‐frequency ground motion observed during the Amatrice earthquake (or, for that matter, during any of the mainshocks of the sequence).

Obtaining a diffusion time in the order of 10 s from equation (4) would require D=1·105  m2/s. We may then assume porosity ϕ2% (a representative estimate for limestones at depths between 5 and 10 km: see Giao and Nguyen, 2022), fluid viscosity η=1.9·104  Pa·s, fluid compressibility βf=5·1010  Pa1, and rock compressibility βr=2·1011  Pa1, all referred to a temperature, or to a depth of about 5 km (Townend and Zoback, 2000). Under these conditions, a diffusion time τ=10  s would correspond to an extremely high permeability value at depth (see Townend and Zoback, 2000), in the order of κ109  m2 (see equation 3).

Noir et al. (1997) calculated an even more extreme value of permeability (κ108  m2) along the faulted area of the 1989 Dobi earthquake sequence in central Afar. However, the basement of that region is primarily composed of basaltic lava flows, which typically form prismatic structures due to cooling; these rocks are subjected to tensile stress, as shown by the numerous fissures seen in the field. Near the surface and down to a depth of several tens of meters, these fissures often attain a width of several meters; this creates a geological environment that differs significantly from that of the central Apennines, in which permeability is expected to be orders of magnitude lower.

According to Manning and Ingebritsen (1999), common values of crustal permeability (i.e., the permeability maintained by critically stressed faults) obey the following power‐law decrease with depth (the variable z in equation 5, expressed in kilometers) as
This inference allowed Townend and Zoback (2000) to conclude that over 1–10 km scales, the permeability of the upper crust falls in the range 10171016  m2. They also pointed out that these values are typically associated with permeabilities in fractured, stressed regions of the crust. Finally, Brace (1980) highlighted that the permeability of the upper crust is generally too high to sustain pore pressures much above hydrostatic levels for extended periods, a key observation in geomechanics.

Based on the values reported earlier for the upper crust, it is also unlikely that fluids in a wedge could be maintained at lithostatic pressure for a significant portion of the seismic cycle. The high permeability in critically stressed faults, as first described by Brace (1980), and later by Manning and Ingebritsen (1999), would cause fluids to escape, resulting in pore‐fluid pressures closer to hydrostatic rather than lithostatic. Geological observations support that sustaining lithostatic fluid pressures requires very low permeability, which is inconsistent with many upper crustal measurements.

Role of fluids in fault zone processes

The Graviquake model posits that a sudden expulsion of overpressurized fluids triggers the gravitational collapse of a dilated wedge of crust. Although our analysis shows that rapid fluid expulsion at the rates implied by the model is inconsistent with observed hydraulic diffusion timescales and lacks seismological evidence, it is important to acknowledge that fluid migration plays a significant role in fault zone processes over longer timescales.

Thus, while rapid coseismic fluid expulsion capable of initiating a large‐scale wedge collapse remains unsupported, longer‐term fluid migration processes are both expected and observed in the aftermath of major earthquakes. These processes contribute to time‐dependent fault healing, permeability evolution, and fluid‐driven aftershock sequences.

In this context, fluid migration represents an important secondary mechanism that modulates fault behavior, but does not fundamentally replace the elastic strain accumulation and release process as the primary driver of normal‐faulting earthquakes. These considerations are consistent with numerous studies emphasizing the role of fluids in fault zone dynamics (e.g., Sibson, 1992; Cappa and Rutqvist, 2011; De Barros et al., 2016) without the need to invoke gravitational collapse as a dominant mechanism.

We conclude that our estimate for the duration of the volume collapse, together with the above discussion on crustal permeability of the central Apennines (based on our estimate from seismological observations and on what can be found in the literature), demonstrates that the hydraulic assumptions underlying the Graviquake model are physically implausible.

Implications of the volume imbalance

An important question to be asked is, what is the volume of fluids that one needs to relocate during a Graviquake? In keeping with the calculations by Bignami et al. (2019), who obtained a volume imbalance of 0.14  km3 between 24 August 2016 and November 2016 as a result of the Amatrice–Visso–Norcia sequence, and hypothesizing that the width of the bottom of the collapsed subvertical wedge was about ∼1 m (as assumed by Doglioni et al., 2015), to meet the volume imbalance we would need a surface of 100  km2, and the volume of water to be expelled in a matter of a few seconds would be approximately equal to the roughly 0.1  km3 imbalance. Again, we note that to sustain the crustal overburden, the dilated wedge needs to be fluid‐saturated at lithostatic pressure.

Doglioni (2024) envisioned an elongated and steeply‐dipping (or even subvertical) dilated wedge that reaches the depth of the brittle–ductile transition (see his fig. 10). We prepared the map and the cross sections shown in Figure 2 using a catalog of double‐difference relocated events (Waldhauser et al., 2021), showing that the seismicity during the 2016–2017 earthquake sequence extends well below 10 km depth, which is the maximum depth reached by the dilated wedge.

Furthermore, the potential constraint of the problem is the already mentioned large amount of missing uplifted volume (0.14  km3, as calculated by Bignami et al. (2019). Because the cumulative length of the entire fault front responsible for the 2016 earthquake sequence is around 40 km, and the average maximum rupture depth is around 10 km, then the dilated wedge can be estimated to extend for about 10×40=400  km2 and to have collapsed coseismically by about 0.25 m, as an average over the entire wedge.

We can now compute the isotropic moment of the elastic energy radiated from the volume collapse, using Bignami et al. (2019) estimates of the relevant parameters. Richards and Kim (2005) discuss two proposed relations for the task,
(originally from Müller, 1973), and
and argue they are the same, but are based on two different definitions of ΔV. Ford et al. (2008) used equation (7) for the Crandall Canyon mine collapse.
Assuming a Poissonian solid (λ=μ=3.0·1010  Pa), equation (6) would give
whereas equation (7) would give
Based on the definition of isotropic moment magnitude
We obtain that the isotropic moment magnitude of the wedge implosion is in the range:
(Recall that the reported magnitude of the 30 October 2016 Norcia earthquake was M 6.5).

From our calculations, we obtain a total deviatoric moment magnitude: MTOT=6.63 (M0TOT=9.97·1018  N·m) is the cumulative moment magnitude obtained by adding the moments of the three mainshocks of the sequence (Amatrice M 6.1, Visso M 5.9, and Norcia, M 6.5). In terms of seismic moment, the collapse represents between 50% and 94% of the cumulative deviatoric moment released by the three mainshocks (that is, between 33% and 48% of the total moment: cumulative isotropic + cumulative deviatoric).

Because the mechanics of the Graviquake model prescribes that the collapse and fault slip occur simultaneously, only the following two scenarios are plausible:

  1. The collapse is “instantaneous”, as claimed by Bignami et al. (2019). More precisely, its duration is short compared to the diffusion time. If this is the case, its radiation could not be neglected (see subsequent section on full MT calculations), yet an insurmountable issue would arise concerning the pore pressure inside the dilated wedge.

  2. The collapse occurs over the entire estimated diffusion time (see equation 4). Because the fault is forced to slip at the same rate, the process would generate a slow earthquake, which does not align with the observed high‐frequency ground motion (see subsequent section on earthquake ground motion).

Finally, we note that Lucente et al. (2010) documented strong fluid diffusion before the 6 April 2009 M 6.1 L’Aquila earthquake. According to these investigators, the fluids migrated from the footwall to the hanging wall of the earthquake‐causative fault; this observation does not confirm the migration of fluids out of the dilatant wedge predicted by the Graviquake model, and it is unclear which conduits the fluids should have adopted to leave the wedge coseismically.

Full MT solutions

Here, we examine the MT behavior of a Graviquake. First, we simulate synthetic waveform data for a composite model containing a double couple and a closing crack, and second, we invert the observed waveforms for the three mainshocks of the sequence. For these analyses, we used the velocity structure calibrated by Herrmann, Malagnini, and Munafò (2011), which is routinely used in the monitoring activities of INGV, and we adopted the complete local/regional waveform inversion method developed at UC Berkeley (e.g., Minson and Dreger, 2008; Dreger, 2018) and in use also at INGV (Scognamiglio et al., 2009).

Bignami et al. (2019) wrote that: “At the coseismic stage, both models generate fault motion, double‐couple mechanism, elastic waves generation, deformation of the surrounding volume, but they differ in several further constraints, such as interseismic fracturing, fluid motion, coseismic shear heating, volume folding, hanging‐wall motion, etc.”; and also that: “The double‐couple visible in the MT is required and implicit with the downward motion of the hanging wall along the main normal fault surface.”. However, there is no mention of the enormous isotropic moment that needs to be released at the time of the collapse of the dilatant wedge, implying that the correct tensorial representation of a Graviquake would be no different from that of an elastic dislocation. Large underground cavity collapses do indeed radiate seismic waves if they occur fast enough (e.g., Ford et al., 2008), though they produce far less static stress change in the crust than double‐couple slip, as modeled by Parsons and Velasco (2009).

Bignami et al. (2019) also hypothesized that the “missing volume occurs at the coseismic phase, hence it can be considered instantaneous…”. But an instantaneous closure of the subvertical dilating wedge described in their figure 7 (∼1 m) requires the detection of a significant isotropic component to the full MT (between 33% and 48% of the total moment) which, in turn, demands that first arrivals be controlled by the isotropic source (for a Graviquake this is a collapse with dilational or down first motions). However, this is not observed in the recorded seismograms.

To illustrate the expected behavior of a full MT inversion for the Graviquake model, we performed an inversion of synthetic data. We computed synthetic data for a source model that has a strike of 0°, dip of 50°, and a pure normal dislocation (rake = −90°). This mechanism was then added to waveforms for an inclined closing tensile crack that is roughly perpendicular to the dipping normal fault, as shown in the cartoon of Figure 1. Displacement Green’s functions were computed with the frequency–wavenumber integration code of Herrmann (2013) and the velocity model in routine use by INGV for a source depth of 5 km (as reported for the Amatrice earthquake). In constructing the synthetic data, 70% of the total moment was assigned to the double couple, and 30% to the closing crack. The synthetic data and Green’s functions were band‐pass filtered with a three‐pole, causal Butterworth filter between 0.02 and 0.05 Hz.

Figure 3 shows the results, in which panel a is the deviatoric inversion, which fails to fit the radial and vertical component records well because it is on those components in which the volumetric aspect of the closing crack is manifest. Figure 3b shows the full MT inversion, which of course fits the data perfectly, and the moment is partitioned as 34% volume decrease in a crack mechanism, and 66% deviatoric, with a dominant north–south‐striking normal‐fault mechanism. Figure 3c,d shows the MT decomposition of the full MT inversion result. This test shows that if present, a large volumetric component of the source would be recoverable from long‐period local/regional waveforms in which the volumetric and shear dislocation components are coseismic, as has been demonstrated in previous studies (e.g., Ford et al., 2008, 2009; Minson and Dreger, 2008; Alvizuri and Tape, 2016, 2018; Nayak and Dreger, 2018). It is also noted that even with the considered 30% closing isotropic moment, the first‐motion polarity over most of the focal sphere is predicted to be down, which is not observed.

Next, we invert the local/regional waveform data from the INGV broadband network for the full MT solution for the 24 August 2016 Amatrice mainshock. Five stations with good three‐component waveform coverage were used in the inversion (P, S, and surface waves). The data were instrument corrected and band‐pass filtered with a three‐pole, causal Butterworth filter between 0.02 and 0.05 Hz. In this frequency passband, moderate‐to‐large earthquakes may be considered as a point source in space and time, and simple layered velocity structures for Green’s functions are appropriate. Figure 4 shows the full MT results, which yield only a small isotropic component representing 4% of the total seismic moment. Finally, note that the orientation of the time‐domain moment tensor solution (see Scognamiglio et al., 2006), not shown, is in excellent agreement with the result of our full MT inversion.

To assess the significance of the volumetric component, an F‐test was used (e.g., Dreger, Tkalcic, and Johnston, 2000). In this test, we compare the variance of a deviatoric MT inversion with five parameters to the full MT with six. Because the data are composed of 5 stations by 3 components by 120 samples (1 sample per second), many samples are correlated, considering the 20 s low‐pass filter that was used. We therefore reduce the number of data in the F‐test calculation by a factor of 20 and take that as uncorrelated data. In this inversion, the ratio of the variances divided by the respective degrees of freedom (F‐factor) is 1.005, corresponding to an F‐test confidence of 47%, whereas for 95% confidence in the volumetric component the F‐factor would need to be 1.436 indicating that the small volumetric component is not statistically distinct from zero.

The sign of the volumetric component does indicate volume reduction; however, it is a statistically insignificant component of the full MT inversion and, in fact, a restrained deviatoric inversion fits the data nearly as well at 91.8% deviatoric versus 91.9% full MT variance reduction. Small isotropic components <∼10% can result in regional distance MT inversions from poor station coverage, the assumed velocity model, and noise (e.g., Dreger, Tkalcic, and Johnston, 2000; Panning et al., 2001). The full MT inversion indicates that the coseismic radiation of the Amatrice mainshock is deviatoric in nature, showing no evidence of a dilational volumetric component expected from a Graviquake.

We also inverted the local/regional waveform data from the INGV broadband network for the full MT solution of the 30 October 2016 Norcia mainshock (five good stations with three‐component seismograms, filtered between 0.02 and 0.05 Hz with a causal three‐pole Butterworth filter: see Fig. 5). The deviatoric solution shows a large calibrated linear vector dipole (CLVD) component (55.4% versus 44.6% double couple). This large CLVD can be attributed to the multiplane rupture identified by Scognamiglio et al. (2018) from the joint inversion of strong motion and GPS for the finite‐fault(s) mechanism. These investigators highlighted the activation of both a pure dip‐slip N155° trending, southwest‐dipping fault and of a second, N210° trending and northwest‐dipping fault, this last having a significant left‐lateral strike‐slip component at its start.

Allowing for a volumetric term does produce a 30% isotropic component, but it is not statistically significant (the variance reduction of the full MT solution is 68.7%, versus a 67.2% variance reduction of the purely deviatoric solution). The statistical F‐test reveals an F‐factor of 1.038 and a statistical significance of 56.7%. The F‐factor for 95% confidence is 1.436, and therefore, this larger isotropic component does not provide a statistically significant improvement in fit. However, more importantly, the full MT solution indicates a volume increase, which is opposite to what would be expected from the collapse of the Graviquake model.

A similar result was obtained for the 2016 M 5.9 Visso earthquake. Table 1 summarizes the MT result for the three earthquakes.

Summarizing the seismic MT analysis, the synthetic test for a 30% isotropic component combined with a normal slip mechanism shows that: (1) the P‐wave radiation pattern is mostly dilational or down vertical‐component first‐motion polarity, and that (2) with the full waveform inversion it would be unlikely to miss a significant volumetric component if it occurs coseismically with the fault shear dislocation. Positive or up first‐motion polarities are observed in these earthquakes (e.g., Ciaccio et al., 2021), not supporting a large volumetric decrease component in the seismic MT expected from the Graviquake model. Furthermore, seismic MT inversions of the Amatrice, Norcia, and Visso mainshocks yield no statistically significant volumetric components. The result for Amatrice shows a volume decrease of only 4% of the total moment and does not statistically improve the fit to the data as stated earlier. The MT results for the Norcia and Visso earthquakes yielded moderate volume increase components to the MT, which violates the Graviquake model expectations. In all three cases, the volumetric terms do not improve the fits to data with an acceptable level of statistical significance. Thus, all three events do not depart from deviatoric MTs consisting of double‐couple or compound double‐couple mechanisms. In addition, if the collapse of the Graviquake model dilatant wedge occurred more slowly before the earthquake to weaken the fault, slow enough to prevent coseismic wave radiation, it would have produced a sizable surface deformation field that would be observable prior to the earthquake, through either continuous GPS or DInSAR observations.

Finally, regarding the 4% volume decrease relative to the Amatrice mainshock, this value falls well within the known uncertainty bounds of full MT inversions for regional earthquakes. Studies by Dreger, Ford, and Walter (2000) and Ford et al. (2008) have shown that uncertainties in MT solutions, particularly for non‐double‐couple and isotropic components, can typically reach ±5%–10% due to factors such as velocity model inaccuracies, limited station azimuthal coverage, and signal‐to‐noise limitations. Therefore, the observed small isotropic contribution in the Amatrice event is not significant and does not provide evidence for a volumetric collapse as assumed by the Graviquake model.

High‐frequency ground motion

Several additional studies were published by Doglioni and coworkers after the 2015 article. The most recent one (Doglioni, 2024) reiterated a concept expressed in Doglioni et al. (2015): “Normal faults show higher seismic energy dissipation as the fault becomes steeper: at a given extensional strain rate, the steeper the fault, the larger the gravitational displacement; the longer the fault displacement and the faster the coseismic slip…”.

The expression “higher seismic energy dissipation” is somewhat vague in a seismological context: we assume this statement means “higher apparent stress (σa),” which, according to Doglioni (2024), characterizes the seismic radiation generated by steeper normal faults. To clarify, apparent stress is defined as the energy radiated per unit fault area per unit slip (σa=μErM0), in which Er is the radiated energy, μ is the rigidity, and M0=μAs is the seismic moment, in which A is the fault area and s is the earthquake coseismic slip. In simple words, because higher apparent stress means higher normalized energy radiation, σa is intimately related to the earthquake stress drop (e.g., the Brune stress drop): the higher the apparent stress of an earthquake of a certain magnitude, the higher its corner frequency.

Choy and Kirby (2004) observed that apparent stress is inversely correlated with fault maturity, and Choy et al. (2006) found that the lowest average apparent stress (0.3 MPa) is associated with thrust‐faulting earthquakes in subduction zones. Malagnini et al. (2010) calculated the apparent stress for earthquakes in sequences involving normal, strike‐slip, and reverse‐thrust faulting; they found that among all the events they analyzed, the 26 September 1997 M 6.0 Colfiorito normal‐faulting sequence in the central Apennines exhibited the highest apparent stress. Large normal faults are generally less mature than major faults in compressional or strike‐slip tectonic environments because their cumulative displacement is limited by their geometry. Moreover, we know that in Italy, active and seismogenic normal faulting is generally much younger than thrust and strike‐slip faulting (e.g., Valensise and Pantosti, 2001), and there is no evidence to support the idea that apparent stress is directly related to fault dip.

Using a spectral ratio technique developed by Mayeda et al. (2007) and applied in various research contexts and settings (e.g., Malagnini and Mayeda, 2008; Malagnini et al., 2008, 2011, 2021; Mayeda and Malagnini, 2009; Malagnini, Mayeda, et al., 2014; Malagnini, Munafò, et al., 2014; Malagnini and Dreger, 2016; Munafò et al., 2016), Malagnini and Munafò (2018) computed the corner frequencies of over 400 events from the 2016–2017 Amatrice–Visso–Norcia sequence using precise seismic moment estimates. They used these results to assess a regional ground‐motion predictive equation (GMPE), which turned out to be identical to that proposed by Malagnini et al. (2011). Figure 6 shows that the seismic waves generated by the events of the 2016 Amatrice–Visso–Norcia–Capitignano sequence were not the result of slow earthquakes. Rather, the mainshocks exhibit relatively high Brune stress‐drop values, well above 10 MPa for the highly damaging Amatrice and Norcia earthquakes.

Let us now focus on the 10–100 s collapsing time estimates proposed by Doglioni et al. (2015) and reiterated by Doglioni (2024). A visual inspection of Figure 6 shows that the Brune stress‐drop values predicted for the 30 October 2016 Norcia earthquake using these estimates are unrealistic: they range between 0.01 MPa (τ=100  s), and less than 1 MPa (τ=10  s).

Figure 7 provides an overview of the recording network available for the three mainshocks and their aftershocks, as well as for the events of the Campotosto Lake–Capitignano subsequence and related aftershocks. The histograms in Figure 8 highlight the accuracy of the regional GMPE calibrated by Malagnini et al. (2011) in predicting the ground shaking caused by the 2016–2017 earthquake sequence. These histograms include all available horizontal peak ground acceleration and peak ground velocity data recorded during the three mainshocks (all having M5.9). The remarkable performance of the GMPE at higher magnitudes also suggests that the L’Aquila mainshock, from which the GMPE was originally derived, does not display the expected seismological characteristics of a Graviquake.

Field evidence: An example from the Basin and Range province

Here, we present observations and models from the Basin and Range province in North America, showing that normal‐faulting earthquakes are not unique from thrust or strike‐slip faulting, as contended by Doglioni et al. (2015) and Petricca et al. (2015); rather, they are simply a response to the magnitudes and directions of the principal stresses in the elastic crust.

Like strike‐slip or thrust faulting, normal‐fault earthquakes focus broadly distributed elastic strain by slipping. For example, the 1906 San Andreas fault earthquake rupture concentrated hundreds of years’ worth of broad transform deformation onto a narrow strike‐slip fault zone. Understanding this process requires little more than Anderson (1905) description of fracture mechanics and faulting (Fig. 9), whereby distributed elastic strain is released during the fracture process. It is the conversion of regional elastic strain to slip onto a fault that led Doglioni et al. (2015) to hypothesize the existence of dilatational zones specifically beneath normal faults in the crust. Here, we provide data and modeling that demonstrate that such a zone is unnecessary.

The Basin and Range province is a broad region in western North America subject to continental extension resulting from northwesterly directed motion of the Pacific plate diverging from stable North America. The present‐day velocity field relative to stable North America is shown in Figure 10 and is derived from campaign and continuous GPS data (Hammond and Thatcher, 2004). The data clearly show a wide (∼1000 km) region under steady‐state elastic extension that is crossed by multiple locked normal faults spaced ∼25–35 km apart (e.g., Parsons, 2006). The region is thus subject to an approximately vertical greatest principal stress and an approximately horizontal least principal stress, leading to normal faulting, as shown in Figure 9.

In response to ongoing elastic extensional strain, a north–south chain of normal‐fault M6 earthquakes struck in 1954 in central Nevada (Fig. 10) (Doser, 1986, 1988). The first two earthquakes happened on the Rainbow Mountain fault, west of Dixie Valley: a 6 July M 6.2 shock, and a 24 August M 6.5 event. On 16 December, the M 7.5 Fairview Peak and M 6.5 Dixie Valley earthquakes ruptured minutes apart.

An unusual opportunity allowed for obtaining a well‐constrained set of geodetic observations associated with the 1954 earthquake series. By chance, the U.S. Coast and Geodetic Survey conducted a triangulation survey in central Nevada during the early summer of 1954, hence before the July–December earthquake series (Whitten, 1957). Those sites were remeasured in 1955 after the earthquakes and all relevelled during the follow‐up campaign. Meister et al. (1968) reported on these results because they pertained to coseismic strain, of which the vertical component is shown in Figure 11. A key observation is that most of the vertical coseismic deformation occurred as subsidence of the hanging wall, whereas the footwall position remained almost unchanged; this was determined from profiles that crossed the Fairview and Dixie Valley faults (Whitten, 1957; Meister et al., 1968; Koseluk and Bischke, 1981; Hodgkinson et al., 1996).

The coseismic subsidence of normal‐fault hanging walls is noted by Doglioni et al. (2015), who explained it by a dilatant wedge beneath the fault plane that expels fluids and causes the collapse. However, modeling studies show that hanging‐wall subsidence can be readily explained without that feature. Finite‐element modeling (Fig. 12) by Thompson and Parsons (2016) replicates observations through basic physics. The model consists of a brittle upper crustal layer that can deform by microfracturing, overlaid upon a ductile mafic lower crust, which in turn lies upon an upper‐mantle layer. The model is subjected to gravity and Basin and Range horizontal extensional strain rates, and has a 60°‐dipping frictional fault surface under a Coulomb criterion. There are no additional structural elements added to the model. The model is allowed to stretch because its western boundary moves at relative Pacific–North American rates, putting the entire elastic crust into a dilatant state until frictional fault failure occurs.

The modeling by Thompson and Parsons (2016) closely replicates observations (Fig. 13). Rather than a localized dilatant wedge creating space for the subsiding hanging wall, it instead collapses into the space made by the receding footwall through post‐failure elastic strain release/relaxation. It is very important to note that, while space is created coseismically, that space is not a permanent, pre‐existing feature. It only exists as elastic relaxation during the earthquake and represents the release of distributed extensional strain. The strain release behavior of a normal fault is the same as a vertical strike‐slip or dipping thrust fault in that all faulting styles act to release built‐up strain in the crust.

The model also predicts observed postseismic uplift that results from isostatic unweighting of the footwall. The full deformation cycle is shown in Figure 14. The occurrence of postseismic uplift after normal‐fault slip is clearly inconsistent with the Graviquake model, which is based on the concept that subsidence of the hanging wall occurs because of a dilatant wedge adjacent to the fault that expels fluid during an earthquake. In that scenario, all the vertical displacements are accommodated within the brittle crust through fluid expulsion. This leaves no explanation for postseismic uplift. However, numerical modeling clearly shows that the footwall displaces the viscoelastic substrate, which causes broad uplift over time (e.g., Stein et al., 1988; Gourmelen and Amelung, 2005; Thompson and Parsons, 2016) (Fig. 15).

This study evaluates the Graviquake conceptual model proposed by Doglioni et al. (2015) and later applied to the 2016 central Apennines seismic sequence by Bignami et al. (2019). Although Segall and Heimisson (2019) already dismissed the significance of the volume imbalance hypothesized by Bignami et al. for the central Apennines sequence, our analysis assumes the physical validity of the Graviquake model and scrutinizes its claims. From this standpoint, we highlight inconsistencies and critical issues, particularly in the application to the 2016–2017, Amatrice–Visso–Norcia–Capitignano earthquake sequence in the central Apennines.

The core of the Graviquake model is the concept that normal (progravity) faulting fundamentally differs from strike‐slip (isogravity) and reverse (countergravity) faulting. According to the model, normal faults are passive frictional structures in which the coseismic gravitational collapse of the hanging wall occurs without significant elastic rebound, and therefore, little or no footwall uplift is expected. The main cornerstone of the Graviquake model is the coseismic collapse and closure of a subvertical dilated wedge (see the cartoon in figs. 1 and 7b in Doglioni, 2024). When this collapse occurs, the crack density (vertical cracks) within the dilated wedge exceeds an unspecified threshold, leading to a total loss of strength. Doglioni (2024) suggests that in strike‐slip and reverse faults, seismic energy is released through elastic rebound across the fault, in marked contrast with normal faulting, producing events termed Elastoquakes.

To begin with, we remark that the idea that normal faults partially or fully violate the principles of elastic dislocation is not supported by over four decades of investigations of normal‐faulting earthquakes using standard elastic dislocation modeling. The methodology was popularized in the early 1980s by studies that used combined geodetic, seismographic strong motion, and surface faulting observations to derive the fault parameters of sizable crustal earthquakes (e.g., Stein and Barrientos, 1985; Crosson et al., 1986). Since then, many normal‐faulting earthquakes worldwide have been investigated using various combinations of leveling, GPS, and DInSAR data, generally processed using the convenient formulations and codes made available by Okada (1985). Elastic dislocation calculations were always seen to reproduce the available geodetic evidence while matching all other coseismic observations.

Petricca et al. (2021) used 32 recent earthquakes from different tectonic settings (M range: 5.7–7.3) to explore the geometric relationships between fault slip at depth, coseismic elevation changes, and damage distribution. Because “… coseismic surface deformation can be obtained by dislocation theory of an elastic medium as shown by Okada (1985)…”, they compared the DInSAR signature of each event with the ground shaking caused by that event, reaching the easily foreseeable conclusion that “…higher macroseismic intensity often corresponds to the area of larger vertical displacement (either downward or upward), apart local site amplification effects…,” regardless of the tectonic environment the earthquake is in. Yet, why should we trust the DInSAR signature of normal‐faulting earthquakes obtained from elastic dislocation modeling, if these events are generated by a process that is mostly gravity‐driven? Possibly aware of this conundrum, Petricca et al. (2021) felt the need to specify that “… The deformation of the ground around the fault activated during an earthquake can be described by a dislocation in an elastic half‐space (e.g., Kanamori, 1973; Okada, 1985), regardless of the source energy, which is elastic for thrust and strike‐slip faults, whereas it is gravitational for normal faults…”. This statement seems to contradict the cornerstones of the articles by Doglioni et al. (2015), Petricca et al. (2015), and Bignami et al. (2019). That is, the elastic dislocation, which is well described by Okada (1985) equations, is independent of the energy that causes it.

In the Graviquake model, the existence of a collapsing dilated wedge is hypothesized based on a volume deficit allegedly revealed by DInSAR data (Bignami et al., 2019). If we accept its existence, before the earthquake this dilated wedge must be filled with fluids at lithostatic pressure; a drained wedge would be unable to support its weight and that of the hanging wall above and behind it (see fig. 7 in Doglioni, 2024). In addition, the overpressurized, fluid‐saturated wedge should remain sealed until a Graviquake triggers rapid fluid diffusion. At that point, pressurized fluids start diffusing, flooding the surrounding crustal volume and directly affecting fault slip and slip velocity. The duration of coseismic slip for such a gravity‐driven fault motion must correspond to the fluid diffusion time, whereas in the Graviquake model it is assumed that the slip velocity implicitly depends on instantaneous fluid flow. Furthermore, a wedge collapse time of 10 s would require a staggering value for the diffusivity parameter: D=1·105  m2/s. Note that the duration of slip on the fault and the fluid‐diffusion time could be different only if elasticity were the driving phenomenon.

In an earlier section, we estimated a minimum diffusion time of several hundred seconds, which conflicts with the high‐frequency ground motion observed during the mainshocks, for which characteristics are summarized by the corner frequencies plotted in Figure 6. An estimate of 10–100 s for the duration of the “pro‐gravity coseismic collapse” is presented by Doglioni et al. (2015) and Doglioni (2024), but no further explanation is given. Even these estimates fail to match the observed high‐frequency ground motion.

To deal with the alternative scenario of an “instantaneous” closure of the subvertical dilated wedge, as proposed by Bignami et al. (2019), we defined a point‐source coherent with the Graviquake model: a double‐couple mechanism combined with a significant isotropic component. This model shows that the volumetric component significantly alters the waveforms in the long‐period band used in regional MT inversion, and would lead to only down first‐motion polarities. Bignami et al. (2019) never mention the need for a substantial isotropic component to make the Graviquake model a viable hypothesis. On the contrary, and quite paradoxically, the very first lines of their article read: “Regardless the tectonic setting, any motion along fault planes determines the shear between two volumes generating the double‐couple mechanism”: a very reasonable sentence, reportedly taken from Stein and Wysession (2009), but one that does not lend any support to the Graviquake model.

In reality, the complex source in the Graviquake model must involve both the collapse of the vertically elongated dilated wedge and the simultaneous double‐couple source along the sliding fault plane. The former represents an isotropic source that accounts for the closure of the wedge, whereas the latter accounts for the frictional sliding of the crustal wedge along the fault. For the central Apennines earthquake sequence, the cumulative moment of the isotropic source (M0ISO=5·1018  N·m is a lower bound for this quantity) corresponds to a moment magnitude of MISO6.43, which represents ∼52% of the cumulative deviatoric moment from the three mainshocks in the sequence (M0TOT=9.56·1018  N·m), equivalent to a cumulative deviatoric moment magnitude of MTOT=6.62.

In summary, let us assume once more that the dilated wedge proposed by Doglioni and colleagues exists and behaves as described in the Graviquake model. Our analyses show that the velocity of wedge closure is a decisive parameter, leading to two opposite scenarios:

  • the wedge closes “instantaneously”: this would cause significant energy radiation from the MISO6.43 cumulative isotropic source associated with the specific mainshock, but would also lead to (1) first‐motion dilatations across most of the entire focal sphere, and (2) massive expulsion of fluids from the fault zone. Both these circumstances are unverified: the observed first motions show distinct compressional and dilatational quadrants, and the observed fluid flow increase is smaller and slower than expected, and it is better justified by coseismic modifications of shallow aquifers than by a large increase in permeability;

  • the wedge closes over a minimum diffusion time of several hundred seconds; this would also limit the rupture time on the fault, preventing the radiation of high‐frequency seismic energy, which contradicts the observed ground motions. The largest events in the 2016 central Apennines sequence show a relatively high Brune stress drop, exceeding 10 MPa (Malagnini and Munafò, 2018).

A third scenario is possible, in which the volume collapse is slow and is decoupled from a fast dislocation. Such scenario, however incompatible with the empirical hydrological evidence (the documented amount of fluid mobilized is 4–5 times that of the alleged volume collapse, and the excess discharge does not take place on the activated fault), is incompatible with the Graviquake model.

We also remark that the purported volume imbalance is partly attributable to the modeling strategy of DinSAR data (Segall and Heimisson, 2019), and partly an expected occurrence that is bound to vanish over the next few decades.

Finally, the finite‐element modeling by Thompson and Parsons (2016) replicated observations of normal faulting through basic physics, showing that the void volume of the envisioned dilated crustal wedge—a fundamental feature of the Graviquake model—is unnecessary. Furthermore, the Graviquake model is inconsistent with the widely observed broad postseismic uplift that occurs over decades following normal‐fault earthquakes. The duration and rate of postseismic uplift are consistent with displacement of viscoelastic lower crust and upper‐mantle rocks beneath the footwall caused by downward hanging‐wall displacement, as shown in Figure 14 (e.g., Stein et al., 1988; Gourmelen and Amelung, 2005).

It is important to recognize that gravitational forces do influence the mechanical behavior of normal faults. Gravity affects the orientation of principal stresses, contributes to the effective normal stress acting on fault planes, and plays a role in controlling the dip angles and structural evolution of fault systems (Dahlen, 1977; Ader et al., 2013). In this sense, gravity acts as a modulating factor within the broader framework of elastic strain accumulation and release. However, gravity alone does not provide the primary driving force for coseismic rupture. Instead, elastic energy stored by tectonic loading remains the fundamental source of seismic moment release during normal faulting earthquakes. Our observations and modeling results are fully consistent with this classical understanding of earthquake physics.

In summary, while gravitational forces play an important role in shaping the tectonic stress regime and modulating the mechanical behavior of normal faults, they do not serve as the primary driving mechanism for coseismic rupture. The cumulative seismological, geodetic, and hydrological evidence examined herein supports the classical interpretation of normal faulting earthquakes as elastic dislocation events governed by tectonic loading. Observations from the 2016 central Apennines earthquake sequence, including surface deformation patterns, ground‐motion characteristics, MT solutions, and hydraulic considerations, consistently align with elastic strain accumulation and release models. No independent evidence supports the occurrence of large‐scale volumetric collapse as proposed by the Graviquake hypothesis. Although continued refinement of observational and modeling techniques is essential, the elastic dislocation framework remains the most robust and comprehensive model for understanding normal‐fault seismicity and for informing seismic hazard assessments.

Given all these circumstances, we conclude that the Graviquake model is not physically viable: as such, it should be critically re‐evaluated to avoid premature incorporation into the broader scientific literature (e.g., Muldashev et al., 2022), in dissemination publications (Pierotti and Prandi, 2018) and in online encyclopedias (see Data and Resources).

Waveforms used in this article were downloaded from the European Integrated Data Archive (EIDA) repository (http://www.orfeus-eu.org/data/eida/). Details of processing for the Field evidence: an example from the Basin and Range Province section is described, and time series and velocities are archived at http://quake.usgs.gov/research/deformation/gps/auto/WesternBRPBO/. The online encyclopedias of Graviquake are available at https://en.wikipedia.org/wiki/Graviquake. All websites were last accessed in May 2025. Ground‐motion parameters used to make the histograms shown in Figure 8, from 91,152 single‐component recordings collected during 872 events of the 2016–2017 seismic sequence of the central Apennines (peak ground acceleration [PGA], peak ground velocity [PGV], pseudospectral acceleration, PSA(f)), are available in csv format in the supplemental material.

The authors acknowledge that there are no conflicts of interest recorded.

The authors express their sincere gratitude to the reviewer, and to Associate Editor Yann Klinger, for their thoughtful and constructive comments, which have significantly contributed to improving the clarity and quality of this article. Their critical insights and careful reading helped refine both the scientific content and the presentation of this paper.

Supplementary data