Numerical modelling of thermodynamics and dynamics of sea ice in the Baltic Sea

Introduction Conclusions References


Introduction
Sea-ice phenomena constitute an important component of the physical and biochemical environment of the Baltic Sea.In the Bothnian Bay as well as in the inner parts of the Gulf of Finland and in the Gulf of Riga, sea ice develops every winter; the maximum ice cover during a typical winter extends over 157 000 km 2 , i.e. over about 40% of the total Baltic Sea area (Granskog et al., 2006).Low salinity, high fresh-water input from the surrounding land and the atmosphere, complex bottom topography and coastlines, coupled with strongly variable weather conditions, render the Baltic Sea quite unique in terms of its sea-ice phenomena -see e.g., Granskog et al. (2006), Leppäranta and Myrberg (2009) or Vihma and Haapala (2009).The specific nature of the Baltic has important consequences for the understanding and numerical modelling of the Baltic Sea ice thermodynamics and dynamics; it is hardly possible to directly apply insights provided by numerous studies in other regions, e.g. in the Arctic, to the Baltic because of substantial differences in temporal and spatial scales involved.The small-scale structure and physical properties of the oceanic and brackish ice are different as well.
Owing to its influence on various human activities (shipping and navigation in particular) and on the regional weather and climate, the Baltic sea ice has been in recent decades the focus of intensive theoretical, numerical and field research.A thorough review of the most important geophysical research on the sea ice in the Baltic, complete with a description of general characteristics of the Baltic ice, its long-term variability, and the state-of-the-art numerical modelling was provided by Vihma and Haapala (2009, see also references therein).Therefore, only those aspects of the sea-ice modelling that are most relevant to the problems dealt with in this paper are highlighted below.
Published by Copernicus Publications on behalf of the European Geosciences Union.
A. Herman et al.: Sea ice modelling in the Baltic Sea In most situations, the sea ice covers only a part of the total water surface and is a mixture of ice types differing in structure and properties -level, rafted and ridged ice, possibly with cracks and leads (Leppäranta and Myrberg, 2009).The numerical sea-ice models operating on the scale of tens up to thousands of kilometers reduce this information to a few parameters, typically to the concentration and mean thickness of the specified ice (and snow) classes in a given grid cell (for studies concerning the Baltic Sea see, e.g., Haapala and Leppäranta, 1996;Meier et al., 1999;Haapala, 2000;Lehmann and Hinrichsen, 2000a;Zhang, 2000).As described below, the model used in this study belongs to that class of models.Modelling studies of the Baltic Sea ice thermodynamics and dynamics can be broadly divided into two categories.One concentrates on problems of climate and climate change, seasonal and interannual sea-ice variability, and the influence of large-scale atmospheric circulation patterns on the sea-ice processes in the Baltic Sea (e.g., Haapala and Leppäranta, 1996;Omstedt and Nyberg, 1996;Lehmann and Hinrichsen, 2000a,b;Schrum et al., 2003).These studies are based on medium-and long-term simulations and mostly involve variables such as the maximum annual ice extent or the length of the ice season, i.e., parameters which can be understood as proxies of winter severity and which are therefore good indicators of a climate change.The other category consists of studies which tackle smaller spatial and temporal scales with the aim of analyzing effects of synoptic-scale weather patterns on the thermodynamics and dynamics of the sea ice.Good examples are provided by Uotila (2001), Brümmer et al. (2002), Rudolph and Lehmann (2006), Wang et al. (2006) or Björk et al. (2008).The need for more extensive observational and numerical research is widely recognized, particularly with respect to short-and medium-term ice dynamics in the Baltic Sea.
Because of the high spatial and temporal variability of the Baltic sea ice properties, observational data representative of larger areas, and thus appropriate for verifying the results of numerical models are difficult to acquire.As opposed to single-point measurements (at stations that are usually located directly at the coast, i.e., within the fast-ice zone) that very often do not meet this requirement, application of remote sensing techniques supplies valuable spatially resolved information.However, the use of satellite or airborne data in combination with numerical sea-ice modelling in the Baltic Sea has been rather limited to date.Moreover, verification of modelling results is usually qualitative only, and/or limited to a few arbitrarily selected situations.For example, Zhang (2000) visually compared results provided by his numerical model with the Finnish Institute of Marine Research (FIMR) ice charts and the Synthetic Aperture Radar (SAR) images.Meier et al. (1999) and Meier (2002) qualitatively compared their modelled ice concentrations and thickness with the Swedish Meteorological and Hydrological Institute (SMHI) ice charts for a number of selected situations.In a similar manner, Rudolph and Lehmann (2006) used ice concentration maps of the US National Ice Center.This approach is tainted by a high uncertainty of the data (although the situation has substantially improved in recent years), and is not suitable for systematic, quantitative comparisons with modelling results over longer time periods.Some models are verified against the measured ice thickness at selected locations (e.g., Meier, 2002).Another frequently used approach involves comparing temporal changes in an observed (e.g.satellite-derived) and simulated total ice-covered area, as in Lehmann and Hinrichsen (2000b), Meier (2002) or Schrum et al. (2003).
In view of the above considerations, this work has two objectives.The first is to provide a systematic, quantitative validation of the sea ice model performance (in terms of sea ice concentrations) against satellite-derived data over the entire model domain and the entire time period analyzed.The second objective is to analyze interactions between synoptic-scale atmospheric forcing and freezing, melting and (re)distribution of the sea ice in the sub-basins of the Baltic Sea.In particular, correlations between meteorological variables and the ice-covered area and the role of ice dynamics in the spatio-temporal evolution of the ice cover are investigated in detail.
The paper is structured as follows: the subsequent section, after a short description of the hydrodynamic model used in the study, presents the theory and implementation of the thermodynamic and dynamic parts of the sea ice model.Section 3 provides the description of the model configuration and the input data, followed by validation of the simulated data against measurements.The modeling results are discussed thoroughly in Sect. 4 in the context of meteorological conditions over the Northern Baltic Sea during the period of study.In particular, changes in the total ice cover and the influence of dynamic processes on the temporal and spatial variability of ice thickness and concentration are analyzed.Finally, Sect. 5 contains the summary and conclusions.

Hydrodynamic model M3D UG
The sea-ice model presented here was developed as an additional module of the hydrodynamic model of the Baltic Sea known as M3D UG.The model, coupled with the ProDeMo ecosystem model, is working pre-operationally at the Institute of Oceanography, University of Gdańsk (see http:// model.ocean.univ.gda.pl/indexeng.html).The hydrodynamic model is based on the Princeton Ocean Model (POM; Blumberg and Mellor, 1987) with modifications concerning calculation of the advective terms in the momentum equations which improve the model's performance in regions with strong salinity and temperature gradients, as described in detail by Kowalewski (1997).M3D UG is a σ -coordinate model formulated on the scattered Arakawa-C grid.As opposed to the original POM, formulated on a "traditional" numerical grid, each mesh of the M3D UG is treated as a one-dimensional entity (water column with a specified number of layers) with pointers to its neighbours that make it possible to calculate spatial derivatives.The model covers the Baltic Sea with a resolution of 5 .For computational efficiency, the equations are split into the fully threedimensional part (an "internal mode" integrated with a time step of 1200 s) and the vertically-integrated part (an "external mode" integrated with a time step of 60 s).The model is driven by the time series of hourly atmospheric pressure, wind, air temperature, humidity, and cloudiness fields from a mesoscale atmospheric model of the Interdisciplinary Centre for Mathematical and Computational Modelling (Warsaw, Poland).The model is a limited-area version of the global Unified Model of the UK Meteorological Office.It covers Central and Northern Europe with a resolution of ∼17 km (Herman-Iżycki et al., 2002).The M3D UG is well calibrated to the Baltic Sea conditions and has been successfully validated against available observational data from a number of stations in the Baltic Sea (Kowalewski, 1997(Kowalewski, , 2002;;Jedrasik, 2005).In this work, the hydrodynamic part of the model was not subject to any modifications.

Sea ice model -general remarks
The sea ice model developed in this study describes the spatial and temporal evolution of the following four variables: ice concentration A i (the part of a grid cell covered with ice); mean ice thickness h i ; mean (vertically averaged) ice temperature T i ; and ice velocity U i = [U i ,V i ] (where U i , V i denote the east-west and south-north components of the velocity vector, respectively).The present version of the model does not include the snow cover.For obvious reasons, a onelayer sea-ice model is not capable of modelling the vertical temperature distribution within the ice; whenever the temperature of the upper surface of the ice, denoted further with T i,up , is required, it is estimated as T i,up = 2T i −T f (where T f denotes the freezing temperature of sea water).
The sea ice model and the hydrodynamic model are coupled at each internal time step by means of air-ice-water heat and momentum flux.The sea ice model itself consists of two largely independent parts: the thermodynamic and the dynamic one (see Sects.2.3-2.5 below), executed one after the other at each internal time step of the hydrodynamic model, similarly as in, e.g., Haapala (2000).h i and A i are the only parameters that are exchanged between the two parts and modified by both of them.Ice temperature is calculated in the thermodynamic part, the dynamic part computing ice velocity.The theoretical background and numerical implementation of the two parts are described in detail below.
The equations of the ice model are discretized on an Arakawa-C grid, which has two important advantages.Firstly, interpolation of various quantities between the grid points of the sea-ice and the hydrodynamic model can be avoided.Secondly, a numerical scheme for the ice rheology developed by Bouillon et al. (2009), more efficient and robust than the available numerical schemes for other grid arrangements, can be used.

Heat flux at the atmosphere-water interface
The total heat flux F tot,aw through the atmosphere-water boundary is a sum of five components: the short-wave solar radiation penetrating through the atmosphere and absorbed by the water surface F s , the long-wave outgoing radiation F Lw , the long-wave incoming radiation F La , the latent heat flux F e , and the sensible heat flux F wt : (1) Negative sign of F tot,aw means a heat loss (cooling) of the water surface.
The short-wave solar radiation penetrating through the atmosphere and reaching the sea surface, F s,in , is calculated as described in Kreżel (1997).Thus: with α w denoting the sea surface albedo (dependent on the sun's position above the horizon).The long-wave outgoing radiation is calculated from the Stefan-Boltzmann law: where T w is the water temperature of the uppermost layer of the hydrodynamic model, and w and C SB denote emissivity of the water surface and the Stephan-Boltzmann constant (C SB = 5.67 × 10 −8 W m −2 K −4 ), respectively.The longwave incoming radiation F La is a function of the cloudiness N ∈ [0,1], air temperature T a and humidity (expressed in terms of water vapour pressure e a ), see e.g., Meier et al. (1999).The present model utilizes an expression adjusted to the Baltic Sea conditions (see Zapadka et al., 2001Zapadka et al., , 2007)): where a denotes emissivity of the atmosphere, e a is the water vapour pressure at T a , and d a , d C are empirical constants, adjusted during calibration of the model.The latent heat flux F e is given by Fischer et al. (1979) as: where L e = (595.9−0.54T w )C w , c e denotes the transfer coefficient, ρ a is the air density, L e is the latent heat of evaporation, C w is the specific heat of water at constant pressure and U a denotes the wind speed 10 m above the water surface.
The absolute humidities q a and q w are calculated according to the formulae derived directly from the equation of state for humid air (see, e.g., Maykut, 1986):  where p is the atmospheric pressure, and E(T w ) is the saturation vapour pressure at T w , calculated from the Clausius-Clappeyron equation: where E 0 = 6.1078 hPa denotes the saturation vapour pressure at temperature T 0 =0 • C, and R is the universal gas constant.The sensible heat flux F wt depends on the temperature difference between the water and the atmosphere.It is positive during stable and neutral conditions T a ≥T w and negative otherwise: Generally, the empirically adjusted transfer coefficients corresponding to stable and unstable conditions, c wts and c wtu , respectively, may have different values (typically c wtu >c wts ).
In the present model, they were set during the model calibration (see further Sect.3.1 and Table 1).

Heat flux at the atmosphere-ice interface
The total heat flux through the air-ice boundary, F tot,ai , is given, by analogy to Eq. (1), by: with F La given by Eq. ( 4) above.The expressions for the remaining four heat flux components are modifications of Eqs.
(2), (3), ( 5) and ( 6): F e,i = c e ρ a L s U a (q a − q i ), (10) where α i denotes the ice surface albedo, L s is the latent heat of sublimation, C i is the specific heat of ice and q i is the saturated absolute humidity at the ice temperature T i,up .The transfer coefficients c e , c wts and c wtu are assumed universal for both air-water and air-ice interface.

Heat transport through the ice
The energy flux through the ice consists of the conduction heat flux F h,i , dependent on conductive properties of the ice, and the short-wave solar radiation flux F s,b , dependent on light extinction within the ice.The heat conduction through the ice is given by: where κ T is the thermal conductivity of ice, and T f denotes the temperature at the bottom ice surface, assumed equal to the freezing temperature, see Eq. ( 16) below.The heat conductivity κ T is assumed constant.The short-wave solar radiation flux through the ice is calculated according to Sahlberg (1988) as: where f sol and κ denote the penetration factor and the bulk extinction coefficient for ice, respectively.

Heat flux at the ice-water interface
The net heat flux at the bottom of the ice is the sum of the two flux components described above, F h,i and F s,b , and the turbulent heat flux F t,b : The turbulent heat flux depends on the ice-water temperature difference and their relative velocity: where c is denotes the bulk heat transfer coefficient from the water into the ice, ρ i is the ice density, and U w is the current velocity in the uppermost layer of the model.

Freezing and melting of sea ice
The freezing point T f for the Baltic Sea water is given by Omstedt (1994) as: where S denotes salinity.Freezing and melting of the sea ice depend on the sign of the heat fluxes at the air-water-ice interfaces and on the temperature of the ice and of the surrounding water.The conditions that must be met for a new ice to develop are: F tot,aw < 0, T w ≤ T f , and A i < 1.The "new" ice forms in the open-water part of the grid cell analyzed (i.e. over an area A i,new = 1 − A i ), its thickness being calculated from: where L i denotes the latent heat of freezing.Freezing occurs only if h i,new = t∂h i,new /∂t ( t is the time step) exceeds a threshold value h i,tres = 10 −4 m; otherwise h i,new is set to zero.Additionally, if h i,new <h i,min = 0.05 m, A i,new is multiplied by h i,new / h i,min and then h i,new is set to h i,min (as in, e.g., Meier et al., 1999).
If A i >0 (ice is present in a given grid cell) and the net heat flux at the lower surface of the ice is different from zero, the ice thickness changes from below according to: Additionally, if F tot,ai −F tot,i >0 and T i ≥T f , the ice melts at its upper surface: If A i >0, F tot,aw >0 and T i ≥T f , the ice is subjected to lateral melting, resulting in reduction of the ice concentration A i .To estimate the lateral melting rate, it is necessary to estimate the lateral surface area of ice floes, a kind of information not available in the type of model used here.To solve this problem, we followed the approach of Steele (1992) (see also Ohshima et al., 1998).If melting occurs uniformly at a rate F lat (assumed equal to F tot,aw ) around the perimeter of the floe, and if the sea ice is adequately described by a single mean caliper diameter L, the rate of change of the ice concentration under melting conditions can be expressed as: where c lat (proportional to L −1 ) can be assumed constant; see Steele (1992) for details.Taking all the above-mentioned processes into account, the changes of A i and h i during a time step t can be expressed as: and: Finally, the rate of change of the ice temperature is calculated from the vertically averaged heat conduction equation:

Ice dynamics
As already stated, the dynamic module of the sea ice model is responsible for simulating changes in the ice concentration A i and thickness h i , produced by the action of wind, currents and stresses within the ice itself.Consequently, in the dynamic module, the total ice mass in the modelled area does not change; it may only be redistributed between the adjacent meshes of the numerical grid.In addition to the modified A i and h i values, the calculations produce components of the ice velocity U i , V i and the surface stress τ , given as weighted average of the stress at the air-water and air-ice interface.The dynamic module is based on two equations: the mass (or, equivalently, volume) conservation equation and the momentum equation.As in the hydrodynamic part, the ice module is formulated on the scattered Arakawa-C grid, with the scalar quantities located in the mid-point of the grid cell, and U i and V i located at its "western" and "southern" boundary, respectively.

Ice volume
Ice conservation law is considered separately for h i (volume conservation): and for A i (area conservation): see e.g., Leppäranta and Myrberg (2009).The index "d" denotes changes due to the dynamics alone, without accounting for thermodynamic effects described in the previous section.Numerically, exchange of the ice volume at each time step of the model is calculated for each pair of neighbouring "wet" grid cells based on the ice velocity component through the common boundary of these cells.This simple approach ensures conservation of the total ice volume at each time step.The ice volume in a given grid box is V = A i h i x y (where x, y denote the grid size in the x-and y-direction, respectively).Analogously, the ice volume in the "western" The volume of ice exchanged between these two boxes in time t is: Thus, the updated ice concentration and thickness are: and: Expressions for A i,W (t+ t) and h i,W (t+ t) are analogous.

Momentum balance equation
The present model assumes the advective acceleration of the ice to be negligible, compared with the inertial acceleration.
The gravitational force produced by the sea surface tilt is neglected as well (see scaling arguments in, e.g., Leppäranta, 2005).Thus, the motion of the ice is governed by the following momentum balance equation: where f is the Coriolis parameter, k = [0,0,1] and σ is the internal ice stress tensor.The motion of the ice is driven by the wind stress acting on the ice upper surface, τ ai , and by the bottom stress at the ice-water interface, τ iw .The air-water stress vector τ aw , together with τ ai and τ wi , are the three stress components acting on the air-ice-water interfaces.They are approximated by the following bulk formulae: where ρ w denotes the water density, C aw , C ai , C iw are nondimensional drag coefficients for the air-water, air-ice, and ice-water boundaries, respectively; and U a is the wind velocity vector at 10 m above the water surface.The total surface stress τ , returned to the hydrodynamic model by the icedynamics module, is calculated as the weighted average of τ iw and τ aw : The last term in Eq. ( 29) is calculated with the Elastic-Viscous-Plastic (EVP) model of Hunke and Dukowicz (1997).Compared with the "traditional" Viscous-Plastic (VP) formulation by Hibler (1979), it has two important advantages.Firstly, it produces more accurate results for shortterm variability of the forcing conditions (the VP model being rather suitable for studies on climatological time scales).
Secondly, and most importantly, EVP is suitable for a computationally effective and numerically stable formulation.This work adopted the numerical scheme for the EVP model on the Arakawa-C grid developed by Bouillon et al. (2009), in which the definition of the elastic parameter E (the Young modulus) is redefined in terms of a damping time scale T d .The rheology module is sub-cycled with a time step t r .The values of T d and t r are chosen so that t r :T d : t = 1:40:120, as suggested by the stability analysis of Hunke and Lipscomb (2006).The ice strength constant P * was set to 2.5×10 4 N m −2 , as in Leppäranta et al. (1998) who calibrated their model (with spatial resolution similar to that used here) to the Bothnian Bay conditions based on the ERS1-SAR-derived ice velocity fields.The spatial resolution of the model used here ( 5) is a compromise between the need to resolve details of the bathymetry and requirements of the continuum approximation underlying the EVP model, which is valid only if there is a large number of ice floes per model grid cell (Hibler, 1979).Typical floes in the Baltic Sea are up to 1-3 km in diameter (Leppäranta et al., 1998;Uotila, 2001).

Model simulations and sea-ice data
The model was run for a period of July 2005-May 2008, i.e., for three winter seasons.The maximum total ice extent S max during those winters was 210×10 3 , 139×10 3 and 49×10 3 km 2 (the smallest area ever recorded), respectively.According to the S max -criterion, the winter 2005/2006 was an average one (the median S max for the period 1980-2008 is ∼140×10 3 km 2 ); it was therefore chosen in this study for a detailed analysis of the sea-ice variability in relation to atmospheric forcing.However, evaluation of the model's performance is based on data from the whole period of simulation.
The adjustable model coefficients (empirical constants), listed in Table 1, were calibrated to the measured temperature data from eight stations (circled dots in Fig. 1) and to the satellite-derived ice concentrations.Initially, a set of realistic values was defined for each of the adjustable coefficients (based on the values reported in the literature).The model was run for each possible combination of values from these sets.The Nash-Sutcliffe effectiveness coefficient was used as a measure of the model performance in order to select the optimal combination of the coefficients.The calibration was performed in two stages.At the first, the coefficients important for the air-water fluxes ( w , a , c e etc.), calibrated previously to the old, hydrodynamics-only model version, were re-calibrated to optimize the water-temperature development during the whole period of study, including the summer seasons (Fig. 2).After those parameters were set, the remaining ones, directly related to the ice development, were determined at the second stage, based on the satellite data.The calibrated values of the coefficients lie within the range of values reported in the literature (e.g., Meier et al., 1999;Wang et al., 2006).
In this study, we used the satellite-derived ice concentrations from the Nimbus-7 Scanning Multichannel Microwave Radiometer (SMMR) and the Defense Meteorological Satellite Program (DMSP) SSM/I passive microwave data, available from the National Snow and Ice Data Center (NSIDC; Cavalieri et al., 1996Cavalieri et al., , updated 2008)).The ice concentration fields are available daily on a regular grid in the polar stereographic projection, with spatial resolution in the Baltic Sea region of about 25 km (black points in Fig. 1).For this study, the data were linearly interpolated onto the model grid, with nearest-neighbor extrapolation along the coasts.
Although the performance of the algorithm used to calculate the ice concentrations from the brightness temperature was assessed in numerous studies (see http://nsidc.org/data/docs/daac/nsidc0051 gsfc seaice.gd.html and references therein), most of those studies dealt with the Arctic region, i.e., areas with consolidated and relatively thick pack ice.The accuracy of the NSIDC data is estimated to be ±5% for such regions and ±15% for regions with thin ice and/or with ice melting and/or very low ice concentrations.Thus, in the case of the Baltic Sea, the data should be interpreted with some caution (Rudolph and Lehmann, 2006, provide an example of applying these data for the Baltic Sea).To assess reliability of the NSIDC data, we performed a qualitative, visual comparisons with the SMHI ice charts covering the time period analyzed.In most cases, the NSIDC data agreed with the SMHI charts very well.Some discrepancies observed during this quality check are described in Sects.3.2 and 4 below, together with discussion of the results.
Instead of comparing the simulated and the satellitederived ice concentrations in a point-by-point manner, we divided the ice concentrations into five classes (0-20%, 20-40%; etc.) and compared the area covered with ice concentrations from a given class, S 0-20 i , S 20-40 i , . . ., S 80-100 i .
The analysis was performed separately for each of the five subbasins shown in Fig. 1 (no sea ice was present in the time period considered in the Southern Baltic, except for shallow coastal locations not resolved by the model).We find this approach to be less controversial in view of the data uncertainties mentioned above.The simulated distribution of the ice thickness was additionally qualitatively compared with the SAR-based ice-thickness charts available at http://haavi.fimr.fi/polarview/charts.php.The selected cases analyzed showed a very close agreement between those data and the SMHI ice charts.

Model performance
Figure 2 shows the simulated and observed water temperature and sea-ice parameters at four selected stations.The simulated water temperatures at the Gulf of Finland's stations Muuga and Sillamäe (Fig. 2a, b; the results for Paldiski were very similar to those at Muuga and are therefore not shown) were slightly underestimated in spring and overestimated in late summer and autumn.The correlation coefficients between the observed and the modelled water temperature at these stations were high and exceeded 88% (Table 2) 2a, b).As shown by the observed data, the ice was present there in all three winter seasons analyzed.The observed and the modelled ice concentrations at Muuga, Sillamäe and Paldiski were scattered over the entire range 0.0-1.0; the respective correlation coefficients were 76%, 45% and 72% (Table 2).
The results for the stations Soru (Fig. 2c) and Lehtma (not shown), both at the entrance to the Gulf of Riga, were very similar, the modelled water temperatures being underestimated during spring and early summer and overestimated in late summer.As far as the sea ice is concerned, the modelling results showed only short episodes of the ice presence, www.ocean-sci.net/7/257/2011/Ocean Sci., 7, 257-276, 2011 as opposed to the observations revealing the ice to have been continuously present at those stations during the first two winter seasons analyzed.The model performed best at the station Pärnu (Fig. 2d), both in terms of the water temperatures and the length of the simulated ice season.The modelled ice thickness was roughly inversely proportional to the ice temperature.Generally, the simulated ice concentrations tended to be higher than the observed ones at stations located in the Gulf of Riga, and lower in the Gulf of Finland (bias q m in Table 2).At the two offshore stations NBPro and SBotSea (not shown), the water temperature variability was simulated very accurately (with correlation coefficients exceeding 98.4%).Apart from the first season in the period analyzed, when the The comparative statistics of the satellite-derived and simulated S i in the five basins analyzed are given in Table 3.The coefficients were calculated for the daily data between 1 January and 15 May during the three winters analyzed, for those data pairs in which at least one value was different from zero (number N p in Table 3).Figures 3 and 4 show the time series of the satellite-derived and modelled areas from different iceconcentration classes, during the winter 2005/2006.As can be seen, the model reproduces the general features of the ice variability in the Baltic Sea quite satisfactorily.The correlation coefficients r between the modelled and observed data are generally high for S 0-20 i and S 80-100 i , and substantially lower for the remaining concentration classes (Table 3).Similarly, the mean difference q m between the observations and the modelling results is negative for areas with the lowest and the highest ice concentrations, and positive otherwise.This feature can be partially attributed to the model deficiencies, and partially to inaccuracies of the satellite data: visual comparison with SMHI ice charts leads to the conclusion that areas with medium ice concentrations tend to be overestimated in the SIDADS data -see Fig. 5 and Supplement Figs. 1  and 2 for some examples.Additionally, in some situations and at some locations, the position of the ice edge, as compared with ice charts, is better reproduced by the model than by the SIDADS data (Fig. 5).Generally, the lowest-quality results were obtained for the Gulf of Finland and the Gulf of Riga, the highest-quality results being those concerning the more open basins: the Bothnian Sea and the Baltic Proper.166 * r (%): correlation coefficient; q m (10 3 km 2 ): mean difference between the observations and the model; σ rms (10 3 km 2 ): standard deviation of differences between the observations and the model; N p : number of data pairs.One of the general shortcomings of the model is its inability to reproduce properly the timing of ice melting in the spring months, especially in the Gulf of Finland and the Bothnian Bay, where the ice may remain at some locations until the second week of May (the SIDADS data agree very well with the SMHI ice charts for that period).The modelling results show the ice to disappear almost completely by mid-April, i.e., up to three weeks earlier than observed.Future versions of the model should definitely get rid of this undesired behaviour.
As already mentioned, and as shown in Figs. 3, 4 and Table 3, the model tends to overestimate/underestimate the high/low ice concentrations.After the onset of freezing in the early winter, which is reproduced with a reasonable accuracy, the simulated ice concentrations grow too fast, especially in smaller basins.For example, in the second half of January 2006, the observed and modelled total ice-covered Only the compact ice cover (mainly along the coasts) could be reproduced quite accurately.The same can be said about the Baltic Proper.It should be mentioned, however, that the SIDADS ice concentrations for mid-March are generally lower than those in the SMHI charts.The modelling results fall in between the two observation-based data sets (Supplement Fig. 2).
In addition to the ice concentrations, the model also produces realistic ice thickness and velocity distributions (Fig. 5 and Supplement Figs. 1 and 2).Typically, the simulated ice www.ocean-sci.net/7/257/2011/Ocean Sci., 7, 257-276, 2011 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 .4 0.5 0.6 0.7 0.8 0.9 1 thickness varies between 10 and 40 cm over most of the icecovered areas, with higher values occurring locally in some coastal locations (temporarily exceeding 1 m in narrow embayments during strong wind events).The modelled ice velocities were typically lower than 0.3 m s −1 ; only in open areas with low ice concentrations were higher values (up to 0.6 m s −1 ) observed.Although it is difficult to quantify this effect, inaccuracies in the meteorological data used to drive the model must have played some part in the emergence of the dis-crepancies between the modelled and observed sea ice distributions.The atmospheric model uses climatic seasurface temperature fields as the lower boundary conditions.Hence, the surface air temperatures are likely to be overestimated/underestimated in exceptionally cold/warm winters, which may -at least in part -explain the patterns observed in the modelled sea-ice cover, i.The meteorological conditions over the northern part of the Baltic Sea during the winter of 2005/2006 were typical of the region at that time of the year.A number of synoptic-scale low and high pressure systems, fronts, ridges and troughs passed over the study area, leading to a high variability of wind speed and direction, which resulted in the advection of air masses with different temperature, humidity, and cloudiness.Atmospheric conditions seldom remained fairly constant for more than one or two days.Considering a typical size of the synoptic-scale pressure systems over the Northern Europe, the fields of the atmospheric parameters can be regarded as time-varying, but spatially uniform within the Baltic Sea subbasins.The time series of the air temperature, atmospheric pressure, wind speed and direction in five locations representative of the basins analyzed here (stars in Fig. 1) are shown in Fig. 6.
A number of specific synoptic situations occurred during the analyzed period (see e.g.weather charts available at http://www.wetter3.de).The first was a cold-air outbreak on 19-21 January, associated with a high pressure system centered to the north-east of the Baltic Sea.Advection of an extremely cold Arctic air mass from the north of Asia led to a sudden temperature drop, especially over the Bothnian Bay (down to −27 • C) and the Gulf of Finland (down to −20 • C).Although the event was relatively short-lived and the air temperatures returned quickly to above −5 • C, it initiated rapid freezing in the Bothnian Bay, Gulf of Finland and the Gulf of Riga (Figs. 3 and 4), where it can be regarded as the onset of the ice season that year.Another cold period occurred in the first week of February, when a large cyclonic meander of the polar front occupied most of the North-Eastern Europe.A few smaller pressure systems moved over the Baltic Sea during that period, leading to variable and relatively strong winds from various directions.Generally, between the second week of January and the first week of February, the winds were strong, particularly in the Baltic Proper and the Bothnian Sea, where about ten events with wind speed over 10 m s −1 were recorded (Fig. 6c).Crucial for the sea-ice development was the period in mid-March.From 1 March, the air temperatures remained below 0 • C over the whole area analyzed, which -in combination with relatively calm weather -led to sea-ice formation not only in small basins, but also in the Bothnian Sea and Northern Baltic Proper (Supplement Fig. 2).Between 13 and 17 March, the Baltic Sea was under the influence of a stable high-pressure system located to the south of the main flow within the polar front.The winds were weak and, due to low cloudiness, a strong diurnal cycle of the air temperatures could be observed, particularly over smaller basins affected by the surrounding land.The total ice-covered area in the Baltic Sea reached its seasonal maximum (210×10 3 km 2 ) at the end of that period, which ended abruptly on 18 March, when a deep low developed over Northern Scandinavia.Within one day, the atmospheric pressure dropped from over 1030 hPa to below 995 hPa, and strong westerly and northwesterly winds (Fig. 6b-d) destroyed the ice cover, relatively uniform until that time, especially in the Bothnian Sea (Fig. 3b).The rest of the ice season was characterized by highly variable weather conditions.As of the second week of April, the air temperatures remained generally above 0 • C (except for the night-time temperatures in the Gulf of Finland), leading to gradual melting of the ice cover.

Sea-ice-meteorology relationships
Not surprisingly, the ice concentrations in the Baltic Sea depend strongly on the atmospheric forcing.In this section, the time-lagged correlations between S 0-20 i , . . ., S 80-100 i and the selected meteorological parameters representative of the five basins (see Fig. 1) are analyzed for all three winters in the simulation period.Prior to the analysis, the atmospheric data were smoothed with a three-day moving-average filter; the smoothed quantities are henceforth denoted with overbars.In terms of the correlation patterns obtained, the basins may be subdivided into two groups: one containing the larger, more open basins, i.e., the Bothnian Sea and the Baltic Proper; the other grouping the remaining, smaller ones.The results for the Bothnian Bay and the Baltic Proper are presented below as representative of those two groups of basins.
In the case of the air temperature T a , its relationship with S i is quite straightforward: during periods with low/high T a , the ice concentrations increase/decrease, which explains the correlation patterns in Fig. 7a, showing linear correlation coefficients between the five S i classes and T a for the Bothnian Bay and the Baltic Proper.The correlation coefficients vary slightly with the time lag, but generally exceed +0.6 for areas with high A i and −0.6 for areas with low A i or ice-free.The difference between the panels in Fig. 7a results from the fact that, whereas the Bothnian Bay is covered with ice throughout most of the winter season, the Baltic Proper is typically ice free.Thus, at low temperatures, the consolidation of the already existing ice cover in the Bothnian Bay and the formation of a new ice in the Baltic Proper are the dominant processes responsible for the patterns illustrated in Fig. 7a.Apart from T a , the smoothed diurnal temperature variations T a (with T a defined as the difference between the maximum and minimum T a on a given day) affect the ice cover, especially in the more open basins (Fig. 7b).Importantly, there is no statistically significant correlation between T a and T a , i.e., the influence of those parameters on S i may be regarded as independent of each other.
The influence of the wind is much more pronounced in the larger basins (the Bothnian Sea and the Baltic Proper) than in the smaller ones.Strong winds and wind-induced currents tend to break and fragment the ice cover, thus expanding/reducing areas with low/high A i .In the large basins, typically only partly covered with ice, the ice broken by winds tends to disperse over large, previously ice-free areas (especially by changing wind directions).This explains the positive correlations for S 0-20 i and the negative correlations for the remaining four S i classes (the right-hand panels in Fig. 7c, d).In smaller basins, with ice covering most of their area, the change from a positive to a negative correlation is shifted towards higher ice concentrations.Moreover, in the smaller basins, the wind direction seems to be more important for the correlations analyzed: in the Bothnian Bay the N-S wind component (i.e.approximately along the basin axis) has an understandably stronger influence than the W-E component (the left-hand panels in Fig. 7c, d).
Generally, the strong correlations (and clear-cut correlation patterns) obtained here suggest an interesting possibility of developing a multiple-regression or neural-network model capable of predicting the basin-wide ice conditions from the time series of atmospheric parameters.It would definitely be an appealing and relevant, practice-oriented, line of further research.Results are for all three winters in the simulation period.Correlations higher than 0.2 are significant at the 95% confidence level.

The role of ice dynamics
The ice in the Baltic Sea is relatively thin and thus can be easily deformed by dynamic processes such as winds, currents, waves, and water-level variations.During strong wind events, even areas with a consolidated ice cover are prone to the influence of the dynamic factors, as clearly shown by Uotila (2001), Wang et al. (2006) and others.Therefore, it is of crucial importance for a Baltic sea-ice model to reproduce the ice dynamics properly.This is succinctly demonstrated by the results of simulations performed with the icedynamics module turned off (Fig. 8).Even though the correlation coefficients between the observations and the results of the thermodynamic-only simulations are in some cases equal to or even higher than those yielded by the full sea-ice model (Tables 3 and 4), the dynamic-thermodynamic model produces a lower bias (q m ) and a lower standard deviation of differences (σ rms ).The positive role of the dynamic module is most pronounced in areas with compact ice (S 80-100 i ), where q m values are up to ten times lower than those obtained with the thermodynamic module alone (Table 4).
The smooth course of the S i curves is not surprising in view of the relationships between the meteorological forcing and the ice concentrations discussed in the previous section.Generally, dynamic processes lead to a higher variability of A i and h i , both spatially (large differences of A i and h i in neighbouring grid points) and temporally.A purely thermodynamic ice model produces h i and A i distributions that are spatially very uniform.This is in agreement with results reported by Meier et al. (1999) who performed a similar numerical experiment.Moreover, although no ice is produced or melted due to the dynamics alone, the dynamic redistribution of ice indirectly affects the temporal evolution of the total ice volume due to the modified surface fluxes and the total stress τ "experienced" by the hydrodynamic model.Each of the two factors leads to differing water temperatures T w in the upper layers of the model: changes of the surface fluxes influence T w directly, whereas changes in τ produce different surface circulation and mixing patterns.Moreover, redistribution of the ice in the h i and A i space modifies F tot,ai and F tot,iw , directly influencing the freezing/melting rates.The role of ice dynamics is well illustrated by the rate of change of the total ice volume V ice in a given basin -Fig.9 presents an example for the Bothnian Bay, obtained with the dynamic part of the ice model turned on and off.Temporal changes of V ice illustrate the problem better than the V ice itself, because they enable to easily identify periods with nearly thermodynamic ice growth/decay from periods during which the dynamics is important.In the case of the Bothnian Bay in the winter analyzed, thermodynamic processes clearly dominated during the cold-air outbreak around 20 January 2006 and during another cold period in early February (Fig. 6a).Between those two periods, the air temperatures were close to the freezing point and, most importantly, the wind direction changed from E-NE to W-NW (Fig. 6d), which was followed by a reversal of the surface currents (Fig. 10), whereby cracks and leads (seen both in the model results and observations) were produced in the ice cover in the northern part of the basin.The modelling results obtained with and without the ice dynamics differed in some important aspects.If the dynamics is ignored, the ice cover acts as a "lid" on the sea surface, cutting off the momentum flux from the atmosphere.As a consequence, the surface currents under the ice are very weak (usually under 10 cm s −1 ) and the influence of the ice edge on the water circulation can be clearly seen (white contours in Fig. 10).If the ice dynamics is included, the surface currents form a more uniform system encompassing the whole basin and the currents in the ice-covered zone are stronger, up to 30-35 cm s −1 in the situation analyzed (Fig. 10b, d).Differences in the circulation led to different sea surface temperatures, particularly after the wind direction reversal on 26 January 2006 when, if the dynamics is turned on, water with temperatures up to 0.3 • C appears in the northern part of the Bothnian Bay.Not surprisingly, without the dynamics, water temperature under the ice hardly changes during the period analyzed.Taken together, the two factors -mechanical fractionation and higher surface water temperatures -led to intensified bottom and lateral melting and, consequently, to the decrease of V ice that could not be modelled with the thermodynamic module alone (Fig. 9).A similar set of factors was responsible for the "destruction" of the ice cover in the second half of March.
The results presented here illustrate complex interrelationships between the thermodynamic and dynamic sea-ice processes.They are definitely sensitive to the details of model formulation (i.a., to the rate of the lateral melting), but the general conclusions are undoubtedly valid.Qualitatively, the results of the full thermodynamic-dynamic ice model show a better agreement with observations than the results produced by the purely thermodynamic model do.

Summary and conclusions
The numerical dynamic-thermodynamic sea-ice model presented in this study was applied to the analysis of short-term variability of the sea ice in the Baltic Sea.The study demonstrated strong relationships between the ice cover variability and the external forcing on the synoptic time scale.As already mentioned, these results suggest the feasibility of constructing a non-physical (e.g.neural network-based) ice model for the Baltic Sea, linking time series of selected meteorological parameters with time series of ice concentrations and volume.Another important aspect of the sea-ice phenomena in the Baltic Sea discussed in this paper concerns the role of ice dynamics.It was shown to influence not only the ice cover itself, but also the thermal regime of the surface water layers and the three-dimensional circulation patterns in the basins affected.
) * θ denotes the sun's position above the horizon (in degrees).

Fig. 1 .
Fig. 1.Bottom topography of the Baltic Sea (in m).Small black dots show the SIDADS-data grid.Thick lines mark the boundaries of the five sub-regions analyzed in the paper (BB: Bothnian Bay, BS: Bothnian Sea, BP: Northern Baltic Proper, GF: Gulf of Finland, GR: Gulf of Riga); asterisks denote the location of sites for which time series of the meteorological data representative of these regions were analyzed.The locations of the measurement stations used for model evaluation are shown with numbered, circled dots.See text for a detailed description.
Figure2shows the simulated and observed water temperature and sea-ice parameters at four selected stations.The simulated water temperatures at the Gulf of Finland's stations Muuga and Sillamäe (Fig.2a, b; the results for Paldiski were very similar to those at Muuga and are therefore not shown) were slightly underestimated in spring and overestimated in late summer and autumn.The correlation coefficients between the observed and the modelled water temperature at these stations were high and exceeded 88%(Table 2).The modelled ice thickness values ranged from 12 to 20 cm and increased eastwards, the highest values being obtained at Sillamäe.The model showed the ice at Paldiski and Muuga to have been present during the 2005/2006 winter only, whereas Sillamäe should have experienced ice during two winter seasons: 2005/2006 and 2006/2007 (Fig.2a, b).As shown by the observed data, the ice was present there in all three winter seasons analyzed.The observed and the modelled ice concentrations at Muuga, Sillamäe and Paldiski were scattered over the entire range 0.0-1.0; the respective correlation coefficients were 76%, 45% and 72% (Table2).The results for the stations Soru (Fig.2c) and Lehtma (not shown), both at the entrance to the Gulf of Riga, were very similar, the modelled water temperatures being underestimated during spring and early summer and overestimated in late summer.As far as the sea ice is concerned, the modelling results showed only short episodes of the ice presence,

Fig. 2 .
Fig. 2. Time series of the selected measured and modelled parameters at stations Muuga (a), Sillamäe (b), Soru (c) and Pärnu (d), see Table 2 and numbered sites in Fig. 1.The measured surface water temperatures were supplied by the BOOS data.The observed ice concentrations were extracted from the NSIDC data.

Fig. 3 .Fig. 4 .
Fig.3.Time series of the area (in 10 2 km 2 ) covered with ice concentrations higher than 0.2 (continuous lines), 0.4 (dashed lines), 0.6 (dotted lines) and 0.8 (dash-dotted lines) for the five analyzed subbasins of the Baltic Sea (a-e), marked in Fig.1.Values calculated from the satellite-derived (NSIDC) daily ice concentration fields interpolated onto the model grid.
e. the underestimated total S i during the rather cold winter of 2005/2006 and the overestimated S i during the extremely mild winter of 2007/2008.www.ocean-sci.net/7

Fig. 6 .
Fig. 6.Time series of the air temperature (a), atmospheric pressure (b), wind speed (c) and direction (d) in winter 2005/2006 at five sites marked with asterisks in Fig. 1.The temporal resolution is 1 h in (a, b) and 12 h in (c, d).
Fig. 7. Correlation between S 0-20 i , . . ., S 80-100 i (area covered with ice concentration from a given range, marked on the vertical axis) and the meteorological parameters in the Bothnian Bay (left) and the Baltic Proper (right) for time lags up to 14 days.The parameters are: the air temperature T a (a), the daily air temperature range T a (b) and the W-E (c) and S-N (d) components of the wind velocity, respectively.Results are for all three winters in the simulation period.Correlations higher than 0.2 are significant at the 95% confidence level.

Fig. 9 .Fig. 10 .
Fig. 9.The rate of change in the total ice volume V ice (in km 3 d −1 ) in the Bothnian Bay in winter 2005/2006, as obtained with the dynamicthermodynamic (continuous line) and purely thermodynamic (dashed line) ice model.

Table 1 .
Adjustable coefficients of the sea-ice model and their calibrated values.

Table 2 .
Measurement stations used for the model calibration and comparative statistics between the observations and modelling results for these stations (r: correlation coefficient in %; q m : mean difference between observations and model; N p : number of data points).Station numbers as in Fig.1.

Table 3 .
Comparative statistics between the satellite-derived and modelled S 0-20

Table 4 .
Comparative statistics between the satellite-derived and modelled S 0-20