On the response of Black Sea level to external forcing: altimeter data and numerical modelling

Abstract

In this work, we address the Black Sea setup of Nucleus for European Modelling of the Ocean (NEMO), and in particular some model enhancements associated with the most important characteristic of ocean dynamics in this semi-enclosed basin, that is the sea-level variability and its relationship with water cycles and wind. Forcing data are presented in detail and compared with previously used coarser-resolution data. One emphasis in this paper is on the statistical analyses of forcing data and outputs from simulations with a focus on the sea level and its change. Numerical simulations are carried out as free run, and alternatively, altimeter data assimilation based on displacement of water properties in the pycnocline is used. Comparisons between the two runs identify the robustness of circulation driven by water balance and winter intensification. Problems in the model to replicate the redistribution of water properties between the two sub-basins in free-run mode are also discussed, which are observed during years with extreme climatic conditions.

Appendices

Appendix 1: Description of surface heat flux components

In our model, we use the same parameterisations for the individual components of heat fluxes as in Stanev et al. (2003). Following Reed (1977), the solar radiation is parametrised as:

$$\label{heat_tot2} Q_{s}=Q_{\rm tot}(1-0.62{}Cl+0.0019\times S_{\rm el})(1-A), $$

(16)

where Q _tot is the solar radiation reaching the sea surface under a clear sky, Cl is the fractional cloud cover (Fig. 3), S _el is the solar elevation at noon, and A is the sea surface albedo. It is noteworthy that mean and variations for clouds used in the present study are consistent with the data presented by Schrum et al. (2001).

The net long wave radiation loss Q _b can be represented as the sum of the long wave radiation emitted by the sea surface and the amount of long wave radiation absorbed from the sea. It was parametrised here following Berliand and Berliand (1952) and modified as from Rosati and Miyakoda (1988):

$$ \begin{array}{lll}\label{Q_b} Q_b &=& {\varepsilon} \sigma {T}^4(0.39 -0.05\sqrt{0.01 r_h e_{sa}(T_a)} (1-0.8Cl)\\ &&+\, 4 {\varepsilon} \sigma {T}^3 (T-T_a), \end{array} $$

(17)

where ε is the emissivity of the ocean, σ is the Stefan–Bolzman constant, and r _h is the relative humidity. The saturation vapour pressure at air temperature e _sa (T _a ) is approximated by the polynomial function from Loewe (1977).

The sensible heat flux H _a and the evaporation E _a are given by:

$$\label{sens} H_a = \rho_a {c_p}^a C_h |W|\times(T-T_a) $$

(18)

$$\label{evap} E_a = \rho_a {c_p}^a C_h |W|\times\left[e_{sa}(T)-r e_{sa}(T_a)\right]0.622/p_a, $$

(19)

where ${c\!_p}^a= \mbox{1.005} \times 10^ 3 \mbox{J kg}^{-1} \mbox{K}^{-1}$ is the specific heat capacity of the air, p _a is the atmospheric pressure at sea surface, and $C_h=1.1\times10^{-3}$ is the climatological value of drag coefficient C _d parametrising the exchange in the boundary layer between air and ocean.

Appendix 2: Reconstruction of river runoff from atmospheric analysis data

We use a straightforward solution to reconstruct the river runoff from independent observations by performing an EOF analysis on historical hydro-meteorological data to statistically link the two variables: precipitation from one side and river runoff from another. We would like to underline again that, in the following, we use precipitation over the sea because long-term data series are available to us only for this area (Altman and Kumish 1986; Simonov and Altman 1991; Stanev and Peneva 2002).

Preliminary analyses of the data (see also Peneva et al. 2001) reveal a clear seasonal and interannual variability. Response times of the Black Sea water balance to meteorological forcing identified in this previous research motivates us to first decompose the hydro-meteorological data in two data series: mean-year (or perpetual-year) for the period 1923–1994, T _py (t _m ), where t _m  = 1,...,12 (Fig. 14a) and annual mean T _am (t _y ), where t _y  = 1923, 1924,...,1994 (Fig. 14b). The notation T above could mean either precipitation or river runoff.

Fig. 14

Climatological values of monthly mean (a) and annual mean (b) water fluxes calculated from long-time observations 1923–1994

The linear correlation between perpetual-year of precipitation over the sea (one modal curve) and perpetual-year of river runoff (two modal curve, like the evaporation) is poor. Therefore, we only analyse the correlation between annual mean precipitation and river runoff. For the EOF analysis, we have a state vector with two entries, that is precipitation and river runoff; thus, only two modes can be derived. Afterwards, we are able to reconstruct the mean year river runoff by:

$$\label{rivrec_1} R^{r}_{am}(t_{y}) = \overline{R}_{am} + \sum\limits_{j=1}^{N} \eta_{j} \psi_{j}(t_{y})\;;\;\;N=2, $$

(20)

where $\overline{R}_{am}$ is the temporal mean of the observed annual mean river runoff and $R^{r}_{am}(t_{y})$ is its reconstruction based on the N EOF${}^{\text{pr}}$s (η _j ) and PC${}^{\text{pr}}$s (ψ _j ). The index ‘${}^{\text{pr}}$’ indicates that the EOFs and PCs refer to EOF analysis of precipitation and river runoff.

In order to obtain a reconstruction of river runoff based on the observed precipitation, a linear regression of annual precipitation (P _am ) and PC${}^{\text{pr}}$-1 and PC${}^{\text{pr}}$-2, respectively, is carried out, which explains 49% of the total variance for PC${}^{\text{pr}}$-1 (Fig. 6a) and 51% for PC${}^{\text{pr}}$-2 (Fig. 6b)). This result demonstrates that ≈ 50% of the observed annual mean river runoff variance at most may be expressed by the variability of precipitation.

Based on the above linear relationship, a reconstruction of river runoff can be given as

$$\label{rivrec_2} R^{r}_{am}(t_{y}) = \overline{R}_{am} + P_{am}(t_{y}) \sum\limits_{j=1}^{N} \eta_{j} b_{\psi_{j}} \;;\;\;N=2, $$

(21)

where $b\!_{\psi_{j}}$ is the correlation coefficient between the yearly amount of precipitation and individual PC${}^{\text{pr}}$s (ψ _j ).

At the end, we superimpose the seasonal signal from the perpetual year

$$ \label{rivrec_3} R^{r}(t) = R^{r}_{am}(t_{y}) \dfrac{R_{py}(t_{m})}{12 \overline{R_{py}}} $$

(22)

in order to obtain the reconstructed monthly river runoff for the whole period R ^r(t).

We calculate the mean error of reconstruction ($\varGamma$) as:

$$ \label{rerr_fqn} \varGamma = 100\% \times \dfrac{1}{N} \sum\limits^{N}_{j=1} \sqrt{\left(\frac{R^{r}(t)-R(t)}{R(t)}\right)^{2} }, $$

(23)

which gives us the mean error with respect to the original river runoff signal. The corresponding reconstruction errors are 18% of the original signal for reconstruction based on Eq. 20, which is the best possible, and 23% for Eq. 21, which is when the correlation method is used. Figure 15 shows the validation of the described reconstructions. The black curve displays the observed river runoff, the green line shows the reconstruction based on Eq. 20, and the red curve shows the reconstruction of runoff data based on Eq. 21.

Fig. 15

Comparison of river runoff from observations (black line) with reconstruction of river runoff based on Eq. 20 (green line) and Eq. 21 (red line)

Appendix 3: Rearrangement of water characteristics by vertical displacement

The implementation of the assimilation method needs some clarification of technical aspects and assumptions. When assimilating altimeter data, we focus on the time change of the dynamical integral, which can be presented as

$$ \label{dyn_int2} \dfrac{\delta \zeta}{\delta t} = - \dfrac{1}{\rho_{\text{0}}} \int_{0}^{H} \dfrac{\delta \rho}{\delta t}\,dz; $$

(24)

for the notations, see Section 3.2.

The general concept when assimilating $\dfrac{\delta \zeta}{\delta t}$ into the model uses Eqs. 13 and 24 and implies a displacement velocity

$$\label{displ_vel} w_{d}=\dfrac{\delta h}{\delta t} $$

(25)

in the area of the pycnocline (see Eq. 15). The change of density in the ‘climatic’ pycnocline due to displacement during unit time would then be

$$\label{cl_ch} \Delta \rho=\Delta h \dfrac{ \delta \overline{\tilde{\rho}_{cl}}}{\delta z}. $$

(26)

In the model, we rather displace temperature and salinity profiles instead of density profiles. Therefore, we introduce the assumption that, in the assimilation process, we can roughly use linear equation of state, thus adding displacement terms in the equation of temperature and salinity. We introduce relaxation terms in the equations of temperature and salinity:

$$ \gamma \left(T_{as} - T\right) \text{ ; } \gamma \left(S_{as} - S\right) \text{ with } $$

(27)

$$ \gamma = \dfrac{10}{\Delta t_{\rm obs}} \text{ and } (T, S)_{as}=(T, S)_{cl}+\Delta h \left(\dfrac{\delta \overline{\tilde{T}_{cl}}}{\delta z}, \dfrac{\delta \overline{\tilde{S}_{cl}}}{\delta z}\right), $$

where Δt _obs is the sampling rate of altimeter observations.

Assuming that climate is correctly considered by the free run, we substitute the needed climatic profiles in the last equation. In this way, we provide enough base to compare assimilation and free run.

To further apply the original method developed for reduced gravity models into the Black Sea z-level model (NEMO), we add some features, which do not change the general idea. We choose a depth range where displacement of temperature and salinity profiles is applied. The prescribed depth range is in the pycnocline and insures that the stratification of deep and surface near ocean layers will remain unaffected by data assimilation. This is a valid assumption for the deep layers because, in the Black Sea, deep profiles are almost homogeneous in the vertical. As for the surface layer, avoiding data assimilation ensures that the upper layer structure is formed mainly as a result of model fluxes at sea surface (mixed layer–seasonal thermocline dynamics). This ensures a performance mode that is physically correct and superior to the one where we disturb surfaces in the upper layer by data assimilation.

Consistent with the above description, we vary the displacement velocity w _d (Eq. 15) in a way which ensures that the major part of displacement takes place only close to the depth of the greatest salinity gradient (h _p ), which, for the Black Sea, is the core of the pycnocline at about 170 m. From the depth h _p , the vertical displacement is forced to decrease upwards and downwards by a factor of ξ(z) until it reaches zero at some prescribed depths (100 m and 500 m):

$$\label{beta} \xi(z) = \frac{\overline{\tilde{S}_{cl}(z)}-\overline{\tilde{S}_{cl}(100m)}} {\overline{\tilde{S}_{cl}(h_{p})}-\overline{\tilde{S}_{cl}(100m)}} $$

(28)

for 100m < z ≤ h _p , and

$$\label{beta_1} \xi(z) = \frac{\overline{\tilde{S}_{cl}(500m)}-\overline{\tilde{S}_{cl}(z)}} {\overline{\tilde{S}_{cl}(500m)}-\overline{\tilde{S}_{cl}(h_{p})}} $$

(29)

for h _p  < z ≤ 500 m, where h _p is 170 m.

To reduce the induced noise from data assimilation, we assimilate the running mean of three weekly data instead of assimilating the data point-by-point.