52 |
\maketitle |
\maketitle |
53 |
|
|
54 |
\begin{abstract} |
\begin{abstract} |
55 |
Some blabla |
|
56 |
|
As part of ongoing efforts to obtain a best possible synthesis of most |
57 |
|
available, global-scale, ocean and sea ice data, a dynamic and thermodynamic |
58 |
|
sea-ice model has been coupled to the Massachusetts Institute of Technology |
59 |
|
general circulation model (MITgcm). Ice mechanics follow a viscous plastic |
60 |
|
rheology and the ice momentum equations are solved numerically using either |
61 |
|
line successive relaxation (LSR) or elastic-viscous-plastic (EVP) dynamic |
62 |
|
models. Ice thermodynamics are represented using either a zero-heat-capacity |
63 |
|
formulation or a two-layer formulation that conserves enthalpy. The model |
64 |
|
includes prognostic variables for snow and for sea-ice salinity. The above |
65 |
|
sea ice model components were borrowed from current-generation climate models |
66 |
|
but they were reformulated on an Arakawa C-grid in order to match the MITgcm |
67 |
|
oceanic grid and they were modified in many ways to permit efficient and |
68 |
|
accurate automatic differentiation. This paper describes the MITgcm sea ice |
69 |
|
model; it presents example Arctic and Antarctic results from a realistic, |
70 |
|
eddy-permitting, global ocean and sea-ice configuration; it compares B-grid |
71 |
|
and C-grid dynamic solvers in a regional Arctic configuration; and it presents |
72 |
|
example results from coupled ocean and sea-ice adjoint-model integrations. |
73 |
|
|
74 |
\end{abstract} |
\end{abstract} |
75 |
|
|
76 |
\section{Introduction} |
\section{Introduction} |
77 |
\label{sec:intro} |
\label{sec:intro} |
78 |
|
|
|
more blabla |
|
|
|
|
79 |
\section{Model} |
\section{Model} |
80 |
\label{sec:model} |
\label{sec:model} |
81 |
|
|
327 |
state variables to be advected by ice velocities, namely enthalphy of |
state variables to be advected by ice velocities, namely enthalphy of |
328 |
the two ice layers and the thickness of the overlying snow layer. |
the two ice layers and the thickness of the overlying snow layer. |
329 |
|
|
|
\section{Funnel Experiments} |
|
|
\label{sec:funnel} |
|
|
|
|
|
For a first/detailed comparison between the different variants of the |
|
|
MIT sea ice model an idealized geometry of a periodic channel, |
|
|
1000\,km long and 500\,m wide on a non-rotating plane, with converging |
|
|
walls forming a symmetric funnel and a narrow strait of 40\,km width |
|
|
is used. The horizontal resolution is 5\,km throughout the domain |
|
|
making the narrow strait 8 grid points wide. The ice model is |
|
|
initialized with a complete ice cover of 50\,cm uniform thickness. The |
|
|
ice model is driven by a constant along channel eastward ocean current |
|
|
of 25\,cm/s that does not see the walls in the domain. All other |
|
|
ice-ocean-atmosphere interactions are turned off, in particular there |
|
|
is no feedback of ice dynamics on the ocean current. All thermodynamic |
|
|
processes are turned off so that ice thickness variations are only |
|
|
caused by convergent or divergent ice flow. Ice volume (effective |
|
|
thickness) and concentration are advected with a third-order scheme |
|
|
with a flux limiter \citep{hundsdorfer94} to avoid undershoots. This |
|
|
scheme is unconditionally stable and does not require additional |
|
|
diffusion. The time step used here is 1\,h. |
|
|
|
|
|
\reffig{funnelf0} compares the dynamic fields ice concentration $c$, |
|
|
effective thickness $h_{eff} = h\cdot{c}$, and velocities $(u,v)$ for |
|
|
five different cases at steady state (after 10\,years of integration): |
|
|
\begin{description} |
|
|
\item[B-LSRns:] LSR solver with no-slip boundary conditions on a B-grid, |
|
|
\item[C-LSRns:] LSR solver with no-slip boundary conditions on a C-grid, |
|
|
\item[C-LSRfs:] LSR solver with free-slip boundary conditions on a C-grid, |
|
|
\item[C-EVPns:] EVP solver with no-slip boundary conditions on a C-grid, |
|
|
\item[C-EVPfs:] EVP solver with free-slip boundary conditions on a C-grid, |
|
|
\end{description} |
|
|
\ml{[We have not implemented the EVP solver on a B-grid.]} |
|
|
\begin{figure*}[htbp] |
|
|
\includegraphics[width=\widefigwidth]{\fpath/all_086280} |
|
|
\caption{Ice concentration, effective thickness [m], and ice |
|
|
velocities [m/s] |
|
|
for 5 different numerical solutions.} |
|
|
\label{fig:funnelf0} |
|
|
\end{figure*} |
|
|
At a first glance, the solutions look similar. This is encouraging as |
|
|
the details of discretization and numerics should not affect the |
|
|
solutions to first order. In all cases the ice-ocean stress pushes the |
|
|
ice cover eastwards, where it converges in the funnel. In the narrow |
|
|
channel the ice moves quickly (nearly free drift) and leaves the |
|
|
channel as narrow band. |
|
|
|
|
|
A close look reveals interesting differences between the B- and C-grid |
|
|
results. The zonal velocity in the narrow channel is nearly the free |
|
|
drift velocity ( = ocean velocity) of 25\,cm/s for the C-grid |
|
|
solutions, regardless of the boundary conditions, while it is just |
|
|
above 20\,cm/s for the B-grid solution. The ice accelerates to |
|
|
25\,cm/s after it exits the channel. Concentrating on the solutions |
|
|
B-LSRns and C-LSRns, the ice volume (effective thickness) along the |
|
|
boundaries in the narrow channel is larger in the B-grid case although |
|
|
the ice concentration is reduces in the C-grid case. The combined |
|
|
effect leads to a larger actual ice thickness at smaller |
|
|
concentrations in the C-grid case. However, since the effective |
|
|
thickness determines the ice strength $P$ in Eq\refeq{icestrength}, |
|
|
the ice strength and thus the bulk and shear viscosities are larger in |
|
|
the B-grid case leading to more horizontal friction. This circumstance |
|
|
might explain why the no-slip boundary conditions in the B-grid case |
|
|
appear to be more effective in reducing the flow within the narrow |
|
|
channel, than in the C-grid case. Further, the viscosities are also |
|
|
sensitive to details of the velocity gradients. Via $\Delta$, these |
|
|
gradients enter the viscosities in the denominator so that large |
|
|
gradients tend to reduce the viscosities. This again favors more flow |
|
|
along the boundaries in the C-grid case: larger velocities |
|
|
(\reffig{funnelf0}) on grid points that are closer to the boundary by |
|
|
a factor $\frac{1}{2}$ than in the B-grid case because of the stagger |
|
|
nature of the C-grid lead numerically to larger tangential gradients |
|
|
across the boundary; these in turn make the viscosities smaller for |
|
|
less tangential friction and allow more tangential flow along the |
|
|
boundaries. |
|
|
|
|
|
The above argument can also be invoked to explain the small |
|
|
differences between the free-slip and no-slip solutions on the C-grid. |
|
|
Because of the non-linearities in the ice viscosities, in particular |
|
|
along the boundaries, the no-slip boundary conditions have only a small |
|
|
impact on the solution. |
|
|
|
|
|
The difference between LSR and EVP solutions is largest in the |
|
|
effective thickness and meridional velocity fields. The EVP velocity |
|
|
fields appears to be a little noisy. This noise has been address by |
|
|
\citet{hunke01}. It can be dealt with by reducing EVP's internal time |
|
|
step (increasing the number of iterations along with the computational |
|
|
cost) or by regularizing the bulk and shear viscosities. We revisit |
|
|
the latter option by reproducing some of the results of |
|
|
\citet{hunke01}, namely the experiment described in her section~4, for |
|
|
our C-grid no-slip cases: in a square domain with a few islands the |
|
|
ice model is initialized with constant ice thickness and linearly |
|
|
increasing ice concentration to the east. The model dynamics are |
|
|
forced with a constant anticyclonic ocean gyre and by variable |
|
|
atmospheric wind whose mean direction is diagnonal to the north-east |
|
|
corner of the domain; ice volume and concentration are held constant |
|
|
(no thermodynamics and no advection by ice velocity). |
|
|
\reffig{hunke01} shows the ice velocity field, its divergence, and the |
|
|
bulk viscosity $\zeta$ for the cases C-LRSns and C-EVPns, and for a |
|
|
C-EVPns case, where \citet{hunke01}'s regularization has been |
|
|
implemented; compare to Fig.\,4 in \citet{hunke01}. The regularization |
|
|
contraint limits ice strength and viscosities as a function of damping |
|
|
time scale, resolution and EVP-time step, effectively allowing the |
|
|
elastic waves to damp out more quickly \citep{hunke01}. |
|
|
\begin{figure*}[htbp] |
|
|
\includegraphics[width=\widefigwidth]{\fpath/hun12days} |
|
|
\caption{Ice flow, divergence and bulk viscosities of three |
|
|
experiments with \citet{hunke01}'s test case: C-LSRns (top), |
|
|
C-EVPns (middle), and C-EVPns with damping described in |
|
|
\citet{hunke01} (bottom).} |
|
|
\label{fig:hunke01} |
|
|
\end{figure*} |
|
|
|
|
|
In the far right (``east'') side of the domain the ice concentration |
|
|
is close to one and the ice should be nearly rigid. The applied wind |
|
|
tends to push ice toward the upper right corner. Because the highly |
|
|
compact ice is confined by the boundary, it resists any further |
|
|
compression and exhibits little motion in the rigid region on the |
|
|
right hand side. The C-LSRns solution (top row) allows high |
|
|
viscosities in the rigid region suppressing nearly all flow. |
|
|
\citet{hunke01}'s regularization for the C-EVPns solution (bottom row) |
|
|
clearly suppresses the noise present in $\nabla\cdot\vek{u}$ and |
|
|
$\log_{10}\zeta$ in the |
|
|
unregularized case (middle row), at the cost of reduced viscosities. |
|
|
These reduced viscosities lead to small but finite ice velocities |
|
|
which in turn can have a strong effect on solutions in the limit of |
|
|
nearly rigid regimes (arching and blocking, not shown). |
|
|
|
|
|
\ml{[Say something about performance? This is tricky, as the |
|
|
perfomance depends strongly on the configuration. A run with slowly |
|
|
changing forcing is favorable for LSR, because then only very few |
|
|
iterations are required for convergences while EVP uses its fixed |
|
|
number of internal timesteps. If the forcing in changing fast, LSR |
|
|
needs far more iterations while EVP still uses the fixed number of |
|
|
internal timesteps. I have produces runs where for slow forcing LSR |
|
|
is much faster than EVP and for fast forcing, LSR is much slower |
|
|
than EVP. EVP is certainly more efficient in terms of vectorization |
|
|
and MFLOPS on our SX8, but is that a criterion?]} |
|
330 |
|
|
331 |
\subsection{C-grid} |
\subsection{C-grid} |
332 |
\begin{itemize} |
\begin{itemize} |
373 |
\subsection{Arctic Domain with Open Boundaries} |
\subsection{Arctic Domain with Open Boundaries} |
374 |
\label{sec:arctic} |
\label{sec:arctic} |
375 |
|
|
376 |
The Arctic domain of integration is illustrated in Fig.~\ref{???}. It |
The Arctic domain of integration is illustrated in Fig.~\ref{fig:arctic1}. It |
377 |
is carved out from, and obtains open boundary conditions from, the |
is carved out from, and obtains open boundary conditions from, the |
378 |
global cubed-sphere configuration of the Estimating the Circulation |
global cubed-sphere configuration of the Estimating the Circulation |
379 |
and Climate of the Ocean, Phase II (ECCO2) project |
and Climate of the Ocean, Phase II (ECCO2) project |
380 |
\citet{menemenlis05}. The domain size is 420 by 384 grid boxes |
\citet{menemenlis05}. The domain size is 420 by 384 grid boxes |
381 |
horizontally with mean horizontal grid spacing of 18 km. |
horizontally with mean horizontal grid spacing of 18 km. |
382 |
|
|
383 |
|
\begin{figure} |
384 |
|
%\centerline{{\includegraphics*[width=0.44\linewidth]{\fpath/arctic1.eps}}} |
385 |
|
\caption{Bathymetry of Arctic Domain.\label{fig:arctic1}} |
386 |
|
\end{figure} |
387 |
|
|
388 |
There are 50 vertical levels ranging in thickness from 10 m near the surface |
There are 50 vertical levels ranging in thickness from 10 m near the surface |
389 |
to approximately 450 m at a maximum model depth of 6150 m. Bathymetry is from |
to approximately 450 m at a maximum model depth of 6150 m. Bathymetry is from |
390 |
the National Geophysical Data Center (NGDC) 2-minute gridded global relief |
the National Geophysical Data Center (NGDC) 2-minute gridded global relief |