The Expected Dynamics for the Extreme Wind and Wave Conditions at the Mouths of the Danube River in Connection with the Navigation Hazards

1. Introduction

The climate dynamics of the last decades indicates that the climate changes affect in a significant way most of the natural systems and this especially concerns the global warming associated with higher temperature variations [1]. Furthermore, various climate scenarios (see for example [2,3]) indicate that a negative impact of climate change is expected in almost all regions of the world inducing a visible enhancement of the intensity and frequency of the extreme events [4]. Coastal environments are more sensitive to such extreme events [5]. They represent the interface between land and water and, on the other hand, the wind is usually stronger in the marine environment [6] and that is why its impact on the land located in the vicinity of the sea is higher.

From this perspective, the objective of the present work is to provide a more comprehensive picture concerning the expected dynamics of the extreme wind and wave conditions at the mouths of the Danube River in the Black Sea in connection with the navigation hazards. The focus is on the entrance in Sulina channel, which represents the main navigation exit from the Danube to the Black Sea [7], being at the same time the southern gate of the seventh Pan-European transportation corridor linking the Black Sea to the North Sea via the Rhine-Maine-Danube navigation system (inaugurated in 1992) [8]. The location of the target area in the Black Sea is illustrated in Figure 1. Since this represents the main southern entrance in the largest European inland navigation system it is subjected to high navigation traffic [9]. Furthermore, for a distance of about 170 kilometers, from Sulina at the Black Sea (considered kilometer zero) and up to the Romanian port of Braila (kilometer 170), the so called Maritime Danube extends and the Danube is also navigable by maritime ships [10], while river ships can navigate upstream to Ulm, in Germany [11]. Moreover, since the beginning of the war in Ukraine in February 2022, this traffic was significantly increased [12], the Danube being an alternative and more secure route for the cereals and other goods exported by Ukraine.

Since the most important water resource for the Black Sea is represented by the Danube [13], the river outflow generates relatively strong currents in the coastal environment neighboring the deltaic zone, which is the largest in Europe [14]. The delta is delimited by the three principal branches of the Danube River, Chilia, Sulina (which represents the main navigation waterway) and Saint George [15], and there are often risks of floods [16,17] or various other hazards [18]. The wind and wave climate are characterized through relatively high energy conditions [19], and extreme storms also occur with a certain frequency in this nearshore [20]. Furthermore, the waves are often increased at the mouths of the Danube by the wave-current interactions that occur between the incoming waves and the coastal currents induced by the river discharge into the sea [21]. Taking into account the complexity of this coastal environment and the processes associated, the relatively high frequency of the severe weather conditions and also the quite high navigation traffic and the proximity of an outgoing war, various studies (see for example [22]) indicate this area as being subjected to high risks from the point of view of the navigation hazards.

In this context, an analysis of the characteristics of the environmental matrix and its expected dynamics in the coastal area at the mouths of the Danube River, and especially at the entrance in Sulina channel (considered as the zero kilometer of the Danube River), is made in the present work. This analysis includes first some in situ measurements for the wind and water level, reanalysis wind data from the European Centre for Medium-Range Weather Forecast (ECMWF), and projections of the future expected wind fields provided by regional climate models (RCMs) under various RCP (Representative Concentration Pathway) and SSP (Shared Socioeconomic Pathway) scenarios. Furthermore, an analysis of the wave conditions is also made by performing simulations over the entire basin of the sea with SWAN (acronym from Simulating WAves Nearshore) spectral model [23] and focusing the modeling system with higher resolution domains on the coastal environment at the mouths of the Danube River.

2. Materials and Methods

2.1. In Situ Measurements

The bathymetric map of the target area considered in the present study is illustrated in Figure 2. Based on the information presented in [24,25,26] the average currents generated by the Danube River outflow in the Black Sea are also represented on the map. Two types of data were available for analysis performed for the 15-year time interval 2009-2023.

The first parameter analyzed is the water level. These are values resulted from the daily measurements [12] at the entrance in Sulina channel. The water level is considered the elevation of the free surface of the Danube relative to a specified vertical datum and refers to the difference in centimeters between the sampling point' vertical datum in reference to a standard Black Sea elevation datum. The wind measurements, data provided by [27], represent the second parameter analyzed. This includes the wind speed provided with one hour time step and also the corresponding maximum speed of the wind gust. The average wind direction is also provided by the data analyzed. Both the water level and wind measurements correspond to the zero kilometer of the Danube, denoted in Figure 2 as RP0 (reference point zero).

2.2. Wind Model Data Considered

The next analysis is based on data coming from wind models. These are carried out in two different locations, reference point 1 (RP1) offshore Sulina channel with the coordinates (30.18°E, 45.12°N) and reference point 2 (RP2) offshore the Saint George arm of the Danube with the coordinates (30.01°E, 44.82°N). Thus, a first analysis of the wind climate relates to the hindcast period of 30 years (the interval 1976-2005). The ERA5 database is considered for this analysis. This is a relatively new reanalysis wind product being an upgrade of the former ERA-Interim [28]. This data starts in 1950 and ends five days behind the current time. It is provided in the public domain by the European Centre for Medium-Range Weather Forecast and is the most recent and one of the most reliable global atmospheric reanalyses databases. The resolution of the ERA5 wind data is 0.25° × 0.25° in the geographical space and the time step is one hour. For each hour it gives wind and atmospheric pressure fields consistent with the previous evolution of the parameters modelled. In the present work, the wind speed at 10m height above the sea level is considered for the analysis with a 3-hour time step.

The second 30-year time period taken into account for analysis considering various climate scenarios (RCP4.5, RCP8.5, and SSP5-8.5) and atmospheric models is (2041-2070). Thus, the first wind field considered is that given by the regional climate atmospheric model version 4 (RCA4) provided by Rossby Centre of the Swedish Meteorological and Hydrological Institute (SMHI) [29]. This regional climate model is driven by the global climate model (GCM), MPI-M-MPI-ESM-LR [30]. MPI-M means the Max Planck Institute for Meteorology, while MPI-ESM indicates the Max-Planck-Institute Earth System Model and LR low resolution. This data is made available in the public domain via the COPERNICUS database [31]. A second RCM considered is ALADIN (versions 6 for RCPs and 6x for SSP5-8.5) [32]. ALADIN stands for Aire Limitée Adaptation dynamique Développement InterNational and is a regional climate model forced by the global atmospheric model GCM CNRM-CM5 (Centre National de Recherches Météorologiques), which is an Earth system model designed to run climate simulations. This GCM comprises ARPEGE Climate Model [33]. ARPEGE is run by Météo-France and is the acronym from Action de Recherche Petite Echelle Grande Echelle, with the meaning research project on small and large scales. As regards the SSP climate scenarios, the data considered come from the fully coupled regional climate system model CNRM-RCSM6/ALADIN [34,35] (where RCSM is the acronym for regional climate system model). Beside the analysis of the data for the 30-year future time period (2041-2070), comparisons have been performed also with the ERA5 data for the past. Taking into account the availability of the data from the climate wind models, which is different for the SSP scenarios than the data for the RCP scenarios, the 30-year time interval (1979-2008) was considered for SSP, while for RCP (either 4.5 and 8.5) the initial 30-year time interval (1976-2005) is kept.

2.3. Wave Model Simulations

The third generation phase averaged model SWAN [23] is considered in the present work to predict the wave conditions. This was implemented and validated against satellite data and in situ measurements in the entire basin of the Black Sea and further on it was focused considering a higher resolution domain on the coastal area at the mouths of the Danube. This wave model is based on the spectrum concept and solves an advection type equation in fifth dimensions (time, geographical and spectral spaces) [36,37]. As in most of the spectral wave models, the spectrum considered is the action density spectrum (N) and not the energy density spectrum, because in the presence of the currents, the action density is conserved while the energy density is not. The action density is equal to the energy density (E) divided by the relative frequency ($\sigma$). Although initially SWAN was designed for the nearshore areas, during the time the model was developed in a significant way and now can be used successful also for sub-oceanic and even oceanic scales, having performances comparable with the standard large scale wave models as WAM (Wave Model) [38] or WW3 (Wave Watch 3) [39]. As regards the sub oceanic scales, SWAN becomes probably the most appropriate model. From this perspective, SWAN represents at this moment the most effective and reliable numerical model for predicting wave conditions and climate in the Black Sea, which is and enclosed sea basin [40].

The advection-type equation considered in SWAN has the general expression:

$$\frac{\partial N}{\partial t}+\nabla \left(N\right)+\frac{\partial }{\partial \sigma }\dot{\sigma }N+\frac{\partial }{\partial \theta }\dot{\theta }N=\frac{S}{\sigma },$$

(1)

For larger scale applications, the spherical coordinates, longitude $\left(\lambda \right)$ and latitude$\left(\phi \right)$, are considered, and the operator $(\nabla )$has the expression:

$${\nabla }_{\mathrm{S}\mathrm{p}\mathrm{h}}\left(N\right)=\frac{\partial }{\partial \lambda }\dot{\lambda }N+\frac{1}{\cos \phi }\frac{\partial }{\partial \phi }\dot{\phi }N,$$

(2)

For coastal applications, the Cartesian coordinates (x) and (y) are more appropriate and the operator $(\nabla )$ becomes in this case:

$${\nabla }_{\mathrm{C}\mathrm{a}\mathrm{r}\mathrm{t}}\left(N\right)=\frac{\partial }{\partial x}\dot{x}N+\frac{\partial }{\partial y}\dot{y}N,$$

(3)

The left side of equation (1) represents the kinematic part, and indicates the propagation of the wave action in time, geographical and spectral spaces. Some relevant phenomena, as wave diffraction or refraction, are also included. On the right hand side is the source (S) expressed as energy density components. Three components are more relevant in deep water. They correspond to the atmospheric input, dissipation by whitecapping and nonlinear quadruplet interactions. In shallow and intermediate water additional terms are included. They correspond especially to the finite depth effects and comprise phenomena as depth induced wave breaking, bottom friction or triad nonlinear wave-wave interactions. Thus, the expression of the source term is:

$$S=S_{in}+S_{nl}+S_{diss}+S_{fd}$$

(4)

Such a wave prediction system SWAN based have been implemented in the basin of the Black Sea and focused on its western coast considering various computational levels [41,42]. Table 1 presents the characteristics of the three computational domains considered in the present work. The first domain, defined in spherical coordinates (longitude and latitude, Δλ and Δφ), corresponds to the entire basin of the Black Sea (illustrated in Figure 1 and Figure 2). The second computational domain, defined also in spherical coordinates, is focused on the coastal area at the mouths of the Danube River and the corresponding geographical space is illustrated in Figure 2. Finally, a higher resolution domain focused only on Sulina channel, was defined in Cartesian coordinates and the corresponding geographical space is illustrated in Figure 1. Table 1 indicates the resolution in the geographical space (Δλ and Δφ for spherical and Δx and Δy for Cartesian coordinates), Δt is the time resolution, nf the number of frequencies considered in the spectral space, nθ number of directions, ngλ (or ngx) the number of grid points in longitude (or x-axis), ngφ (or ngy) the number of grid points in latitude (or y-axis, respectively) and np total number of grid points.

The physical processes activated in the SWAN simulations, corresponding to each of the three computational domains considered in the present work are presented in Table 2. In this table: Wave indicates the wave forcing, Wind the wind forcing, Td the tide forcing, Cr the current field input, Gen generation by wind, Wc the whitecapping process, Qd –the quadruplet nonlinear interactions (interactions between four waves that usually occur in deep water), Tri the triad nonlinear interactions (reflecting the interactions between three waves, which are usually characteristic to intermediate and shallow water), Dif the diffraction process, Bfr the bottom friction, Set up – the wave induced set up and Br the depth induced wave breaking.

Intensive validations have been performed to assess the performances of this wave prediction system in the Black Sea, indicating that the SWAN results are in general reliable for both average energy and storm conditions. For example, in [43] comparisons against satellite data have been performed for the 20-year period 1997-2016. In terms of significant wave height (Hs) the following values resulted for the main statistical parameters: Bias= 0.05m, RMSE=0.38m, SI=0.36, R=0.87 and S=0.99. The parameters above presented are: Bias, Root Mean Square Error (RMSE), scatter index (SI), correlation coefficient (R), and the regression slope (S), all of them being computed according to their standard definitions, as given below, where X_i represents the measured value, Y_i the simulated value and n the number of observations:

$$X_{med}=\tilde{X}=\frac{\sum _{i=1}^n{X}_i}{n},Bias=\frac{\sum _{i=1}^n\left(X_i-Y_i\right)}{n},RMSE=\sqrt{\frac{\sum _{i=1}^n{\left(X_i-Y_i\right)}^2}{n}}$$

(5)

$$SI=\frac{RMSE}{\tilde{X}},r=\frac{\sum _{i=1}^n\left(X_i-\tilde{X}\right)\left(Y_i-\tilde{Y}\right)}{{\left(\sum _{i=1}^n{\left(X_i-\tilde{X}\right)}^2\sum _{i=1}^n{\left(Y_i-\tilde{Y}\right)}^2\right)}^{\frac{1}{2}}},S=\sqrt{\frac{\sum _{i=1}^n{Y_i}^2}{\sum _{i=1}^n{X_i}^2}}$$

(6)

It has to be highlighted at this point that the values of the Bias, RMSE and SI are better when they are smaller, while the correlation coefficient and the regression slope are better when they are closer to the unity. For the case of the waves with significant wave heights higher than 3m (associated with the storm conditions), for the 20-year time interval considered the following values resulted for the main statistical parameters: Bias= -0.07m, RMSE=0.47m, SI=0.21, R=0.78 and S=1.05. The simulation results show also that, in the Black Sea basin, the storm waves represent about 2% of the total. These statistical results can be considered as reasonable accurate and they give a good degree of credibility to the analyses related to the past and future storm conditions in the Black Sea, which will be presented next. Furthermore, an assimilation scheme of the satellite data, based on an optimal interpolation approach, has been implemented [44] improving the statistical results in relative terms with (Bias 57%, RMSE 17%, SI 20%, R 3% and S 1%). It has to be highlighted also that the above results were performed when the wave modeling system was forced with ERA5 reanalysis data, while for the future wave projections RCM wind data have been used [45] considering the SWAN version 41.31AB.

3. Results

3.1. Analysis of the In Situ Measurements

An analysis of some in situ measurements performed at the entrance in Sulina channel (the zero kilometer of the Danube, location denoted in Figure 1 as RP0) is performed next for the 15-year time interval 2009-2023. The first parameter considered is the water level and Figure 3 illustrates the variation of the water level at the entrance in Sulina channel (RP0). The data processed represent the results of the daily measurements carried out for the 15-year time interval 2009-2023. The results show that the maximum value of the water level measured in this time interval is 125cm, the minimum is 13cm and the average value is 71.6cm.

The next analysis relates the wind speed measured in the same reference point (RP0). The data processed provides hourly values of the wind speed (Uw), maximum value of the wind gust (Uwg) registered in the respective hour, and the mean wind direction. Thus, Figure 4 presents the wind roses for the two parameters, Uw and Uwg. Figure 5 illustrates the annual maximum series (bars) for both wind parameters considered (Uw and Uwg), corresponding to the 15-year time interval considered. In this figure, the linear regression (indicating the trend) is also illustrated. A trend line, or the line of best fit [46], is given by the relationship:

y=sx+m,

(7)

where x is the independent variable, y is the dependent variable, s is the slope of the line, and m is the y-intercept. The regression parameters (s and m) are computed with the relationships:

$$s=\frac{\sum _{i=n}^n\left(x_i-\overline{x}\right)\left(y_i-\overline{y}\right)}{\sum _{i=1}^n{\left(x_i-\overline{x}\right)}^2}$$

(8)

$$m=\overline{y}-s\overline{x}$$

(9)

whith n the number of years for which the maximum annual values were calculated (in the present case 15 years), x is the independent variable represented by each year, y is the dependent variable represented by the maximum value corresponding to each year.

In order to complete the picture, Table 3 presents the maximum values of the wind speed and wind gust, corresponding to each month for the 15-year period 2009-2023. The values of the ratio between Uwg and Uw are also presented in this table. As it can be seen in Table 3, the maximum value of this ratio is 1.44. However, it can be noticed that such value is mainly characteristic to the high winds, because the average value of the same ratio for the entire data set is 1.73. This means that in normal conditions we can very often expect at the entrance in Sulina channel wind gusts with intensity almost double than the hourly registered value of the wind speed.

3.2. Analysis of Wind Data Offshore the Danube’s Mouths, Recent Past against near Future

The next analysis is made considering model wind data in the reference points RP1 (offshore the entrance in Sulina channel) and RP2 (offshore Saint George arm of the Danube). The analysis is made first for the past processing ERA5 data. The regional climate wind models, considered to evaluate the future expected wind conditions for the near future period 2041-2071, are RCA4 and ALADIN6, both for RCP4.5 and RCP8.5, and ALADIN6x for SSP5-8.5. Since the historical data provided by the climate models are available for the 30-year time interval 1976-2005, this period is considered for comparison with the ERA5 reanalysis to assess the convergence, but also the differences, between hindcast and historical data of the wind models. On the other hand, the SSP data are provided starting with 1979 and in this case, the comparison with ERA5 is made for the 30-year period 1979-2008. The wind roses corresponding to the reference points RP1 and RP2 (Sulina and Saint George, respectively) are illustrated in Figure 6 and Figure 7. In order to cover the entire period 1976-2008, Figure 6a and Figure 7a were designed for this 33-year period, Figure 6b, Figure 6c, Figure 7b and Figure 7c for the 30-year period 1976-2005, while Figure 6d and Figure 7d for the 30-year time interval 1979-2008. First, it can be noticed from these Figures the fact that the wind conditions in the two locations (RP1 and RP2) are in general similar. As regards the values of the maximum wind speed, for Sulina they are U_10ERA=20.94m/s, U_10RCA=27.49m/s, U_10AL6=26.16m/s and U_10AL6x=24.2m/s, while for Saint George U_10ERA=21.55m/s, U_10RCA=22.01m/s, U_10AL6=25.77m/s and U_10AL6x=25.4m/s. According to the results presented, the highest value of the maximum wind speed (27.49 m/s), between the two reference points considered, coresponds to RP1 (Sulina) and it was provided by the RCA4 climate wind model.

The annual maximum series corresponding to the same time intervals form the recent past considered before are presented in Figure 8 (for RP1) and Figure 9 (for RP2). At this point, it should be also highlighted that the time step of all wind model data considered and processed is of 3 hours, so considerably higher wind speed values can occur in this time window. It should be also noticed that there is an acceptable concordance between the models, the RCMs (and especialy ALADIN) providing in general higher wind speed values than ERA5.

The wind roses expected for the near future period (2041-2070) are illustrated in Figure 10 and Figure 11, for RP1 and RP2, respectively. The results presented indicate similarity with the recent past, with the observation that, according to the results, a higher percentage of the wind coming from the Northeast direction is expected in the future.

The corresponding annual maximum series and linear trends are presented in Figure 12 and Figure 13, respectively. For Sulina the maximum values of the wind speed evaluated by the climatic models are: U_10RCA4.5=26.11m/s, U_10RCA8.5=23.96m/s, U_10AL4.5=24.65 m/s, U_10AL8.5=24.48 m/s and U_10AL5-8.5=29.85m/s, while for Saint George U_10RCA4.5=25.35m/s, U_10RCA8.5=22.16m/s, U_10AL4.5=25.02 m/s, U_10AL8.5=22.83 m/s and U_10AL5-8.5=28.32 m/s.

As for the recent past period, the results show that in the near future time interval considered the highest wind speed is expected also in Sulina (29.85m/s). However, this time the highest wind conditions are provided by the SSP5-8.5 climate scenario.

It can be also noticed that, for the near future period analysed, the SSP5-8.5 data provide systematically higher wind speed values than all the other climate scenarios. On the other hand, the SSP data indicate higher wind speeds than in the past, a relative increase with about 22% in RP1 and 11.5% in RP2. From this perspective, the results indicate a slight tendency of enhancement of the extreme winds in the coastal environment at the mouths of the Danube River, a tendency also noticed in other works (as for example [47]).

3.3. Analysis of Wave Data Offshore the Danube’s Mouths, Recent Past against near Future

Using the wind fields provided by ERA5 for the recent past and those coming from the regional climate models (RCA4 and ALADIN) under the scenarios previously considered (RCP4.5, RCP8.5 and SSP5-8.5), simulations with the wave modeling system SWAN based have been carried out. The focusing of the system towards the target area, for a case of a typical storm in the Black Sea, corresponding to the time frame 21.12.2004h03, is illustrated in Figure 14 in terms of significant wave height fields and wave vectors. According to the results presented in this figure, for the case of the storms propagating from east to west (which are typical for this region), due to the effects of the wave-current interactions at the mouths of the Danube River the significant wave height values are very high, very close to the higher value in the entire sea basin. For example, in Figure 14 the maximum Hs value in the Black Sea is 7.26m, while the maximum Hs value at the mouths of the Danube is 7.24m.

Based on the results of the wave modeling system in the two reference points defined (RP1 and RP2) for the two 30-year time intervals (recent past and near future), and considering as forcing factors the wind fields provided by the models mentioned in the previous section, an analysis of the past and future expected wave climate in the target area is next presented.

The Hs wave roses corresponding to the reference points RP1 and RP2 (Sulina and Saint George, respectively) are illustrated in Figure 15 and Figure 16. In order to cover the entire period 1976-2008, Figure 15a and Figure 16a were designed for this 33-year period, Figure 15b, Figure 15c, Figure 16b and Figure 16c for the 30-year period 1976-2005, while Figure 16d and Figure 17d for the 30-year time interval 1979-2008. As in the case of the wind, it can be noticed from these figures the fact that the wave conditions in the two locations (RP1 and RP2) are in general similar in terms of both Hs and wave directions. As regards the values of the maximum Hs values, for Sulina they are H_sERA=4.5m, H_sRCA=5.65m, H_sAL6=6.27m/s and H_sAL6x=7.04m, while for Saint George the results are H_1s0ERA=4.91m, H_sRCA=6.061m, H_sAL6=6.96m and H_sAL6x=7.33m/s. According to the results presented the highest value of the significant wave height (7.33 m/s), between the two reference points considered, corresponds to RP2 and was provided by the ALADIN6x wind model as forcing system for SWAN.

The annual maximum series corresponding to the same time intervals from the recent past considered before are presented in Figure 17 (for RP1) and Figure 18 (for RP2). It should be mentioned that, as in the case of the wind, the time step considered for the wave model output is 3 hours. It should be also noticed that there is an acceptable concordance between the models. However, it can be observed at this point that when forcing the wave modeling system, the RCMs (and especially ALADIN) are providing in general sensible higher Hs values than ERA5.

The Hs wave roses expected for the near future period (2041-2070) are illustrated in Figure 19 and Figure 20, for RP1 and RP2, respectively. The results presented indicate similarity with those from the recent past. The corresponding Hs annual maximum series and linear trends are presented in Figure 21 and Figure 22, respectively. For Sulina the maximum Hs values evaluated are: H_sRCA4.5=6.19m, H_sRCA8.5=6.13m, H_sAL4.5=7.13m, H_1AL8.5=5.7m and H_sAL5-8.5=7.05m, while for Saint George H_sRCA4.5=6.86m, H_sRCA8.5=6.51m, H_sAL4.5=7.29m, H_1AL8.5=6.19m and H_sAL5-8.5=7.31m.

At this point, it can be noticed that, for the near future period analysed, a similar tendency as in the case of the wind fields occurs, in the sense that ALADIN data provides higher values for the significant wave heights than when the RCA4 model is considered as forcing wind field for the wave model. However, unlike in the case of the wind analysis, ALADIN6 for the RCP4.5 climate scenario provides higher Hs values in the target area than SSP5-8.5 data. Finaly, it can be also underlined that the results provided by the RCMs for the recent past in terms of the maximum values of the significant wave height, which may occur offshore the coastal area at the mouths of the Danube River, appear more realistic than the ones provided when ERA5 is considered as the wind driver for the wave model. This is probably because ERA5 data have a tendency to underestimate systematically the maximum wind speed values. This observation is in line with the findings of various previuos works (for example [43,44,45] or [47,48]), but also with some in situ observations.

4. Discussion

The analysis of the wind data shows that there is in general a good match between the measurements and the model results, especially as regards the intensity of the maximum wind speeds. However, by comparing Figure 4 and Figure 6 we can notice that in RP0 the dominant wind direction is from North-Northwest while in RP1 is from North-Northeast. Furthermore, unlike in RP0, in RP1 significant winds are coming also from the Southeast, while in both cases significant winds are coming from the South. The main explanation for the differences is that while RP1 is located 30km offshore, RP0 represents the zero kilometer of the Danube, where both the influence of the coast and also of the local wind currents propagating along the river are consistent.

As regards the waves, the results show that we can expect significant wave heights even higher than 7m in the coastal environment offshore the mouths of the Danube River. Furthermore, as Figure 14 clearly illustrates, the local processes, and especially the wave-current interactions and the shallow water effects, induce in the nearshore considerable enhancements of the significant wave height. From this perspective and in order to assess better the influence of the local effects at the entrance in Sulina channel, wave model simulations have been carried out considering the most common wave patterns from the point of view of the significant wave height (Hso) and mean wave direction of the incoming waves on the offshore boundary (Wdir), as indicated by the Figure 15, Figure 16, Figure 19 and Figure 20. Furthermore, previous results in terms of wave modeling presented in [48,49] have been also considered. Thus, in the high resolution computational domain described in Table 1, SWAN simulations have been performed considering significant wave heights and wave directions in the ranges [1m-5m] and [30⁰-150⁰]. Based on the results of these simulations, an index giving the relative enhancement due to the wave-current interactions of the significant wave height at the entrance in Sulina channel (REHs) has been evaluated. The expression of this index is given by the equation below, while the corresponding values are provided in Table 4.

$${RE}_{Hs}=\left({Hs}_{max}-Hso\right)/Hso$$

(10)

On the other hand, while the spectral wave models, such as SWAN, are used to provide predictions in terms of the significant wave height, which represents four times the area under spectral density, from the point of view of the navigation hazards we are more interested in the prediction of the larger wave heights, such as the maximum wave height. Usually, this is estimated based on some theoretical probabilistic distributions [50] among which the most common is the Rayleigh distribution [51], according to which waves with maximum wave height almost double the significant wave height can be often expected in a wave group. There are however situations when the ratio between the maximum and the significant wave height can be higher, or even much higher, than 2. Such kind of waves are called freak or rogue waves and they are very dangerous inducing a high risk of hazards in both navigation and any other marine activities. Thus, passing now from a statistical approach to a probabilistic one, we can estimate higher risks of occurrences for such kind of waves via the Benjamin-Feir Index (BFI). BFI, or the steepness-over-randomness ratio, has been introduced by Jansen [52] and is defined as:

$$BFI=\sqrt{2\pi }St·Q_p$$

(11)

St represents the integral wave steepness and represents the ration between the significant wave height and the wavelength while Qp is the peakedness of the wave spectrum defined as:

$$Q_p=2\frac{\iint \sigma E^2\left(\sigma ,\theta \right)d\sigma d\theta }{{\left(\iint \sigma E\left(\sigma ,\theta \right)d\sigma d\theta \right)}^2}$$

(12)

In fact, BFI is a spectral shape parameter that can be related to the kurtosis of the wave height distribution, the kurtosis depending on the square of BFI. Experimental results indicate that for BFI = 0.2 the maximum wave heights are very well described by the Rayleigh distribution while for values of BFI greater than 0.9 the ratio Hmax/Hs is substantially underestimated. The values of the BFI index corresponding to the most relevant wave propagation patterns are presented also in Table 4. The results indicate that BFI values higher than 1.5 occur in this area for waves with significant wave heights between 3 and 4 meters and when the mean wave direction on the offshore boundary is between 90 and 120 degrees. From this perspective, two additional simulations have performed considering two situations presenting the highest risk form the point of view of the rogue waves’ occurrence, and the results are illustrated in Figure 23. Thus, Figure 23a presents the bathymetric map of the computational domain and the average current field. Figure 23b illustrates in terms of significant wave height scalar fields and wave vectors the results of the SWAN simulations considering the offshore wave conditions Hs=3.5m and mean wave direction of 90 degrees, while in Figure 23c the same results considering the offshore wave conditions Hs=3.5m and mean wave direction 120 degrees. The values of the BFI index are indeed higher than those presented in Table 4, being 1.95 (with RE_Hs=28.6%) for the simulations from Figure 23b and 1.9 (with RE_Hs=31.4%) for that illustrated in Figure 23c. The above results indicate that, besides the enhancement in terms of Hs expressed by the values of the RE_Hs index, there is also a high risk of occurrence of abnormal waves for which the ratio between the maximum and the significant wave height is higher than 2.

5. Conclusions

The objective of this work was to study the past and future dynamics of the environmental matrix, especially wind and waves, at the entrance in Sulina channel. This is one of the three arms defining the Danube Delta in the Black Sea. At the same time, Sulina channel represents the main southern entrance in the navigation system Rhine-Maine-Danube, which represents the largest inland waterway in Europe linking the Black and the North seas. Hence, this sector is subjected to high navigation traffic. On the other hand, strong environmental conditions are very often encountered in this area inducing a high risk of navigation hazards and many accidents and incidents occurred in the past in this sector.

The analysis of 15 years of recent in situ measurements at the zero kilometer of the Danube showed that the maximum value of the wind speed (one hour step) registered in this point was almost 31m/s, while the maximum speed of the wind gust was almost 40m/s being with about 30% higher than the hourly value. The results also indicate that, for the entire data set analyzed, the average value of ratio between the maximum speed of the wind gust and the hourly wind speed value is of 1.7. However, in the case of the high wind speeds the intensity of the wind gust is relatively lower and the values of the ratio between the maximum speed of the wind gust and the hourly wind speed value is in general in the range [1.3-1.4].

A second analysis was performed in two locations (one offshore Sulina and the other offshore Saint George) and it is based on data provided by two regional wind climate models (RCA4 and ALADIN) under several scenarios (RCP4.5, RCP8.5 and SSP5-8.5). The analysis was performed for the future 30-year time interval (2041-2070). In order to compare the evolution of wind climate in relationship with the past conditions, a 30-year time interval from the past was also considered when ERA5 reanalysis data have been processed and analyzed in parallel with the RCM data. The climate model results indicate for the past similar intensities with those provided by the measurements for the maximum wind speed, all the models indicating maximum values of the wind speed higher than 25m/s. It has to be highlighted also that the time step considered for all models was 3 hours. The results show that slightly higher values than in the past are expected in the future for the wind speed, especially in the case of SSP5-8.5.

The last analysis relates the waves. In this connection, a wave modeling system that was previously validated against both in situ measurements and satellite data, was forced successively with each wind data field considered and simulations have been performed for the two 30-year time periods previously considered in the case of the wind analysis. The results show that Hs values of about 7 meters have been encountered in the past offshore close to the mouths of the Danube and such extreme values might be also expected more often in the future. The wave model simulations performed in the high resolution computational domain focused on Sulina channel also showed that, due to the nearshore effects, an enhancement between 30% and 40% in terms of significant wave height may be expected at the entrance in Sulina channel. Furthermore, not only the significant wave heights are higher, but also, due to the strong wave-current interaction processes, the probability of rogue waves occurrence is very high. This indicates that waves with heights of 20 meters, or even more, and also wind speeds higher than 30m/s, associated with wind gusts that might exceed 40m/s, can be expected in extreme storm conditions in this coastal environment.

Finally, the general conclusion of this work is that the navigation conditions in this sector present high risk of hazards, especially in the winter time. In fact, for this reason very often, even in the case of the moderate storms, the navigation is closed at the entrance in Sulina channel.