Numerical Estimation of the Black Sea Circulation near the Continental Slope Using SKIRON and ERA5 Atmospheric Forcing

: Assessments of the state of sea waters and complex studies of the marine environment in various ocean basins are often based on hydrophysical ﬁelds (currents, temperature, salinity, etc.) obtained through the use of numerical modeling. The regular ﬁelds of currents are of particular importance for assessing the transport of impurities in sea waters at different depths, including pollutants of various origins. The results of hydrophysical ﬁeld modeling, in turn, depend on the conditions set at the boundaries of the basin. Therefore, the correct setting of rapidly changing atmospheric conditions is extremely important for the reconstruction of marine dynamics. This paper presents model estimates of the Black Sea circulation obtained using two different datasets, SKIRON and ERA5, as atmospheric forcing. Numerical experiments for 2016 are carried out based on the eddy-resolving MHI-model. ARGO ﬂoats and R/V Cruises data are used to validate the simulation results. It was discovered that temperature and salinity RMSE between the model and measurement data are decreased under ERA5 forcing. Near the northeastern continental slope, a change in the direction of the alongshore subpycnocline current, which is detected in the ARGO ﬂoat trajectory, is modeled using ERA5 rather than SKIRON. Therefore, for a more accurate reconstruction of the Black Sea circulation, ERA5 atmospheric forcing is recommended.


Introduction
The circulation in the upper layer of each ocean is in direct contact with the atmosphere and is related to the distribution of meteorological parameters [1][2][3]. At the same time, the influence of atmospheric forcing on the structure of deepwater circulation is not so clear. For the Black Sea, this problem is complicated by the presence of a strong vertical density gradient (permanent pycnocline) at 50-100 m horizons, which blocks vertical seawater exchange [3].
Regional features of density stratification often arise near the Black Sea continental slope due to the mixing and lowering processes of surface waters along the slope into the deep sea layers [3,4]. The formation of density anomalies here can be caused by external forcing at the boundaries of the basin (including wind, river runoff, etc.), the sinking of denser waters down the continental slope, and the transfer of water with thermohaline characteristics that differ from the ones in eddies [4,5]. These processes are especially important in the northeastern part of the basin due to the narrow and steep continental slope in this region. The seawater density anomalies formed near the slope can lead to the transformation of the velocity field at deep horizons [5]. Thus, the generation of unsteady deepwater undercurrents is found out there [6].
Below the permanent pycnocline, the Black Sea waters become warmer and more saline accumulates towards the bottom [2]. At the same time, anticyclones can form near the shelf edge and then can move along the slope [4,5,7]. In the centers of the anticyclones, subpycnocline waters that are colder and contain less saline deepen, and their movement contributes to the transfer of thermohaline anomalies and the corresponding perturbations of dynamical fields. Such complicated dynamics near the continental slope require a detailed and accurate reconstruction of all hydrophysical characteristics, which is only possible if boundary conditions are correctly specified.

Materials and Methods
The Black Sea circulation was reconstructed by an eddy-resolving model from the Marine Hydrophysical Institute (MHI-model) [8]. The model was based on the Navier-Stokes equations in Boussinesq and hydrostatic approximations. Vertical turbulent mixing was described by the Mellor-Yamada closure model 2.5, and horizontal mixing was described using a bilaplacian operator with constant coefficients. The model circulation was driven by atmospheric forcing, including wind stress, heat fluxes, precipitation, and evaporation, on the sea surface. The climatological Black Sea rivers runoff and exchange through the straits were considered. Data assimilation (except for the satellite sea surface temperature data) was not used in the discussed numerical experiments. The MHI-model was implemented on a C-grid with a resolution of (1/48) • longitude, (1/66) • latitude, and 27 z-levels vertically. The detailed model description is presented in [8].
Basin bathymetry was built from EMODnet data [9]. The initial data were obtained from the Black Sea Physical Reanalysis CMEMS [10]. All initial and input fields were linearly interpolated in the MHI-model grid nodes.
In this work, two numerical simulations with identical model setups but different atmospheric forcing were carried out for the year 2016. In the first simulation (SKIRONexperiment), the forcing included 2 h of data on wind velocity, thermal, latent, sensible, and solar heat fluxes, evaporation, and precipitation provided by the SKIRON/Dust modeling system (Greece), with a spatial resolution of 0.1 • [11]. In the second simulation (ERA5-experiment), the forcing was based on the freely available hourly data of reanalysis supported by the European Centre for Medium-Range Weather Forecasts for the global climate, with a resolution of 0.25 • [12].
Comparative analysis of the SKIRON and ERA5 data showed a significant difference in wind forcing in the Black Sea region. As can be seen in Figure 1, the ERA5 wind stress is stronger than SKIRON one by about 25-30%, and the repeatability of NN-E and N-E wind directions (forming surface cyclonic circulation of the Black Sea) is higher. The remaining fluxes in ERA5 and SKIRON are close to each other, with there being some excess (15-20%) ERA5 data on total heat flux during the year and mass flux (precipitation minus evaporation) in autumn and winter.

Results
Daily data on sea surface height and three-dimensional fields of seawater temperature, salinity, and current velocity for the year 2016 were obtained in the two numerical experiments described above. The next stages of the study compare the simulation results with observational data (validation) and analyze deepwater circulation with a focus on the continental slope region, where the most interesting features of the currents are observed. In the northeastern part of the sea, so-called undercurrents (opposite to the Black

Results
Daily data on sea surface height and three-dimensional fields of seawater temperature, salinity, and current velocity for the year 2016 were obtained in the two numerical experiments described above. The next stages of the study compare the simulation results with observational data (validation) and analyze deepwater circulation with a focus on the continental slope region, where the most interesting features of the currents are observed. In the northeastern part of the sea, so-called undercurrents (opposite to the Black Sea surface basin-scale cyclonic gyre-the Rim Current [1][2][3]) are detected at a depth of 200 m for the period 9 June 2016-14 October 2016, according to ARGO data. Validation of the model fields was performed based on temperature and salinity measurement data obtained by ARGO profiling floats [13] and R/V «Professor Vodyanitsky» Cruises 87, 89, and 91 [14] in 2016. Our validation methodology is described in [8] (Section 2.2). Root mean square errors RMSE between the model and in situ data for both experiments are presented in Table 1. The temperature RMSE in the upper layer (0-300 m) decrease in the ERA5-experiment compared to the SKIRON-experiment. The highest decrease in temperature error was observed at a depth of 0-30 m. The model salinity in the ERA5-experiment correlated more at a depth of 30-300 m. Therefore, the permanent pycnocline and seasonal thermocline layers in the ERA5-experiment are closer to the measurement data. A difference between the simulation results was primarily found in the velocity fields due to the strong influence of the wind on Black Sea dynamics [1,2]. The increasing wind velocity in ERA5 (Figure 1b) led to a more typical structure of the Rim Current at the end of 2016, when the basin-scale cyclonic gyre was propagated above the continental slope ( Figure 2b). The Rim Current was not regenerated in winter, and mesoscale eddies were developed in the central sea part in the SKIRON-experiment (Figure 2a) due to insufficient kinetic energy inflow from the wind [8].
The model circulation in the upper layer was generally cyclonic for both experiments. At the same time, the most significant difference of the current velocity fields was detected below the permanent pycnocline core. Thus, at deepwater horizons, in the ERA5-experiment (Figure 2d), the current field was more intense, and maximal velocity was higher than in the SKIRON one (Figure 2c).
Analysis of ARGO float ID6901833 trajectory data [13] revealed a change in the direction of the alongshore subpycnocline current from the northwestern (cyclonic) to the northeastern (anticyclonic) near the northeastern continental slope. Thus, from 6 September to 14 October, 2016 the float drifted anticyclonicaly at its parking depth of 200 m (Figure 3a, red arrows). Such behavior of the alongshore current was not modeled in the SKIRONexperiment (Figure 3b), but was clearly reconstructed in the ERA5-experiment (Figure 3c near the northeastern slope of the Black Sea seem to be of an anticyclonic nature, similar to the undercurrents formed by anticyclones in the western part of the Bay of Bengal [15]. ( Figure 2b). The Rim Current was not regenerated in winter, and mesoscale eddies were developed in the central sea part in the SKIRON-experiment (Figure 2a) due to insufficient kinetic energy inflow from the wind [8].
The model circulation in the upper layer was generally cyclonic for both experiments. At the same time, the most significant difference of the current velocity fields was detected below the permanent pycnocline core. Thus, at deepwater horizons, in the ERA5-experiment (Figure 2d), the current field was more intense, and maximal velocity was higher than in the SKIRON one (Figure 2c). Analysis of ARGO float ID6901833 trajectory data [13] revealed a change in the direction of the alongshore subpycnocline current from the northwestern (cyclonic) to the northeastern (anticyclonic) near the northeastern continental slope. Thus, from 6 September to 14 October, 2016 the float drifted anticyclonicaly at its parking depth of 200 m (Figure 3a, red arrows). Such behavior of the alongshore current was not modeled in the SKIRON-experiment (Figure 3b), but was clearly reconstructed in the ERA5-experiment (Figure 3c). Averaged over the period of anticyclonical movement of the float, the model velocity of the undercurrent riches 0.03-0.05 m/s with instant value up to 0.08 m/s. The undercurrent generation near the Black Sea continental slope is probably associated with the intense mesoscale variability under the permanent pycnocline in the ERA5-experiment (Figure 3c). Here, some eddies were observed along the continental slope. The undercurrents that form near the northeastern slope of the Black Sea seem to be of an anticyclonic nature, similar to the undercurrents formed by anticyclones in the western part of the Bay of Bengal [15]. The structure of the circulation is inextricably linked with the spatiotemporal variability of seawater thermohaline characteristics [3,7]. The model temperature and salinity fields on the zonal cross-section along 44°N averaged over the period of existing undercurrent are shown in Figure 4. As seen in temperature fields (Figure 4a,b), the upper mixed layer reached a depth of 20-25 m in both experiments, but in the ERA5-experiment its thickness was larger near the eastern coast (up to 25-30 m), and its temperature was higher here as well. The mesoscale anticyclones shown in Figure 3c led to the deepening of isotherms and isohalines near the eastern coast and the formation of an undercurrent along the slope. There is a downward deflection of the isotherms at zone of 38.3-39.0° N in Figure 4b that corresponds to the anticyclonic current. A similar deflection is also visible in the salinity field (Figure 4d). Thus, the distribution of temperature and salinity in the ERA5-experiment is consistent with the anticyclonic current near the continental slope The structure of the circulation is inextricably linked with the spatiotemporal variability of seawater thermohaline characteristics [3,7]. The model temperature and salinity fields on the zonal cross-section along 44 • N averaged over the period of existing undercurrent are shown in Figure 4. As seen in temperature fields (Figure 4a,b), the upper mixed layer reached a depth of 20-25 m in both experiments, but in the ERA5-experiment its thickness was larger near the eastern coast (up to 25-30 m), and its temperature was higher here as well. The mesoscale anticyclones shown in Figure 3c led to the deepening of isotherms and isohalines near the eastern coast and the formation of an undercurrent along the slope. There is a downward deflection of the isotherms at zone of 38.3-39.0 • N in Figure 4b that corresponds to the anticyclonic current. A similar deflection is also visible in the salinity field ( Figure 4d). Thus, the distribution of temperature and salinity in the ERA5-experiment is consistent with the anticyclonic current near the continental slope detected in the ARGO float ID6901833 data [13].
bility of seawater thermohaline characteristics [3,7]. The model temperature and salinity fields on the zonal cross-section along 44°N averaged over the period of existing undercurrent are shown in Figure 4. As seen in temperature fields (Figure 4a,b), the upper mixed layer reached a depth of 20-25 m in both experiments, but in the ERA5-experiment its thickness was larger near the eastern coast (up to 25-30 m), and its temperature was higher here as well. The mesoscale anticyclones shown in Figure 3c led to the deepening of isotherms and isohalines near the eastern coast and the formation of an undercurrent along the slope. There is a downward deflection of the isotherms at zone of 38.3-39.0° N in Figure 4b that corresponds to the anticyclonic current. A similar deflection is also visible in the salinity field (Figure 4d). Thus, the distribution of temperature and salinity in the ERA5-experiment is consistent with the anticyclonic current near the continental slope detected in the ARGO float ID6901833 data [13].

Conclusions
The important outcome of the study is that atmospheric fluxes can affect the circulation of both the surface and deepwater layers of the Black Sea, and the choice of atmospheric forcing data can be decisive for the correct modeling of hydrophysical fields in the entire basin. As was determined through numerical analysis, with a significant influence of wind forcing, both the upper layer circulation and the deepwater dynamics in the Black Sea depend on the characteristics of the atmosphere. Despite strong density stratification and difficult vertical exchange with deep layers, the atmospheric forcing also affects the circulation at a horizon of 200 m and deeper. Thus, mesoscale features of the model dynamics near the continental slope, such as subpycnocline undercurrents detected from the ARGO observations in the northeastern part of the sea, appear only when using ERA5 forcing. Additionally, the Black Sea thermohaline structure is more accurately reconstructed under ERA5 forcing, and this was confirmed by the TS-measurement data. Atmospheric fluxes in ERA5 were more likely to be intense compared to SKIRON. Thus, for more accurate modeling of Black Sea circulation and its subsequent application for complex studies, from the two widely used meteorological datasets, the use of ERA5 atmospheric forcing data is recommended rather than SKIRON.  Data Availability Statement: The datasets generated and analyzed during the current study are available from the MHI administration upon reasonable and/or special request.