/[MITgcm]/manual/s_phys_pkgs/text/seaice.tex
ViewVC logotype

Diff of /manual/s_phys_pkgs/text/seaice.tex

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph | View Patch Patch

revision 1.21 by mlosch, Tue Apr 1 07:30:32 2014 UTC revision 1.25 by mlosch, Tue Mar 29 14:50:54 2016 UTC
# Line 63  no additional CPP options are required. Line 63  no additional CPP options are required.
63  Parts of the SEAICE code can be enabled or disabled at compile time  Parts of the SEAICE code can be enabled or disabled at compile time
64  via CPP preprocessor flags. These options are set in  via CPP preprocessor flags. These options are set in
65  \code{SEAICE\_OPTIONS.h}.  \code{SEAICE\_OPTIONS.h}.
66  Table \ref{tab:pkg:seaice:cpp} summarizes the most important ones.  Table \ref{tab:pkg:seaice:cpp} summarizes the most important ones. For
67    more options see the default \code{pkg/seaice/SEAICE\_OPTIONS.h}.
68    
69  \begin{table}[!ht]  \begin{table}[!ht]
70  \centering  \centering
# Line 527  non-linear iterations $\code{SEAICEnewto Line 528  non-linear iterations $\code{SEAICEnewto
528  number of Krylov iterations $\code{SEAICEkrylovIterMax} = 50$, because  number of Krylov iterations $\code{SEAICEkrylovIterMax} = 50$, because
529  the Krylov subspace has a fixed dimension of 50.  the Krylov subspace has a fixed dimension of 50.
530    
531    Setting \code{SEAICEuseStrImpCpl = .TRUE.,} turns on ``strength
532    implicit coupling'' \citep{hutchings04} in the LSR-solver and in the
533    LSR-preconditioner for the JFNK-solver. In this mode, the different
534    contributions of the stress divergence terms are re-ordered in order
535    to increase the diagonal dominance of the system
536    matrix. Unfortunately, the convergence rate of the LSR solver is
537    increased only slightly, while the JFNK-convergence appears to be
538    unaffected.
539    
540  \paragraph{Elastic-Viscous-Plastic (EVP) Dynamics\label{sec:pkg:seaice:EVPdynamics}}~\\  \paragraph{Elastic-Viscous-Plastic (EVP) Dynamics\label{sec:pkg:seaice:EVPdynamics}}~\\
541  %  %
542  \citet{hun97}'s introduced an elastic contribution to the strain  \citet{hun97}'s introduced an elastic contribution to the strain
# Line 592  the damping time scale $T$ accordingly, Line 602  the damping time scale $T$ accordingly,
602  $E_{0}\Delta{t}=\mbox{forcing time scale}$, or directly  $E_{0}\Delta{t}=\mbox{forcing time scale}$, or directly
603  \code{SEAICE\_evpTauRelax} ($T$) to the forcing time scale.  \code{SEAICE\_evpTauRelax} ($T$) to the forcing time scale.
604    
605    \paragraph{More stable variants of Elastic-Viscous-Plastic Dynamics:
606      EVP* , mEVP, and aEVP \label{sec:pkg:seaice:EVPstar}}~\\
607    %
608    The genuine EVP schemes appears to give noisy solutions \citep{hun01,
609      lemieux12, bouillon13}. This has lead to a modified EVP or EVP*
610    \citep{lemieux12, bouillon13, kimmritz15}; here, we refer to these
611    variants by modified EVP (mEVP) and adaptive EVP (aEVP)
612    \citep{kimmritz16}. The main idea is to modify the ``natural''
613    time-discretization of the momentum equations:
614    \begin{equation}
615      \label{eq:evpstar}
616      m\frac{D\vec{u}}{Dt} \approx m\frac{u^{p+1}-u^{n}}{\Delta{t}}
617      + \beta^{*}\frac{u^{p+1}-u^{p}}{\Delta{t}_{\mathrm{EVP}}}
618    \end{equation}
619    where $n$ is the previous time step index, and $p$ is the previous
620    sub-cycling index. The extra ``intertial'' term
621    $m\,(u^{p+1}-u^{n})/\Delta{t})$ allows the definition of a residual
622    $|u^{p+1}-u^{p}|$ that, as $u^{p+1} \rightarrow u^{n+1}$, converges to
623    $0$. In this way EVP can be re-interpreted as a pure iterative solver
624    where the sub-cycling has no association with time-relation (through
625    $\Delta{t}_{\mathrm{EVP}}$) \citep{bouillon13, kimmritz15}. Using the
626    terminology of \citet{kimmritz15}, the evolution equations of stress
627    $\sigma_{ij}$ and momentum $\vec{u}$ can be written as:
628    \begin{align}
629      \label{eq:evpstarsigma}
630      \sigma_{ij}^{p+1}&=\sigma_{ij}^p+\frac{1}{\alpha}
631      \Big(\sigma_{ij}(\vec{u}^p)-\sigma_{ij}^p\Big),
632      \phantom{\int}\\
633      \label{eq:evpstarmom}
634      \vec{u}^{p+1}&=\vec{u}^p+\frac{1}{\beta}
635      \Big(\frac{\Delta t}{m}\nabla \cdot{\bf \sigma}^{p+1}+
636      \frac{\Delta t}{m}\vec{R}^{p}+\vec{u}_n-\vec{u}^p\Big).
637    \end{align}
638    $\vec{R}$ contains all terms in the momentum equations except for the
639    rheology terms and the time derivative; $\alpha$ and $\beta$ are free
640    parameters (\code{SEAICE\_evpAlpha}, \code{SEAICE\_evpBeta}) that
641    replace the time stepping parameters \code{SEAICE\_deltaTevp}
642    ($\Delta{T}_{\mathrm{EVP}}$), \code{SEAICE\_elasticParm} ($E_{0}$), or
643    \code{SEAICE\_evpTauRelax} ($T$). $\alpha$ and $\beta$ determine the
644    speed of convergence and the stability. Usually, it makes sense to use
645    $\alpha = \beta$, and \code{SEAICEnEVPstarSteps} $\gg
646    (\alpha,\,\beta)$ \citep{kimmritz15}. Currently, there is no
647    termination criterion and the number of mEVP iterations is fixed to
648    \code{SEAICEnEVPstarSteps}.
649    
650    In order to use mEVP in the MITgcm, set \code{SEAICEuseEVPstar =
651      .TRUE.,} in \code{data.seaice}. If \code{SEAICEuseEVPrev =.TRUE.,}
652    the actual form of equations (\ref{eq:evpstarsigma}) and
653    (\ref{eq:evpstarmom}) is used with fewer implicit terms and the factor
654    of $e^{2}$ dropped in the stress equations (\ref{eq:evpstresstensor2})
655    and (\ref{eq:evpstresstensor12}). Although this modifies the original
656    EVP-equations, it turns out to improve convergence \citep{bouillon13}.
657    
658    Another variant is the aEVP scheme \citep{kimmritz16}, where the value
659    of $\alpha$ is set dynamically based on the stability criterion
660    \begin{equation}
661      \label{eq:aevpalpha}
662      \alpha = \beta = \max\left( \tilde{c}\pi\sqrt{c \frac{\zeta}{A_{c}}
663        \frac{\Delta{t}}{\max(m,10^{-4}\text{\,kg})}},\alpha_{\min} \right)
664    \end{equation}
665    with the grid cell area $A_c$ and the ice and snow mass $m$.  This
666    choice sacrifices speed of convergence for stability with the result
667    that aEVP converges quickly to VP where $\alpha$ can be small and more
668    slowly in areas where the equations are stiff. In practice, aEVP leads
669    to an overall better convergence than mEVP \citep{kimmritz16}.
670    %
671    To use aEVP in the MITgcm set \code{SEAICEaEVPcoeff} $= \tilde{c}$;
672    this also sets the default values of \code{SEAICEaEVPcStar} ($c=4$)
673    and \code{SEAICEaEVPalphaMin} ($\alpha_{\min}=5$). Good convergence
674    has been obtained with setting these values \citep{kimmritz16}:
675    \code{SEAICEaEVPcoeff = 0.5, SEAICEnEVPstarSteps = 500,
676      SEAICEuseEVPstar = .TRUE., SEAICEuseEVPrev = .TRUE.}
677    
678    Note, that probably because of the C-grid staggering of velocities and
679    stresses, mEVP may not converge as successfully as in
680    \citet{kimmritz15}, and that convergence at very high resolution
681    (order 5\,km) has not been studied yet.
682    
683  \paragraph{Truncated ellipse method (TEM) for yield curve \label{sec:pkg:seaice:TEM}}~\\  \paragraph{Truncated ellipse method (TEM) for yield curve \label{sec:pkg:seaice:TEM}}~\\
684  %  %
685  In the so-called truncated ellipse method the shear viscosity $\eta$  In the so-called truncated ellipse method the shear viscosity $\eta$
# Line 600  is capped to suppress any tensile stress Line 688  is capped to suppress any tensile stress
688    \label{eq:etatem}    \label{eq:etatem}
689    \eta = \min\left(\frac{\zeta}{e^2},    \eta = \min\left(\frac{\zeta}{e^2},
690    \frac{\frac{P}{2}-\zeta(\dot{\epsilon}_{11}+\dot{\epsilon}_{22})}    \frac{\frac{P}{2}-\zeta(\dot{\epsilon}_{11}+\dot{\epsilon}_{22})}
691    {\sqrt{(\dot{\epsilon}_{11}+\dot{\epsilon}_{22})^2    {\sqrt{\max(\Delta_{\min}^{2},(\dot{\epsilon}_{11}-\dot{\epsilon}_{22})^2
692        +4\dot{\epsilon}_{12}^2}}\right).        +4\dot{\epsilon}_{12}^2})}\right).
693  \end{equation}  \end{equation}
694  To enable this method, set \code{\#define SEAICE\_ALLOW\_TEM} in  To enable this method, set \code{\#define SEAICE\_ALLOW\_TEM} in
695  \code{SEAICE\_OPTIONS.h} and turn it on with  \code{SEAICE\_OPTIONS.h} and turn it on with
# Line 803  $\sigma_{21}^{Z}=0$, or equivalently $\e Line 891  $\sigma_{21}^{Z}=0$, or equivalently $\e
891    
892  \paragraph{Thermodynamics\label{sec:pkg:seaice:thermodynamics}}~\\  \paragraph{Thermodynamics\label{sec:pkg:seaice:thermodynamics}}~\\
893  %  %
894    \noindent\textbf{NOTE: THIS SECTION IS TERRIBLY OUT OF DATE}\\
895  In its original formulation the sea ice model \citep{menemenlis05}  In its original formulation the sea ice model \citep{menemenlis05}
896  uses simple thermodynamics following the appendix of  uses simple thermodynamics following the appendix of
897  \citet{sem76}. This formulation does not allow storage of heat,  \citet{sem76}. This formulation does not allow storage of heat,
# Line 818  way to that of \citet{parkinson79} and \ Line 907  way to that of \citet{parkinson79} and \
907  The conductive heat flux depends strongly on the ice thickness $h$.  The conductive heat flux depends strongly on the ice thickness $h$.
908  However, the ice thickness in the model represents a mean over a  However, the ice thickness in the model represents a mean over a
909  potentially very heterogeneous thickness distribution.  In order to  potentially very heterogeneous thickness distribution.  In order to
910  parameterize a sub-grid scale distribution for heat flux  parameterize a sub-grid scale distribution for heat flux computations,
911  computations, the mean ice thickness $h$ is split into seven thickness  the mean ice thickness $h$ is split into $N$ thickness categories
912  categories $H_{n}$ that are equally distributed between $2h$ and a  $H_{n}$ that are equally distributed between $2h$ and a minimum
913  minimum imposed ice thickness of $5\text{\,cm}$ by $H_n=  imposed ice thickness of $5\text{\,cm}$ by $H_n= \frac{2n-1}{7}\,h$
914  \frac{2n-1}{7}\,h$ for $n\in[1,7]$. The heat fluxes computed for each  for $n\in[1,N]$. The heat fluxes computed for each thickness category
915  thickness category is area-averaged to give the total heat flux  is area-averaged to give the total heat flux \citep{hibler84}. To use
916  \citep{hibler84}. To use this thickness category parameterization set  this thickness category parameterization set \code{SEAICE\_multDim} to
917  \code{\#define SEAICE\_MULTICATEGORY}; note that this requires  the number of desired categories (7 is a good guess, for anything
918  different restart files and switching this flag on in the middle of an  larger than 7 modify \code{SEAICE\_SIZE.h}) in
919  integration is not possible.  \code{data.seaice}; note that this requires different restart files
920    and switching this flag on in the middle of an integration is not
921    advised. In order to include the same distribution for snow, set
922    \code{SEAICE\_useMultDimSnow = .TRUE.}; only then, the
923    parameterization of always having a fraction of thin ice is efficient
924    and generally thicker ice is produced \citep{castro-morales14}.
925    
926    
927  The atmospheric heat flux is balanced by an oceanic heat flux from  The atmospheric heat flux is balanced by an oceanic heat flux from
928  below.  The oceanic flux is proportional to  below.  The oceanic flux is proportional to
Legend:
Removed from v.1.21  
changed lines
  Added in v.1.25
  ViewVC Help
Powered by ViewVC 1.1.22