Sea-Ice Model

Sea-Ice Model

The sea-ice component of ROMS is a combination of the elastic-viscous-plastic (EVP) rheology ( Hunke and Dukowicz, 1997, Hunke, 2001) and simple one-layer ice and snow thermodynamics with a molecular sublayer under the ice ( Mellor and Kantha, 1989). It is tightly coupled, having the same grid (Arakawa-C) and timestep as the ocean and sharing the same parallel coding structure for use with MPI or OpenMP.

Dynamics

The momentum equations describe the change in ice/snow velocity due to the combined effects of the Coriolis force, surface ocean tilt, air and water stress, and internal ice stress:

$M{du \over dt}=Mfv-Mg{\partial \zeta _{w} \over \partial x}+\tau _{a}^{x}+\tau _{w}^{x}+{F}_{x}$

$M{dv \over dt}=-Mfu-Mg{\partial \zeta _{w} \over \partial y}+\tau _{a}^{y}+\tau _{w}^{y}+{F}_{y}$

In this model, we neglect the nonlinear advection terms as well as the curvilinear terms in the internal ice stress. The force due to the internal ice stress is given by the divergence of the stress tensor $\sigma \!\,$. The rheology is given by the stress-strain relation of the medium. We would like to emulate the viscous-plastic rheology of Hibler (1979):

$\sigma _{ij}=2\eta {\dot {\epsilon }}_{ij}+(\zeta -\eta ){\dot {\epsilon }}_{kk}\delta _{ij}-{P \over 2}\delta _{ij}$

${\dot {\epsilon }}_{ij}\equiv {1 \over 2}\left({\partial u_{i} \over \partial x_{j}}+{\partial u_{j} \over \partial x_{i}}\right)$

$P=P^{*}Ah_{i}e^{-C(1-A)}$

while having an explicit model that can be solved efficiently on parallel computers. The EVP rheology has a tunable coefficient E (the Young's modulus) which can be chosen to make the elastic term small compared to the other terms. We rearrange the VP rheology:

${1 \over 2\eta }\sigma _{ij}+{\eta -\zeta \over 4\eta \zeta }\sigma _{kk}\delta _{ij}+{P \over 4\zeta }\delta _{ij}={\dot {\epsilon }}_{ij}$

then add the elastic term:

${1 \over E}{\partial \sigma _{ij} \over \partial t}+{1 \over 2\eta }\sigma _{ij}+{\eta -\zeta \over 4\eta \zeta }\sigma _{kk}\delta _{ij}+{P \over 4\zeta }\delta _{ij}={\dot {\epsilon }}_{ij}$

Much like the ocean model, the ice model has a split timestep. The internal ice stress term is updated on a shorter timestep so as to allow the elastic wave velocity to be resolved.

Once the new ice velocities are computed, the ice tracers can be advected using the MPDATA scheme ( Smolarkiewicz). The tracers in this case are the ice thickness, ice concentration, snow thickness, internal ice temperature, and surface melt ponds.

The ice model variables:

Name	Description
$u,v\!\,$	Horizontal ice velocity components
$M\!\,$	Ice mass, $\rho _{i}Ah_{i}\!\,$
$A\!\,$	Ice concentration, $0.0\leq A\leq 1.0\!\,$
$h_{i}\!\,$	Ice thickness
$f\!\,$	Coriolis parameter
$g\!\,$	Gravity
$\zeta _{w}\!\,$	Surface elevation of the underlying water
$\tau \!\,$	Air and water stresses
$F_{x},F_{y}\!\,$	Forces due to internal ice stress
$\sigma _{ij}\!\,$	Internal ice stress tensor
$\epsilon _{ij}\!\,$	Strainrate tensor
$\delta _{ij}\!\,$	Kronecker delta function
$\zeta ,\eta \!\,$	Nonlinear ice viscosities
$P\!\,$	Ice strength
$P^{*},C\!\,$	Ice strength parameters
$E\!\,$	Young's modulus

Thermodynamics

The thermodynamics is based on calculating how much ice grows and melts on each of the surface, bottom, and sides of the ice floes, as well as frazil ice formation:

 [Figure with W_xx]

Once the ice tracers are advected, the ice concentration and thickness are timestepped according to the terms on the right:

${DAh_{i} \over Dt}={\rho _{o} \over \rho _{i}}\left[A(W_{io}-W_{ai})+(1-A)W_{ao}+W_{fr}\right]$

${DA \over Dt}={\rho _{o}A \over \rho _{i}h_{i}}\left[\Phi (1-A)w_{ao}+(1-A)W_{fr}\right]\qquad \qquad 0\leq A\leq 1$

The term $Ah_{i}\!\,$ is the "effective thickness", a measure of the ice volume. Its evolution equation is simply quantifying the change in the amount of ice. The ice concentration equation is more interesting in that it provides the partitioning between ice melt/growth on the sides vs. on the top and bottom. The parameter $\Phi \!\,$ controls this and has differing values for ice melt and retreat. In principle, most of the ice growth is assumed to happen at the base of the ice while rather more of the melt happens on the sides of the ice due to warming of the water in the leads.

The heat fluxes through the ice are based on a simple one layer Semtner (1976) type model with snow on top. The temperature is assumed to be linear within the snow and within the ice. The ice contains brine pockets for a total ice salinity of 5. The surface ocean temperature and salinity is half a dz below the surface. The water right below the surface is assumed to be at the freezing temperature; a logarithmic boundary layer is computed having the temperature and salinity matched at freezing.

[Note: I have gobs more here...]