1 Introduction

River mouth deviation is not a simple geometric deflection of the terminal part of a river, involving its hydrodynamics. Indeed, it could induce important effects within the environmental, ecological, geomorphological, and social frameworks.

Several causes may lead to a relocation of a river at its natural mouth. Prevention of the flood is one of these, especially when heavily populated areas are concerned. This is the case of Guadalhorce River in southern Spain. Nieto-López et al. (2020) described the hydrological and ecological responses to Guadalhorce River mouth channelization, which was concluded in 2003 to mitigate the flood risk in Malaga. Years later, the authors observed how the modification of the natural river configuration caused alterations in the hydrology (i.e. surface water and groundwater) and, consequently, in the biodiversity (i.e. flora and fauna) of the area.

Even land exploitation can lead to river mouth relocation. In 1956, the Kaituna River in New Zealand was diverted at its mouth to make the estuarine areas drainable and prone to be farmed. Nevertheless, the loss of freshwater flows, together with the increase in sedimentation rate, undermined the estuarine ecological health (Mawer 2012). Locals began to complain, motivated by the will to restore the natural estuarine ecosystem. Their claims were heard: the authorities developed a plan to increase the freshwater flow into the estuary by 2018, maximizing the ecological and cultural benefits, while mitigating the adverse environmental effects. The project was completed in 2020 (Barrett et al. 2021).

Another important issue to be considered is the interference that may occur between river mouths and nearby harbors, which materializes in the siltation phenomenon. It is a frequent problem for several ports worldwide (e.g. Winterwerp 2005; Kuijper et al. 2005; Di Risio et al. 2017b; Lisi et al. 2019), indeed, “silt in harbor basins is a problem that exists as long as harbors exist” (Van Rijn 2005). In this context, if harbors are deployed close to river mouths, the siltation phenomenon could be enhanced, as it can be seen as a source of sediment supply (e.g. Truong and Tanaka 2007). According to the environmental conditions and the access channel configuration, the amount of siltation may be such that navigation can be compromised, with the consequent risk of downtime. Such a problem may become even more severe in the case of channel harbors, i.e. for harbors placed along river channels. This is the case of the port of Pescara, a city in the Abruzzo region, central Italy. It lies along the Adriatic Sea at the mouth of the Pescara River, along which the harbor is placed (see Fig. 1).

Fig. 1
figure 1

Location of Pescara municipality. Source Google Earth

Full size image

In the 1990s, a rubble mound breakwater was deployed in front of the river mouth to counteract wave action (e.g. Celli et al. 2018, 2020; Fischione et al. 2022), facilitating the harbor entrance. Nevertheless, the breakwater turned out to be an obstacle for the river outflow (e.g. Lalli et al. 2001; Gallerano et al. 2020), making the water quality significantly worse due to the river polluting load (see Fig. 2).

Fig. 2
figure 2

Copernicus Sentinel 2 data (2nd December 2016) retrieved from Copernicus Scihub: detail of the suspended sediments at the Pescara River mouth

Full size image

To preserve port operativity, dredging activities were carried out, a commonly employed procedure, especially for harbors situated near river mouths (e.g. Van Schijndel and Kranenburg 1998). However, these operations can be costly and may have adverse environmental effects due to the potential re-suspension of contaminated sediments (e.g. Di Risio et al. 2017b; Lisi et al. 2019).

Basically, dredging should be seen as a maintenance operation. However, sustainable solutions can be identified. The modification of the harbor configuration is one of these. For instance, Romdani et al. (2022) showed how the extension of the North Jetty at Mailiao Port (Taiwan) can mitigate the siltation problems related to the sediments originating from the nearby Zhuoshui River.

Alternatively, river mouth deviations in combination (or not) with new harbor configurations, can be carried out. This is what is going to be realized at the port of Pescara (see Fig. 3). On the one side, the re-organization of the internal layout, with the relocation of the port to a new site, is supposed to be a suitable solution to significantly reduce the siltation issues. On the other, with the Pescara River mouth deviation, the water quality is expected to significantly improve.

Fig. 3
figure 3

Comparison between the old and the future Pescara port layouts. Adaptation from the Pescara port masterplan 2008

Full size image

Deviating a river mouth is however a highly demanding challenge. Different from a river diversion, in which the river flow is redirected into a new channel before returning to the original channel downstream (possibly to the mouth, e.g. Flatley et al. 2018), the river mouth deviation involves the nearshore hydrodynamics effects too. In this context, work activities can involve variations of river banks or coastal structures. Any layout modification may affect river hydrodynamics and the propagation of waves, thereby influencing nearshore hydrodynamics. The phenomena are intimately related to each other: wave-driven nearshore hydrodynamics affects the sea level that represents the downstream boundary condition of the river flow.

Furthermore, due to the large-scale of these projects, river mouth deviations often require a significant amount of time to be completed, typically spanning several years. Hence, the intermediate configurations may persist over time.

Then, the main objective of this paper is to highlight the importance of simulating potential flood scenarios that may occur during the stages of river mouth deviation. Additionally, it presents a purpose-built modeling framework, the adoption of which can assist engineers and regulatory authorities in planning and managing construction activities. The paper is organized as follows: Sect. 2 details the modeling framework. Then, its application is described in Sect. 3: the case study of Pescara River mouth deviation is addressed in its implementation stages, to catch the effects of the intermediate configurations on the river hydrodynamics, of which the flood hazard is a function of. Section 4 presents the conclusions.

2 The proposed framework

Since a river mouth deviation typically spans a significant duration (e.g. Flatley et al. 2018), it is crucial to analyze both the final configuration and the construction stages, as they can have implications in various areas (even if only the effects on the river hydrodynamics are considered herein).

The realization steps may influence the river flow in such a way as to cause flooding (e.g. Yanmaz 2000). Then, the flood hazard evaluation and risk assessment (e.g. Di Risio et al. 2017a; Li et al. 2020) have to be carried out (beyond the scope of this manuscript).

Therefore, how the water discharge flows within intermediate configurations should be assessed.

The sea level evolution is affected by the nearshore hydrodynamics, which in turn depends on wave propagation. Acting on both phenomena, even coastal structures (emerged or submerged) could influence the evolution of sea levels (e.g. Marini et al. 2022). Basically, all the phenomena considered herein are physically-coupled. The wave motion forces the nearshore circulation (sea levels and water currents), which, however, influences the wave propagation. The sea levels affect the river flow, which, however, modifies the levels at the outlet.

Although the interaction between coastal and fluvial hydrodynamics at the mouth would require the analysis of the sea-river interplay (Melito et al. 2020), it is worth noting that in this context, the focus is on the possible flooding areas and not on the hydrodynamics at the mouth. Then, as usual for the estimation of floodable areas along coastal river reaches (e.g. Christian et al. 2015; Feng and Brubaker 2016; Kumbier et al. 2018; Sopelana et al. 2018; Santiago-Collazo et al. 2019), the involved phenomena can be analyzed by adopting a one-way coupled approach. Numerical models should be employed, resulting in a modeling chain (see Fig. 4). It should be stressed that its definition does not rely on the selection of specific numerical models, whose choice is arbitrary.

Fig. 4
figure 4

The purpose-built modeling framework

Full size image

In the first phase of the chain, the aim is to define the wave features at the nearshore boundary. Then, the use of a large-scale wave propagation model is suggested, by neglecting the interaction between the wave motion and the nearshore circulation. The one-way approach may be assumed as reliable due to the different spatial scales of large-scale wave propagation and the nearshore circulation.

Within the second phase, the goal is to assess the sea level variation due to the nearshore circulation (e.g. Pasquali et al. 2015, 2019), to be superposed to the astronomic tide components. Then, the use of a small-scale wave propagation and forced nearshore circulation model (up to the mouth) is suggested, by neglecting the interaction with the river flow.

In the last phase, the sea levels should be considered as the downstream boundary condition for detecting the hydrodynamics features of the river flow and the river flood evolution, by neglecting the nearshore circulation coupling but the given downstream boundary conditions. The one-way coupled approach is based on the definition of the downstream boundary condition related to the nearshore circulation. Even if approximated, it may be considered reliable due to the different velocity intensities, with the river flow currents larger than the wave-driven kinematic field in the nearshore.

The obtained results can be used to evaluate the effects of the mouth deviation on the river flow during its stages of realization (i.e. for different configurations of the river mouth), as well as the possible flooding events. As usual, flood maps are generally obtained by neglecting sediment transport and morphological changes (e.g. Manfreda et al. 2015; Chen et al. 2021; Khojeh et al. 2022). This is mainly due to the shortage of sediment transport data and the lack of standardized approaches for its integration within river flood analysis (Vázquez-Tarrío et al. 2023). Even in this work, the sediment transport phenomena have not been considered within the river flood analysis. Despite this, it is worth mentioning that sediment transport may influence the channel morphology (e.g. Baldoni et al. 2021, 2022) during floods and should therefore be considered in flood hazard analysis (e.g. Vázquez-Tarrío et al. 2023; Nones and Guo 2023), whenever possible.

The wave features at the small-scale area boundary can be evaluated by using a large-scale wave propagation model, once defined an offshore boundary. Just as an example, the SWAN (Simulating Waves Nearshore) numerical model (e.g. Booij et al. 1997; Zijlema 2010) can be employed. It solves the wave action balance equation (SWAN team 2021) for each angular frequency and for each direction.

Once defined the characteristics of the waves at the small-scale area boundary, a circulation model useful to evaluate the features of the wave-induced nearshore hydrodynamics (e.g. water levels) inside the area affecting the river flow, should be employed. For instance, the XBeach numerical model (Roelvink et al. 2010) can be used. This model is commonly used in coastal engineering problems involving wave propagation processes within the nearshore area (e.g. Vousdoukas et al. 2012; Jamal et al. 2014; Saponieri et al. 2018). The steady state wave-induced hydrodynamics is reproduced by solving the non-linear shallow water equations. To account for the wave-induced mass-flux and the subsequent (return) flow, these are cast into a depth-averaged Generalized Lagrangian Mean (GLM) formulation (e.g. Andrews and McIntyre 1978).

The final stage of the modeling chain aims to reach the goal of the analysis: to assess the hydrodynamic characteristics of river flow and the evolution of river floods based on predefined upstream conditions (i.e., water discharge) and downstream conditions (i.e., water level). For this purpose, the numerical model BASEMENT (Vetsch et al. 2017) can be used.

In a 2DH domain, this model solves the non-linear shallow water equations, i.e. by neglecting non-hydrostatic pressure terms, and hence vertical velocities and accelerations. The role of turbulence is taken into account via the eddy viscosity model, introducing further flow resistance due to turbulent fluctuations. In order to numerically solve the system of equations, it is necessary to define the closure hypothesis for the bottom shear stress as well as the depth-averaged viscous and turbulent stresses.

3 The case study of Pescara River mouth deviation

3.1 Problem statement

To decrease the intensity of the siltation phenomenon and to significantly improve the water quality at the Pescara channel port, the Pescara River mouth is going to be deviated according to the Pescara port masterplan.

The design configuration relocates the river mouth at the rubble mound breakwater location (see Fig. 3), where the water depth is about 10 m. The new configuration connects with the river channel, where the port channel is now located, with a river channel width equal to about 60 m and a water depth of 4 m. Proceeding downstream, the cross-sections widen, with a maximum width equal to about 65 m and a water depth equal to 5 m. To date, the achieved steps toward the completion of the Pescara River mouth deviation, involved (see Fig. 5):

  • the opening of a breach in the rubble mound breakwater located in front of the actual mouth;

  • the realization of a short groin at the opening to improve the outflow and shelter from wave action;

  • the realization of a temporary groin (initially submerged, see upper right panel of Fig. 5) that will be the starting point for the new left bank.

During this intermediate configuration, which involves partial dismantling of the breakwater and temporary groin construction, the nearshore circulation is likely to experience variations. Hence, it could induce effects on the sea level at the actual river outflow.

To evaluate how these effects could influence river hydrodynamics, the modeling chain illustrated in Fig. 4 has been implemented. In this regard, three different configurations have been considered (see Fig. 5). The first one, to be intended as a reference (configuration A0), is the initial mouth configuration in 2018, before the start of the works. The second one is the configuration in 2020 (configuration A1) characterized by the partial opening of the breakwater, the realization of both the small groin at the opening, and the temporarily submerged groin on the left bank (crest elevation equal to \(-\)0.50 m asl). In the end, the configuration in 2023 (configuration A2) which differs from the previous one for the temporary groin become emerged (crest elevation equal to + 1.40 m asl) and for the completion of the groin at the opening.

Fig. 5
figure 5

Achieved steps toward the completion of the Pescara River mouth deviation. Comparison over the years. Source Google Earth

Full size image

3.2 Computational domains

SWAN numerical model (Booij et al. 1997) has been employed to evaluate the wave features within the large-scale domain by using the finite element method (FEM). The whole domain has been discretized in triangular elements. Each element is characterized by the position of the nodes at its vertices. According to the extension of the investigated area, a variable mesh resolution has been implemented. The triangles area has ranged between \(5\cdot 10^5\) \(\hbox {m}^2\) (offshore) and 100 \(\hbox {m}^2\) (near the breakwater). Then, 47377 triangular elements and 24097 nodes have been adopted. The numerical grid is 60 km long alongshore, whilst it extends 27 km in the cross-shore direction. The whole grid covers an area of 1270 \(\hbox {km}^2\) with the offshore water depth reaching a value greater than 100 m. Figure 6 shows the implemented numerical grid. It should be noted that only one numerical grid (i.e. for all the analyzed configurations) has been implemented for the large-scale wave propagation.

Fig. 6
figure 6

Computational grid implemented in SWAN within the large-scale wave propagation. The white dot represents the location where the wave features have been extracted to force the nearshore circulation

Full size image

Further computational grids have been implemented for XBeach simulations, according to the finite difference method, one for each of the tested configurations (i.e. A0, A1, and A2). To catch detailed results at the river mouth area, a variable grid resolution has been employed. Toward offshore and at lateral boundaries, the spatial resolution is equal to 100 m. It reduces toward the coast, reaching 5 m at the river mouth area. Then, 128700 computational nodes have been adopted. The grid extends for 1700 m cross-shore. Along the coast, the grid is 4200 m long. The whole grid covers an area of 7 \(\hbox {km}^2\) and the water depth at the offshore boundary is about 10 m. As a reference, Fig. 7 shows the implemented numerical grid concerning configuration A0 (similar to the other grids, but the configuration of the structures).

Fig. 7
figure 7

Computational grid implemented in XBeach for the configuration A0 (similar to the other grids, but the configuration of the structures). Within the dashed green area, the mean water level has been computed to be applied as downstream boundary condition for the river flow modeling

Full size image

As for the computation of the large-scale wave propagation, the river flow hydrodynamics has been computed (through BASEMENT) by implementing computational grids (one for each of the analyzed configurations) with triangular elements. The computational domain extends 2 km upstream from the river mouth and it includes part of the Pescara downtown and a part of the open sea.

Concerning the probable flooding zone extension, a maximum area of 100 \(\hbox {m}^2\) has been adopted for the mesh elements. Within the main channel and at the river mouth, a maximum area of 50 \(\hbox {m}^2\) has been considered. Then, 63467 triangular elements and 32259 computational nodes have been implemented. With the aim of catching preferential flow paths, the hydraulic characteristics of the main city roads have been also reproduced. The whole grid covers an area of 4.8 \(\hbox {km}^2\). About 2.3 \(\hbox {km}^2\) refers to the river mouth area, whilst about 2.5 \(\hbox {km}^2\) refers to the city area. It should be emphasized that the computational domain extension has been selected large enough to prevent the flow from reaching the boundaries.

As a reference, Fig. 8 shows the implemented numerical grid for configuration A2.

Fig. 8
figure 8

Computational grid implemented in BASEMENT for the configuration A2. Colors represent different bottom roughness. The downstream boundary condition (in terms of water levels) has been applied along the external offshore arc. The section where the upstream boundary condition has been applied is also depicted

Full size image

3.3 Boundary conditions

3.3.1 Large-scale wave modeling

At the offshore boundary of the large-scale domain, JONSWAP spectra have been imposed, according to the selected return period along with the wind forcing. Anemometric and wave data have been taken from the ERA5 database, recently developed by ECMWF (European Center for Medium-range Weather Forecasts). The data, available since 1950, have been extracted at a single point, whose location is illustrated in Fig. 9. For the anemometric characterization of the zone, the information related to the (average) speed of the wind and its (average) direction at an elevation of 10 m asl has been taken into account. Regarding the wave data, the significant wave height \(H_{s}\) (equal to \(H_{m0}\), offshore spectral wave height), the peak period \(T_p\) and the wave direction \(\theta _m\) have been analyzed.

Fig. 9
figure 9

Offshore point location from which anemometric and wave data have been extracted from ERA5 database, and the Ortona wave-buoy location

Full size image

In particular, the wave heights have been firstly calibrated through a comparison with the data recorded by the wave buoy (belonging to the National Wave Network, managed by ISPRA), installed offshore Ortona (see Fig. 9), whose data are available from 1989 to 2014. Extreme Value Analysis (EVA) has been performed for both wind speeds and significant wave heights. To ensure statistical independence, the extreme events have been identified by the POT (Peak Over Threshold) method. The wind speed and significant wave height excesses over the threshold have been considered following the GPD (Generalized Pareto Distribution, e.g. Coles et al. 2001) probability distribution function. A GEV (Generalized Extreme Value Distribution, e.g. Coles et al. 2001; Celli et al. 2021) probability distribution has been also considered. The statistical inference on both the probability distribution functions (i.e. for wind speed and significant wave height) allowed identifying the values of the stochastic variables as a function of the return period \(T_R\). The peak periods associated with the calculated significant wave heights have been determined according to the relationship \(T_p = aH_s^b\) (e.g. Goda 2010), where the parameters a and b have been obtained evaluating the statistical correlation between the involved variables (\(a=5.44\) and \(b=0.34\)). Similar wave parameters (\(H_s\), \(T_p\)) have been obtained calculating the environmental contours (e.g. Eckert-Gallup et al. 2016; Vanem 2019), according to selected return periods (see Fig. 10).

The wave directions have been selected in order to consider the worst condition in terms of wave penetration through the breakwater breach (existing in configurations A1 and A2). Table 1 shows the wave parameters which the JONSWAP spectrum is a function of (a peak enhancement factor equal to 3.3 has been considered) together with the wind speed U\(_W\).

Fig. 10
figure 10

100 years environmental contour and the adopted 100 years sea state condition

Full size image
Table 1 Wave features and wind speed at the offshore boundary of the large-scale domain, according to selected return periods (\(T_R\))
Full size table

Along with the boundary conditions, even the (static) water levels have to be specified. In this regard, the effects of long-term mean sea level associated to climate change may be taken into account if the expected duration of the intermediate river mouth configuration is significantly high. In this case, since the temporary nature of the considered short-term river mouth configurations, long-term mean sea-level changes have not been considered.

Then, the tide components have been estimated on the basis of the data collected by the tide gauge installed in Ortona (20 km south of Pescara). The astronomical components are illustrated in Table 2. The meteorological components (i.e. residual levels, see Table 3) have been subjected to EVA. To ensure their statistical independence, the POT (Peak Over Threshold) method has been used. The residual levels excesses over the threshold have been considered following the GPD (e.g. Coles et al. 2001). Then, the statistical inference on the probability distribution function allowed identifying the values of the residual levels as a function of the return period \(T_R\), as synthesized in Table 3.

Table 2 Characteristics astronomical tide levels
Full size table
Table 3 Meteorological tides and the associated return periods
Full size table

The superposition of the Mean High Water Springs (MHWS, + 0.19 m asl) and the meteorological tides with the associated return periods have then been imposed as static water levels.

3.3.2 Nearshore hydrodynamics modeling

Within the nearshore modeling, at the offshore boundary of the small-scale domain, the results achieved by using SWAN can be applied. Basically, the wave features in the nearshore area of the large-scale domain represent the boundary conditions (i.e. forcing term with a given return period) for the nearshore circulation model. In this context, the XBeach model has been implemented by choosing a stationary wave model, including wave propagation, directional spreading, shoaling, refraction, bottom dissipation and wave breaking. To achieve the wave field, at the lateral boundaries, the gradient in the wave energy along the wave crest can be set to zero as a boundary condition. As far as the hydrodynamics is concerned, at the offshore boundary, the radiation condition has been applied, hence no net water flux can come into the model domain. To avoid the influence of lateral boundaries, the Neumann condition has been applied so that no local change in surface elevation and velocity takes place (i.e. the spatial derivative along the normal to the boundaries is set to zero).

As for the large-scale wave propagation model, even for the nearshore hydrodynamics modelling, the same static water levels have been imposed (i.e. MHWS, +0.19 m asl, and the meteorological tides with the associated \(T_R\)).

3.3.3 River modeling

In BASEMENT, the water levels obtained from the nearshore simulations serve as the downstream boundary condition (i.e., outflow) for the river flow modeling at the downstream boundary. It should be underlined that the sea level obtained by the nearshore numerical model has not been directly applied at the river mouth to correctly reproduce the influence of the river flow on the local water levels. In particular, for each simulation, the water level has been uniformly applied along the offshore arc of Fig. 8. It represents the mean value of the water levels computed within the green dashed area of Fig. 7 for the steady solution (area of river outflow). At the upstream boundary, the temporal variations of the water discharge (i.e. hydrographs) have been imposed for the river flow investigation. The results provided by the regional planning have been then used. The peak water discharges \(Q(T_R)\), listed in Table 4 as a function of the return period \(T_R\), have been estimated according to the VAPI procedure, a regional methodology based on two-component extreme value (TCEV) probability distribution (e.g. Rossi and Villani 1994). A two-parameter Gamma distribution has been used to calculate the synthetic hydrographs, depicted in Fig. 11. Once known the peak water discharges \(Q(T_R)\) and the time of concentration of the Pescara River watershed basin, the generic synthetic hydrograph q(t) is given by (e.g. Croley II 1980):

$$\begin{aligned} q(t) = \frac{Q(T_R)}{(\alpha -1)^{\alpha -1}\exp \left[ {-(\alpha -1)}\right] }\left( \frac{t}{\beta }\right) ^{(\alpha -1)}\exp {\left( \frac{-t}{\beta }\right) } \end{aligned}$$
(1)

where \(\alpha\) and \(\beta\) are the gamma function parameters, iteratively computable as reported in Croley II (1980).

Fig. 11
figure 11

River upstream boundary conditions according to the selected return period

Full size image
Table 4 Upstream water discharges and the associated return periods
Full size table

Concerning the hydraulic roughness, the following Manning’s coefficient values have been adopted (e.g. Chow 1959; Jens et al. 1979):

  • 0.035 s/\(\hbox {m}^{1/3}\) for the river channel and the seabed;

  • 0.060 s/\(\hbox {m}^{1/3}\) for the floodplains and the roads;

  • 0.50 s/\(\hbox {m}^{1/3}\) for the urban areas (buildings).

It should be stressed that the effects of the buildings on the hydrodynamics have been taken into account by means of high Manning’s coefficient value (Beretta et al. 2018) and the role of moored boats on the river flow (Sammarco and Di Risio 2017) has been neglected. Lateral boundaries has been modeled as impervious walls.

3.3.4 Selection of the return period

To sum up, the wave propagation at large-scale and the nearshore circulation require:

  • the definition of the wave features at the offshore boundary of the large-scale domain for given return periods;

  • the definition of the astronomical component of tide levels;

  • the definition of the meteorological component of tide levels for given return periods.

On the other hand, the evaluation of the river flow hydrodynamics requires:

  • the definition of flow discharge at the upstream boundary, according to given return periods;

  • the definition of the downstream boundary condition.

The selection of an appropriate return period is the last issue to be addressed. Within a flood hazard evaluation, the scenarios to be investigated depend on the selected probability of occurrence, i.e. on the return period, which the boundary conditions are function of. In this specific case, both hydrological and met-ocean variables are involved (e.g. Orton et al. 2020). Then, evaluating if the sea levels at the mouth and the river discharge are two statistically correlated variables, is paramount.

For the problem at hand, the nearshore circulation is influenced by the meteorological perturbation occurring in the Adriatic Sea basin. On the other hand, the river flow is induced by perturbations acting on the wide Pescara River watershed basin. For their different extents, the considered perturbations are likely to be not physically correlated, as associated to different meteorological genesis. This would mean a low probability of simultaneous occurrence of river discharge and sea level (meteorological component) with high return periods.

To support this hypothesis, a joint probability analysis of extreme events should be carried out (e.g. Ward et al. 2018; Sopelana et al. 2018; Couasnon et al. 2020. Such analysis requires an extensive set of high-quality observations.

However, the available river discharge data for the current case are limited.

Then, by following the approach of Pasquier et al. (2019), a comparison of the available data on past storm surges and discharge has been carried out.

In particular, the maximum daily river discharge values (from 1999 to 2005) have been compared to the sea level data. The last have been collected by the tide gauge deployed in Ortona (20 km south of Pescara) belonging to the National Mareographic Network. The tidal time series have been analyzed employing standard harmonic analysis (e.g. Codiga 2011). The residual tide has been extracted by assuming the linear superposition of deterministic and stochastic components (e.g. Pasquali et al. 2019). The maximum residual tide values within a time window spanning from the day before to the day after the peak river discharge events have been subsequently taken into account.

Figure 12 compares the observed daily water discharge along the Pescara River with the simultaneous storm surges evaluated from the sea levels collected at Ortona. Such a comparison relies on a limited available daily discharge dataset. Nevertheless, large values of water discharge are still observable (roughly equal to 500 \(\hbox {m}^3\)/s), to which small sea levels correspond (less than \(z_w\) = +0.75 m, characterized by a return period of 2 years). Basically, intense river discharges and intense storm surge events are not synchronous, at least for almost all data. Then, the selection of the return periods has been carried out assuming the nearshore hydrodynamics and the rainfall-runoff transformation process as uncorrelated phenomena. This reasonable assumption, although not supported by a joint probability analysis, has already been adopted in other cases (e.g. Klerk et al. 2015; Webster et al. 2014; Couasnon et al. 2020)

Then, the offshore wave features and the nearshore hydrodynamics are evaluated according to return periods ranging from 2 years (high probability of occurrence) up to 100 years (design condition). Concerning the river flow, consistent with the regional law for flood hazard assessment, the following return periods have been considered: 50, 100, and 200 years.

Fig. 12
figure 12

Observed daily water discharge along the Pescara River and the simultaneous storm surges evaluated from the sea levels collected at Ortona

Full size image

3.4 Scenarios

For the sake of synthesis, Table 5 synthesizes the parameters adopted within the large-scale wave propagation model and the nearshore circulation simulations. In the Table, the value \(z_w\) represents the static water level, equal to the superposition of the Mean High Water Springs (MHWS, equal to +0.19 m asl) and meteorological tides with the associated return periods; \(U_W\) is the wind speed, \(\Theta\) is the wave and wind direction. In this context, it should be specified that the considered short-term river mouth configurations may persist over time, but they still remain temporary configurations. For this reason, long-term mean sea-level changes have not been taken into account.

Table 5 Scenarios for the large-scale wave propagation model and the nearshore circulation, according to selected return periods (\(T_R\))
Full size table

In the same way, Table 6 illustrates the parameters adopted for the river hydrodynamics numerical simulations.

Since the low probability of simultaneous occurrence of high return period river discharge and storm surge events (as discussed in Sect. 3.3), a return period ranging from 50 to 200 years for the river discharge has been considered. For the storm surge, a return period of 2 years (\(z_w\) = +0.75 m) has been considered.

Table 6 Scenarios for the river flow modeling according to the selected return period
Full size table

3.5 Numerical results

The main aim of the large-scale wave propagation simulations is to estimate the wave features just offshore the breakwater, representing the forcing of the nearshore circulation. Table 7 shows the calculated values for a point where the water depth is 10 m (see the white dot in Fig. 6).

To evaluate the effects of the different configurations (i.e. A0, A1, and A2) on the nearshore circulation at the mouth, Figs. 13 and 14 illustrate the wave-induced currents. By looking at Fig. 13, the sea states coming from N-E are depicted (see scenario M002.1 in Table 5), inducing currents mainly directed from North-East to South-East. In the upper panel, where configuration A0 is analyzed, the confinement effect provided by the breakwater could be appreciated. It materializes with an increase of the velocity magnitude (about 1.4 m/s) at the southern section of the breakwater, i.e. at the port channel entrance. In the middle panel, configuration A1 is depicted. An increased velocity at the submerged groin could be observed due to the locally limited water depth (about 0.83 m). In contrast with configuration A0, the velocity field at the southern section of the breakwater decreases due to the breach which reduces the confinement effect (the water current is partially diverted northward). In the lower panel, configuration A2 is depicted. The intensity and the extension of the water current are influenced by its interaction with the emerged groin. High velocities could be observed in the gap between the groin and the breakwater. The velocities are very low near the mouth. Even in this configuration, the water current is partially diverted through the breach. For return periods (of sea states) starting from 100 years (not shown in the Figure), the currents return to being directed southward: the water levels are high enough to make the groin submerged. In Fig. 14, the currents coming from the East are depicted (scenario M010.2 in Table 5). For configurations A0 and A1 (upper and middle panel) the velocity fields approach the mouth entering at the southern section of the breakwater. Progressively, they both tend to spread. Also in this case, the effect of the submerged groin in increasing the velocity field within configuration A1 is clearly reproduced (middle panel). The emerged groin again influences the velocity field in configuration A2 (lower panel), by forcing the current to flow through the gap between the groin and the breakwater. Basically, the emerged groin limits the northward drainage, increasing the sea levels at the mouth. For return periods greater than 50 years (not shown in the Figure), the water levels are such to make the groin submerged and the effects on the sea level are similar to the ones in configurations A0 and A1.

It is worth recalling that the nearshore circulation simulations have been carried out to define the downstream boundary conditions for the river hydrodynamics simulations. Table 8 shows the mean sea levels at the terminal section of the river, for each of the considered configurations. Since the sea levels and the water discharges are two statistically independent variables for the problem at hand, within the river hydrodynamics simulations, the downstream boundary conditions are represented hereinafter by mean sea levels with \(T_R\) equal to 2 years (i.e. M002.1 and M002.2 for each configuration–see Sect. 3.4).

Table 7 Wave features at a point located just offshore the breakwater, where the water depth is 10 m
Full size table
Table 8 Nearshore circulation: mean sea levels at the terminal section of the river for each of the considered configurations
Full size table
Fig. 13
figure 13

Nearshore circulation velocity fields: comparison between the different configurations for the scenario M002.1 (\(T_R\) = 2 years - waves coming from N-E)

Full size image
Fig. 14
figure 14

Nearshore circulation velocity fields. Comparison between the different configurations for the scenario M010.2 (\(T_R\) = 10 years - waves coming from East)

Full size image

The river hydrodynamic simulations have been carried out to reach the goal of the analysis: to evaluate the effects of the river mouth execution steps on the potentially flooded area within the domain of interest. This forms the basis for any possible evaluation of the flood hazard and risk. According to the duration of the selected hydrographs (see Fig. 11), 36 h of river flow hydrodynamics have been simulated for each configuration. Then, the temporary evolution of the flooded areas has been examined. In Fig. 15, the water depth values after 16 h are depicted, for the configuration A1 in Scenario F200 - M002.1 (see Table 6). Similar results are achieved for all the configurations. This is confirmed by the maximum flooded area values, reported in Table 9 for the three configurations. Therefore, the implementation of the river mouth execution steps do not significantly impact the extent of flooded areas. As expected, the extension of the flooded area increases as the hydrograph return period grows. Nevertheless, for the same scenario, similar values are achieved for the three configurations, with a negligible decrease concerning configuration A2 (i.e. last step). This can be related to the slight variations in the downstream boundary conditions. For configurations A0 and A1, the mean sea levels are very similar. For configuration A2, the mean sea levels show a difference that does not exceed the value of + 0.15 m. As a result, the flooded areas do not change between the selected configurations, except for a small zone located near the emerged groin (see Fig. 16). Although marginally, the emerged groin influences the flooding on the left bank in configuration A2. However, the slight differences show up only on the beach. No difference can be appreciated in terms of flow velocity as well.

Table 9 Maximum flooded areas according to the simulated scenario, for each of the considered configurations
Full size table
Fig. 15
figure 15

River hydrodynamics: water depth values after 16 h. Scenario F200 - M002.1 (see Table 6) for the configuration A1

Full size image
Fig. 16
figure 16

River hydrodynamics: river flow velocity values after 16 h. Comparison between the different configurations for the scenario F200 - M002.1 (see Table 6)

Full size image

4 Concluding remarks

The present work aims to illustrate the importance of simulating potential scenarios of river flood that may occur during the execution of a river mouth deviation. This complex process is generally completed over a relatively long time frame, characterized by stepwise time-limited river mouth layouts. Then, evaluating how such short-term (i.e. intermediate) configurations could influence the river hydraulics, could be of interest, especially in terms of potential variation of the flood hazard in coastal urban areas. Since the short-term nature of the considered configurations, the influence of the morphological evolution of the site on the river hydraulics has not been considered.

To achieve the goal, a purpose-built modeling framework has been adopted and applied to the Pescara River mouth deviation, considering three different consecutive mouth configurations.

At first, the wave features at the nearshore boundary have been estimated by means of a large-scale wave propagation model, by neglecting the interaction between the wave motion and the nearshore circulation. These represent the forcing of the nearshore circulation and do not vary between the different configurations of the river mouth.

Within the second phase, the effects of the different mouth configurations on the nearshore circulation have been analyzed by a 2DH numerical model. The definition of the boundary conditions is based on the results of the first phase. The interaction with the river flow is considered negligible.

In the last phase, the sea levels estimated in the previous phase have been imposed as the downstream boundary condition for detecting the river flood evolution. Hence, the coupling with the nearshore circulation has been neglected.

For the case at hand, the sea levels at the mouth, i.e. the downstream boundary condition for the river flow, do not significantly change between the considered configurations. Basically, it comes out that the flooded areas do not vary between the considered configurations. Changes might still take place during the following stage of realization, eventually evaluable with the same rationale.