Storm Surges and Extreme Wind Waves in the Caspian Sea in the Present and Future Climate

The Caspian Sea is of particular interest. Against the background of long-term sea level changes, low-lying coastal areas in the northern part are subject to constant flooding as a result of storm surges. The elongation of the sea in the meridional direction allows the development of strong waves in the middle and southern parts. A comprehensive understanding of the characteristics of storm surges and storm waves is especially important in the context of ongoing climate change. This study is devoted to the analysis of storm surges for the time period from 1979 up to 2017 and wind waves from 1979 to 2020 in the Caspian Sea region. The circulation model ADCIRC and the wave model WAVEWATCH III with wind and pressure forcing from the NCEP/CFSR reanalysis were used. The modeling is performed on different unstructured grids with spacings of 500–900 m in the coastal zone. Mean and extreme values of surges, wave parameters and storm activity are provided in the research. The maximum significant wave height for the whole period was 8.2 m. The average long-term SWH did not exceed 1.1 m. No significant trends in the storm activity were found. The maximum surge height was 2.7 m. The analysis of the interannual variability of the surges' occurrence showed that 7–10 surges with a height of more than 1 meter were detected every year. The total duration of these surges was 20–30 days per year. Assessment of the risks of coastal flooding was carried out by calculating the extreme values of the sea for different return periods: 5, 10, 25, 50, and 100 years. The extreme sea level values in the northern part of the Caspian Sea for the 100-year return period are close to 3 m, and the areas with big surges are located along the eastern and western coasts. A forecast is made for the recurrence of storm wind waves in the 21st century based on climatic scenarios in CMIP5. A statistically significant increase in the recurrence of storm waves is to be expected in the near future, but that increase is not severe.


Introduction
The investigation of specific features of surges and wind waves in oceans and seas is an issue for modern society. The data collected for mean and extreme characteristics of sea level and wind waves are needed for naval navigation, marine works, offshore and coastal facility construction works, exploration and extraction of minerals, coastal protection strategies, etc. Extreme wind waves and surges get an important place among the factors risk and hazards assessing [1,2]. Currently, the generally accepted method for evaluating level and wave parameters in oceans and seas is mathematical modeling. This is due to the deficiency of instrumental observations data with the high spatial resolution and rather a long series often absent.
Some articles present the results of wind wave studies for the Caspian Sea [18][19][20][21][22][23]. However, the wave climate has not been properly researched. Lopatuhin et al. (2003) provided a long series of wave hindcasts [19]. There, some detailed tables of the repeatability of wave height depending on wind direction and the season are constructed using the WAVEWATCH III model and NCEP/NCAR reanalysis data. In this case, the grid spacing corresponded to the input data regarding the wind (~ 1.8°). Another article [21] demonstrated wave climate changes for the second half of the 20th century-up to the beginning of the 21st century and also used the NCEP/NCAR (National Centers for Environmental Prediction/National Center for Atmospheric Research) reanalysis with the resolution of 2.5° and the SWAN (Simulating WAves Nearshore) wave model to simulate the wave climate. The authors of [24] used the MIKE 21 SW model and ERA-40 reanalysis with a spacing of ~1.25°. They modified the coefficients of whitecapping and wave breaking based on the comparison of simulation and observational data, which reduced modeling errors for specific points and for a certain time period. Kudryavtseva et al. (2016) showed the main features of the wave climate in the Caspian Sea for the 2002-2013 time period based on altimetry data [25]. It is not surprising that the results of the analysis and forecasting of waves highly depend on the quality of the wind data input. Unfortunately, wind errors from the reanalysis occur in space and time factors, so that does not allow using the fixed coefficients to increase the accuracy of wave simulation [26]. The data of the NCEP/NCAR reanalysis considerably underestimate wind speed, which is proved by the results presented in [26][27][28]. The low correlation coefficient (0.7) obtained during the assessment of the skill of wave simulation based on this reanalysis was also confirmed for the Caspian Sea [21]. Some articles [26,29,30] demonstrate that the use of the data from the new-generation reanalysis (for example, NCEP/CFSR or ERA-Interim) leads to a considerable improvement in wave simulation quality. Thus, the earlier data on the Caspian Sea wave climate based on the low-quality reanalysis should be specified. The most modern analysis of waves in the Caspian Sea is presented in Bruneau & Toumi (2016) [31], where the authors used an ensemble of models (ROMS, WRF, and SWAN). However, the authors of Bruneau & Toumi (2016) [31] provided a hindcast only for three years. Such length of the series is insufficient for the analysis of mean and extreme characteristics. Some publications deal with the hindcast and climate projections of the frequency of synoptic conditions' occurrence, which cause severe hydrometeorological events, including storm waves [32,33]. However, these articles also used wind data from the NCEP/NCAR reanalysis. The article [34] provides the wave energy analysis for the Caspian Sea. Ivkina and Galaeva (2017), Strukov et al. (2013), and Zamani et al. (2009) [35][36][37] describe the systems of operational analysis and forecasting of waves in the Caspian Sea. These articles present simulation results; however, the analysis of mean wave characteristics for a long time period is absent. The task of forecasting storm waves is extremely important since their fluctuations lead to changes in marine ecosystems, the operating conditions of coastal facilities, and shipping. Not only the short-term forecast for the upcoming hours and days is important, but so is the climate one. This is the basis for developing a climate change adaptation and mitigation strategy [38,39].
Atmospheric circulation is one of the main factors forming the structure of the wind field and its intensity. In this work, patterns of atmospheric circulation are revealed and classified for storm cases in the Caspian Sea. The sea level pressure data for the present climate and climate projections are used.
Priorly, preliminary results were obtained on the modelling of storm surges and wind waves [40,41]. However, this article provides a more detailed and extended analysis. The motivation of this research is to simultaneously assess the impact of climate change on storm surge variability, storm activity over the last few decades, and the forecast of storm wind waves in the 21st century in the Caspian Sea. This issue is quite new and specific to the Caspian Sea. State-of-theart models and reanalyzes were used to calculate natural hazards for the Caspian Sea and coastal zones, which is an important regional contribution for ocean engineering in general.
The main objective of the present article is to determine the mean and extreme parameters of wind waves and storm surges in the Caspian Sea. Calculations were performed using the circulation model ADCIRC and the WAVEWATCH III spectral wave model and data of the NCEP/CFSR/CFSv2 reanalysis

ADCIRC Circulation Model Setup
The ADvanced CIRCulation model was used in order to estimate interannual and seasonal changes of wind-surges fluctuations in the Caspian Sea [42,43]. ADCIRC is a numerical model used for calculating water circulation and sea level fluctuations. It solves the complete equations of motion for moving fluid on a rotating earth.
The model equations were formulated using the hydrostatic pressure and Boussinesq approximations. It uses a finite element method for discretization by spatial variables, which allows the use of the unstructured grids. The approximation in time is carried out by the finite difference method. It takes into account such parameters as Coriolis force, atmospheric pressure, wind stress, tidal potential, and bottom friction. In the model, one can set the properties of the underlying surface. ADCIRC also includes flooding and drying of low-lying areas, as well as river flow.
There are two options for using the ADCIRC model: as a two-dimensional depth-integrated model (2DDI), and as a three-dimensional model (3D). In this work, a two-dimensional model was used. Level exceedance was obtained by solving the depth-integrated continuity equation in the Generalized Wave-Continuity Equation form (GWCE) [44]. Velocity was obtained from the solutions of momentum equations. All non-linear terms are preserved. The ADCIRC architecture successfully allows to use this model complex in parallel computing. ADCIRC can be run with either a Cartesian or a spherical coordinate system. The initial equations of the 2DDI ADCIRC model are: where x, y and t are horizontal grid points and time; H = h + η is total water depth; η is surface elevation; h is bathymetric depth; U and V are depth-integration currents in the x-and y-directions, respectively; Qx = UH and Qy = VH are fluxes per unit width; f is the Coriolis parameter; g is gravitational acceleration; Pa is atmospheric pressure at the surface; ρ0 is the reference density of water; (τsx,wind, τsy,wind) and (τsx,wave, τsy,wave) are surface stresses due to winds and waves, respectively; τbx and τby are bottom stresses; Mx and My are horizontal eddy viscosity; Dx and Dy are momentum dispersion terms; τ0 is a numerical parameter that optimizes the phase propagation properties.
Calculations were performed on the unstructured grid based on the bathymetric data of detailed navigation maps of the Department of Navigation and Oceanography (https://navysoft.ru/). On these maps, the depth was presented relative to the sea level of -28 m (Baltic Level System), taking into account flooding and drying of the coastal land. For the ADCIRC model, the grid consisted of 71523 points ( Figure 1). The grid spacing varied from 10 km in the open sea to 500 m in the coastal zone ( Figure 1). The computational grid was created using the Aquaveo Surface-water Modeling System (SMS 11.0) application. The similar computational grids used for the wave and surge modeling in the other seas within the Russian Federation proved their efficiency [30,45,46].

Figure 1. The computational unstructured grid for ADCIRC model and the map of the Caspian Sea depth
The input data were used for the fields of surface wind (at a height of 10 meters) and the atmospheric pressure of the NCEP (National Centers for Environmental Prediction) CFSR (Climate Forecast System Reanalysis) reanalysis ( Figure  2 and Table 1). The NCEP/CFSR reanalysis is a modern technology of the National Centers for Environmental Prediction implemented in 2010. The CFSR is a global atmosphere-ocean-land-sea-ice system with high resolution that provides the best estimate of the state of these interconnected systems. Data covering the period from 1979 to 2010 have a 1-hour time interval and a spatial resolution of ∼ 0.3° × 0.3° [47].  For numerical calculations from 2011 to 2017, the NCEP/CFSv2 (Climate Forecast System Version 2) reanalysis was used with a spatial resolution of ∼ 0.2° × 0.2° [48]. Sea ice concentration was also set as input (OSI-450), which reprocesses the brightness temperature data based on passive microwave radiometer SMMR, SSM/I, and SSMIS and ECMWF ERA-Interim reanalysis [49]. Ice accounting is very important in the study of the Caspian Sea [50].
A separate experiment was carried out using coupled model version ADCIRC+SWAN. The following configuration was used for the SWAN model: GEN3, KOMEN (cds2 = 2.36 − 5, stpm = 3.02 − 3), Quadrupl, Triad, Breaking constant (alfa = 1.0, gamma = 0.73) and Friction JONSWAP constant (cf = 0.067). The spectral resolution of the model is 36 directions (Δθ = 10°), the frequency range includes 36 intervals (from 0.03 to 1.1 Hz). The general time step for the integration of the full wave equation was 15 minutes.
Water level fields were obtained for every one hour from 1979 to 2017 (total 39 years) for each point of calculation grid as one of the ADCIRC model output.

Extreme Sea Level Calculation
In this work, the extreme sea level values were calculated to assess the risk of coastal flooding for different return periods: 5, 10, 25, 50 and 100 years. The extreme sea level analysis based on ADCIRC model output data with one hour time step. In accordance with the National Standard of the Russian Federation statement [51], the regime distributions of the total sea level heights belong to the type of exponential distributions, and, therefore, the Gumbel distribution is used to determine the extreme characteristics: Parameters A, B can be determined from a ranked sample of N annual maximum heights of the total level using the least squares method: after determining the parameters A and B, VT values possible once in T years, are determined as quantiles (1 − 1 ) • 100% probability distribution defined by F(V), V. Langbein established the relationship between the period of the frequency and number of values used in the sample that exceeds a certain value of the statistical variable: where N -is the total number of data.

WAVEWATCH III Wave Model Setup
The WAVEWATCH III (WW3) version 6.07 third-generation spectral wave model was used to calculate wind wave parameters for the Caspian Sea [52]. It is known that the energy inflow to waves is provided by wind energy and its dissipation is caused by the breaking of wave crests due to the bottom friction or wave breaking at the critical depth.
This model is based on a numerical solution of the equation of the wave action density spectrum: where, N(k, θ) ≡ F(k, θ)/σ, k -wavenumber, θ -propagation direction, σ -relative frequency, D/Dt represents the total derivative (moving with a wave component) and S represents the net effect of sources and sinks for the spectrum [52]. S function describes the transfer of energy from the wind to the waves, nonlinear wave interactions, and dissipation of the energy through the collapse of the crests at a great depth and in the coastal zone, friction against the bottom and ice, wave scattering by ground relief forms and reflection from the coastline and floating objects. The energy balance equation is integrated using finite-difference schemes by the geographic grid and the spectrum of wave parameters.
The ST6 scheme was used to generate waves during modelling. The nonlinear interactions were calculated using the DIA (Discrete Interaction Approximation) parameterization being a standard approximation for the computation of nonlinear interactions in all modern wave models.
The influence of the sea ice on the wave development was considered by the IC0 scheme, where a grid point is considered as ice-covered if ice concentration was > 0.5.
The increase in wave height as waves approach the shore and the related wave breaking after waves reach the critical value of the steepness is taken into consideration in the shallow water coastal zone, excluding the wave breaking due to the long wind effect on the sea surface (it is taken into account in the ST6 scheme). The standard JONSWAP scheme was used to take the bottom friction into account. The spectral resolution of the model is 36 directions (Δθ = 10°), the frequency range includes 36 intervals (from 0.03 to 0.83 Hz). The default wave model setup provided by ST6 scheme was used. Calculations were performed on the unstructured grid that consisted of 15792 points. The grid spacing varied from 10 km in the open sea to 900 m in the coastal zone ( Figure 3).

Figure 3. The computational unstructured grid for the WAVEWATCH III model and the map of the Caspian Sea depth [41]
The general time step for the integration of the full wave equation is 30 minutes, the time step for the integration of functions of sources and sinks of wave energy is 30 seconds, the time step for the spectral energy transfer and for satisfying the Courant-Friedrichs-Lewy condition is 900 seconds. The wind and ice data were the same as used for ADCIRC model (NCEP/CFSR/CFSv2 and OSI-450,430) ( Figure 2 and Table 1). The Caspian Sea level has considerably varied for the latest 40 years [5]. Therefore, different levels were specified in the model for every year to provide more accurate calculations. Year average of Caspian Sea level was calculated by using the data from several gauging stations (Makhachkala, Tyulenii, Peshnoi, and Baku) which was collected from "http://www.caspcom.com/". A similar implementation of the WAVEWATCH III model was successfully used for studying wave parameters in the other Russian Seas [53][54][55][56]. Wind wave fields were obtained every three hours from 1979 to 2020 as a model output (total 42 years). The data with time step of 1 hour and 3 hours were tested and did not reveal a significant change in the extremes (no more than 0.1 m). The model results include the SWH (4√ 0 , where m0 is the zero-order moment of the wave spectrum, approximately SWH is the mean value from 1/3 of the highest waves), the wave propagation direction, the mean wave period (WP) Tm02=(2 √ 2 ̅̅̅ ), and mean wavelength (WL)= (2 −1 ̅̅̅̅̅ ). The data were used to compute maximum and average long-term values as well as extreme characteristics.
The storm activity analysis was held according to the Peak Over Threshold (POT) method, which is widely used [57]. The essence of the POT method is to find extreme values of those samples that exceed a certain threshold value. Previously, POT was used for the Barents and Kara Seas wave analysis [53,56]. The number of storm waves with different SWH from 2 to 5 m was calculated for each year in the whole Caspian Sea region. The calculation procedure included the following steps: if at least one node in the investigated sea area had the SWH exceeding for example a 2 m (or a different threshold from 2 to 5 m), then such an event was attributed to the storm case with SWH threshold of 2 m. This event continued until the SWH was not less than the threshold at all nodes of the investigated area. Then, if the SWH threshold exceeded in one of the nodes again, then this event was added to the following case. A period of at least 9 hours passed between two storm cases in order to eliminate any possible errors. This technique has some inaccuracies associated with storms running in a row or from different directions at the same time. However, such cases are rare. The proposed algorithm works correctly; it was validated by visual analysis that has been conducted for several years.

Future Climate Data and Methods
Mathematically constructed models are effective tools for predicting the state of the climate system. However, even if one uses the most advanced climate models, only those indicators that are based on average values of climatic variables related to large territories or accumulated amounts can be predicted with confidence (with sufficient accuracy for practical use). Specific problems arise when predicting extreme events. In this article, one of those is consideredthe storm waves event. A necessary (but not sufficient) condition for the development of the storm is high wind speed, but its direct model forecast has a lower quality than the forecast of the atmospheric pressure field that forms the surface wind field.
Climate projection of weather hazards needs their classification. There are two fundamental approaches that allow to investigate the link between the large-scale circulation and environmental variables [58][59][60]. In the framework of the first approach, so called "circulation -to environment" approach [58]. There, the arrangement of the circulation data of interest (e.g. sea level pressure, geopotential height, etc.) is carried out to group them into circulation types (CTs) according to a selected methodology (clustering, principal component analysis (PCA), regression, etc.). Then, one searches for relations of CTs with the local-scale environmental variable (e.g. storm waves). On the contrary, the "environment -to circulation" approach is based on the classification of circulation data for certain criteria of the environmental variable, so that composite maps of the circulation variable can be derived for a specific environmental condition. Both approaches are widely used in various fields of atmospheric sciences. The first one was successfully implemented in the framework of the COST733 action (http://cost733.met.no) entitled "The harmonization and application of weather types classifications for European regions", an extensive review of existing classification methods and those used in COST733 is performed by Huth et al. (2008) [61]. In this study, the second "environmentto circulation" approach was used for storm events in the Caspian Sea region.
The classifications of atmospheric circulation are very useful for climate change researches; they are used for reconstruction of the past climate, analysis of variability of the present climate, and to estimate the future climate. The practical usage -to produce climate projection and for a more specific purposes, for example, for storm waves frequency in the future, as it was done in this study.
The wave height of 3 m and more is chosen as a criterion of a storm day. According to the data, the storm calendar was created. For these days sea level pressure (SLP) fields are used to classify circulation types. Sea level pressure daily data was taken from NCEP/CFSR reanalysis [47,48] for years 1979-2017. Then, basing on SLP fields, atmospheric circulation patterns were obtained by cluster analysis (k-means approach (e.g., [62]) pre-processed by Empirical Orthogonal Function (EOF) analysis (e.g., [63]) to reveal few leading modes determining the most part of the variance. These techniques of EOF decomposition and k-means cluster analysis together or in combinations with other techniques are widely used in circulation types classifications [59,[64][65][66][67][68][69]. At the first stage of classification, a dataset consisting of 30 daily Sea Level Pressure (SLP) grids was prepared for every storm case from the calendar -15 days before and after a storm day. After EOF decomposition of daily SLP grids, the first three eigenvectors explaining more than 70% of the variability were retained [33], thus filtering high-frequency perturbations (SLP-EOF fields).
SLP-EOF fields for storm days (according to the storm calendar of the sea) were used as input variables to classify circulation patterns. The definition of circulation types was carried out using the k-mean cluster analysis. The k-means algorithm starts with a preset number of k clusters and then it moves objects between clusters with the intention to, firstly, minimize the variance within clusters and, secondly, maximize the variance between clusters. Cluster centroids (ensemble mean of cluster members) were constructed for each circulation type by averaging the SLP grids of all days that belonged to the same circulation type. The dependence of assigning SLP field to a certain cluster within the area was checked: starting from initial territory 0-90E, 30-80N and gradually reducing it until it was comparable with the sea size. It barely influenced the result of sorting. Finally, the area 0-90 E, 30-80 N was left for centroids construction to have a full large-scale synoptic view of the circulation pattern.
After that, the daily SLP data of CSFR reanalysis for the time period of 1979-2017 was sorted considering previously derived circulation patterns. As one of the distance measures [61,70,71], space correlation was used between model data and reanalysis SLP-EOF fields for the storm calendar days. To eliminate 'noise' on the boundaries of a rather large initial domain (0-90E, 30-80N) and yet to save individual features of the SLP field, a correlation was estimated for 30-70E, 30-70N for the Caspian Sea. The spatial scale of the smaller domains is comparable to the size of such typical midlatitude synoptic structures as cyclones and anticyclones governing surface winds and therefore -storm waves.
Next, the spatial correlation was calculated for the storm days and CMIP5 models for Historical climate experiments and for the RCP8.5 scenario of the future climate [72]. The data of the following models were used: ACCESS1.

Storm Surges
To assess the quality of the ADCIRC model, the   The main task of this work was to analyze the interannual and spatial distribution of storm surges for the time period of 1979-2017. The distribution of the maximum surge heights for the Caspian Sea for the modeling period (1979-2017) is shown in Figure 6. Spatial distribution analysis of the maximum surges height showed that in the central part of the Caspian Sea maximal surges height did not exceed the 0.3 m. In the south part of the Sea, the maximal surges height did not exceed 0.8-0.9 m and it was located in the south-eastern part. In the northern part of the Caspian Sea, there were two areas of maxima. The first area was in the north-western part along the coast near the Volga River, and the second was in the east, where surges reached 2.5 and 2.7 m respectively ( Figure 6).

Figure 6. Distribution of the maximum surge heights for period 1979 -2017
The depth in the northern part of the Caspian Sea does not exceed 20 m. That part of the sea is semi-closed, thus it is probably the main reason for the strong surges to occur in this area. The next step of our research was to analyze the main synoptic situations which lead to the formation of storm surges. Storm surge is understood as a deviation of the sea level from the average annual level by a certain amount. This analysis is based on a model output sea level data from 1979 to 2017 in all points of the computational grid. The storm surge was calculated using POT method. All cases of surges exceeding 1 m were found, so synoptic situations that lead to these surges were analysed.
Three main synoptic situations that lead to the formation of surges of more than 1 meter have been identified:  For the formation of surge on the west coast, a strong anticyclone is formed north of the Caspian Sea, which comes either from the Asian maximum or from the north. At the same time, a cyclone is formed west of the Caspian Sea (in most cases on the European territory).
 The opposite situation is formed for the formation of the sea level decline. A strong cyclone is formed to the north of the Caspian Sea, as well as an anticyclone -to the west.
 The formation of a cyclone directly above the Caspian Sea, which moves to its northern part and forms a setdown of the sea level on the west coast.
One of the major tasks in the study of storm surges is to identify the contribution of the main factors to the surge development. To estimate the contribution of wind and atmospheric pressure to the formation of setup and setdown sea level, numerical experiments were carried out where the circulation model is forced only by the wind or atmospheric pressure. In both cases, ice fields were taken into consideration in calculations. As a result of comparing sea level according to different experiments, it was revealed that wind forcing is responsible for the largest changes in water level due to surge 92 -100% of total surge height and the atmospheric pressure induced 0 -8% in different synoptic situations. We understand that there are nonlinear and resonant effects under the simultaneous influence of pressure and wind; however, even the analysis showed that the contribution of wind is certainly much greater.
The number of surges events per year was calculated for the Caspian Sea according to the POT method (the technique is described in Section 2.2 and it was the same for surges and wind wave events). The storm surges events were calculated with different thresholds of 0.5 and 1 m. Setup and setdown events were calculated and analyzed separately. Analysis of the results showed that an average number of surges are 7-10 surges per year, with a height of more than 1 meter and the total duration of all surges -20-30 days per every year. Figure 7 shows an example of the distribution of the number of surges for year 1981.  The number of surges (setup and setdown events) for each year for 6 different points is shown on Figure 9. From 1979 up to 2005-2008, there was a long-term significant trend of decreasing the number of surges in all points. For example, the number of surges with a threshold 0.5 m reduced from ~25 to 15 cases in point 1 ( Figure  9). The significance of trends was assessed by the F-test. The F-statistic is the standard test statistic for testing statistical significance test of any linear model. The F-test statistics of the analysis of variance (ANOVA) approach was applied, which based on the null-hypothesis that the means of a given set sample of normally distributed populations having the same standard deviation are equal. Trends for all points except point 2 are significant at the level p = 0.05.

Figure 9. Number of surges per year at points 1-6
After years 2005-2008, the rise of the number of surges was detected in all points (Figure 9). If we consider the changes in the long-term average annual Caspian Sea level fluctuations, we can see the sea level rise from 1978 up to 1995 and a sustained decrease from 2005 to 2017 ( Figure 10). These fluctuations are associated with climate change [5]. A slight negative correlation can be observed between sea level and the number of surges. Since the surges are observed mainly in the northern part of the Caspian Sea, they are formed by winds that blow southwards. On the other hand, there is usually less precipitation during southern wind blowing, also less humidity and more evaporation. Thus, with an increase in the recurrence of southern winds and storm surges, the reduction of sea level is observed due to climate changes in the water balance. Then the extreme sea level values were calculated for different return periods: 5, 10, 25, 50 and 100 years based on probability analysis described in section 2.2. Figure 11 shows a graph of extreme sea level values for different return periods 5, 10, 25, 50 and 100 years. The minimum extreme values are observed at points 3 and 6 in the north and southeast, the maximum -at points 4 and 5 in the east. This is due to the prevailing winds causing the largest surges on the western and eastern coasts while the wave of level rise passes along the northern and south-eastern coasts. On the east coast, the level with a return period 100 years corresponds to 2.9 m, on the west -2.6 m. Figure 11. Extreme sea level values at points 1-6 for the return periods 5, 10, 25, 50 and 100 years Figure 12 shows the extreme sea level values in the northern part of the Caspian Sea for the return periods 10 and 100 years. The maximum value for a return period 100 years is close to 3 m and the areas with big surges are located along the eastern and western coasts.

Extreme Wind Waves
At the first stage of the wind wave climate studying, the quality of wave model results was assessed (due to the absence of digital data of direct measurements, observational data from   [74] were used). The results of the comparison of model and observational data for the wave heights of 3% probability of exceedance for the point located in the Central Caspian Sea are presented in Figure 13-a. The simulation quality may be visually assessed as satisfactory (we do not have any digital data, only images). The model simulates the phase of the storm beginning quite well; the model does not simulate storm peaks very accurately, but no systematic underestimation or overestimation of data is observed. Hindcasts were also compared with AltiKa altimeter data (rads.tudelft.nl). The satellite data of significant wave height (SWH) contain 34 990 points for the time period of 2013 to 2016 (Figure 13-b). The correlation coefficient was 0.918, the root-mean-square error, 0.28 m and the BIAS 0.07 m Scatter Index 0.29. Such accuracy of simulation generally corresponds to the modern realizations of wave models [26,36,53]. According to section 3.1, strong storm surges are observed in the northern part of the Caspian Sea, so that the sea level can vary by 0.5-1.5 m during several days. The water depth changes, which may affect the wave parameters. To estimate the effect of strong surges on the wave simulation for the northern part of the Caspian Sea, a special numerical experiment was performed. For October 1984, when the strong storm surge was registered, waves were simulated using data on the average annual sea level; the sea level increased by 1.5 m. The analysis of obtained data revealed that in the northern part of the Caspian Sea, where the depth is more than 7 m and waves with the maximum height (>2 m) are observed, the depth increase does not affect the wave height. Evidently, the main factor limiting the wave growth in this area of the Caspian Sea is the short fetch. In the areas where the sea depth is less than 2-3 m, its increase led to the wave height growth by 5-10%. Therefore, the effect of surges on mean and extreme characteristics of waves is not considered in the present article. The results of model calculations for the period of 1979 to 2020 were statistically processed. Figures 14-a and 14-b present the long-term mean and maximum SWH (13% probability of exceedance) over the whole simulation period. The multiyear mean wave height in the central part of the Caspian Sea reaches 1.1 m. The maximum value of wave height is 8.2 m, which is observed in the central part of the Caspian Sea. The local maximum SWH is equal to 6.9 m near the coast of Iran in the south part of the Caspian Sea. The wave development in the northern part of the Caspian Sea is essentially limited due to the short fetch, small depth and ice being present; therefore, the maximum wave height is 2. According to the data of wave climate book [19], the maximum height of 3% probability waves (which can be observed in the Caspian Sea once in 50 years) is 7 m. The results show a maximum value of SWH 8.2 m, which can be transformed to 10.8 m of 3% probability of exceedance (the wave heights of 3% probability of exceedance calculated as 1.32*SWH [75]. The use of NCEP/NCAR reanalysis data on the wind by the authors of (2003) for calculations probably considerably underestimated the wave height. Previous results had a maximum SWH about 8.9 m [41], but there an outdated ST1 scheme for WW3 model with the same wind data was used. The number of storm events per year was calculated for the whole Caspian Sea region according to the POT method (the technique is described in section 2.2). The events have different SWH thresholds, from 2 to 5 m. For the rest of the work, these storm events with a different wave height will be called storms. First, the number of storms for each year was analyzed (Figure 16  It is worth noting that there is some high interannual variability in the number of storms. The average variance is about 8% from year to year for storms with an SWH threshold 2 m and 60% for SWH thresholds 5 m. The most interesting feature is the significant linear trend for the storms with SWH threshold ≥ 3 m from 2003 to 2016 ( Figure 16). An increase for this period was from 30 up to 65. It correlated well with an increase in the number of surges of more than 0.5 meters which was also obtained for the period from 2005-2008 to 2016.

Future Climate Projections
For the situations of storm waves events, the typization of the surface atmospheric pressure fields was performed using the cluster analysis method. Three main types were identified ( Figure 17). These types are consistent with the five synoptic types developed earlier by genetically and circulation features [76] and five types of wind fields [77]. However, three types identified in this work are more generalized. The first type of circulation pattern is characterized by the development of a powerful anticyclone in the east of the European part of Russia, the Middle Urals and the Western Siberia. Cyclonic activity is forced out to the north and to the Mediterranean Sea. The southern periphery of the anticyclone determines the wind regime over the Caspian region. Almost half of the storm cases are attributed to the second type. Its distinguishing features are: the extensive cyclonic region over the Arctic, cyclones over the Middle East and Iran, two large anticyclones with centers over south-eastern Europe and Mongolia. Between the anticyclones in Central Asia, there is a low-pressure jumper that extends to the Caspian. The third type has some similarities with the second type, but it is distinguished by a less pronounced anticyclone over Europe and a pronounced region of reduced pressure over Central Asia. In this case, increased cyclonic activity contributes to the development of high wind speeds above the sea. In percentage numbers, type II and type III are dominant (Figure 17). Figure 18 shows the long-term dynamics of the anomalies for the general (accumulated) decades of repeatability for synoptic situations, which can be compared with events that include storm waves in the nowadays climate (with a SWH threshold ≥ 3 m). The results obtained from the ensemble of climate models indicate that, when the most unfavourable climatic scenario RCP8.5 is implemented during the 21st century, a gradual increase in the frequency of occurrence of such situations is possible. Most of the scenarios give an increase in synoptic situations leading to storms. Regarding the present climate situation, the number of storms with the SWH threshold ≥ 3 m is 40-60 cases and if a positive anomaly of 10-12 cases per 10 years is detected, but it does not grow significantly.

Discussions
However, there are several discussion points in the obtained results. The first point is a quality assessment of the circulation and wave models for the extremely high surges or waves. Unfortunately, the authors do not have full-scale direct measurement of sea level and wave height data for the Caspian Sea. The maximum value of sea level was 1.4 m, according to direct measurements. The calculated RMSE was about 0.1 m for sea level from the ADCIRC model. What will be the quality of the model for a surge height of more than 2-2.5 m? Unfortunately, we do not know the answer to this. The solution to this problem could be long-term sea level monitoring, which is currently not very active in the northern part of the Caspian Sea, or the data is not open to public.
In the case of the quality assessment of the wave model, the maximum SWH of 5.5 m was noted from direct measurements and the maximum SWH of 4.8 m from satellite data. An additional problem arises there -that the presence in the measurement data of 1-2 cases of storms with SWH 5-6 m height or a surge with a height of 1.2-1.4 m is completely insufficient for statistically acceptable quality assessments for the range.
Prior results had a maximum SWH of about 8.9 m [41], but an outdated ST1 scheme for the WW3 model was used. The quality assessment of the previous version of calculations based on the same satellite data showed exactly the same quality of the model for the whole range of SWH. However, it turned out that the ST1 scheme overestimated the SWH for more than 4-5 m (Figure 19), so it was decided to use the ST6 scheme, which provides the same quality over the whole range and a slightly better quality for high waves ( Figure 13). However, we are not sure that the modeling of extreme waves has been successful since there is very little measurement data for the high wave range.

Figure 19. Scatter plot of SWH derived from satellite and model data based on previous model configuration [41]
The same goes for the calculation of surge heights, especially for extreme events with return periods of 50 or 100 years. If you slightly change the maximum surge height (which is provided by long-term modeling), then the distribution function will change because the distribution function is highly dependent on the values of extreme surges.

Conclusions
This article presents new information about storm surges, wind waves, and their recurrence in the Caspian Sea region based on the results of numerical modelling. Long-term calculations were performed using the ADCIRC and WAVEWATCH III models.
The storm surge's maximum was 2.7 m. It was observed in the northern part of the sea for the modelling period . The northern part of the Caspian Sea is shallow water and semi-closed, thus it is probably the main reason for the strong surges in this area. The contribution of the wind stress to sea level fluctuations due to the surge is 92-100%, and the contribution of the atmospheric pressure is 0-8% in different synoptic situations.
The number of surge events per year was calculated according to the POT method for the whole sea or for the different points in the northern part of the Caspian Sea. There are 7-10 surges per year (on average), with a height of more than 1 meter for the whole sea. Based on analysis of the number of surges in the different points, it was found that from 1979 to 2005-2008 there was a long-term significant trend to reduce the number of surges in all points. After years 2005-2008, the rise of the number of surges was obtained in all points. A slight negative correlation between the average annual sea level and the number of surges was found. The number of surges on the east coast is greater than on the west for the northern part of the sea.

Data Availability Statement
The data presented in this study are available on request from the corresponding author.

Funding
The authors received no financial support for the research, authorship, and/or publication of this article.

Acknowledgements
The surge simulation was carried out by V. S. Arkhipkin and A. V. Pavlova and was supported by the Russian Foundation for Basic Research (grant 18-05-80088). The work of S.A. Myslenkov and G.V. Surkova was supported by the Interdisciplinary Scientific and Educational School of Moscow State University "The Future of the Planet and Global Environmental Changes."

Conflicts of Interest
The authors declare no conflict of interest.