/[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.22 by mlosch, Wed Jan 21 18:18:22 2015 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 variant of Elastic-Viscous-Plastic Dynamics:  EVP*\label{sec:pkg:seaice:EVPstar}}~\\  \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,  The genuine EVP schemes appears to give noisy solutions \citep{hun01,
609    lemieux12, bouillon13}. This has lead to a modified EVP or EVP*    lemieux12, bouillon13}. This has lead to a modified EVP or EVP*
610  \citep{lemieux12, bouillon13, kimmritz15}; here, refer to these  \citep{lemieux12, bouillon13, kimmritz15}; here, we refer to these
611  variants by EVP*. The main idea is to modify the ``natural''  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:  time-discretization of the momentum equations:
614  \begin{equation}  \begin{equation}
615    \label{eq:evpstar}    \label{eq:evpstar}
# Line 605  time-discretization of the momentum equa Line 617  time-discretization of the momentum equa
617    + \beta^{*}\frac{u^{p+1}-u^{p}}{\Delta{t}_{\mathrm{EVP}}}    + \beta^{*}\frac{u^{p+1}-u^{p}}{\Delta{t}_{\mathrm{EVP}}}
618  \end{equation}  \end{equation}
619  where $n$ is the previous time step index, and $p$ is the previous  where $n$ is the previous time step index, and $p$ is the previous
620  sub-cycling index. The term allows the definition of a residual  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  $|u^{p+1}-u^{p}|$ that, as $u^{p+1} \rightarrow u^{n+1}$, converges to
623  $0$ and a re-interpretation of EVP as a pure iterative solver where  $0$. In this way EVP can be re-interpreted as a pure iterative solver
624  the sub-cycling has lost all time-relation \citep{bouillon13,  where the sub-cycling has no association with time-relation (through
625    kimmritz15}. Using the terminology of \citet{kimmritz15}, the  $\Delta{t}_{\mathrm{EVP}}$) \citep{bouillon13, kimmritz15}. Using the
626  evolution equations of stress $\sigma_{ij}$ and momentum $\vec{u}$ can  terminology of \citet{kimmritz15}, the evolution equations of stress
627  be written as:  $\sigma_{ij}$ and momentum $\vec{u}$ can be written as:
628  \begin{align}  \begin{align}
629    \label{eq:evpstarsigma}    \label{eq:evpstarsigma}
630    \sigma_{ij}^{p+1}&=\sigma_{ij}^p+\frac{1}{\alpha}    \sigma_{ij}^{p+1}&=\sigma_{ij}^p+\frac{1}{\alpha}
# Line 620  be written as: Line 633  be written as:
633    \label{eq:evpstarmom}    \label{eq:evpstarmom}
634    \vec{u}^{p+1}&=\vec{u}^p+\frac{1}{\beta}    \vec{u}^{p+1}&=\vec{u}^p+\frac{1}{\beta}
635    \Big(\frac{\Delta t}{m}\nabla \cdot{\bf \sigma}^{p+1}+    \Big(\frac{\Delta t}{m}\nabla \cdot{\bf \sigma}^{p+1}+
636    \frac{\Delta t}{m}\vec{R}^{p+1/2}+\vec{u}_n-\vec{u}^p\Big).    \frac{\Delta t}{m}\vec{R}^{p}+\vec{u}_n-\vec{u}^p\Big).
637  \end{align}  \end{align}
638  $\vec{R}$ contains all terms in the momentum equations except for the  $\vec{R}$ contains all terms in the momentum equations except for the
639  rheology terms and the time derivative, $\alpha$ and $\beta$ are free  rheology terms and the time derivative; $\alpha$ and $\beta$ are free
640  parameters (\code{SEAICE\_evpAlpha}, \code{SEAICE\_evpBeta}) that  parameters (\code{SEAICE\_evpAlpha}, \code{SEAICE\_evpBeta}) that
641  replace the time stepping parameters \code{SEAICE\_deltaTevp}  replace the time stepping parameters \code{SEAICE\_deltaTevp}
642  ($\Delta{T}_{\mathrm{EVP}}$), \code{SEAICE\_elasticParm} ($E_{0}$), or  ($\Delta{T}_{\mathrm{EVP}}$), \code{SEAICE\_elasticParm} ($E_{0}$), or
643  \code{SEAICE\_evpTauRelax} ($T$). $\alpha$ and $\beta$ determine the  \code{SEAICE\_evpTauRelax} ($T$). $\alpha$ and $\beta$ determine the
644  speed of convergence and the stability. Usually, it makes sense to use  speed of convergence and the stability. Usually, it makes sense to use
645  $\alpha = \beta$, and \code{SEAICEnEVPstarSteps} $>> \alpha = \beta$  $\alpha = \beta$, and \code{SEAICEnEVPstarSteps} $\gg
646  \citep{kimmritz15}.  (\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 EVP* in the MITgcm, set \code{SEAICEuseEVPstar =  In order to use mEVP in the MITgcm, set \code{SEAICEuseEVPstar =
651    .TRUE.,} in \code{data.seaice}. \code{SEAICEuseEVPrev =.TRUE.,} uses    .TRUE.,} in \code{data.seaice}. If \code{SEAICEuseEVPrev =.TRUE.,}
652  the actual form of equations (\ref{eq:evpstarsigma}) and  the actual form of equations (\ref{eq:evpstarsigma}) and
653  (\ref{eq:evpstarmom}) with fewer implicit terms and the factor of  (\ref{eq:evpstarmom}) is used with fewer implicit terms and the factor
654  $e^{2}$ dropped in the stress equations (\ref{eq:evpstresstensor2})  of $e^{2}$ dropped in the stress equations (\ref{eq:evpstresstensor2})
655  and (\ref{eq:evpstresstensor12}).  This turns out to improve  and (\ref{eq:evpstresstensor12}). Although this modifies the original
656  convergence \citep{bouillon13}.  EVP-equations, it turns out to improve convergence \citep{bouillon13}.
657    
658  Note, that for historical reasons, \code{SEAICE\_deltaTevp} needs to  Another variant is the aEVP scheme \citep{kimmritz16}, where the value
659  be set to some value in order to use also EVP*. Also note, that  of $\alpha$ is set dynamically based on the stability criterion
660  probably because of the C-grid staggering of velocities and stresses,  \begin{equation}
661  EVP* does not converge as successfully as in \citet{kimmritz15}.    \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  %  %
# Line 654  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 857  $\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 872  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.22  
changed lines
  Added in v.1.25
  ViewVC Help
Powered by ViewVC 1.1.22