Chapter 2 Model Description
2.2. LIM model
LIM is a dynamic-thermodynamic sea ice model designed specially for climate studies. Its thermodynamic part is based on Semtner’s three-layer model (1976), taking into account snow-ice formation and assuming brine pockets in the ice, and its dynamic part basically follows the viscous-plastic model of Hibler (1979). LIM has been used for sea ice studies in both polar regions through coupling with a one-dimensional upper ocean model (Fichefet and Morales Maqudea, 1997), for investigation of the importance of sea ice-ocean interaction for the global ocean circulation (Goosse and Fichefet, 1999), and for process studies in high latitudes (Timmermann et al., 2005) through coupling with a global, free surface OGCM.
2.2.1. Sea ice dynamics
The motion of sea ice can be described by a two-dimensional momentum equation. The advection term is usually neglected in the global climate studies according to the scaling arguments stated in the last part of this section.
u k u F water, Ca and Cw are the drag coefficients for air and water, θ is the turning angle in the ocean boundary layer beneath sea ice.
30 The internal force F is given by,
F=∇⋅σ (2.32)
in which σ is a two-dimensional internal ice stress tensor, depending on the relationship between stress and strain-rates. Hibler’s viscous-plastic rheology (1979) treats sea ice as linear viscous for very small strain rates and as plastic for large strain rates. This rheology has been widely used so far and is described as,
I
in which ec is the ratio of the principal axes of the elliptical yield curve, and,
where ε&0is a parameter termed “creep limit”, and D(ε&) is the determinant of the strain rate tensor. The compressive strength of ice P is given by,
P=P*Ahie−C(1−A) (2.35)
where P* and C are empirical constants, and A is the ice concentration representing the grid cell area covered by ice. The formulation above illustrates that the ice strength is strongly dependent on the amount of ice. When the ice becomes thicker, sea ice is strengthened.
Table 2.1. Scale analysis Internal friction Ph/L (P taking 10 kpa) 10-1
The air and water stresses are independent on the study scale and the ice thickness, and the other terms are all dependent on the ice thickness. In addition, the ice advection, sea surface tilt and internal friction are also related to the study scale (Table 2.1). To make scale analysis for the ice momentum equation, the following typical scales are selected: sea ice thickness hi=1 m, density ρi=900 kg/m3, grid spacing L=100 km, wind speed Ua=10 m/s, current velocity U=0.1 m/s, and f=10-4, water density ρw=1000 kg/m3, air density ρa=1.3 kg/m3. The air and water stresses are dominant forces together with the internal friction, in the order of 10-1. Coriolis acceleration and sea surface tilt are an order smaller (10-2), and ice advection is even smaller (10-4). Therefore, ice advection can be neglected.
The magnitudes of ice advection, Coriolis acceleration, sea surface tilt and internal friction are proportional to the ice thickness (Table 2.1). In the case of thin-ice cover (0.1 m), these terms drop down by one order of magnitude, while the air and water stresses keep unchanged.
This results in a governing balance between the air and water stresses. In the case of thick-ice
31
cover (10 m), the Coriolis acceleration, sea surface tilt, internal friction and ice advection increase by one order of magnitude. As a result, ice advection is still negligible (10-3). The internal friction is dominant in compact ice, while the Coriolis acceleration and sea surface tilt rank to the same order of magnitude as the air and water stresses. More discussion of the significance of the terms can be found in Leppäranta (2005).
2.2.2 Sea ice thermodynamics
Sea ice is thermodynamically considered to be a homogenous slab over which snow can accumulate. Because the temperature gradient is much smaller in the horizontal direction than in the vertical direction, the just vertical flow of heat is needed in the heat diffusion equation,
( )
st (including snow-ice), Tc is the ice temperature, and G(he) is a correction factor representing the effect of subgrid-scale snow and ice thickness distributions on the average heat conduction.In case of a uniform thickness distribution over the ice-covered part of the grid cell,
( )
thickness for heat conduction, defined as,
qst in Eq. (2.36) is an internal source term accounting for solar radiation penetrating into the ice when ice is free of snow. The penetrated solar radiation has an important effect on ice thickness (Maykut and Untersteiner, 1971). It is usually not for immediate surface melting, but for warming the subsurface interior of ice and initiating the internal melting (Semtner, 1976). The trapped brine pockets behave as a thermal reservoir to retard the heating or cooling of the ice (Maykut and Untersteiner, 1971). The effects of penetrated solar radiation and brine pockets in the ice are parameterized as (Fichefet and Morales Maqueda, 1997),
( )
in case of snow-covered ice (Maykut and Unstersteiner, 1971), and
i0=0.18
(
1−c)
+0.35c (2.41)in case of free of snow-covered ice (Ebert and Curry, 1993), where c is cloud cover. Flbp is the heat flux used to prevent the upper ice layer from cooling below the melting point of the ice.
The variation of ice thickness due to thermal effects is determined by the heat balance at the top surface and at the base of the ice. At the ice base, the variation of ice thickness is
where Li is the volumetric latent heat of ice freezing, Qcb and Qw are the conductive heat flux within the ice and the oceanic heat flux, respectively.
32
At the top surface, the heat balance Qsu is determined by,
Qsu =(1−i0)(1−αsu)Fsw+εsu
(
Flw −σTsu4)
+Fs+Fl+Fc (2.43) where the subscripts su is the surface, sw, lw, l and s have the same meaning as in Eq. (2.21).Thus, the second term on the right hand is the net long-wave radiation (Flwnet). ε is the emissivity, σ is the Stefan-Boltzmann constant, and Fc is the conductive heat flux from below the surface. The heat flux is assumed positive toward the surface. The energy imbalance reflects phase changes (freezing or melting) or warming or cooling (storage or release of heat) of the ice. When Qsu>0, the excessive heat is used for melting the top surface while keeping its temperature at the freezing point,
where Ls is the volumetric latent heat diffusion of snow. The heat balance at the top and the base of sea ice is summarized as in Fig. 2.2.
Fig. 2.2. Heat balance of sea ice.
Snow accumulated over the ice has two effects, as an insulator between ice and air and as a contributor to the ice growth through the snow-ice formation. When the mass of snow is heavy enough, the lower part of the snow layer submerges into the water. Slush develops and eventually forms snow-ice. Snow ice has remarkable contribution to the ice cover in the Baltic Sea (Leppäranta, 1983) and in the Southern Ocean (e.g., Massom et al., 2001 for review;
Worby et al., 2008). Neglecting the density difference between snow-ice and sea ice, the variation of snow and ice thickness after snow-ice formation is given by,
where the subscript si represents snow-ice.
The variation of A depends on the heat budget of the open water Bl. It is parameterized as, Snow Ts=Ti;
33 denoting as (Fichefet and Morales Maqueda, 1997),
Φ(A)=
(
1−A2)
1/2 (2.47)h0 is the thickness of newly formed ice in the leads.
2.2.3. Sea ice thickness categories and transportation
There has been a variety of approaches to describe ice thickness distribution in sea ice models. The classical one is the two-level approach (Doronin, 1971), which includes open water (ice concentration) and mean ice thickness.
The conservation law for ice/snow is,
while the conservation law for ice concentration is, iA SthA Di A SdyA
t
A =−∇⋅ + + ∇ +
∂
∂ (u ) 2 (2.49)
where Di is the horizontal diffusivity, used for numerical stability (Hibler, 1979; Fichefet and Morales Maqueda, 1997), S represents the source term due to thermodynamic (th) and shear deformation (dy) processes. The dynamic source Sdy in Eq. (2.48-2.49) is generally taken as zero in the standard two-level ice models, although the shear deformation contributes to the formation of open water and ice ridges.
2.2.4. Numerical methods and algorithms
Differing from the ocean model ORCA, Arakawa B grid is used for arrangement of the variables in space (Fig. 2.1b) in the sea ice model LIM, where the vectors are located at the corner of the grid cell and the scalars are situated at the center of the grid cell. Because of the development of numerical instabilities when snow or ice thickness becomes small (Fichefet and Morales Maqueda 1997), a fully implicit time stepping scheme is applied to solve the heat diffusion equation (2.36). In order to strictly satisfy the heat balance at the top surface of the snow-ice system, a Newton-Raphson iterative procedure is employed to solve Eq. (2.45).
Basically similar to that of Hibler (1979), a semi-implicit predictor-corrector procedure is adopted to solve the ice momentum balance Eq. (2.30) with a simultaneous under-relaxation technique. While for the ice transportation equations (2.48)–(2.49), a forward time-stepping scheme is used. To conserve the second-order moment of the spatial distribution of the advected quantities within each grid cell, the Prather technique (1986) is adopted to solve the advection term. More details can be found in Fichefet and Morales Maqueda (1997).