# Evidence of supershear during the 2018 magnitude 7.5 Palu earthquake from space geodesy

Anne Socquet, James Hollingsworth, Erwan Pathier, Michel Bouchon

► To cite this version:

Anne Socquet, James Hollingsworth, Erwan Pathier, Michel Bouchon. Evidence of supershear during the 2018 magnitude 7.5 Palu earthquake from space geodesy. Nature Geoscience, 2019, 12 (3), pp.192-199. 10.1038/s41561-018-0296-0. hal-02416371

HAL Id: hal-02416371

https://hal.science/hal-02416371v1

Submitted on 2 Oct 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Evidence of supershear during the 2018 magnitude 7.5 Palu earthquake from space geodesy

Anne Socquet*, James Hollingsworth, Erwan Pathier and Michel Bouchon

A magnitude 7.5 earthquake hit the city of Palu in Sulawesi, Indonesia on 28 September 2018 at 10:02:43 (coordinated universal time). It was followed a few minutes later by a 4-7-m-high tsunami. Palu is situated in a narrow pull-apart basin surrounded by high mountains of up to 2,000 m altitude. This morphology has been created by a releasing bend in the Palu-Koro fault, a rapidly moving left-lateral strike-slip fault. Here we present observations derived from optical and radar satellite imagery that constrain the ground surface displacements associated with the earthquake in great detail. Mapping of the main rupture and associated secondary structures shows that the slip initiated on a structurally complex and previously unknown fault to the north, extended southwards over 180 km and passed through two major releasing bends. The 30 km section of the rupture south of Palu city is extremely linear, and slightly offset from the mapped geological fault at the surface. This part of the rupture accommodates a large and smooth surface slip of 4-7 m, with no shallow slip deficit. Almost no aftershock seismicity was recorded from this section of the fault. As these characteristics are similar to those from known supershear segments, we conclude that the Palu earthquake probably ruptured this segment at supershear velocities.

Indonesia is a country sadly known for its exposure to extreme earthquake and tsunami hazard¹,². Eastern Indonesia is located in a complicated plate tectonic context; the Sunda, Philippine and Australian plates meet at a triple junction, and regional deformation is characterized by rapid microblock rotations and intense faulting across an area of ~500 km² (refs. ³⁻⁷).

The 2018 earthquake occurred on the active strike-slip fault of Palu-Koro⁸,⁹, which is the main plate boundary structure that accommodates the relative motion between the Makassar block to the west and the North Sula block to the east⁷,¹⁰,¹¹ (Fig. 1). The strike-slip fault connects to the Minahassa subduction zone to the North, and to the Matano strike-slip fault to the south⁷,⁸,¹⁰,¹², and is proposed to have accommodated a geological offset of 120–250 km (ref. ¹³). Three moment magnitude ($M_w$) 6.8–8.0 earthquakes have been identified during the last 2,000 years by palaeoseismological investigations⁸,¹⁴. Global positioning system (GPS) measurements made between 1992 and 2003 constrained the slip rate of the Palu-Koro fault to 40 mm yr⁻¹ (consistent with the Holocene strike-slip rate of 35 ± 8 mm yr⁻¹ (ref. ⁸)), and showed that the displacement profile across the fault was compatible with a locking depth of ~12 km (refs. ⁷,¹⁵,¹⁶) (Fig. 1). A large earthquake was therefore expected on this fault.

The significant seismic and tsunami hazard posed to the area has long been known¹,². Over the past century, six tsunamigenic earthquakes have hit the region² associated with the Minahassa subduction zone, the offshore portion of the Palu-Koro fault or offshore faults in the Makassar straight to the west. An earthquake of comparable magnitude (7.4) occurred on the Palu-Koro fault in 1968, and also resulted in a tsunami with wave heights up to 10 m (ref. ²). The most recent significant events to occur in the region involved a sequence of earthquakes (1996–1998) that broke both the Minahassa trench and the Palu-Koro fault, and have been shown to interact¹⁰.

The epicentre of the 2018 earthquake was located ~75 km north of the city of Palu (Fig. 1). The moment tensor indicates a dominant left-lateral strike-slip, with a small normal component on an

east-dipping fault. The initial analyses of the teleseismic data indicate that most of the moment release occurred in less than 30 seconds, with propagation of the rupture at an average velocity of 4–4.1 km s⁻¹ (refs. ¹⁷⁻¹⁹). The mainshock was preceded by many smaller magnitude events that occurred close to the epicentre the days before (Fig. 1). The largest foreshock ($M_w$ 6.1) occurred 3 hours before the mainshock, and ~20 km from the mainshock epicentre. The $M_w$ 7.5 earthquake was followed by a series of aftershocks with a maximum magnitude of $M_w$ 5.8 (Fig. 1), which is low compared to the magnitude of the mainshock, as expected for supershear earthquakes²⁰. Aftershocks are mostly located north of Palu city along the so-called Sulawesi neck, and 40–50 km south of Palu city, whereas few have occurred so far in the Palu basin area (Fig. 1).

## Intermediate displacement field as seen from InSAR

Wide-swath synthetic aperture radar interferograms computed using pre- and post-earthquake images (21 August 2018 and 2 October 2018) acquired by the Advanced Land Observing Satellite-2 (ALOS-2) satellite in ScanSAR mode constrain the coseismic deformation field in the satellite line-of-sight (Fig. 1 and Supplementary Fig. 1). Significant deformation is visible up to ~50 km away from the epicentre. This medium- and far-field deformation pattern is particularly important for constraining the distribution and magnitude of the slip on the fault plane at depth below the surface (typically, >3 km). Close to the rupture, large strains result in multiple phase cycles between neighbouring pixels, which causes, together with intense shaking and damage, decorrelation of the interferometric synthetic aperture radar (InSAR) signal. Surprisingly, the earthquake did not appear to initiate on the mapped geological fault (located offshore at the latitude of the epicentre⁵,⁸,⁹), but rather nucleated inland in a complicated setting within the Sulawesi neck, close to the Balaesang peninsula, which corresponds to a releasing bend. The asymmetry of fringes (Fig. 1 and Supplementary Fig. 1), with more fringes on the eastern side of the fault, indicates that the fault plane dips eastward, in agreement with the focal mechanism. Four main lobes of deformation are clearly visible. The orientation
<figure><img src="image_2.png" /></figure>

**Fig. 1** | **Setting and measured surface displacements associated with M<sub>w</sub> 7.5 Palu earthquake.** The epicentre is indicated by a black star. The focal mechanism is from the USGS. The coloured dots represent one month of foreshock (blue) and aftershock (yellow and red) seismicity (different temporal colour scales are provided on each panel to show the detailed evolution of the seismic sequence). The topography is from the Shuttle Radar Topography Mission. **a**, The unwrapped ALOS-2 interferogram (21 August 2018 to 2 October 2018) showing the coseismic displacement in the line of sight (© JAXA). A relative increase between two pixels means a relative displacement towards the satellite. The surface trace of the ruptured fault is shown as a red line. Grey circles represent the background seismicity. The geological traces of the main faults are shown in black. **b**, Map of the horizontal surface displacement computed by the correlation of Sentinel-2 images. The arrows show the horizontal displacement, whereas the colours correspond to the N-S component of the displacements. **c,d**, Distribution along the rupture of the surface coseismic offsets extracted at 0.01° of the fault (red, north-south; blue, east-west) (c) and of the azimuth (green, azimuth of the slip vector; black, azimuth of the fault trace) (d). LOS, line-of-sight.

of the northern lobe is towards the northwest, which suggests that the rupture turns towards the northwest at latitude 0.1°S to continue offshore and ultimately connect up with the Minahassa trench further north. East of the fault, south of Palu city, the phase gradient indicates a range increase compatible with a left-lateral strike-slip and possible subsidence. North of 0.7°S, the phase gradient shows a range decrease, which suggests that uplift occurred to the east of the fault, leading for local transpression in the area of the Sulawesi neck.

## Rupture characteristics and evidence for supershear

Subpixel optical image correlation²¹ from Landsat-8, Sentinel-2 and WorldView satellite imagery (Figs. 1 and 2 and Supplementary Fig. 1) allows us to map the trace of the rupture in great detail and to quantify the horizontal component of the coseismic displacement at the surface. Critically, the optical correlation data provide valuable constraints on the near-field displacement pattern (in the region where the InSAR decorrelates), which thus reveals details of how the fault broke the surface and gives resolution to the shallow part of our fault slip model. South of Palu city (0.9°S), the north-south

(N-S) displacement field (Fig. 2) shows that the trace of the 2018 rupture reproduces, to the first order, the shape of the Palu-Koro fault, as documented from tectonic geomorphology and geological investigation⁵,⁸,⁹. South of the Palu coastline (that is, south of 0.9°S), the rupture is linear and strikes ~N172°. Near 1.187°S (~33km south of the coastline), the rupture bends sharply to the southeast for ~9km, where it forms a major releasing bend before recovering its initial azimuth and continuing for another ~20km. North of Palu city, the rupture disappears offshore within the Palu bay, and reappears 21km further north within the Sulawesi neck, where a much smoother displacement gradient can be followed northwards for 60km (probably due to a buried slip that does not come up to the surface), until it reaches the Balaesang Peninsula releasing bend at latitude 0.1°S. The azimuth of this northern segment (~182°) is rotated ~10° clockwise from that of the Palu segment to the south (~172°); the change occurs somewhere in the bay of Palu.

The high-resolution displacement field (40 m resolution) allows a structural examination of the rupture south of Palu, which is made up of six segments (A–F, Figs. 2 and 3) (detailed description in Supplementary Information). Some segments appear to be
<figure><img src="image_3.png" /></figure>

**Fig. 2 | Detailed features of the Palu rupture.** **a,** Slip vectors determined from the near-field (primary rupture, blue arrows) and medium field (primary + secondary, red arrows) ruptures. **b,c,** N-S (**b**) and E-W (**c**) displacement (m) fields from the correlation of pre- and post-earthquake Sentinel-2 images. The fault rupture is shown in black (bold is the primary rupture and dashed is the secondary one). **d,** Topography with the different segment lengths (m) labelled as in Fig. 3. **e,** Azimuth of the ground displacement (°). **f,** High-resolution E-W displacement map for the northern Palu section of the rupture, from correlation of WorldView satellite images (©2018 DigitalGlobe, a Maxar company; pre-image, 20 February 2018; post-image, 2 October 2018). The white decorrelated patch corresponds to a large mudslide, visible in the post-event WorldView imagery.

exceptionally sharp, expressed as a narrow surface rupture, with little coseismic deformation taken up by distributed shear off the main fault trace. Significantly, segments A and B do not follow the long-term fault morphology, lying 2–3 km from the bedrock range front within the basin (Fig. 2d). The coseismic displacement (4–7 m) is remarkably smooth and almost pure left-lateral strike-slip, whereas the fault-normal component is accommodated on secondary structures located off the main fault (Figs. 2 and 3). This is precisely what was observed along the supershear segments of the 1999 Izmit<sup>22-24</sup> and 2001 Kunlun<sup>25</sup> earthquakes. In both cases, the strike-slip rupture trace was located a few kilometres from the obvious bedrock fault structure, whereas the fault perpendicular component was accommodated on secondary structures through a slip partitioning mechanism<sup>26</sup>.

Fault segments that have hosted supershear ruptures have been shown to share specific structural characteristics (straight and

without major structural complexity) and very smooth coseismic slip in pure mode II, which reaches the surface (even though some faults may be dipping or accommodating oblique movement over the long term<sup>27,28</sup>). As pure strike-slip does not produce topography, ancient traces of such ruptures are challenging to identify, in particular in a humid climate setting. However, subtle evidence for laterally confined river channels in Palu basin were found<sup>9</sup>, suggesting that a fault branch crossed straight through the sediments, which in turn led the authors to conclude that the Palu-Koro fault may be capable of generating supershear ruptures.

Our observations of the M<sub>w</sub> 7.5 surface rupture south of Palu city are remarkably consistent, and together suggest that a rupture at a supershear velocity on some segments is highly probable. Segments A, B and E are all characterized by minor normal slip, a very straight azimuth, no differential slope, nearly uniform slip (especially when considering the primary and secondary slip) and
<figure><img src="image_4.png" /><figcaption>**Fig. 3** | Characteristics of surface slip for the Palu rupture. From left-to-right. Left-lateral (red) and fault-normal (blue) slip profiles plotted along strike by latitude. Orange lines show the slip along the secondary fault in the releasing bend. Solid lines are the total displacement over a 4 km swath (accounting for the primary and secondary displacement); dashed lines are the near-field displacement on the primary rupture. The OFD for the strike-slip component is computed from the difference between the primary and primary and secondary displacement estimates. Fault-strike azimuth, topography, local slope and change in slope across the fault are shown, respectively. The far right shows the location of the primary fault rupture, which is made up of several straight segments generally less than 15 km. Orange bars show the divisions between different fault segments; the grey bar highlights the releasing bend.</figcaption></figure>

a slip vector parallel to rupture trace (Figs. 2 and 3). Along segments A and B, apparent normal faulting occurs west of the rupture trace (although some liquefaction and mudsliding may contribute to this short wavelength signal (Fig. 2f)). Nevertheless, they are indicative that some limited but clear partitioning along A and B occurs.

The slip vector orientation suddenly changes at the start of segment C when the rupture meets the bedrock and begins to run along the bottom of the hills, along the geological trace of the Palu-Koro fault. Clearly, strike-slip and some normal component of slip are now on the same fault on segment C (at least at the surface)—with some additional minor thrust motion east of the primary rupture. This is also consistent with the low component of off-fault deformation (OFD) on this segment (7% (Fig. 3)). Enhanced slip localization may result from the strong contrast in material rigidity (that is, bimaterial interface) along this segment<sup>29,30</sup>. This may suggest that the rupture velocity decreases on segment C, which would be consistent with the stepover the rupture will undertake at the end of segment C (at supershear speed the rupture would try to go as straight as possible). Nevertheless, longer-wavelength slip vectors (red arrows in Fig. 2a) indicate a pure strike-slip motion at rela-

tively shallow depths along segment C, which is also consistent with supershear rupture.

Near 1.187°S, the rupture bends sharply to the southeast for ~9 km along segment D, where it forms a major releasing bend. The slip vectors are relatively parallel along segment D, with normal slip partitioned onto a neighbouring structure 5 km to the northeast (Figs. 2a and 3). Further south, the fault recovers its initial azimuth and the slip diminishes before terminating in the bedrock mountains near 1.41°S. Topographic amplification of the horizontal shear waves has been proposed as a mechanism for inducing off-fault damage in bedrock mountain settings<sup>31,32</sup> and may contribute to the increase in OFD for segments E (25%) and F (33%) (Fig. 3).

**Slip distribution.** To obtain a macroscopic view of the rupture process, the surface displacements from satellite geodesy were used to invert for the slip distribution on a buried fault plane (Fig. 4 and Methods). The fault geometry was constrained by the trace of the rupture at the surface and implied a varying strike along the fault that reproduces the observed releasing and transpressive bends. The inclusion of surface displacements near the fault (from optical
<figure><img src="image_5.png" /></figure>

<figure><img src="image_6.png" /></figure>

<figure><img src="image_7.png" /><figcaption>**Fig. 4** | Model of slip distribution on the fault. **a,b,** Model of slip distribution on the fault for the dip-slip (**a**) and strike-slip (**b**) components. **c-e,** Curves (dip-slip in blue and strike-slip in red) show the surface slip distribution along latitude (**d**) and the moment distributions along latitude (**c**) and along depth (**km**) (**e**). Continuous lines, the whole rupture; dashed lines, the southern part of the rupture in Palu Bay and south of it; crosses, the northern part of the rupture in Sulawesi neck and Balaesang Peninsula.</figcaption></figure>

image correlation) and at distance from the fault (InSAR) allowed us to constrain well the slip at both deep and shallow depths, which provided valuable resolution throughout the entire upper crust³³⁻³⁵ (Supplementary Fig. 4). The Palu earthquake therefore represents

one of the first instances in which both the shallow and deeper slip can be very well resolved for a mature continental plate-boundary strike-slip fault. The obtained geodetic magnitude is $M_w$ 7.62, with a maximum strike-slip of 5.25 m (obtained at shallow depths
below Palu basin, segment B in Fig. 3). Both reverse and normal slip are required on the fault plane to explain the measured InSAR field. The maximum normal dip-slip reached ~3 m (in the releasing bend area, segment D in Fig. 3), and the maximum reverse slip is ~2 m.

The surface slip distribution obtained from the model nicely complements the detailed one obtained directly from the offsets south of Palu city. A striking feature is the complementarity between the strike-slip and dip-slip at the surface. The dip-slip component is generally very small at the surface, except in three areas, where the slip is oblique with comparable strike-slip and dip-slip components (Fig. 4): (1) in the releasing bend area at 1.2° S, (2) in the Palu bay at 0.8° S where transpression is observed and (3) in the releasing bend at 0.2° S. When integrated over the entire seismogenic depth interval, the moment release in dip-slip evolves smoothly along the fault, except for a peak in the releasing bend at 1.2° S.

The moment distribution with depth indicates that most of the moment is released in the upper crust, between 0 and 15 km depth, the maximum being released between 0 and 10 km. This is compatible with the 12 km depth of interseismic locking previously estimated on this fault.⁷ It is noticeable that this depth is of the same order as the segment lengths deduced from rupture mapping from offset fields (typically less than 15 km (Fig. 3)), consistent with observations from other earthquakes³⁶.

If we consider that the earthquake actually broke two different faults (the mature Palu-Koro fault segment and a less mature unknown fault to the northeast), we see that the mature Palu-Koro fault experienced the peak moment release in the shallowest crust, but the northeast fault slipped at greater depth (down to 14–20 km), with a significant dip-slip component and a moderate deficit of slip near the surface (Fig. 4).

These results on the mature Palu-Koro fault contrast with previous findings on less structurally mature strike-slip faults, for which a systematic shallow slip deficit (SSD) has been observed in the upper 3 km (ref. ³⁷). The Palu earthquake represents a rare opportunity to examine the slip distribution with depth on a mature plate boundary fault, for which we also have high resolution in the shallowest part of our model. The distinct lack of SSD on the Palu fault segment is significant, and has important implications for both the seismic hazard of the region, as well as for our understanding of how the seismic cycle is accommodated on major tectonic structures globally. Structural maturity probably plays a key role in how deformation is distributed both on and off the fault plane, which in turn takes energy out of the system and eventually leads to rupture termination. Mature strike-slip faults are known to simplify geometrically (that is, smooth) with time³⁸,³⁹, and eventually result in smooth fault surfaces with few structural complexities and with well-developed gouge zones. Such faults are expected to produce similarly smooth slip distributions in earthquakes, both along-strike and up-dip; indeed, this is a requirement for the fault to rupture at supershear velocities²⁷. Therefore, SSD is not necessarily expected for earthquakes that occur on mature strike-slip faults⁴⁰.

## Tsunami

The tsunami generated by this earthquake took the population by surprise and is puzzling for the scientific community given the strike-slip mechanism of the earthquake, which produces limited uplift of the seafloor. Nevertheless, an earthquake of comparable magnitude ($M_w$ 7.4) on the Palu-Koro fault in 1968 also created a tsunami with wave heights up to 10 m (ref. ²). Strike-slip ruptures may generate some vertical land motion if the rake is not purely horizontal, notably in fault jog areas, or if they laterally displace submarine slopes and thereby displace the water column¹¹. Although our inversion results are less well constrained in the offshore region, the section north of Palu hosts a significant dip-slip displacement, notably in the Palu bay between 0.7° S and 0.8° S, and close to the epicentre

between 0.1° S and 0.2° S. In addition, landslides triggered by strong ground motions probably also contributed to the increase the local tsunami run-up⁴². In the case of Palu City, the shape of the narrow bay coupled with the rapid shallowing near the coastline most probably contributed to amplify the size of the tsunami. Initial numeric simulations of the tsunami waves indicate that the wave heights increase towards the south, when they reach the termination of the bay⁴³. However, to elucidate what were the critical factors that resulted in the Palu tsunami, further data acquisition is needed, notably offshore surveys, as well as modelling, which includes the effect of local bathymetry.

## Seismic cycle balance

Similar events are likely to occur on a regular basis on the Palu fault. Given a fully locked fault, a slip deficit of 2 m has accumulated since the previous tsunamigenic event in 1968, which corresponds to a moment deficit of $2.4 \times 10^{20}$ N m (~$M_w$ 7.51), assuming a down-dip width of 12 km and an along-strike length of 150 km. Simple seismic cycle balance calculations indicate a potential for a $M_w$ 7.5 earthquake about every 50 years on this fault, or a $M_w$ 7.25 event every 20 years. In 1996, a $M_w$ 7.9 event on the Minahassa trench triggered a $M_w$ 6.6 earthquake on the Palu fault two years later¹⁰. Interactions with nearby structures, such as the Minahassa trench to the north, or the Matano strike-slip fault to the south must therefore be taken into account in multisegment rupture scenarios.

## Conclusion

The rupture and slip measurements of the $M_w$ 7.5 Palu earthquake display the typical characteristics of supershear earthquakes. Like the 1999 Izmit (Turkey), the 2001 Kunlun (Tibet) and the 2002 Denali (Alaska) earthquakes, a long segment of the fault exhibits remarkably simple faulting characteristics: a sharp, straight, pure mode II rupture with a large smoothly varying slip, not directly associated with the morphology or with significant aftershocks²²,²⁴,²⁵,²⁷,²⁸,⁴⁴⁻⁴⁷. This finding is consistent with the short duration of the moment release evidenced by the source time function inferred from teleseismic records, which does not exceed 25–30 seconds¹⁷,¹⁸. The simple ratio between the distance from the epicentre to the southern termination of the slip inverted here (130 km) leads to an average rupture velocity of 4.3–5.2 km s⁻¹. The distribution of slip for this supershear segment shows no surface slip deficit, and contrasts with other known strike-slip events on relatively immature fault systems. The lack of observed SSD is probably due to a lack of geometric complexity on the fault plane, which is expected to be much smoother for structurally mature plate boundary faults.

Like many earthquakes, notably Kunlun and Denali²⁵,⁴⁵,⁴⁷,⁴⁸, the Palu earthquake nucleated in a complex area. The rupture propagated bilaterally at first. To the north, the rupture went in the direction of the Minahassa through a large releasing bend, a structurally complex area that contributed to arrest the rupture. Geometric complexities are known to play an important role in the initiation and termination of ruptures (especially supershear events) because they do not permit significant stress release⁴⁹. The rupture, therefore, mostly propagated southward. Its velocity probably increased dramatically to supershear speed when it reached the straight segment in the Palu basin. The rupture propagation at a shallow depth at supershear velocity within weakly consolidated alluvial sediments in the basin resulted in widespread liquefaction and severe damage from shock-wave-induced ground shaking. Finally, the rupture terminated after passing the jog at the southern end of the pull apart, again in a structurally complex area.

## Online content

Any methods, additional references, Nature Research reporting summaries, source data, statements of data availability and
associated accession codes are available at https://doi.org/10.1038/
s41561-018-0296-0.

Received: 25 October 2018; Accepted: 21 December 2018;
Published online: 4 February 2019

References

1. Cipta, A. et al. in Geohazards in Indonesia: Earth Science for Disaster Risk Reduction (eds Cummins, P. R. & Meilano, I.) 133–152 (Geological Society Special Publications Vol. 441, Geological Society, London, 2017).

2. Prasetya, G. S., de Lange, W. P. & Healy, T. R. The Makassar Strait tsunamigenic region, Indonesia. *Nat. Hazards* **24**, 295–307 (2001).

3. Puntodewo, S. S. O. et al. GPS measurements of crustal deformation within the Pacific–Australia Plate boundary zone in Irian Jaya, Indonesia. *Tectonophysics* **237**, 141–153 (1994).

4. Rangin, C. et al. Plate convergence measured by GPS across the Sundaland/Philippine Sea Plate deformed boundary; the Philippines and eastern Indonesia. *Geophys. J. Int.* **139**, 296–316 (1999).

5. Hall, R. Cenozoic geological and plate tectonic evolution of SE Asia and the SW Pacific: computer-based reconstructions, model and animations. *J. Asian Earth Sci.* **20**, 353–431 (2002).

6. Wallace, L. M., McCaffrey, R., Beavan, J. & Ellis, S. Rapid microplate rotations and backarc rifting at the transition between collision and subduction. *Geology* **33**, 857–860 (2005).

7. Socquet, A. et al. Microblock rotations and fault coupling in SE Asia triple junction (Sulawesi, Indonesia) from GPS and earthquake slip vector data. *J. Geophys. Res.* **111**, B08409 (2006).

8. Bellier, O. et al. High slip rate for a low seismicity along the Palu-Koro active fault in central Sulawesi (Indonesia). *Terra Nova* **13**, 463–470 (2001).

9. Watkinson, I. M. & Hall, R. in Geohazards in Indonesia: Earth Science for Disaster Risk Reduction (eds Cummins, P. R. & Meilano, I.) 71–120 (Geological Society Special Publications, Vol. 441, Geological Society, London, 2017).

10. Vigny, C. et al. Migration of seismicity and earthquake interactions monitored by GPS in SE Asia triple junction: Sulawesi, Indonesia. *J. Geophys. Res.* **107**, 2231 (2002).

11. Simons, W. et al. A decade of GPS in Southeast Asia: resolving Sundaland motion and boundaries. *J. Geophys. Res.* **112**, B06420 (2007).

12. Walpersdorf, A., Vigny, C., Manurung, P., Subarya, C. & Sutisna, S. Determining the Sula Block kinematics in the triple junction area in Indonesia by GPS. *Geophys. J. Int.* **135**, 351–361 (1998).

13. Silver, E. A., McCaffrey, R. & Smith, R. B. Collision, rotation, and the initiation of subduction in the evolution of Sulawesi, Indonesia. *J. Geophys. Res. Solid Earth* **88**, 9407–9418 (1983).

14. Bellier O. et al. in GEODYSSEA Final Report GFZ Scientific Technical Report 98/14 (eds Wilson, P. & Michel, G. W.) 276–312 (GeoForshingsZentrum, 1988).

15. Walpersdorf, A., Vigny, C., Subarya, C. & Manurung, P. Monitoring of the Palu-Koro Fault (Sulawesi) by GPS. *Geophys. Res. Lett.* **25**, 2313–2316 (1998).

16. Stevens, C. et al. Rapid rotations about a vertical axis in a collisional setting revealed by the Palu Fault, Sulawesi, Indonesia. *Geophys. Res. Lett.* **26**, 2677–2680 (1999).

17. $M_w$ 7.5 Palu earthquake, Indonesia. *US Geological Survey* https://earthquake.usgs.gov/earthquakes/eventpage/us1000h3p4/executive#executive (2018).

18. $M_w$ 7.5 earthquake, Sulawesi 2018/09/28 10:02:43 UTC. *Geoscope Observatory* http://geoscope.igpp.fr/index.php/en/catalog/earthquake-description?seis=us1000h3p4 (2018).

19. Bao, H. et al. Early and persistent supershear rupture of the 2018 magnitude 7.5 Palu earthquake. *Nat. Geosci.* https://doi.org/10.1038/s41561-018-0297-z (2019).

20. Bouchon, M. & Karabulut, H. The aftershock signature of supershear earthquakes. *Science* **320**, 1323–1325 (2008).

21. Leprince, S., Barbot, S., Ayoub, F. & Avouac, J. Automatic and precise orthorectification, coregistration, and subpixel correlation of satellite images, application to ground deformation measurements. *IEEE Trans. Geosci. Remote Sens.* **45**, 1529–1558 (2007).

22. Emre, Ö., Awata, Y., & Duman, T. Y. (eds) Surface Rupture Associated with the 17 August 1999 Izmit Earthquake (General Directorate of Mineral Research and Exploration, Ankara, 2003).

23. Barka, A. et al. The surface rupture and slip distribution of the 17 August 1999 Izmit earthquake (M 7.4), North Anatolian fault. *Bull. Seismol. Soc. Am.* **92**, 43–60 (2002).

24. Michel, R. & Avouac, J. P. Deformation due to the 17 August 1999 Izmit, Turkey, earthquake measured from SPOT images. *J. Geophys. Res. Solid Earth* **107**, 2062 (2002).

25. Klinger, Y. et al. High-resolution satellite imagery mapping of the surface rupture and slip distribution of the $M_w$ 7.8, 14 November 2001 Kokoxili earthquake, Kunlun fault, northern Tibet, China. *Bull. Seismol. Soc. Am.* **95**, 1970–1987 (2005).

26. King, G., Klinger, Y., Bowman, D. & Tapponnier, P. Slip-partitioned surface breaks for the $M_w$ 7.8 2001 Kokoxili earthquake, China. *Bull. Seismol. Soc. Am.* **95**, 731–738 (2005).

27. Bouchon, M. et al. Faulting characteristics of supershear earthquakes. *Tectonophysics* **493**, 244–253 (2010).

28. Robinson, D. P., Das, S. & Searle, M. P. Earthquake fault superhighways. *Tectonophysics* **493**, 236–243 (2010).

29. Andrews, D. J. & Ben-Zion, Y. Wrinkle-like slip pulse on a fault between different materials. *J. Geophys. Res.* **102**, 553–571 (1997).

30. Milliner, C. W. et al. Quantifying near-field and off-fault deformation patterns of the 1992 $M_w$ 7.3 Landers earthquake. *Geochem. Geophys. Geosyst.* **16**, 1577–1598 (2015).

31. Boore, D. M. A note on the effect of simple topography on seismic SH waves. *Bull. Seismol. Soc. Am.* **62**, 275–284 (1972).

32. Ma, S. A physical model for widespread near-surface and fault zone damage induced by earthquakes. *Geochem. Geophys. Geosyst.* **9**, Q11009 (2008).

33. Xu, X. et al. Refining the shallow slip deficit. *Geophys. J. Int.* **204**, 1867–1886 (2016).

34. Fialko, Y. Probing the mechanical properties of seismically active crust with space geodesy: study of the coseismic deformation due to the 1992 $M_w$ 7.3 Landers (southern California) earthquake. *J. Geophys. Res.* **109**, B03307 (2004).

35. Gombert, B. et al. Revisiting the 1992 Landers earthquake: a Bayesian exploration of co-seismic slip and off-fault damage. *Geophys. J. Int.* **212**, 839–852 (2018).

36. Klinger, Y. Relation between continental strike-slip earthquake segmentation and thickness of the crust. *J. Geophys. Res.* **115**, B07306 (2010).

37. Fialko, Y., Sandwell, D., Simons, M. & Rosen, P. Three-dimensional deformation caused by the Bam, Iran, earthquake and the origin of shallow slip deficit. *Nature* **435**, 295–299 (2005).

38. Stirling, M. W., Wesnousky, S. G. & Shimazaki, K. Fault trace complexity, cumulative slip, and the shape of the magnitude-frequency distribution for strike-slip faults: a global survey. *Geophys. J. Int.* **124**, 833–868 (1996).

39. Brodsky, E., Gilchrist, J., Sagy, A. & Collettini, C. Faults smooth gradually as a function of slip. *Earth Planet. Sci. Lett.* **302**, 185–193 (2011).

40. Dolan, J. F. & Haravitch, B. D. How well do surface slip measurements track slip at depth in large strike-slip earthquakes? The importance of fault structural maturity in controlling on-fault slip versus off-fault surface deformation. *Earth Planet. Sci. Lett.* **388**, 38–47 (2014).

41. Hooper, A. et al. Importance of horizontal seafloor motion on tsunami height for the 2011 $M_w$=9.0 Tohoku-Oki earthquake. *Earth Planet. Sci. Lett.* **361**, 469–479 (2013).

42. Carvajal, M. et al. The extraordinary tsunami of Palu 2018 and its implications for Chile. In XV Congreso Geologico Chileno (2018).

43. Schäfer, A. M. & Wenzel, F. Tsupy: computational robustness in tsunami hazard modelling. *Comp. Geosci.* **102**, 148–157 (2017).

44. Bouchon, M. & Vallée, M. Observation of long supershear rupture during the magnitude 8.1 Kunlunshan earthquake. *Science* **301**, 824–826 (2003).

45. Dunham, E. M. & Archuleta, R. J. Evidence for a supershear transient during the 2002 Denali fault earthquake. *Bull. Seismol. Soc. Am.* **94**, S256–S268 (2004).

46. Konca, A. O., Leprince, S., Avouac, J. P. & Helmberger, D. V. Rupture process of the 1999 $M_w$ 7.1 Duzce earthquake from joint analysis of SPOT, GPS, InSAR, strong-motion, and teleseismic data: a supershear rupture with variable rupture velocity. *Bull. Seismol. Soc. Am.* **100**, 267–288 (2010).

47. Lassere, C. et al. Coseismic deformation of the 2001 $M_w$=7.8 Kokoxili earthquake in Tibet, measured by synthetic aperture radar interferometry. *J. Geophys. Res. Solid Earth* **110**, B12408 (2005).

48. Ozacar, A. A. & Beck, S. L. The 2002 Denali fault and 2001 Kunlun fault earthquakes: complex rupture processes of two large strike-slip events. *Bull. Seismol. Soc. Am.* **94**, S278–S292 (2004).

49. King, G. & Nábělek, J. Role of fault bends in the initiation and termination of earthquake rupture. *Science* **228**, 984–987 (1985).

Acknowledgements

Our first thoughts go to the population of Sulawesi and to the victims of this earthquake. Satellite optical imagery has been made freely available by ESA (Sentinel-2, https://schub.copernicus.eu/dhus/#/home), USGS (Landsat-8, https://earthexplorer.usgs.gov/) and DigitalGlobe (WorldView, www.digitalglobe.com/opendata/indonesia-earthquake-tsunami/). We thank JAXA for providing ALOS-2 data under PI investigation no. 3328 (ALOS 6th RA) and the AW3D30 digital elevation model. We also thank the contributors to GMT and GMTSAR opensource software.
**Author contributions**

J.H. processed and analysed the image correlation data. E.P. processed and analysed ALOS-2 interferograms. A.S. coordinated the study, analysed the displacement maps and inverted for coseismic slip distributions. M.B. cross-examined the observations and results. All the authors discussed the content of the paper and shared the writing.

**Competing interests**

The authors declare no competing interests.

**Additional information**

**Supplementary information** is available for this paper at https://doi.org/10.1038/s41561-018-0296-0.

**Reprints and permissions information** is available at www.nature.com/reprints.

**Correspondence and requests for materials** should be addressed to A.S.

**Publisher's note:** Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

© The Author(s), under exclusive licence to Springer Nature Limited 2019
## Methods

**InSAR processing.** We used data from the ALOS-2 satellite, operated by the Japan Aerospace Exploration Agency (JAXA), acquired by the PALSAR2 instrument in the L band (wavelength = 24 cm, right-side looking) in descending observation mode WD1 (ScanSAR, Beam W2, Path 26, Frame 3650 with the AUIG frame shift = -5). Only the first scan (over five) is used in this study as the other scans are offshore to the west. The ALOS interferogram was processed using GMTSAR (open-source InSAR processing code<sup>[20,21]</sup>) using images acquired on the 21 August 2018 and on the 2 October 2018, separated by a perpendicular baseline of 60 m, and the ALOS 30 m resolution digital elevation model from the Global Digital Surface Model<sup>[52]</sup> (AW3D30) computed from JAXA stereo images of the optical PRISM instrument (www.eorc.jaxa.jp/ALOS/en/aw3d30/data/index.htm). Unwrapping was performed with the SNAPHU Algorithm for Phase Unwrapping<sup>[23]</sup>, and carefully controlled by visual checks. Decorrelating areas, as well as zones in which fringes were too narrow or discontinuous, were masked to prevent unwrapping errors.

**Optical images processing.** We correlated Sentinel-2 and Landsat-8 optical satellite images that covered the Palu region acquired before (S2, 17 September 2018; L8, 16 September 2018) and after (S2 and L8, 2 October 2018) the earthquake. Due to the wider spatial coverage and higher resolution (10 m) available in four different wavelengths (bands 2, 3, 4 and 8), we focused our optical analysis on Sentinel-2 imagery, which were mosaicked using the Geospatial Data Abstraction Library<sup>[54]</sup>, and then correlated using the COSI-Corr software package<sup>[21,55]</sup>, which features a phase correlator with subpixel accuracy. We used a sliding multiscale window of 64 × 32 pixels, a step size of four pixels and four robustness iterations (which optimizes the masking of noisy frequencies). We then median stack the four correlations (using an Ames Stereo Pipeline) to help reduce noise<sup>[56]</sup>. Similar parameters were used for Landsat-8 imagery, albeit with just one band (B8, 15 m resolution). The Landsat-8 correlation allowed a consistency check with the Sentinel-2 data and yielded identical results.

Several postprocessing steps were undertaken to reduce the noise further: (1) destriping 1 deals with the misalignments between the CCD (charge-coupled device) arrays on the Sentinel-2 sensor that produce striping in the along-track direction of the correlation map, which is removed by subtracting the average along-track value from each column in the along-track geometry; (2) destriping 2 deals with the additional striping that occurs in the satellite across-track direction due attitude oscillations of the satellite ('jitter'), which is removed in a similar manner to the CCD misalignments (we determined the striping pattern using stable areas of the correlations); (3) for outliers 1, we identified and discarded outlier pixels if >50% of its neighbours (within a 32×32 window) deviated from the central pixel by 0.8 m (that is, the s.d. of stable regions in the correlation maps); (4) for ramps, we detrended the correlation maps by removing a linear ramp from the image, determined from stable areas far from the fault rupture and (5) for outliers 2, we removed any values that fell outside the probable range expected from the earthquake (estimated from the histogram of displacement values).

Displacements across the fault were measured directly from the correlation maps (E–W and N–S components) using fault-perpendicular profiles. We interactively fitted a linear trend to each side of the fault, and retrieved the on-fault displacement projected into the fault-parallel and fault-normal directions (COSI-Corr features a 'stacking' tool for this purpose).

**Fault-rupture mapping and analysis.** The rupture extends southwards from the Palu coastline for ~33 km, breaking three fault segments (A, B and C (Figs. 2 and 3)), which both display very linear (albeit slightly different) strikes as they cut across the sedimentary deposits of the valley floor (segments A and B), before reaching the range front and the geological trace of the Palu-Koro fault (segment C). At the termination of segment C, the fault turns sharply to the southeast, where it crosses the valley in a releasing bend (segment D), before continuing south (segments E and F) into the mountains along the same fault azimuth as A, B and C. Each segment displays a primary rupture core, which accommodates the vast majority of deformation (that is, strike-slip) over a narrow zone. Minor secondary faulting either side of the primary rupture appears to accommodate largely fault-normal (that is, dip-slip) motion (Fig. 2b,f) over a wider zone. Horizontal slip vectors measured across the primary rupture (blue arrows (in Fig. 2a) indicate small variations in the fault-normal slip component. However, slip vectors measured across the rupture (using profiles that extend 4 km from the rupture trace), with only the longer-wavelength (that is, deeper) displacement considered (red arrows in Fig. 2a), indicate a remarkably smooth displacement profile (Fig. 3a,b), which is almost pure strike-slip (that is, pure Mode II). Minor variations in the slip vector on the primary rupture at the surface may also reflect additional ground deformation effects, such as landslides and mudflows (Fig. 2f), which are documented from the northern end of the Palu fault segment (near Palu city, www.youtube.com/watch?v=2TvDBGkfBgs).

Partitioning of the slip onto off-fault structures is also highlighted by the difference between the primary slip (short wavelength), and the cumulative primary-plus-secondary slip (long wavelength) (Fig. 3). Segments A and B both show ~16% OFD, which is probably related to the weak water-saturated sediment through which the rupture passes.

In general, OFD is larger and more fluctuating when the rupture is in the sediments than when it follows the topography. This suggests that it may be due to the way the sediment cover is deforming: either the slip amplitude may decrease at a shallow depth or the sudden decomposition of the shallow sediments during the passage of the shock wave expected from a supershear rupture may produce it, such as the extended cracked surface observed in Kunlun, for example<sup>[25]</sup>

**Slip distribution inversion and resolution.** The surface deformation fields associated with the earthquake were modelled as a series of dislocations buried in an elastic half space<sup>[57]</sup>.

Surface displacement fields were subsampled to manage a reasonable amount of data in the inversion. A subsampling scheme that depends on the distance to the fault<sup>[58]</sup> was applied to the ALOS-2 interferogram and to the correlation offset fields from Sentinel-2 (Supplementary Fig. 2). For distances less than 5 km from the fault plane, one displacement value every 1 km was kept when the images were coherent. The subsampling was then increased to one point every 2 km for distances between 5 and 20 km from the fault plane, then one point every 4 km from 20 to 45 km from the fault plane and finally to one point every 8 km for distances farther than 45 km from the fault plane.

The fault geometry was constrained by the trace of the rupture at the surface, derived from offset fields. North of 0.2° S, the offsets failed to constrain the trace of the fault, but Insar fringes infer that the fault is running offshore to connect the Minahassa trench to the North. The offshore trace of the fault was constrained from multibeam bathymetry<sup>[59]</sup>. We assumed a uniform dip of 60° eastward, close to the W-phase moment tensor from the United States Geological Survey (USGS). This dip is also compatible with the morphological trace of the fault that generates a clear scarp and triangular facets compatible with a significant amount of normal faulting south of Palu city<sup>[3]</sup>. To the north, in the Palu bay and Sulawesi neck areas, the fault orientation is compatible with a restraining bend, which implies a transpressive movement that is compatible with the change in sign of the InSAR phase gradient east of the fault at those latitudes. Finally, at 0.2° S, the fault forms a major releasing bend when its trace rotates towards the northwest through the Balaesang peninsula.

The fault was discretized into an array of 42×7 elements, their size varying locally because the fault follows the rupture trace (Fig. 4 and Supplementary Fig. 4). To solve for the slip distribution along the 294 fault patches, we used a least squares minimization with a non-negativity constraint on the slip to stabilize the solution. In this model, the dip-slip sign is forced to positive values (that is, normal faulting) south of Palu (segments A to F) and north of 0.2° S (the sea and peninsula segments), whereas the dip-slip is negative (that is, reverse faulting) in the Palu bay and Sulawesi neck segments. The rake is allowed to vary. Alternative models with uniform dip-slip sign have been performed but fail to explain the InSAR data. Inversions without the non-negativity constraint have also been conducted and show compatible results, but with some instabilities in the slip distribution where it is less well constrained.

We estimated the sensitivity of our data set to unit displacements on each node of the grid by summing the horizontal deformation on the whole network following Loveless and Meade<sup>[60]</sup>. The power of data to constrain the fault interface behaviour depends on the size of the patches, on the number of data available and on the distance between the patches and the data points. The power of our data to constrain the slip on the fault interface is very high at the surface and decreases with depth. The strike-slip component is extremely well constrained by our data set, but the dip-slip component is slightly less well resolved. Expectedly, the resolution at the surface nodes is less good where the fault trace is offshore, notably in the Palu bay and the slight passage at sea south of the Peninsula. The northern extremity of the fault, where it runs offshore to connect the Minahassa trench, is not resolved.

We performed a joint inversion, which included the subsampled ALOS descending InSAR track and N–S and E–W surface offsets derived from subpixel correlation of Sentinel-2 images. For each offset field an ambiguity value was also estimated in the inversion scheme, whereas a ramp was estimated for the interferogram. Values of tilts and ambiguities estimated during inversion are indicated in Supplementary Fig. 3. We experimented with different weightings. We searched for a compromise between the model-data root mean square (r.m.s.) misfit and the spatial density of each type of data (Supplementary Fig. 5), and chose a weight of 50% for InSAR data and 50% for Sentinel offsets. To limit oscillations of the solution, we applied smoothing by minimizing the second-order derivative of the fault slip. Slip was forced to zero at the bottom and edges of the fault, but left free at the surface. We determined the optimal solution roughness<sup>[61]</sup> that was used in our final models, searching for a compromise between the roughness and misfit of the solution. The roughness of the preferred coseismic distribution is 1.9 cm km<sup>-1</sup> for a total r.m.s (L2-norm misfit) of 59 cm. Supplementary Fig. 3 shows the fit of the model to the data (r.m.s.<sub>inSAR</sub> = 8.8 cm, r.m.s.<sub>offset_NS</sub> = 37 cm and r.m.s.<sub>foeget_EW</sub> = 32 cm). The seismic moment was 3.4 × 10<sup>20</sup> N m, and corresponds to a magnitude of 7.618.

## Code availability

Synthetic aperture radar data were processed using GMTSAR software, freely available from https://topex.ucsd.edu/gmtsar/. Optical satellite images were
mosaicked using the Geospatial Data Abstraction Library (http://gdal.org) and then correlated using the COSI-Corr software package available at www.tectonics.caltech.edu/slip_history/spot_coseis/download_software.html. Deformation generated by a static earthquake source was modelled using Okada's equations implemented in the disloc program (www.physics.hmc.edu/GL/disloc/disloc.c).

## Data availability

The data sets generated during the current study (displacement fields from Landsat-8, Sentinel-2 and WorldView image correlation and from the ALOS-2 interferogram, as well as the static slip distribution) are available from the corresponding author upon request. Raw satellite optical imagery was made freely available by ESA (Sentinel-2, https://scihub.copernicus.eu/dhus/#/home), USGS (Landsat-8, https://earthexplorer.usgs.gov/) and DigitalGlobe (WorldView, www.digitalglobe.com/opendata/indonesia-earthquake-tsunami/). Raw ALOS-2 data availability is restricted to PI investigation at www.eorc.jaxa.jp/ALOS/en/aw3d30/data/index.htm.

## References

50. Sandwell, D., Mellors, R., Tong, X., Wei, M. & Wessel, P. Open radar interferometry software for mapping surface deformation. *Eos Trans. AGU 92*, 28 (2011).

51. Sandwell, D., Mellors, R., Tong, X., Wei, M., & Wessel, P. GMTSAR: An *InsAR Processing System Based on Generic Mapping Tools* Scripps Institution of Oceanography Technical Report (2011).

52. Tadono, T. et al. Precise global DEM generation by ALOS PRISM. *ISPRS Ann. Photogram Remote Sensing Spatial Inf. Sci. II-4*, 71–76 (2014).

53. Chen, C. W. & Zebker, H. A. Phase unwrapping for large SAR interferograms: statistical segmentation and generalized network models. *IEEE Trans. Geosci. Remote Sensing 40*, 1709–1719 (2002).

54. GDAL/OGR contributors. GDAL/OGR Geospatial Data Abstraction (Open Source Geospatial Foundation).

55. Ayoub, F., Leprince, L. & Avouac, J.-P. *User's Guide to COSI-CORR: Co-registration of Optically Sensed Images and Correlation* (California Institute of Technology, 2015).

56. Beyer, R. A., Alexandrov, O. & McMichael, S. The Ames Stereo Pipeline: NASA's open source software for deriving and processing terrain data. *Earth Space Sci.* 5, 537–548 (2018).

57. Okada, Y. Surface deformation due to shear and tensile faults in a half-space. *Bull. Seismol. Soc. Am.* 75, 1135–1154 (1985).

58. Grandin, R. et al. The September 2005 Wal’is-Dabbahu rifting event, Afar (Ethiopia): constraints provided by geodetic data. *J. Geophys. Res.* 114, B08404 (2009).

59. Tiranda, H. *Tectonostratigraphic Evolution of Offshore NW Sulawesi with Implications for the Hydrocarbon Prospectivity*. MSc thesis, Univ. London (2017).

60. Loveless, J. P. & Meade, B. J. Spatial correlation of interseismic coupling and coseismic rupture extent of the 2011 Mw=9.0 Tohoku Oki earthquake. *Geophys. Res. Lett.* 38, L17306 (2011).

61. Jonsson, S., Zebker, H., Segall, P. & Amelung, F. Fault slip distribution of the 1999 Mw 7.1 Hector Mine, California, earthquake, estimated from satellite radar and GPS measurements. *Bull. Seismol. Soc. Am.* 92, 1377–1389 (2002).