# Simulation of the Mediterranean tsunami generated by the Mw 6.0 event offshore Bejaia (Algeria) on 18 March 2021

P Heinrich, A Dupont, M Menager, A Trilla, A Gailler, B Delouis, H Hébert

► To cite this version:

P Heinrich, A Dupont, M Menager, A Trilla, A Gailler, et al. Simulation of the Mediterranean tsunami generated by the Mw 6.0 event offshore Bejaia (Algeria) on 18 March 2021. Geophysical Journal International, 2024, 237 (3), pp.1400-1413. 10.1093/gji/ggae121. hal-04595626

HAL Id: hal-04595626

https://u-paris.hal.science/hal-04595626

Submitted on 31 May 2024

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
# Simulation of the Mediterranean tsunami generated by the Mw 6.0 event offshore Bejaia (Algeria) on 18 March 2021

P Heinrich, A Dupont, M Menager, A Trilla, B. Delouis, B Delouis, H Hébert

► To cite this version:

P Heinrich, A Dupont, M Menager, A Trilla, B. Delouis, et al.. Simulation of the Mediterranean tsunami generated by the Mw 6.0 event offshore Bejaia (Algeria) on 18 March 2021. Geophysical Journal International, 2024, 237 (3), pp.1400-1413. 10.1093/gji/ggae121. hal-04595626

HAL Id: hal-04595626

https://u-paris.hal.science/hal-04595626v1

Submitted on 31 May 2024

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
# Simulation of the Mediterranean tsunami generated by the $M_w$ 6.0 event offshore Bejaia (Algeria) on 18 March 2021

P. Heinrich,<fnref n="1" /> A. Dupont,<fnref n="1" /> M. Menager,<fnref n="2" /> A. Trilla,<fnref n="1" /> A. Gailler,<fnref n="1" /> B. Delouis<fnref n="3" /> and H. Hébert<fnref n="1" />

Accepted 2024 March 25. Received 2024 March 21; in original form 2024 January 16

## SUMMARY

On 18 March 2021 an earthquake of magnitude $M_w$ = 6.0 occurred offshore the Algerian coasts and generated a tsunami with offshore amplitudes smaller than a few millimetres crossing the western Mediterranean Sea. The objective of this study is threefold: first, to determine whether seismic sources calculated in the context of tsunami early warning are relevant; secondly, to determine whether tsunami simulations are able to reproduce tide-gauge observations and thirdly, to define the sensitivity of simulations to the grid resolutions and tsunami parameters. In the Mediterranean Sea, a very small number of available coastal tide gauges recorded the tsunami. Among them, a few French tide gauge stations recorded water waves with amplitudes smaller than a few centimetres and with periods ranging from 5 to 20 min associated to harbour or bay resonances. Numerical simulations of the tsunami are performed by the operational code Taitoko for seven different source fault models. Three of them allow for a rapid source detection and characterization in the framework of tsunami warning at CENALT (Centre National d'Alerte aux Tsunamis, France). The integrated code Taitoko uses a system of multiple nested grids. Standard Boussinesq equations are solved in the Mediterranean grid, whereas non-linear shallow water equations are solved in coastal and harbour grids with 25 and 5 m resolutions, respectively. Whatever the fault model, the observed time-series of water heights are reproduced satisfactorily both in phase and amplitude by the model at Nice and Monaco but poorly at Port Mahon (Minorca) and Toulon.

**Key words:** Numerical modelling; Tsunamis; Earthquake early warning.

## 1 INTRODUCTION

On 18 March 2021, at 00:04:06 UTC, an offshore earthquake occurred at about 20 km away from the town of Bejaia in northeastern Algeria (Yelles-Chaouche et al. 2023). This earthquake ($M_w$ = 6.0) was felt at distances up to 900 km from the epicentral area along the Italian and French coasts (De Caro et al. 2023). The earthquake occurred on the Eurasian–African plate boundary, also called the Algerian margin, which is the location of an important seismic activity (Delouis et al. 2004; Braunmiller & Bernardi 2005; Stich et al. 2005, 2006; Buforn et al. 2015, 2017). The margin extends from the Strait of Gibraltar to Tunisia, passing through southern Spain, northern Morocco and northern Algeria (Fig. 1). It is characterized by a continuous occurrence of low to moderate magnitude earthquakes ($M$ < 5), and by numerous events with a magnitude $M$ > 5 (Buforn et al. 1995). More precisely, 59 seismic events with $M$ > 5 have been listed since 1960 in northern Algeria according to ISC (2024). The convergence between the African and Eurasian plates

causes a NNW–SSE regional compression and creates multiple reverse faults, which are completed by strike-slip faults (Buforn et al. 1995, 2004; Alasset et al. 2006; Stich et al. 2006). The convergence is estimated at about 6 mm yr⁻<fnref n="1" /> in northeastern Algeria, and at about 7.6 mm yr⁻<fnref n="1" /> in northwestern Algeria (Buforn et al. 2004).

While many large earthquakes occurring in northern Algeria are located onshore, some others are located offshore (Fig. 1). Historical tsunamis are reported within the region after seismic events in 1365, 1802, 1856, 1954, 1980, 2003 and 2016 (Roger & Hebert 2008; Sahal et al. 2009; Gailler et al. 2013; Masina et al. 2020). The regional tsunami occurrence is likely related to active submarine thrust-and-fold structures that can be connected to continental faults (Alasset et al. 2006). Along the margin of north Africa, the 2021 event is the second strongest seismic event recorded since the devastating earthquake ($M_w$ = 6.8) in Boumerdes in 2003. The 2003 event generated tsunami heights of several tens of centimetres along the coasts of the western Mediterranean Sea (Alasset et al. 2006; Sahal et al. 2009; Heidarzadeh & Satake 2013). Farther to

<footnote n="1">¹CEA, DAM, DIF, F-91297 Arpajon, France. E-mail: philippe.heinrich@cea.fr</footnote>

<footnote n="2">²Institut de Physique du Globe de Paris, Observatoire volcanologique du Piton de la Fournaise, 97418 La Plaine des Cafres, La Réunion</footnote>

<footnote n="3">³Université Côte d'Azur, CNRS, IRD, Observatoire de la Côte d'Azur, Géoazur, 250 rue Albert Einstein, 06560 Valbonne, France</footnote>
# Simulation of the Mediterranean tsunami generated by the $M_w$ 6.0 event offshore Bejaia (Algeria) on 18 March 2021

P. Heinrich,<fnref n="1" /> A. Dupont,<fnref n="1" /> M. Menager,<fnref n="2" /> A. Trilla,<fnref n="1" /> A. Gailler,<fnref n="1" /> B. Delouis<fnref n="3" /> and H. Hébert<fnref n="1" />

Accepted 2024 March 25. Received 2024 March 21; in original form 2024 January 16

## SUMMARY

On 18 March 2021 an earthquake of magnitude $M_w$ = 6.0 occurred offshore the Algerian coasts and generated a tsunami with offshore amplitudes smaller than a few millimetres crossing the western Mediterranean Sea. The objective of this study is threefold: first, to determine whether seismic sources calculated in the context of tsunami early warning are relevant; secondly, to determine whether tsunami simulations are able to reproduce tide-gauge observations and thirdly, to define the sensitivity of simulations to the grid resolutions and tsunami parameters. In the Mediterranean Sea, a very small number of available coastal tide gauges recorded the tsunami. Among them, a few French tide gauge stations recorded water waves with amplitudes smaller than a few centimetres and with periods ranging from 5 to 20 min associated to harbour or bay resonances. Numerical simulations of the tsunami are performed by the operational code Taitoko for seven different source fault models. Three of them allow for a rapid source detection and characterization in the framework of tsunami warning at CENALT (Centre National d'Alerte aux Tsunamis, France). The integrated code Taitoko uses a system of multiple nested grids. Standard Boussinesq equations are solved in the Mediterranean grid, whereas non-linear shallow water equations are solved in coastal and harbour grids with 25 and 5 m resolutions, respectively. Whatever the fault model, the observed time-series of water heights are reproduced satisfactorily both in phase and amplitude by the model at Nice and Monaco but poorly at Port Mahon (Minorca) and Toulon.

**Key words:** Numerical modelling; Tsunamis; Earthquake early warning.

## 1 INTRODUCTION

On 18 March 2021, at 00:04:06 UTC, an offshore earthquake occurred at about 20 km away from the town of Bejaia in northeastern Algeria (Yelles-Chaouche et al. 2023). This earthquake ($M_w$ = 6.0) was felt at distances up to 900 km from the epicentral area along the Italian and French coasts (De Caro et al. 2023). The earthquake occurred on the Eurasian–African plate boundary, also called the Algerian margin, which is the location of an important seismic activity (Delouis et al. 2004; Braunmiller & Bernardi 2005; Stich et al. 2005, 2006; Buforn et al. 2015, 2017). The margin extends from the Strait of Gibraltar to Tunisia, passing through southern Spain, northern Morocco and northern Algeria (Fig. 1). It is characterized by a continuous occurrence of low to moderate magnitude earthquakes ($M$ < 5), and by numerous events with a magnitude $M$ > 5 (Buforn et al. 1995). More precisely, 59 seismic events with $M$ > 5 have been listed since 1960 in northern Algeria according to ISC (2024). The convergence between the African and Eurasian plates

causes a NNW–SSE regional compression and creates multiple reverse faults, which are completed by strike-slip faults (Buforn et al. 1995, 2004; Alasset et al. 2006; Stich et al. 2006). The convergence is estimated at about 6 mm yr⁻<fnref n="1" /> in northeastern Algeria, and at about 7.6 mm yr⁻<fnref n="1" /> in northwestern Algeria (Buforn et al. 2004).

While many large earthquakes occurring in northern Algeria are located onshore, some others are located offshore (Fig. 1). Historical tsunamis are reported within the region after seismic events in 1365, 1802, 1856, 1954, 1980, 2003 and 2016 (Roger & Hebert 2008; Sahal et al. 2009; Gailler et al. 2013; Masina et al. 2020). The regional tsunami occurrence is likely related to active submarine thrust-and-fold structures that can be connected to continental faults (Alasset et al. 2006). Along the margin of north Africa, the 2021 event is the second strongest seismic event recorded since the devastating earthquake ($M_w$ = 6.8) in Boumerdes in 2003. The 2003 event generated tsunami heights of several tens of centimetres along the coasts of the western Mediterranean Sea (Alasset et al. 2006; Sahal et al. 2009; Heidarzadeh & Satake 2013). Farther to

<footnote n="1">¹CEA, DAM, DIF, F-91297 Arpajon, France. E-mail: philippe.heinrich@cea.fr</footnote>

<footnote n="2">²Institut de Physique du Globe de Paris, Observatoire volcanologique du Piton de la Fournaise, 97418 La Plaine des Cafres, La Réunion</footnote>

<footnote n="3">³Université Côte d'Azur, CNRS, IRD, Observatoire de la Côte d'Azur, Géoazur, 250 rue Albert Einstein, 06560 Valbonne, France</footnote>
<figure><img src="image_2.png" /><figcaption>**Figure 1.** Map of the seismic activity recorded between January 2012 and March 2022 with magnitudes larger than 2 (white dots). Location of stations used for the GRiD MT and FMNEAR inversions is given by yellow and violet triangles, respectively. Tide-gauge locations are indicated by blue squares. The red star is the epicentre of the Bejaia earthquake published by Yelles-Chaouche et al. (2023). On the bottom-right frame, published locations of the earthquake are presented: Yelles-Chaouche (Y-C, red), GCMT (orange), FMNEAR (violet), USGS (green) and GRiD MT (yellow). A summary of published sources is also available online: https://seismicportal.eu/eventdetails.html?unid=20210318.0000001.</figcaption></figure>

the west, the 2016 Al Hoceima earthquake ($M_w$ = 6.4) occurred off the Moroccan coasts (Stich et al. 2020) and generated a small tsunami (with amplitudes of a few centimetres) recorded by a few tide gauges of the Spanish coast. Despite its moderate moment magnitude ($M_w$ = 6.0), the 2021 earthquake triggered a small-amplitude tsunami travelling across the Mediterranean Sea toward Spain and France. This tsunami was recorded by a few tide gauges in the Balearic Islands, Sardinia, Spain and France.

Given its potential tsunamigenic hazards, the northern Algerian margin is a region of interest for the French Tsunami Warning Center (Gailler et al. 2013; Roudil et al. 2013; Schindelé et al. 2015), also known as ‘Centre National d’Alerte aux Tsunamis’ (CENALT) since 2012. Its main task is to issue tsunami bulletins and warnings within 15 min after an event of interest to the French authorities. As a Tsunami Service Provider (TSP) recognized by IOC/Unesco (as Portugal, Italy, Greece and Turkey) it also sends information to any other subscribing country or institute willing to receive alert messages. CENALT is mainly focused on early seismological characterization, followed by estimation of a possible tsunami trigger based on a decision matrix defined from past known events, taking into account magnitude, location and depth basically (Roudil et al. 2013; Schindelé et al. 2015). These parameters have to be obtained very quickly after the main shock, and be the most reliable. Beyond the message issued within 15 min after the main shock, the tsunami simulation is essential to estimate possible coastal effects, but the results are not considered in the evacuation procedures that

rely on rapid response only so far. Usually focal mechanisms are not available within the first minutes, therefore the first simulations rely either on pre-defined fault rupture patterns, or they are carried out once a focal mechanism is available, generally after 20–30 min (Gailler et al. 2013). However tsunami warning procedures can be expected to incorporate more and more simulation results in the future, and CENALT is also exploring such capacities, in addition to the offline assessment of tsunami coastal impact. In order to validate modelling methods and seismological analyses, the comparison on actual events is therefore essential. Since the birth of CENALT in 2012, the 2021 Bejaia event triggered the largest tsunami of tectonic origin recorded along the French coasts, and it is thus a significant event for CENALT. We present here the results of two methods, GRiD MT and FMNEAR, that are able to rapidly detect or characterize an event (and to provide complete information to CENALT) if any. The GRiD MT (Grid-based Real-time Determination of Moment Tensors) method is based on a seismic deviatoric moment tensor inversion scheme performing a continuous long-period seismic waveform scanning over a region of interest (Kawakatsu 1998; Guilhem & Dreger 2011; Menager et al. 2023). The FMNEAR method, developed by the OCA (Observatoire de la Côte d’Azur) institute, is based on both the waveform inversion of near-source seismic records and on a linear finite-source model (Delouis 2014). The primary source parameters that are determined are the moment magnitude ($M_w$), the double-couple focal mechanism (strike, dip and rake), the source depth and the distribution of seismic moment
<figure><img src="image_2.png" /><figcaption>**Figure 1.** Map of the seismic activity recorded between January 2012 and March 2022 with magnitudes larger than 2 (white dots). Location of stations used for the GRiD MT and FMNEAR inversions is given by yellow and violet triangles, respectively. Tide-gauge locations are indicated by blue squares. The red star is the epicentre of the Bejaia earthquake published by Yelles-Chaouche et al. (2023). On the bottom-right frame, published locations of the earthquake are presented: Yelles-Chaouche (Y-C, red), GCMT (orange), FMNEAR (violet), USGS (green) and GRiD MT (yellow). A summary of published sources is also available online: https://seismicportal.eu/eventdetails.html?unid=20210318.0000001.</figcaption></figure>

the west, the 2016 Al Hoceima earthquake ($M_w$ = 6.4) occurred off the Moroccan coasts (Stich et al. 2020) and generated a small tsunami (with amplitudes of a few centimetres) recorded by a few tide gauges of the Spanish coast. Despite its moderate moment magnitude ($M_w$ = 6.0), the 2021 earthquake triggered a small-amplitude tsunami travelling across the Mediterranean Sea toward Spain and France. This tsunami was recorded by a few tide gauges in the Balearic Islands, Sardinia, Spain and France.

Given its potential tsunamigenic hazards, the northern Algerian margin is a region of interest for the French Tsunami Warning Center (Gailler et al. 2013; Roudil et al. 2013; Schindelé et al. 2015), also known as ‘Centre National d’Alerte aux Tsunamis’ (CENALT) since 2012. Its main task is to issue tsunami bulletins and warnings within 15 min after an event of interest to the French authorities. As a Tsunami Service Provider (TSP) recognized by IOC/Unesco (as Portugal, Italy, Greece and Turkey) it also sends information to any other subscribing country or institute willing to receive alert messages. CENALT is mainly focused on early seismological characterization, followed by estimation of a possible tsunami trigger based on a decision matrix defined from past known events, taking into account magnitude, location and depth basically (Roudil et al. 2013; Schindelé et al. 2015). These parameters have to be obtained very quickly after the main shock, and be the most reliable. Beyond the message issued within 15 min after the main shock, the tsunami simulation is essential to estimate possible coastal effects, but the results are not considered in the evacuation procedures that

rely on rapid response only so far. Usually focal mechanisms are not available within the first minutes, therefore the first simulations rely either on pre-defined fault rupture patterns, or they are carried out once a focal mechanism is available, generally after 20–30 min (Gailler et al. 2013). However tsunami warning procedures can be expected to incorporate more and more simulation results in the future, and CENALT is also exploring such capacities, in addition to the offline assessment of tsunami coastal impact. In order to validate modelling methods and seismological analyses, the comparison on actual events is therefore essential. Since the birth of CENALT in 2012, the 2021 Bejaia event triggered the largest tsunami of tectonic origin recorded along the French coasts, and it is thus a significant event for CENALT. We present here the results of two methods, GRiD MT and FMNEAR, that are able to rapidly detect or characterize an event (and to provide complete information to CENALT) if any. The GRiD MT (Grid-based Real-time Determination of Moment Tensors) method is based on a seismic deviatoric moment tensor inversion scheme performing a continuous long-period seismic waveform scanning over a region of interest (Kawakatsu 1998; Guilhem & Dreger 2011; Menager et al. 2023). The FMNEAR method, developed by the OCA (Observatoire de la Côte d’Azur) institute, is based on both the waveform inversion of near-source seismic records and on a linear finite-source model (Delouis 2014). The primary source parameters that are determined are the moment magnitude ($M_w$), the double-couple focal mechanism (strike, dip and rake), the source depth and the distribution of seismic moment
along strike, which can provide a first-order estimate of rupture length and directivity.

This paper presents results from a short seismological analysis of the 2021 Bejaia main shock and it then focuses on tsunami simulations. Its first part describes the earthquake sources available in the literature and more particularly the ones calculated by the GRiD MT and FMNEAR methods. Both methods are under testing and implementation at CENALT. Tsunami simulations are performed for seven seismic sources and results are compared to observed tide-gauge records. Sensitivity studies to grid resolutions and some relevant tsunami parameters are also analysed. The objective of this study is threefold: first, to determine whether seismic sources calculated in the context of tsunami early warning are relevant; secondly, to determine whether tsunami simulations are able to reproduce tide-gauge observations and thirdly, to define the sensitivity of simulations to the grid resolutions and tsunami parameters.

## 2 FOCAL MECHANISMS

This section presents published focal mechanisms and three methods implemented or under implementation at CENALT: the ‘Polarities’, GRiD MT, OCA/FMNEAR methods. Only the first two are described in the section.

### 2.1. Published focal mechanisms solutions

In their paper, Yelles-Chaouche et al. (2023) summarizes epicentres and focal mechanisms of eight seismic institutes for the 2021 event. The early CMT focal solutions calculated by these institutes show dominant dip-slip mechanisms with an assumed homogeneous slip distribution. According to Yelles-Chaouche et al. (2023), the Bejaia main shock and its aftershocks sequence in 2021 and 2022 constitute an episode of a long-lasting crisis, which started in the Bejaia region since 2012. The calculation by Yelles-Chaouche et al. (2023) of the 2021 source parameters is based on an extensive network of Algerian stations and their solution is considered in this paper as the reference solution. This latter is characterized as a reverse faulting earthquake on an almost east-west offshore fault (Fig. 2, Table 1).

Six other focal mechanisms are also analysed in this study: the solutions published the day of the earthquake, USGS, GCMT, OCA/FMNEARb and those not published or under test at CENALT, the ‘Polarities’ method (see the Annex), OCA/FMNEARa and GRiD MT. Parameters of the seven tested seismic sources are given in Table 1. The GCMT and Yelles-Chaouche sources are located offshore within the Bejaia Bay (Fig. 1) whereas the other epicentres are located further north, off the Bejaia Bay.

The USGS, OCA/FMNEARb and GCMT solutions were synthesized by EMSC (European Mediterranean Seismological Centre) within a few hours after the origin time (2h40min, 6h39min and 12h57min, respectively). The ‘Polarities’ method, GRiD MT and OCA/FMNEARa solutions would be potentially available earlier at CENALT in the context of early tsunami warning (12, 15 and 52 min, respectively).

### 2.2. GRiD MT approach for an automated source characterization

The GRiD MT solution can be obtained within the 15 min timeframe required for the CENALT’s mission of first alert. Indeed, it can be used to determine the event’s origin time, location (latitude,

longitude and depth), moment magnitude as well as its focal mechanism. The approach makes use of automated time-domain moment tensor inversions performed at different potential source hypocentres over a 3-D grid. Several seismological centres around the world use it in their seismic monitoring algorithms (Tsuruoka et al. 2009; Guilhem & Dreger 2011; Guilhem et al. 2013; MacPherson et al. 2013; Lee et al. 2013). Guilhem et al. (2013) also showed its interest for the detection of large and potentially tsunamigenic earthquakes, and in particular for the 2011 M<sub>w</sub> 9.1 Tohoku-oki earthquake in Japan.

The approach GRiD MT has recently been tested for earthquakes occurring along the Algerian margin (Menager et al. 2023). It works with a limited number of five broad-band stations that are selected over a pool of regional stations, called CENALT’s backbone network, providing reliable realtime data to CENALT in the western Mediterranean region: IM.MDT, PM.PESTR, IU.PAB, RD.MTT2 and RD.ORT2 (here, network code followed by station name). The three-component seismic waveforms from these five selected stations are deconvolved with the instrument response, integrated to displacement, and filtered between periods of 25 and 50 s. The grid covering the Algerian margin ranges from 34.50° to 38.00°N in latitude, and from -5.75° to 7°E in longitude (Figs 1 and 2). Because CENALT is interested in potentially tsunamigenic earthquakes, only four source depths are considered: 10, 20, 30 and 40 km.

The GRiD MT solution is defined by calculating the variance reduction value (VR) (eq. 1). It assesses the fit between the observed and the predicted data during the inversion. The higher the VR value is, the better is the fit, and so is the solution’s confidence.

$$VR = \left[ 1 - \frac{\sum_{i} (\text{observed}_i - \text{predicted})^2}{\sum_{i} \text{observed}_i^2} \right] \times 100. \quad (1)$$

Fig. 2 shows the results of the rapid GRiD MT approach defined for the region of interest between 100 s before and after the event origin time (OT). The maximum peak of VR (68.4 per cent) obtained for the Bejaia earthquake is detected 5 s after the OT, at the grid point having the coordinates 37.00°N/5.25°E (Fig. 2). This optimal position is in good agreement with locations of other sources (Fig. 1 and Table 1) and is quickly calculated (within 15 min). The earthquake focal mechanism is reverse, West-East oriented, with a high double-couple component, and is very similar to published solutions (Table 1 and Fig. 2).

The best GRiD MT solution (i.e. the one having the maximal VR value) is found around 10–20 km depth, with a magnitude M<sub>w</sub> 5.9 (Fig. 2). It is worth noting that the depth is difficult to constrain at the frequency band used in this approach.

### 2.3. The FMNEAR method

An instance of FMNEAR (Delouis 2014) is running at the Observatoire de la Côte d'Azur (OCA), automated for the near-real time but also with the possibility to be relaunched manually with a modified parametrization. A first automatic solution was obtained at T0 + 52 min the day of the earthquake, but with quality criteria just below the thresholds defined for transmission to the EMSC (cf. OCA/FMNEARa solution, Table 1, Fig. 3). Obtained with five stations only, this solution is nonetheless meaningful and would have been of interest for the purpose of tsunami warning. A solution with a higher degree of confidence was obtained at T0 + 6h39, by allowing the use of stations at larger distances in order to compensate the absence of available records in the south. This later
along strike, which can provide a first-order estimate of rupture length and directivity.

This paper presents results from a short seismological analysis of the 2021 Bejaia main shock and it then focuses on tsunami simulations. Its first part describes the earthquake sources available in the literature and more particularly the ones calculated by the GRiD MT and FMNEAR methods. Both methods are under testing and implementation at CENALT. Tsunami simulations are performed for seven seismic sources and results are compared to observed tide-gauge records. Sensitivity studies to grid resolutions and some relevant tsunami parameters are also analysed. The objective of this study is threefold: first, to determine whether seismic sources calculated in the context of tsunami early warning are relevant; secondly, to determine whether tsunami simulations are able to reproduce tide-gauge observations and thirdly, to define the sensitivity of simulations to the grid resolutions and tsunami parameters.

## 2 FOCAL MECHANISMS

This section presents published focal mechanisms and three methods implemented or under implementation at CENALT: the ‘Polarities’, GRiD MT, OCA/FMNEAR methods. Only the first two are described in the section.

### 2.1. Published focal mechanisms solutions

In their paper, Yelles-Chaouche et al. (2023) summarizes epicentres and focal mechanisms of eight seismic institutes for the 2021 event. The early CMT focal solutions calculated by these institutes show dominant dip-slip mechanisms with an assumed homogeneous slip distribution. According to Yelles-Chaouche et al. (2023), the Bejaia main shock and its aftershocks sequence in 2021 and 2022 constitute an episode of a long-lasting crisis, which started in the Bejaia region since 2012. The calculation by Yelles-Chaouche et al. (2023) of the 2021 source parameters is based on an extensive network of Algerian stations and their solution is considered in this paper as the reference solution. This latter is characterized as a reverse faulting earthquake on an almost east-west offshore fault (Fig. 2, Table 1).

Six other focal mechanisms are also analysed in this study: the solutions published the day of the earthquake, USGS, GCMT, OCA/FMNEARb and those not published or under test at CENALT, the ‘Polarities’ method (see the Annex), OCA/FMNEARa and GRiD MT. Parameters of the seven tested seismic sources are given in Table 1. The GCMT and Yelles-Chaouche sources are located offshore within the Bejaia Bay (Fig. 1) whereas the other epicentres are located further north, off the Bejaia Bay.

The USGS, OCA/FMNEARb and GCMT solutions were synthesized by EMSC (European Mediterranean Seismological Centre) within a few hours after the origin time (2h40min, 6h39min and 12h57min, respectively). The ‘Polarities’ method, GRiD MT and OCA/FMNEARa solutions would be potentially available earlier at CENALT in the context of early tsunami warning (12, 15 and 52 min, respectively).

### 2.2. GRiD MT approach for an automated source characterization

The GRiD MT solution can be obtained within the 15 min timeframe required for the CENALT’s mission of first alert. Indeed, it can be used to determine the event’s origin time, location (latitude,

longitude and depth), moment magnitude as well as its focal mechanism. The approach makes use of automated time-domain moment tensor inversions performed at different potential source hypocentres over a 3-D grid. Several seismological centres around the world use it in their seismic monitoring algorithms (Tsuruoka et al. 2009; Guilhem & Dreger 2011; Guilhem et al. 2013; MacPherson et al. 2013; Lee et al. 2013). Guilhem et al. (2013) also showed its interest for the detection of large and potentially tsunamigenic earthquakes, and in particular for the 2011 M<sub>w</sub> 9.1 Tohoku-oki earthquake in Japan.

The approach GRiD MT has recently been tested for earthquakes occurring along the Algerian margin (Menager et al. 2023). It works with a limited number of five broad-band stations that are selected over a pool of regional stations, called CENALT’s backbone network, providing reliable realtime data to CENALT in the western Mediterranean region: IM.MDT, PM.PESTR, IU.PAB, RD.MTT2 and RD.ORT2 (here, network code followed by station name). The three-component seismic waveforms from these five selected stations are deconvolved with the instrument response, integrated to displacement, and filtered between periods of 25 and 50 s. The grid covering the Algerian margin ranges from 34.50° to 38.00°N in latitude, and from -5.75° to 7°E in longitude (Figs 1 and 2). Because CENALT is interested in potentially tsunamigenic earthquakes, only four source depths are considered: 10, 20, 30 and 40 km.

The GRiD MT solution is defined by calculating the variance reduction value (VR) (eq. 1). It assesses the fit between the observed and the predicted data during the inversion. The higher the VR value is, the better is the fit, and so is the solution’s confidence.

$$VR = \left[ 1 - \frac{\sum_{i} (\text{observed}_i - \text{predicted})^2}{\sum_{i} \text{observed}_i^2} \right] \times 100. \quad (1)$$

Fig. 2 shows the results of the rapid GRiD MT approach defined for the region of interest between 100 s before and after the event origin time (OT). The maximum peak of VR (68.4 per cent) obtained for the Bejaia earthquake is detected 5 s after the OT, at the grid point having the coordinates 37.00°N/5.25°E (Fig. 2). This optimal position is in good agreement with locations of other sources (Fig. 1 and Table 1) and is quickly calculated (within 15 min). The earthquake focal mechanism is reverse, West-East oriented, with a high double-couple component, and is very similar to published solutions (Table 1 and Fig. 2).

The best GRiD MT solution (i.e. the one having the maximal VR value) is found around 10–20 km depth, with a magnitude M<sub>w</sub> 5.9 (Fig. 2). It is worth noting that the depth is difficult to constrain at the frequency band used in this approach.

### 2.3. The FMNEAR method

An instance of FMNEAR (Delouis 2014) is running at the Observatoire de la Côte d'Azur (OCA), automated for the near-real time but also with the possibility to be relaunched manually with a modified parametrization. A first automatic solution was obtained at T0 + 52 min the day of the earthquake, but with quality criteria just below the thresholds defined for transmission to the EMSC (cf. OCA/FMNEARa solution, Table 1, Fig. 3). Obtained with five stations only, this solution is nonetheless meaningful and would have been of interest for the purpose of tsunami warning. A solution with a higher degree of confidence was obtained at T0 + 6h39, by allowing the use of stations at larger distances in order to compensate the absence of available records in the south. This later
<figure><img src="image_3.png" /><figcaption>**Figure 2.** Map of the GRiD MT results for the 18 March 2021 Bejaia earthquake showing the moment tensor solutions obtained for each grid point at the optimal depth (10–20 km) and origin time (00:04:06 UTC), colour-coded and sized by their goodness of fit (VR). The red star is the epicentre location published by Yelles-Chaouche *et al.* (2023). Focal mechanisms published by EMSC (European-Mediterranean Seismological Centre) from several contributors are shown and compared with the GRiD MT solution. The bottom-left map shows the 5 stations (yellow triangles) and the grid extension (black box) used for the GRiD MT method. The best GRID MT solution is represented by the biggest focal mechanism in black within the grid, with the centre of the MT corresponding to the optimal centroid position (grid point) found. The GRiD MT focal mechanism is compared with FMNEAR (Observatoire Côte d’Azur), USGS (United States Geological Survey), GCMT (Global Centroid Moment Tensor Project, Columbia University) and Y.-C. (Yelles-Chaouche *et al.* 2023) solutions.</figcaption></figure>

solution, labelled OCA/FMNEARb (Table 1) reached the thresholds to be published at the EMSC. In the present day, an instance of FMNEAR is under test at the CENALT, working with a list of pre-selection stations established on the basis of waveform inversions carried out for reference events ($M > 5.5$) in the Mediterranean Sea. Particular attention is paid to incorporate a significant number of stations with accelerometric sensors in order to compensate for the potential saturation of nearby stations in case of an event of large magnitude. By limiting the maximum of stations included in the inversion to six, well enough constrained solution can be obtained in the time frame $T0 + 30$ min.

# 3 TSUNAMI OBSERVATIONS AND SIMULATIONS

## 3.1. Tsunami observations

The 2021 Bejaia tsunami was recorded by several tide-gauges along the French coasts but observed amplitudes are generally smaller than a few centimetres and observed time-series (with a sampling interval of 1 or 2 min) are strongly affected by ambient noise. By default, data are bandpass filtered with cut-off periods between 2 min and 3 hr in order to eliminate oceanic tides and ambient noise. For illustration, Fig. 4 shows the noise and filtered signals recorded over a period of 12 hr for seven stations, 6 hr before and after the

earthquake origin time. Amplitudes of signals with a bandpass of 2–20 min are strongly reduced for signals recorded in the Toulon Bay (stations Ifremer and Toulon harbour), or at Cassis, Sète, Port-Mahon which indicates resonances around 20 min at these locations (Fig. 4).

For most tide-gauges, the tsunami signature is hardly readable regardless of the filter chosen. Finally, only four stations with an acceptable signal-to-noise ratio have been selected: three tide-gauges on the French coast (Nice, Monaco, Ifremer) and one tide-gauge at Port-Mahon (Minorca, Balearic Islands). Bathymetries of the harbours of Nice, Monaco, Port-Mahon and of Toulon Bay are shown in Fig. 5.

## 3.2. Description of the tsunami sources

Seven different focal mechanisms among those proposed in the literature are tested in this paper. The tsunami source is calculated from Okada's formulas (1985) based on the elastic dislocation theory. Tsunami propagation is simulated by the code Taitoko (Heinrich *et al.* 2021) which has been extensively used (Poupardin *et al.* 2018; Jamelot *et al.* 2019; Heinrich *et al.* 2023) and is described in Sections 3.3.1 and 3.3.2.

Parameters of the seven fault models introduced in Okada's formulas are described in Table 1. Only results of the nodal planes south
<figure><img src="image_3.png" /><figcaption>**Figure 2.** Map of the GRiD MT results for the 18 March 2021 Bejaia earthquake showing the moment tensor solutions obtained for each grid point at the optimal depth (10–20 km) and origin time (00:04:06 UTC), colour-coded and sized by their goodness of fit (VR). The red star is the epicentre location published by Yelles-Chaouche *et al.* (2023). Focal mechanisms published by EMSC (European-Mediterranean Seismological Centre) from several contributors are shown and compared with the GRiD MT solution. The bottom-left map shows the 5 stations (yellow triangles) and the grid extension (black box) used for the GRiD MT method. The best GRID MT solution is represented by the biggest focal mechanism in black within the grid, with the centre of the MT corresponding to the optimal centroid position (grid point) found. The GRiD MT focal mechanism is compared with FMNEAR (Observatoire Côte d’Azur), USGS (United States Geological Survey), GCMT (Global Centroid Moment Tensor Project, Columbia University) and Y.-C. (Yelles-Chaouche *et al.* 2023) solutions.</figcaption></figure>

solution, labelled OCA/FMNEARb (Table 1) reached the thresholds to be published at the EMSC. In the present day, an instance of FMNEAR is under test at the CENALT, working with a list of pre-selection stations established on the basis of waveform inversions carried out for reference events ($M > 5.5$) in the Mediterranean Sea. Particular attention is paid to incorporate a significant number of stations with accelerometric sensors in order to compensate for the potential saturation of nearby stations in case of an event of large magnitude. By limiting the maximum of stations included in the inversion to six, well enough constrained solution can be obtained in the time frame $T0 + 30$ min.

# 3 TSUNAMI OBSERVATIONS AND SIMULATIONS

## 3.1. Tsunami observations

The 2021 Bejaia tsunami was recorded by several tide-gauges along the French coasts but observed amplitudes are generally smaller than a few centimetres and observed time-series (with a sampling interval of 1 or 2 min) are strongly affected by ambient noise. By default, data are bandpass filtered with cut-off periods between 2 min and 3 hr in order to eliminate oceanic tides and ambient noise. For illustration, Fig. 4 shows the noise and filtered signals recorded over a period of 12 hr for seven stations, 6 hr before and after the

earthquake origin time. Amplitudes of signals with a bandpass of 2–20 min are strongly reduced for signals recorded in the Toulon Bay (stations Ifremer and Toulon harbour), or at Cassis, Sète, Port-Mahon which indicates resonances around 20 min at these locations (Fig. 4).

For most tide-gauges, the tsunami signature is hardly readable regardless of the filter chosen. Finally, only four stations with an acceptable signal-to-noise ratio have been selected: three tide-gauges on the French coast (Nice, Monaco, Ifremer) and one tide-gauge at Port-Mahon (Minorca, Balearic Islands). Bathymetries of the harbours of Nice, Monaco, Port-Mahon and of Toulon Bay are shown in Fig. 5.

## 3.2. Description of the tsunami sources

Seven different focal mechanisms among those proposed in the literature are tested in this paper. The tsunami source is calculated from Okada’s formulas (1985) based on the elastic dislocation theory. Tsunami propagation is simulated by the code Taitoko (Heinrich *et al.* 2021) which has been extensively used (Poupardin *et al.* 2018; Jamelot *et al.* 2019; Heinrich *et al.* 2023) and is described in Sections 3.3.1 and 3.3.2.

Parameters of the seven fault models introduced in Okada’s formulas are described in Table 1. Only results of the nodal planes south
<table><caption>**Table 1.** Location and focal mechanism solution for the two dipping nodal planes of the Bejaia earthquake. NP1 and NP2 parameters are given in the first and second lines, respectively. The blue MT solutions (2–4) are obtained the day of the seismic event at CENALT and OCA (except solution 3 which has been recomputed). The red MT solutions are summarized from Yelles-Chaouche *et al.* 2023 (solution 1) or issued from EMSC Moment Tensors web service (https://seismicportal.eu/eventdetails.html?unid=20210318.0000001) for solutions 5–7. The publication time is relative to the main shock occurrence time.</caption><thead><tr><th>Source</th><th>Epicentral location</th><th>Depth (km)</th><th>Strike (°)</th><th>NP1/NP2 Dip (°)</th><th>Rake (°)</th><th>M<sub>w</sub></th><th>Publication (relative to OT)</th><th>MT solution</th></tr></thead><tbody><tr><td>(1) CRAAG/Y-C (2023)</td><td>36.87°N; 5.28°E</td><td>9</td><td>71<br>270</td><td>46<br>45</td><td>76<br>104</td><td>6.0</td><td>n.c.</td><td></td></tr><tr><td>(2) CENALT/ Polarities</td><td>36.97°N; 5.20°E</td><td>10</td><td>77<br>257</td><td>11<br>79</td><td>91<br>90</td><td>5.9</td><td>12 min</td><td></td></tr><tr><td>(3) CENALT/ GRiDMT</td><td>37.00°N; 5.25°E</td><td>16</td><td>126<br>266</td><td>15<br>79</td><td>129<br>81</td><td>5.9</td><td>15 min</td><td></td></tr><tr><td>(4) OCA/ FMNEARa</td><td>36.92°N; 5.20°E</td><td>11</td><td>90<br>270</td><td>30<br>60</td><td>90<br>90</td><td>5.9</td><td>52 min</td><td></td></tr><tr><td>(5) USGS</td><td>36.92°N; 5.20°E</td><td>8</td><td>81<br>283</td><td>22<br>69</td><td>69<br>98</td><td>6.0</td><td>2h40</td><td></td></tr><tr><td>(6) OCA/ FMNEARb</td><td>36.90°N; 5.20°E</td><td>11</td><td>95<br>275</td><td>35<br>55</td><td>90<br>90</td><td>5.9</td><td>6h39 min</td><td></td></tr><tr><td>(7) GCMT</td><td>36.84°N; 5.28°E</td><td>12</td><td>85<br>282</td><td>28<br>63</td><td>75<br>98</td><td>6.0</td><td>12h57 min</td><td></td></tr></tbody></table>

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

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

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

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

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

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

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

dipping are presented. Preliminary tests have shown that sources associated with the second nodal plane produce similar initial water surfaces and thus similar water time-series at coastal stations. The fault dimensions have been defined following the empirical relationships of Wells & Coppersmith (1994). The length (L) and width (W) of the focal plane have been set at 15 and 7.5 km, respectively. Considering a standard rigidity of $\mu = 4 \times 10^{10}$ Pa, and an average slip $U = 0.3$ m, the seismic moment $M_0 = \mu ULW$ is $1.2 \times 10^{18}$ N.m which corresponds to a moment magnitude $M_w = 6.0$. The resulting initial water surface elevations are shown in Fig. 6. For each source, a test with $L = 12$ km, $W = 6$ km and $U = 0.5$ m has been performed corresponding to reduced-size sources but with the same seismic moment. Along the French coasts, the corresponding simulations show that time-series at tide-gauge stations are identical

in terms of periods and display amplitudes larger by approximately 10–15 per cent. Finally, tests with a larger fault width ($W = 10$ km, $U = 0.3$ m) show that times-series at Nice and Monaco are similar to signals in Fig. 6 provided that the seismic moment is the same.

The water surface displacement in Taitoko is either copied from the vertical bottom deformation or calculated through a transfer function taking into account the water depth (Heinrich *et al.* 2021). Sensitivity tests show that the transfer function is not required for this earthquake.

Seismic parameters of the sources OCA/FMNEARb and OCA/FMNEARa are very close to each other (Table 1). As expected, numerical results for both sources are then almost identical. In the following, only the results with OCA/FMNEARb are presented, hereafter labelled simply OCA.
<table><caption>**Table 1.** Location and focal mechanism solution for the two dipping nodal planes of the Bejaia earthquake. NP1 and NP2 parameters are given in the first and second lines, respectively. The blue MT solutions (2–4) are obtained the day of the seismic event at CENALT and OCA (except solution 3 which has been recomputed). The red MT solutions are summarized from Yelles-Chaouche *et al.* 2023 (solution 1) or issued from EMSC Moment Tensors web service (https://seismicportal.eu/eventdetails.html?unid=20210318.0000001) for solutions 5–7. The publication time is relative to the main shock occurrence time.</caption><thead><tr><th>Source</th><th>Epicentral location</th><th>Depth (km)</th><th>Strike (°)</th><th>NP1/NP2 Dip (°)</th><th>Rake (°)</th><th>M<sub>w</sub></th><th>Publication (relative to OT)</th><th>MT solution</th></tr></thead><tbody><tr><td>(1) CRAAG/Y-C (2023)</td><td>36.87°N; 5.28°E</td><td>9</td><td>71<br>270</td><td>46<br>45</td><td>76<br>104</td><td>6.0</td><td>n.c.</td><td></td></tr><tr><td>(2) CENALT/ Polarities</td><td>36.97°N; 5.20°E</td><td>10</td><td>77<br>257</td><td>11<br>79</td><td>91<br>90</td><td>5.9</td><td>12 min</td><td></td></tr><tr><td>(3) CENALT/ GRiDMT</td><td>37.00°N; 5.25°E</td><td>16</td><td>126<br>266</td><td>15<br>79</td><td>129<br>81</td><td>5.9</td><td>15 min</td><td></td></tr><tr><td>(4) OCA/ FMNEARa</td><td>36.92°N; 5.20°E</td><td>11</td><td>90<br>270</td><td>30<br>60</td><td>90<br>90</td><td>5.9</td><td>52 min</td><td></td></tr><tr><td>(5) USGS</td><td>36.92°N; 5.20°E</td><td>8</td><td>81<br>283</td><td>22<br>69</td><td>69<br>98</td><td>6.0</td><td>2h40</td><td></td></tr><tr><td>(6) OCA/ FMNEARb</td><td>36.90°N; 5.20°E</td><td>11</td><td>95<br>275</td><td>35<br>55</td><td>90<br>90</td><td>5.9</td><td>6h39 min</td><td></td></tr><tr><td>(7) GCMT</td><td>36.84°N; 5.28°E</td><td>12</td><td>85<br>282</td><td>28<br>63</td><td>75<br>98</td><td>6.0</td><td>12h57 min</td><td></td></tr></tbody></table>

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

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

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

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

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

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

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

dipping are presented. Preliminary tests have shown that sources associated with the second nodal plane produce similar initial water surfaces and thus similar water time-series at coastal stations. The fault dimensions have been defined following the empirical relationships of Wells & Coppersmith (1994). The length (L) and width (W) of the focal plane have been set at 15 and 7.5 km, respectively. Considering a standard rigidity of $\mu = 4 \times 10^{10}$ Pa, and an average slip $U = 0.3$ m, the seismic moment $M_0 = \mu ULW$ is $1.2 \times 10^{18}$ N.m which corresponds to a moment magnitude $M_w = 6.0$. The resulting initial water surface elevations are shown in Fig. 6. For each source, a test with $L = 12$ km, $W = 6$ km and $U = 0.5$ m has been performed corresponding to reduced-size sources but with the same seismic moment. Along the French coasts, the corresponding simulations show that time-series at tide-gauge stations are identical

in terms of periods and display amplitudes larger by approximately 10–15 per cent. Finally, tests with a larger fault width ($W = 10$ km, $U = 0.3$ m) show that times-series at Nice and Monaco are similar to signals in Fig. 6 provided that the seismic moment is the same.

The water surface displacement in Taitoko is either copied from the vertical bottom deformation or calculated through a transfer function taking into account the water depth (Heinrich *et al.* 2021). Sensitivity tests show that the transfer function is not required for this earthquake.

Seismic parameters of the sources OCA/FMNEARb and OCA/FMNEARa are very close to each other (Table 1). As expected, numerical results for both sources are then almost identical. In the following, only the results with OCA/FMNEARb are presented, hereafter labelled simply OCA.
<figure><img src="image_11.png" /><figcaption>**Figure 3.** Moment tensor solution and best-fitted waveforms (observed data in black, synthetics in red) obtained by the FMNEAR method, obtained at T0 + 52 min (cf. OCA/FMNEARa solution). The top panel provides information on the moment tensor parameters. Components displaying flat lines were discarded during the inversion.</figcaption></figure>

## 3.3. Tsunami simulation

Simulations are performed using the nested multigrid code Taitoko developed at CEA (Heinrich *et al. 2021*), which is shortly described in this section.

### 3.3.1. Equations and numerical methods

The code Taitoko solves either the non-linear shallow water or the weakly non-linear Boussinesq equations in a spherical coordinates system. More precisely, the so-called ‘standard Boussinesq equations’ (Peregrine 1967) are solved only in the ≪ mother ≫ grid covering the western Mediterranean Sea. Alternative formulations with improved dispersion properties allow for the propagation of shorter waves (Nwogu 1993) but are not yet implemented in this code.

As regards Boussinesq equations, the implicit nature of dispersive terms in the momentum equations requires solving a linear system of equations. For a given direction, the dispersion terms in the other

direction are discretized explicitly following a standard ADI (Alternating Direction Implicit) method. The two implicit momentum equations are then associated to two tridiagonal equations that are solved by an iterative Gauss-Seidel method. It has been checked than non-linear terms in the momentum equations are negligible in this study.

As regards shallow water propagation, the non-linear long wave equations are solved by means of a staggered-grid finite difference method. The finite-difference scheme is implicit and uses an iterative procedure based on a second-order Strong Stability-Preserving Runge-Kutta method (Heinrich *et al. 2021*).

### 3.3.2. Nested multigrid code

The system of multiple nested grids allows to model multiscale processes with increasing resolution: low-resolution grids cover deep water regions and include coastal areas covered by high-resolution grids.
<figure><img src="image_11.png" /><figcaption>**Figure 3.** Moment tensor solution and best-fitted waveforms (observed data in black, synthetics in red) obtained by the FMNEAR method, obtained at T0 + 52 min (cf. OCA/FMNEARa solution). The top panel provides information on the moment tensor parameters. Components displaying flat lines were discarded during the inversion.</figcaption></figure>

## 3.3. Tsunami simulation

Simulations are performed using the nested multigrid code Taitoko developed at CEA (Heinrich *et al. 2021*), which is shortly described in this section.

### 3.3.1. Equations and numerical methods

The code Taitoko solves either the non-linear shallow water or the weakly non-linear Boussinesq equations in a spherical coordinates system. More precisely, the so-called ‘standard Boussinesq equations’ (Peregrine 1967) are solved only in the ≪ mother ≫ grid covering the western Mediterranean Sea. Alternative formulations with improved dispersion properties allow for the propagation of shorter waves (Nwogu 1993) but are not yet implemented in this code.

As regards Boussinesq equations, the implicit nature of dispersive terms in the momentum equations requires solving a linear system of equations. For a given direction, the dispersion terms in the other

direction are discretized explicitly following a standard ADI (Alternating Direction Implicit) method. The two implicit momentum equations are then associated to two tridiagonal equations that are solved by an iterative Gauss-Seidel method. It has been checked than non-linear terms in the momentum equations are negligible in this study.

As regards shallow water propagation, the non-linear long wave equations are solved by means of a staggered-grid finite difference method. The finite-difference scheme is implicit and uses an iterative procedure based on a second-order Strong Stability-Preserving Runge-Kutta method (Heinrich *et al. 2021*).

### 3.3.2. Nested multigrid code

The system of multiple nested grids allows to model multiscale processes with increasing resolution: low-resolution grids cover deep water regions and include coastal areas covered by high-resolution grids.
<figure><img src="image_12.png" /><figcaption>**Figure 4.** Tide-gauge time-series recorded by the stations of Nice, Monaco, Ifremer (Toulon Bay), Port Mahon, Cassis, Toulon harbour and Sete. The theoretical first-arrival travel times of the tsunami are indicated by the vertical grey lines. Bandpass filters with a passband of 2–180 min (black), 2–40 min (red) and 2–20 min (green) are applied. The tsunami arrival on the French coasts occurs around *t* = 1h–1h30min after the seismic event (*t* = 0). Data at Nice, Monaco, Mahon, Toulon and Sete are sampled at 1 min. The stations Ifremer and Cassis belong to the network of 18 HTM-NET stations located on the French coast and their data are sampled at 2 min.</figcaption></figure>

In this study, the mother grid covers the western Mediterranean Sea using a resolution of 250 m × 250 m and a time step of 0.4 s. The mesh is composed of 4800 × 3400 cells with depths ranging from a minimum of 10 m at the coast (to cope with the insufficient resolution along the coast) to 3530 m offshore. Grid convergence has been checked by performing numerical tests with a resolution of 150 m × 150 m. For Nice, Monaco and Toulon, two additional coastal grids are used with resolutions of 25 and 5 m, respectively. For Minorca, only one coastal grid with a resolution of 30 meters is considered due to the lack of bathymetric data. Propagation across the Mediterranean Sea is simulated by solving either Boussinesq or shallow water equations. The CPU times for a 5-hr shallow water or Boussinesq simulations are about 10 min and 8 hr, respectively, on a HPC machine with 512 cores.

One-way and two-way communications between grids are imple-
mented in the code. The one-way communication is standard and
consists in interpolating velocities and water surface elevations at

the boundaries of the child nested grid from the solution calculated
in the parent grid. The two-way communication is optional. The ad-
ditional step performed in two-way simulations consists in updating
the solution of the parent grid from that of the child grid in regions
where the two grids coincide. More precisely, the values at a coarse
grid point are calculated from a combination of the four points of
the fine domain surrounding the coarse grid point.

3.3.3. *Other tsunami parameters*

Several other parameters have been tested in this study. A bottom
friction force (based on the Manning's formulation) has been intro-
duced in coastal grids but its effect is negligible. A sponge layer at
the open boundaries of the mother grid has been also computed to
limit spurious reflection of water waves but its effect is also minor.
<figure><img src="image_12.png" /><figcaption>**Figure 4.** Tide-gauge time-series recorded by the stations of Nice, Monaco, Ifremer (Toulon Bay), Port Mahon, Cassis, Toulon harbour and Sete. The theoretical first-arrival travel times of the tsunami are indicated by the vertical grey lines. Bandpass filters with a passband of 2–180 min (black), 2–40 min (red) and 2–20 min (green) are applied. The tsunami arrival on the French coasts occurs around *t* = 1h–1h30min after the seismic event (*t* = 0). Data at Nice, Monaco, Mahon, Toulon and Sete are sampled at 1 min. The stations Ifremer and Cassis belong to the network of 18 HTM-NET stations located on the French coast and their data are sampled at 2 min.</figcaption></figure>

In this study, the mother grid covers the western Mediterranean Sea using a resolution of 250 m × 250 m and a time step of 0.4 s. The mesh is composed of 4800 × 3400 cells with depths ranging from a minimum of 10 m at the coast (to cope with the insufficient resolution along the coast) to 3530 m offshore. Grid convergence has been checked by performing numerical tests with a resolution of 150 m × 150 m. For Nice, Monaco and Toulon, two additional coastal grids are used with resolutions of 25 and 5 m, respectively. For Minorca, only one coastal grid with a resolution of 30 meters is considered due to the lack of bathymetric data. Propagation across the Mediterranean Sea is simulated by solving either Boussinesq or shallow water equations. The CPU times for a 5-hr shallow water or Boussinesq simulations are about 10 min and 8 hr, respectively, on a HPC machine with 512 cores.

One-way and two-way communications between grids are imple-
mented in the code. The one-way communication is standard and
consists in interpolating velocities and water surface elevations at

the boundaries of the child nested grid from the solution calculated
in the parent grid. The two-way communication is optional. The ad-
ditional step performed in two-way simulations consists in updating
the solution of the parent grid from that of the child grid in regions
where the two grids coincide. More precisely, the values at a coarse
grid point are calculated from a combination of the four points of
the fine domain surrounding the coarse grid point.

3.3.3. *Other tsunami parameters*

Several other parameters have been tested in this study. A bottom
friction force (based on the Manning's formulation) has been intro-
duced in coastal grids but its effect is negligible. A sponge layer at
the open boundaries of the mother grid has been also computed to
limit spurious reflection of water waves but its effect is also minor.
<figure><img src="image_13.png" /><figcaption>**Figure 5.** Bathymetry of the Nice harbour, Monaco harbour, Toulon Bay and Port-Mahon inlet (Minorca). The location of tide-gauges is indicated by a red dot. Their geographic coordinates are (7.285°, 43.695°) at Nice, (7.4215°, 43.7286°) at Monaco, (5.8856°, 43.1052°) at Ifremer (Toulon Bay) and (4.2705°, 39.8930°) at Port Mahon. Isobaths in solid lines (0, 2, 5 and 10 m).</figcaption></figure>

Finally, the water compressibility has been taken into account following Wang (2015). The associated depth correction delays the arrival of water waves on the French coast by 30 s approximately.

## 3.4. Propagation of the tsunami across the Mediterranean Sea

The tsunami propagation towards French coasts is highly dispersive and is then simulated by solving standard Boussinesq equations. Time-series along the tsunami path from the source to the French coast (Fig. 7) show that the shallow water propagation is characterized by a first wave with the largest water height whereas Boussinesq simulations produce a train of water waves. Maps of water wave elevations show that offshore wavelengths range from 10 to 50 km in the middle of the Mediterranean Sea (with a maximum ocean depth of 3500 m). In harbours (for instance at Nice), differences between simulations are strongly attenuated due to resonances phenomena of harbours which dominate the propagation in the harbour. At Nice, signals calculated by Boussinesq or shallow water models are nearly identical except for the first wave (Fig. 7).

Recorded and simulated time-series are compared at the four selected tide-gauges (Fig. 8) over a period of 5 hr considering the reference source. The computed water waves are calculated using either the one-way or the two-way communication. Comparison between the two simulations show that the two-way mode is required at the harbours of Nice and Monaco located on the Côte d'Azur. These small harbours (their width is around 100 m) are affected by

multiple reflections which results in resonances of harbours with a period of about 5–7 min. The one-way communication seems to amplify these resonance phenomena after a few hours which does not fit with the attenuation of recorded signals. Conversely, few differences between 1-way and 2-way simulations are observed inside the Mahon inlet or in the Toulon Bay. At Nice and Monaco, amplitudes of recorded signals are approximately reproduced. Smaller source dimensions (12 km × 6 km) with a larger slip (U = 0.5 m) have been tested and result in larger amplitudes close to the maximum recorded ones.

At the Ifremer station, the recorded periods (about 15 min) are reproduced but amplitudes are about twice smaller than the recorded ones. At Port-Mahon, the signal is approximately reproduced provided that the 45 min period is filtered with a bandpass of 2–40 min.

## 3.5. Sensitivity to the seismic source

So far, time-series have been calculated using the Yelles-Chaouche et al. (2023) focal mechanism of the Bejaia earthquake (Table 1). Synthetic tsunami time-series for the five other seismic sources are presented only at Nice and Monaco in Figs 9 and 10, respectively. Results at Ifremer station (Toulon Bay) or at Port-Mahon are poor for the six sources and the time-series calculated for the Yelles-Chaouche source (presented in Fig. 8) are representative. At Ifremer station, none of the tested sources produces water wave amplitudes exceeding 2 cm and the effect of source dimensions is not sufficient to account for the maximum recorded amplitudes of about 5 cm.
<figure><img src="image_13.png" /><figcaption>**Figure 5.** Bathymetry of the Nice harbour, Monaco harbour, Toulon Bay and Port-Mahon inlet (Minorca). The location of tide-gauges is indicated by a red dot. Their geographic coordinates are (7.285°, 43.695°) at Nice, (7.4215°, 43.7286°) at Monaco, (5.8856°, 43.1052°) at Ifremer (Toulon Bay) and (4.2705°, 39.8930°) at Port Mahon. Isobaths in solid lines (0, 2, 5 and 10 m).</figcaption></figure>

Finally, the water compressibility has been taken into account following Wang (2015). The associated depth correction delays the arrival of water waves on the French coast by 30 s approximately.

## 3.4. Propagation of the tsunami across the Mediterranean Sea

The tsunami propagation towards French coasts is highly dispersive and is then simulated by solving standard Boussinesq equations. Time-series along the tsunami path from the source to the French coast (Fig. 7) show that the shallow water propagation is characterized by a first wave with the largest water height whereas Boussinesq simulations produce a train of water waves. Maps of water wave elevations show that offshore wavelengths range from 10 to 50 km in the middle of the Mediterranean Sea (with a maximum ocean depth of 3500 m). In harbours (for instance at Nice), differences between simulations are strongly attenuated due to resonances phenomena of harbours which dominate the propagation in the harbour. At Nice, signals calculated by Boussinesq or shallow water models are nearly identical except for the first wave (Fig. 7).

Recorded and simulated time-series are compared at the four selected tide-gauges (Fig. 8) over a period of 5 hr considering the reference source. The computed water waves are calculated using either the one-way or the two-way communication. Comparison between the two simulations show that the two-way mode is required at the harbours of Nice and Monaco located on the Côte d'Azur. These small harbours (their width is around 100 m) are affected by

multiple reflections which results in resonances of harbours with a period of about 5–7 min. The one-way communication seems to amplify these resonance phenomena after a few hours which does not fit with the attenuation of recorded signals. Conversely, few differences between 1-way and 2-way simulations are observed inside the Mahon inlet or in the Toulon Bay. At Nice and Monaco, amplitudes of recorded signals are approximately reproduced. Smaller source dimensions (12 km × 6 km) with a larger slip (U = 0.5 m) have been tested and result in larger amplitudes close to the maximum recorded ones.

At the Ifremer station, the recorded periods (about 15 min) are reproduced but amplitudes are about twice smaller than the recorded ones. At Port-Mahon, the signal is approximately reproduced provided that the 45 min period is filtered with a bandpass of 2–40 min.

## 3.5. Sensitivity to the seismic source

So far, time-series have been calculated using the Yelles-Chaouche et al. (2023) focal mechanism of the Bejaia earthquake (Table 1). Synthetic tsunami time-series for the five other seismic sources are presented only at Nice and Monaco in Figs 9 and 10, respectively. Results at Ifremer station (Toulon Bay) or at Port-Mahon are poor for the six sources and the time-series calculated for the Yelles-Chaouche source (presented in Fig. 8) are representative. At Ifremer station, none of the tested sources produces water wave amplitudes exceeding 2 cm and the effect of source dimensions is not sufficient to account for the maximum recorded amplitudes of about 5 cm.
<figure><img src="image_14.png" /><figcaption>**Figure 6.** Initial surface displacement for the Yelles-Chaouche, USGS, GRiD MT, Polarity, OCA/FMNEARb and GCMT source models (Table 1). Isobaths are shown in dashed line (–1000 m) and in solid lines (–2500 m, –100 m, 0 m).</figcaption></figure>

At Port-Mahon, none of the mechanisms reproduces the 30–45 min period observed wave. These discrepancies are attributed first to the poor quality of bathymetric data at these two locations. The tide-gauge recorded signals could be also questioned since tide-gauges are affected by ambient noise with periods similar to the tsunami ones (Fig. 4).

At Nice, time-series calculated by the GCMT, OCA or Yelles-Chaouche et al. (2023) sources are similar (Fig. 9) since the initial surface displacements of these sources are close to each other (Fig. 6). The POLARITY, GRiD MT and USGS sources are located further north and produce computed signals in nearly perfect phase with data of the first recorded packet of water waves (between $t = 1h30$ and $t = 3h00$). By contrast, the second calculated packet (between $t = 3$ hr and $t = 4h30$) exhibits a phase shift of about 3–4

min whatever the source. This phase shift could be attributed to the bathymetry imprecisions in the Nice area. As regards amplitudes, the largest ones are produced by the GRiD MT source for the two packets of water waves and they compare well with the absolute maximum amplitudes recorded at Nice.

For a given source, sensitivity tests on seismic parameters show that the most sensitive parameter is the source latitude. As a consequence, results calculated for the GRiD MT and POLARITY sources are almost identical. As regards these two solutions, wave amplitudes are strongly reduced after $t = 2h00$ for a source relocation 15 km to the south. This shift over a shallower region produces water waves with longer periods. In the Nice area, water waves are affected by multiple coastal reflections and longer waves finally result in smaller wave amplitudes in the Nice harbour. Conversely, a
<figure><img src="image_14.png" /><figcaption>**Figure 6.** Initial surface displacement for the Yelles-Chaouche, USGS, GRiD MT, Polarity, OCA/FMNEARb and GCMT source models (Table 1). Isobaths are shown in dashed line (–1000 m) and in solid lines (–2500 m, –100 m, 0 m).</figcaption></figure>

At Port-Mahon, none of the mechanisms reproduces the 30–45 min period observed wave. These discrepancies are attributed first to the poor quality of bathymetric data at these two locations. The tide-gauge recorded signals could be also questioned since tide-gauges are affected by ambient noise with periods similar to the tsunami ones (Fig. 4).

At Nice, time-series calculated by the GCMT, OCA or Yelles-Chaouche et al. (2023) sources are similar (Fig. 9) since the initial surface displacements of these sources are close to each other (Fig. 6). The POLARITY, GRiD MT and USGS sources are located further north and produce computed signals in nearly perfect phase with data of the first recorded packet of water waves (between $t = 1h30$ and $t = 3h00$). By contrast, the second calculated packet (between $t = 3$ hr and $t = 4h30$) exhibits a phase shift of about 3–4

min whatever the source. This phase shift could be attributed to the bathymetry imprecisions in the Nice area. As regards amplitudes, the largest ones are produced by the GRiD MT source for the two packets of water waves and they compare well with the absolute maximum amplitudes recorded at Nice.

For a given source, sensitivity tests on seismic parameters show that the most sensitive parameter is the source latitude. As a consequence, results calculated for the GRiD MT and POLARITY sources are almost identical. As regards these two solutions, wave amplitudes are strongly reduced after $t = 2h00$ for a source relocation 15 km to the south. This shift over a shallower region produces water waves with longer periods. In the Nice area, water waves are affected by multiple coastal reflections and longer waves finally result in smaller wave amplitudes in the Nice harbour. Conversely, a
<figure><img src="image_15.png" /><figcaption>**Figure 7.** Time-series along the path from the epicentre to Nice harbour considering the Yelles-Chaouche *et al.* (2023) source. Numerical simulations are performed by solving either Boussinesq equations (black) or shallow water equations (red). The first three time-series are calculated in the mother grid at offshore points (6°E, 39°N), (7°E, 41.5°N) and 10 km south of Nice at (7.28°E, 43.6°N). The 4th time-series is calculated in the child grid at the Nice tide-gauge location (7.285°E, 43.695°N) with a grid resolution of 5 m.</figcaption></figure>

<figure><img src="image_16.png" /><figcaption>**Figure 8.** Comparisons between recorded time-series (dashed black line) and those simulated using the one-way (green) or two-ways (red) communication at Nice, Monaco, Ifremer (Toulon Bay) and Port Mahon stations, considering the Yelles-Chaouche *et al.* (2023) source. Recorded time-series are bandpass filtered. The resolutions of the child grids covering Nice, Monaco and Toulon are 5 m, the resolution of the Minorca grid is 30 m.</figcaption></figure>
<figure><img src="image_15.png" /><figcaption>**Figure 7.** Time-series along the path from the epicentre to Nice harbour considering the Yelles-Chaouche *et al.* (2023) source. Numerical simulations are performed by solving either Boussinesq equations (black) or shallow water equations (red). The first three time-series are calculated in the mother grid at offshore points (6°E, 39°N), (7°E, 41.5°N) and 10 km south of Nice at (7.28°E, 43.6°N). The 4th time-series is calculated in the child grid at the Nice tide-gauge location (7.285°E, 43.695°N) with a grid resolution of 5 m.</figcaption></figure>

<figure><img src="image_16.png" /><figcaption>**Figure 8.** Comparisons between recorded time-series (dashed black line) and those simulated using the one-way (green) or two-ways (red) communication at Nice, Monaco, Ifremer (Toulon Bay) and Port Mahon stations, considering the Yelles-Chaouche *et al.* (2023) source. Recorded time-series are bandpass filtered. The resolutions of the child grids covering Nice, Monaco and Toulon are 5 m, the resolution of the Minorca grid is 30 m.</figcaption></figure>
<figure><img src="image_17.png" /><figcaption>**Figure 9.** Comparison at the tide-gauge of Nice between recorded time-series of surface elevation (dashed black line) and those simulated using six seismic sources (red). The resolution of the child grid is 5 m.</figcaption></figure>

shift 10–15 km to the north of the other sources (Yelles-Chaouche, OCA and GCMT), leads to offshore shorter waves and then larger wave amplitudes in the Nice harbour.

At Monaco, time-series calculated for the six sources display similar amplitudes within the first hours. As regards phases, time-series between $t = 1h30$ and $t = 2h00$ do not present a phase shift for the OCA, GCMT and Yelles-Chaouche sources (Fig. 10) whereas GRiD MT, POLARITY and USGS sources (located further north) are phase shifted by a few minutes. Significant synthetic resonance for a second packet (after $t = 3h00$) is calculated for the USGS and OCA sources (Fig. 10) but is not recorded by the tide-gauge. As for Nice, the location of epicentres is the most sensitive parameter. Results show that resonance after $t = 3$ hr for the USGS and OCA sources is eliminated by moving their epicentres 10 km to the south.

## 4 DISCUSSION AND CONCLUSION

Seismic events in the northern Algerian margin may trigger tsunamis in the western Mediterranean Sea. Compared to the Pacific Ocean, the Mediterranean Sea is a small basin and tsunamis

generated near one coast can reach the opposing shore of the basin within one hour (Heidarzadeh et al. 2017; Heidarzadeh & Gusman 2021; Necmioglu et al. 2023). Observations and/or simulations confirm that the 2003 (Alasset et al. 2006; Masina et al. 2020) and 2021 tsunamis triggered off Algeria took 70–80 min to reach the coasts of southern France. The tsunami alert off the French coasts strongly depends on the rapid determination of a reliable seismic source.

In this paper, seven focal mechanisms of the 2021 Algerian main shock have been tested. For each source, propagation of these water waves is simulated towards French coasts and more precisely towards the few harbours, which detected the tsunami. Numerical results at the four selected tide-gauges depend on the seismic source but do not differ significantly. For a given source, the most sensitive parameter is the distance of the source from Algerian coasts. Resonance in the harbours of Nice and Monaco is activated or not depending on the wavelengths of offshore waves.

Among the seven seismic sources, those of the GRiD MT and FMNEAR methods, potentially obtained within 15 min for the first and effectively calculated in 52 min after origin time for the second, provide synthetic time-series which are close to those calculated
<figure><img src="image_17.png" /><figcaption>**Figure 9.** Comparison at the tide-gauge of Nice between recorded time-series of surface elevation (dashed black line) and those simulated using six seismic sources (red). The resolution of the child grid is 5 m.</figcaption></figure>

shift 10–15 km to the north of the other sources (Yelles-Chaouche, OCA and GCMT), leads to offshore shorter waves and then larger wave amplitudes in the Nice harbour.

At Monaco, time-series calculated for the six sources display similar amplitudes within the first hours. As regards phases, time-series between $t = 1h30$ and $t = 2h00$ do not present a phase shift for the OCA, GCMT and Yelles-Chaouche sources (Fig. 10) whereas GRiD MT, POLARITY and USGS sources (located further north) are phase shifted by a few minutes. Significant synthetic resonance for a second packet (after $t = 3h00$) is calculated for the USGS and OCA sources (Fig. 10) but is not recorded by the tide-gauge. As for Nice, the location of epicentres is the most sensitive parameter. Results show that resonance after $t = 3$ hr for the USGS and OCA sources is eliminated by moving their epicentres 10 km to the south.

## 4 DISCUSSION AND CONCLUSION

Seismic events in the northern Algerian margin may trigger tsunamis in the western Mediterranean Sea. Compared to the Pacific Ocean, the Mediterranean Sea is a small basin and tsunamis

generated near one coast can reach the opposing shore of the basin within one hour (Heidarzadeh et al. 2017; Heidarzadeh & Gusman 2021; Necmioglu et al. 2023). Observations and/or simulations confirm that the 2003 (Alasset et al. 2006; Masina et al. 2020) and 2021 tsunamis triggered off Algeria took 70–80 min to reach the coasts of southern France. The tsunami alert off the French coasts strongly depends on the rapid determination of a reliable seismic source.

In this paper, seven focal mechanisms of the 2021 Algerian main shock have been tested. For each source, propagation of these water waves is simulated towards French coasts and more precisely towards the few harbours, which detected the tsunami. Numerical results at the four selected tide-gauges depend on the seismic source but do not differ significantly. For a given source, the most sensitive parameter is the distance of the source from Algerian coasts. Resonance in the harbours of Nice and Monaco is activated or not depending on the wavelengths of offshore waves.

Among the seven seismic sources, those of the GRiD MT and FMNEAR methods, potentially obtained within 15 min for the first and effectively calculated in 52 min after origin time for the second, provide synthetic time-series which are close to those calculated
<figure><img src="image_18.png" /><figcaption>**Figure 10.** Comparison at the tide-gauge of Monaco between recorded time-series of surface elevation (dashed black line) and those simulated using six seismic sources (red).</figcaption></figure>

using the three sources published a few hours later. Work is underway to recover more rapidly seismic data so that the FMNEAR solution could be obtained within 30 min. An operational implementation of both methods is planned in the framework of CENALT.

The complex North Algeria seismotectonics is expressed by numerous faults systems accommodating convergence between Africa and Europa, inherited from probable inversion of former passive margins (Leffondré et al. 2021). Although tsunami models are not extremely different depending on the source parameters, the latitudinal dependence recalls the need to accurately characterize active faulting off north Algeria. The 2021 rupture may have occurred close to, and west to the 1856 earthquake area, possibly as a result of stress transfer (Leffondré et al. 2021) induced by this major seismic sequence of the last 200 yr associated with a significant tsunami (Yelles-Chaouche et al. 2009). This very active thrusting area can therefore produce again very important tsunamis in the western Mediterranean Sea, posing a high risk level for Algerian coastlines as much as for distant regions (Balearic, Sardinia, French coastlines). The 2021 event also offers the opportunity to strengthen

research studies in this region, and to improve tsunami warning procedures.

So far, maximum water heights along the French coast are inferred from tables decision matrix established for a large number of seismic scenarios. The next step will be to simulate the propagation of tsunamis for both sources in the context of tsunami warning. As regards the 2021 tsunami, shallow water simulations with two levels (a mother grid with a resolution of 500 and 25 m child grids) produce time-series in Nice and Monaco close to results presented in this paper. These simulations are very low CPU-time (around 2 min) and should be used in the future to forecast tsunami amplitudes in French harbours. In any case, most of the methods used in real time have to benefit from the up-to-date knowledge and characterization of the seismically active faulting patterns.

## ANNEX

Solution obtained by the method of ‘polarities’ at CENALT within 12 min after the origin time.
<figure><img src="image_18.png" /><figcaption>**Figure 10.** Comparison at the tide-gauge of Monaco between recorded time-series of surface elevation (dashed black line) and those simulated using six seismic sources (red).</figcaption></figure>

using the three sources published a few hours later. Work is underway to recover more rapidly seismic data so that the FMNEAR solution could be obtained within 30 min. An operational implementation of both methods is planned in the framework of CENALT.

The complex North Algeria seismotectonics is expressed by numerous faults systems accommodating convergence between Africa and Europa, inherited from probable inversion of former passive margins (Leffondré et al. 2021). Although tsunami models are not extremely different depending on the source parameters, the latitudinal dependence recalls the need to accurately characterize active faulting off north Algeria. The 2021 rupture may have occurred close to, and west to the 1856 earthquake area, possibly as a result of stress transfer (Leffondré et al. 2021) induced by this major seismic sequence of the last 200 yr associated with a significant tsunami (Yelles-Chaouche et al. 2009). This very active thrusting area can therefore produce again very important tsunamis in the western Mediterranean Sea, posing a high risk level for Algerian coastlines as much as for distant regions (Balearic, Sardinia, French coastlines). The 2021 event also offers the opportunity to strengthen

research studies in this region, and to improve tsunami warning procedures.

So far, maximum water heights along the French coast are inferred from tables decision matrix established for a large number of seismic scenarios. The next step will be to simulate the propagation of tsunamis for both sources in the context of tsunami warning. As regards the 2021 tsunami, shallow water simulations with two levels (a mother grid with a resolution of 500 and 25 m child grids) produce time-series in Nice and Monaco close to results presented in this paper. These simulations are very low CPU-time (around 2 min) and should be used in the future to forecast tsunami amplitudes in French harbours. In any case, most of the methods used in real time have to benefit from the up-to-date knowledge and characterization of the seismically active faulting patterns.

## ANNEX

Solution obtained by the method of ‘polarities’ at CENALT within 12 min after the origin time.
CENALT - Polarisations

Bejaia
18/03/2021
00:04 UTC

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

Lat / Lon = (36.97°, 5.20°)

Depth = 10 km

M<sub>W</sub> = 6.0

N₁: (Strike, Dip, Rake) = (77, 11, 91)

N₂: (Strike, Dip, Rake) = (257, 79, 90)

REFERENCES

Alasset, P.-J., Hebert, H., Maouche, S., Calbini, V. & Meghraoui, M., 2006. The tsunami induced by the 2003 Zemmouri earthquake (Mw = 6.9, Algeria): modelling and results, *Geophys. J. Int.*, **166**, 213–226.

Braunmiller, J. & Bernardi, F., 2005. The 2003 Boumerdes, Algeria earthquake: regional moment tensor analysis, *Geophys. Res. Lett.*, **32**(6), doi:10.1029/2004GL022038.

Buforn, E., Sanz de Galdeano, C. & Udias, A., 1995. Seismotectonics of the Ibero-Maghrebian region, *Tectonophysics*, **248**, 247–261.

Buforn, E., Bezzeghoud, M., Udias, A. & Pro, C., 2004. Seismic sources on the Iberia-African plate boundary and their tectonic implications, *Pure appl. Geophys.*, **161**, 623–646.

Buforn, E., Udias, A. & Pro, C., 2015. Large earthquakes at the Ibero-Maghrebian region: basis for an EEWS, *Pure appl. Geophys.*, **172**, 2387–2396.

Buforn, E., Pro, C., Sanz de Galdeano, C., Cantavella, J.V., Cesca, S., Caldeira, B., Udias, A. & Mattesini, M., 2017. The 2016 south Alboran earthquake (Mw = 6.4): a reactivation of the Ibero-Maghrebian region?, *Tectonophysics*, **712-713**, 704–715.

De Caro, M., Montuori, C., Frugoni, F., Monna, S. & Giuntini, A., 2023. Seismic T phases in the Western-Central Mediterranean: source of Seismic hazard?, *Seismol. Res. Lett.*, **95**, 859–869.

Delouis, B., 2014. Fmnear: determination of focal mechanism and first estimate of rupture directivity using near-source records and a linear distribution of point sources, *Bull. seism. Soc. Am.*, **104**, 1479–1500.

Delouis, B. et al. 2004. Slip distribution of the 2003 Boumerdes-Zemmouri earthquake, Algeria, from teleseismic, GPS, and coastal uplift data, *Geophys. Res. Lett.*, **31**(18), doi:10.1029/2004GL020687.

Gailler, A., Hebert, H., Loevenbruck, A. & Hernandez, B., 2013. Simulation systems for tsunami wave propagation forecasting within the French tsunami warning center, *Nat. Haz. Earth Syst. Sci.*, **13**, 2465–2482.

Guilhem, A. & Dreger, D., 2011. Rapid detection and characterization of large earthquakes using quasi-finite-source Green's functions in continuous moment tensor inversion, *Geophys. Res. Lett.*, **38**(13), doi:10.1029/2011GL047550.

Guilhem, A., Dreger, D., Tsuruoka, H. & Kawakatsu, H., 2013. Moment tensors for rapid characterization of megathrust earthquakes: the example of the 2011 M 9 Tohoku-oki, Japan earthquake, *Geophys. J. Int.*, **192**, 759–772.

Heidarzadeh, M. & Gusman, A.R., 2021. Source modeling and spectral analysis of the Crete tsunami of 2nd May 2020 along the Hellenic Subduction Zone, offshore Greece, *Earth Planets Space*, **73**, doi:10.1186/s40623-021-01394-4.

Heidarzadeh, M. & Satake, K., 2013. The 21 May 2003 tsunami in the Western Mediterranean Sea: statistical and wavelet analyses, *Pure appl. Geophys.*, **170**, 1449–1462.

Heidarzadeh, M. et al., 2017. Bodrum–Kos (Turkey–Greece) Mw 6.6 earthquake and tsunami of 20 July 2017: a test for the Mediterranean tsunami warning system, *Geosci. Lett.*, **4**, doi:10.1186/s40562-017-0097-0.

Heinrich, P., Jamelot, A., Cauquis, A. & Gailler, A., 2021. Taitoko, an advanced code for tsunami propagation, developed at the French Tsunami Warning Centers, *Eur. J. Mech.—B/Fluids*, **88**, 272–288.

Heinrich, P., Gailler, A., Dupont, A., Rey, V., Hébert, H. & Listowski, C., 2023. Observations and simulations of the meteotsunami generated by the

Tonga eruption on 15 January 2022 in the Mediterranean Sea, *Geophys. J. Int.*, **234**, 2, 903–914,

ISC, 2024. International Seismological Centre. On-line Bulletin. doi.org/10.31905/D808B830.

Jamelot, A., Gailler, A., Heinrich, P., Vallage, A. & Champenois, J., 2019. Tsunami simulations of the Sulawesi Mw7.5 event: comparison of seismic sources issued from a Tsunami warning context versus post-event finite source, *Pure appl. Geophys.*, **176**, 3351–3376.

Kawakatsu, H., 1998. On the realtime monitoring of the long-period seismic wavefield, *Bull. Earthq. Res. Inst. Univ. Tokyo*, **73**, 267–274.

Lee, S.-J. et al., 2013. Towards real-time regional earthquake simulation I: real-time moment tensor monitoring (RMT) for regional events in Taiwan, *Geophys. J. Int.*, **196**(1), 432–446.

Leffondré, P., Déverchère, J., Medaouri, M., Klingelhoefer, F., Graindorge, D. & Arab, M., 2021. Ongoing inversion of a passive margin: spatial variability of strain markers along the Algerian margin and basin (Mediterranean Sea) and seismotectonic implications, *Front. Earth Sci.*, **9**, doi:10.3389/feart.2021.674584.

Macpherson, K.A., Ruppert, N.A., Dreger, D.S., Lombard, P.N., Freymueller, J.T., Nicolsky, D. & Guilhem, A., 2013. Towards implementation of the grid MT algorithm for near real-time calculation of moment tensors at the Alaska earthquake information center, in *American Geophysical Union, Fall Meeting 2013*, abstract id. S51A-2314. https://ui.adsabs.harvard.edu/abs/2013AGUFM.S51A2314M.

Masina, M., Archetti, R. & Lamberti, A., 2020. 21 May 2003 Boumerdès earthquake: numerical investigations of the rupture mechanism effects on the induced tsunami and its impact in harbors, *J. Mar. Sci. Eng.*, **8**, 933.

Menager, M., Aurélie Guilhem Trilla, A. & Delouis, B., 2023. An automated earthquake detection and characterization tool for rapid earthquake and tsunami response in Western Mediterranean, *Geophys. J. Int.*, **234**(2), 998–1011.

Necmioglu, O. et al., 2023. Landslide induced tsunami hazard at volcanoes: the case of Santorini, *Pure appl. Geophys.*, **180**, 1811–1834.

Nwogu, O., 1993. Alternative form of Boussinesq equations for nearshore wave propagation, *J. Waterw. Port, Coastal Ocean Eng.*, **119**(6), 618–638.

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

Peregrine, D.H., 1967. Long waves on a beach, *J. Fluid Mech*, **77**, 417–431.

Poupardin, A., Heinrich, P., Hébert, H., Schindelé, F., Jamelot, A. & Reymond, D., 2018. Traveltime delay relative to the maximum energy of the wave train for dispersive tsunamis propagating across the Pacific Ocean: the case of 2010 and 2015 Chilean tsunamis, *Geophys. J. Int.*, **214**, 1538–1555.

Roger, J. & Hebert, H., 2008. The 1856 Djijelli (Algeria) earthquake and tsunami: source parameters and implications for tsunami hazard in the Balearic Islands, *Nat. Haz. Earth Syst. Sci.*, **8**, 721–731.

Roudil, P. et al. 2013. The French tsunami warning center for the mediterranean and Northeast Atlantic: CENALT, *Sci. Tsunami Haz.*, **32**, 1–7.

Sahal, A., Roger, J., Lemaire, B., Hebert, H., Schindelé, F. & Lavigne, F., 2009. The tsunami triggered by the 21 May 2003 Boumerdès-Zemmouri (Algeria) earthquake: field investigations on the French Mediterranean coast and tsunami modelling, *Nat. Haz. Earth Syst. Sci.*, **9**, 1823–1834.

Schindelé, F. et al. 2015. Implementation and challenges of the Tsunami Warning System in the Western Mediterranean, *Pure appl. Geophys.*, **172**, 821–833.

Stich, D., Mancilla, F., Baumont, D. & Morales, J., 2005. Source analysis of the Mw 6.3 2004 Al Hoceima earthquake (Morocco) using regional apparent source time functions, *J. geophys. Res.*, **110**, doi:10.1029/2004JB003366.

Stich, D., Serpelloni, E., Mancilla, F. & Morales, J., 2006. Kinematics of the Iberia-Maghreb plate contact from seismic moment tensors and GPS observations, *Tectonophysics*, **426**, 295–317.

Stich, D., Martín, R., Morales, J., López-Comino, J.Á. & Mancilla, F.L., 2020. Slip partitioning in the 2016 Alboran Sea earthquake sequence (Western Mediterranean), *Front. Earth Sci.*, **8**, doi:10.3389/feart.2020.587356.
CENALT - Polarisations

Bejaia
18/03/2021
00:04 UTC

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

Lat / Lon = (36.97°, 5.20°)

Depth = 10 km

M<sub>W</sub> = 6.0

N₁: (Strike, Dip, Rake) = (77, 11, 91)

N₂: (Strike, Dip, Rake) = (257, 79, 90)

REFERENCES

Alasset, P.-J., Hebert, H., Maouche, S., Calbini, V. & Meghraoui, M., 2006. The tsunami induced by the 2003 Zemmouri earthquake (Mw = 6.9, Algeria): modelling and results, *Geophys. J. Int.*, **166**, 213–226.

Braunmiller, J. & Bernardi, F., 2005. The 2003 Boumerdes, Algeria earthquake: regional moment tensor analysis, *Geophys. Res. Lett.*, **32**(6), doi:10.1029/2004GL022038.

Buforn, E., Sanz de Galdeano, C. & Udias, A., 1995. Seismotectonics of the Ibero-Maghrebian region, *Tectonophysics*, **248**, 247–261.

Buforn, E., Bezzeghoud, M., Udias, A. & Pro, C., 2004. Seismic sources on the Iberia-African plate boundary and their tectonic implications, *Pure appl. Geophys.*, **161**, 623–646.

Buforn, E., Udias, A. & Pro, C., 2015. Large earthquakes at the Ibero-Maghrebian region: basis for an EEWS, *Pure appl. Geophys.*, **172**, 2387–2396.

Buforn, E., Pro, C., Sanz de Galdeano, C., Cantavella, J.V., Cesca, S., Caldeira, B., Udias, A. & Mattesini, M., 2017. The 2016 south Alboran earthquake (Mw = 6.4): a reactivation of the Ibero-Maghrebian region?, *Tectonophysics*, **712-713**, 704–715.

De Caro, M., Montuori, C., Frugoni, F., Monna, S. & Giuntini, A., 2023. Seismic T phases in the Western-Central Mediterranean: source of Seismic hazard?, *Seismol. Res. Lett.*, **95**, 859–869.

Delouis, B., 2014. Fmnear: determination of focal mechanism and first estimate of rupture directivity using near-source records and a linear distribution of point sources, *Bull. seism. Soc. Am.*, **104**, 1479–1500.

Delouis, B. et al. 2004. Slip distribution of the 2003 Boumerdes-Zemmouri earthquake, Algeria, from teleseismic, GPS, and coastal uplift data, *Geophys. Res. Lett.*, **31**(18), doi:10.1029/2004GL020687.

Gailler, A., Hebert, H., Loevenbruck, A. & Hernandez, B., 2013. Simulation systems for tsunami wave propagation forecasting within the French tsunami warning center, *Nat. Haz. Earth Syst. Sci.*, **13**, 2465–2482.

Guilhem, A. & Dreger, D., 2011. Rapid detection and characterization of large earthquakes using quasi-finite-source Green's functions in continuous moment tensor inversion, *Geophys. Res. Lett.*, **38**(13), doi:10.1029/2011GL047550.

Guilhem, A., Dreger, D., Tsuruoka, H. & Kawakatsu, H., 2013. Moment tensors for rapid characterization of megathrust earthquakes: the example of the 2011 M 9 Tohoku-oki, Japan earthquake, *Geophys. J. Int.*, **192**, 759–772.

Heidarzadeh, M. & Gusman, A.R., 2021. Source modeling and spectral analysis of the Crete tsunami of 2nd May 2020 along the Hellenic Subduction Zone, offshore Greece, *Earth Planets Space*, **73**, doi:10.1186/s40623-021-01394-4.

Heidarzadeh, M. & Satake, K., 2013. The 21 May 2003 tsunami in the Western Mediterranean Sea: statistical and wavelet analyses, *Pure appl. Geophys.*, **170**, 1449–1462.

Heidarzadeh, M. et al., 2017. Bodrum–Kos (Turkey–Greece) Mw 6.6 earthquake and tsunami of 20 July 2017: a test for the Mediterranean tsunami warning system, *Geosci. Lett.*, **4**, doi:10.1186/s40562-017-0097-0.

Heinrich, P., Jamelot, A., Cauquis, A. & Gailler, A., 2021. Taitoko, an advanced code for tsunami propagation, developed at the French Tsunami Warning Centers, *Eur. J. Mech.—B/Fluids*, **88**, 272–288.

Heinrich, P., Gailler, A., Dupont, A., Rey, V., Hébert, H. & Listowski, C., 2023. Observations and simulations of the meteotsunami generated by the

Tonga eruption on 15 January 2022 in the Mediterranean Sea, *Geophys. J. Int.*, **234**, 2, 903–914,

ISC, 2024. International Seismological Centre. On-line Bulletin. doi.org/10.31905/D808B830.

Jamelot, A., Gailler, A., Heinrich, P., Vallage, A. & Champenois, J., 2019. Tsunami simulations of the Sulawesi Mw7.5 event: comparison of seismic sources issued from a Tsunami warning context versus post-event finite source, *Pure appl. Geophys.*, **176**, 3351–3376.

Kawakatsu, H., 1998. On the realtime monitoring of the long-period seismic wavefield, *Bull. Earthq. Res. Inst. Univ. Tokyo*, **73**, 267–274.

Lee, S.-J. et al., 2013. Towards real-time regional earthquake simulation I: real-time moment tensor monitoring (RMT) for regional events in Taiwan, *Geophys. J. Int.*, **196**(1), 432–446.

Leffondré, P., Déverchère, J., Medaouri, M., Klingelhoefer, F., Graindorge, D. & Arab, M., 2021. Ongoing inversion of a passive margin: spatial variability of strain markers along the Algerian margin and basin (Mediterranean Sea) and seismotectonic implications, *Front. Earth Sci.*, **9**, doi:10.3389/feart.2021.674584.

Macpherson, K.A., Ruppert, N.A., Dreger, D.S., Lombard, P.N., Freymueller, J.T., Nicolsky, D. & Guilhem, A., 2013. Towards implementation of the grid MT algorithm for near real-time calculation of moment tensors at the Alaska earthquake information center, in *American Geophysical Union, Fall Meeting 2013*, abstract id. S51A-2314. https://ui.adsabs.harvard.edu/abs/2013AGUFM.S51A2314M.

Masina, M., Archetti, R. & Lamberti, A., 2020. 21 May 2003 Boumerdès earthquake: numerical investigations of the rupture mechanism effects on the induced tsunami and its impact in harbors, *J. Mar. Sci. Eng.*, **8**, 933.

Menager, M., Aurélie Guilhem Trilla, A. & Delouis, B., 2023. An automated earthquake detection and characterization tool for rapid earthquake and tsunami response in Western Mediterranean, *Geophys. J. Int.*, **234**(2), 998–1011.

Necmioglu, O. et al., 2023. Landslide induced tsunami hazard at volcanoes: the case of Santorini, *Pure appl. Geophys.*, **180**, 1811–1834.

Nwogu, O., 1993. Alternative form of Boussinesq equations for nearshore wave propagation, *J. Waterw. Port, Coastal Ocean Eng.*, **119**(6), 618–638.

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

Peregrine, D.H., 1967. Long waves on a beach, *J. Fluid Mech*, **77**, 417–431.

Poupardin, A., Heinrich, P., Hébert, H., Schindelé, F., Jamelot, A. & Reymond, D., 2018. Traveltime delay relative to the maximum energy of the wave train for dispersive tsunamis propagating across the Pacific Ocean: the case of 2010 and 2015 Chilean tsunamis, *Geophys. J. Int.*, **214**, 1538–1555.

Roger, J. & Hebert, H., 2008. The 1856 Djijelli (Algeria) earthquake and tsunami: source parameters and implications for tsunami hazard in the Balearic Islands, *Nat. Haz. Earth Syst. Sci.*, **8**, 721–731.

Roudil, P. et al. 2013. The French tsunami warning center for the mediterranean and Northeast Atlantic: CENALT, *Sci. Tsunami Haz.*, **32**, 1–7.

Sahal, A., Roger, J., Lemaire, B., Hebert, H., Schindelé, F. & Lavigne, F., 2009. The tsunami triggered by the 21 May 2003 Boumerdès-Zemmouri (Algeria) earthquake: field investigations on the French Mediterranean coast and tsunami modelling, *Nat. Haz. Earth Syst. Sci.*, **9**, 1823–1834.

Schindelé, F. et al. 2015. Implementation and challenges of the Tsunami Warning System in the Western Mediterranean, *Pure appl. Geophys.*, **172**, 821–833.

Stich, D., Mancilla, F., Baumont, D. & Morales, J., 2005. Source analysis of the Mw 6.3 2004 Al Hoceima earthquake (Morocco) using regional apparent source time functions, *J. geophys. Res.*, **110**, doi:10.1029/2004JB003366.

Stich, D., Serpelloni, E., Mancilla, F. & Morales, J., 2006. Kinematics of the Iberia-Maghreb plate contact from seismic moment tensors and GPS observations, *Tectonophysics*, **426**, 295–317.

Stich, D., Martín, R., Morales, J., López-Comino, J.Á. & Mancilla, F.L., 2020. Slip partitioning in the 2016 Alboran Sea earthquake sequence (Western Mediterranean), *Front. Earth Sci.*, **8**, doi:10.3389/feart.2020.587356.
Tsuruoka, H., Kawakatsu, H. & Urabe, T., 2009. Grid mt (grid-based real-time determination of moment tensors) monitoring the long-period seismic wavefield, *Phys. Earth planet. Inter.*, **175**, 8–16.

Wang, D., 2015. An ocean depth-correction method for reducing model errors in tsunami travel time: application to the 2010 Chile and 2011 Tohoku tsunamis. *Sci. Tsunami Haz.*, **34**, 1–22.

Wells, D. & Coppersmith, K., 1994. New empirical relationships among magnitude, rupture length, rupture width, rupture area, and surface displacement. *Bull. seism. Soc. Am.*, **84**, 4, 974–1002.

Yelles-Chaouche, A. et al., 2009. The 1856 tsunami of Djidjelli (Eastern Algeria): seismotectonics, modelling and hazard implications for the Algerian coast, *Pure appl. Geophys.*, **166**, 283–300.

Yelles-Chaouche, A. et al., 2023. The 2021–2022 Mw 6.0 Bejaia Bay, NE Algeria, earthquake sequence: tectonic implications at the Algerian margin between lesser and greater Kabylian blocks, *Acta Geophys.*, doi:10.1007/s11600-023-01171-9.
Tsuruoka, H., Kawakatsu, H. & Urabe, T., 2009. Grid mt (grid-based real-time determination of moment tensors) monitoring the long-period seismic wavefield, *Phys. Earth planet. Inter.*, **175**, 8–16.

Wang, D., 2015. An ocean depth-correction method for reducing model errors in tsunami travel time: application to the 2010 Chile and 2011 Tohoku tsunamis. *Sci. Tsunami Haz.*, **34**, 1–22.

Wells, D. & Coppersmith, K., 1994. New empirical relationships among magnitude, rupture length, rupture width, rupture area, and surface displacement. *Bull. seism. Soc. Am.*, **84**, 4, 974–1002.

Yelles-Chaouche, A. et al., 2009. The 1856 tsunami of Djidjelli (Eastern Algeria): seismotectonics, modelling and hazard implications for the Algerian coast, *Pure appl. Geophys.*, **166**, 283–300.

Yelles-Chaouche, A. et al., 2023. The 2021–2022 Mw 6.0 Bejaia Bay, NE Algeria, earthquake sequence: tectonic implications at the Algerian margin between lesser and greater Kabylian blocks, *Acta Geophys.*, doi:10.1007/s11600-023-01171-9.