Changeset 11388 for NEMO/trunk/doc/latex/NEMO

Timestamp:

2019-08-01T18:09:07+02:00 (5 years ago)

Author:

acc

Message:

Add adaptive-implicit vertical advection section to ZDF chapter. May still need a few tweaks next week and refers to an updated traadv_fct.F90 which is not yet submitted. Finals tests are in progress

File:

• NEMO/trunk/doc/latex/NEMO/subfiles/chap_ZDF.tex (modified) (1 diff)

NEMO/trunk/doc/latex/NEMO/subfiles/chap_ZDF.tex

@@ -1269 +1269 @@
 \label{subsec:ZDF_aimp}
 
-This refers to \citep{shchepetkin_OM15}.
-
-TBC ...
-
+The adaptive-implicit vertical advection option in NEMO is based on the work of
+\citep{shchepetkin_OM15}.  In common with most ocean models, the timestep used with NEMO
+needs to satisfy multiple criteria associated with different physical processes in order
+to maintain numerical stability. \citep{shchepetkin_OM15} pointed out that the vertical
+CFL criterion is commonly the most limiting. \citep{lemarie.debreu.ea_OM15} examined the
+constraints for a range of time and space discretizations and provide the CFL stability
+criteria for a range of advection schemes. The values for the Leap-Frog with Robert
+asselin filter time-stepping (as used in NEMO) are reproduced in
+\autoref{tab:zad_Aimp_CFLcrit}. Treating the vertical advection implicitly can avoid these
+restrictions but at the cost of large dispersive errors and, possibly, large numerical
+viscosity. The adaptive-implicit vertical advection option provides a targetted use of the
+implicit scheme only when and where potential breaches of the vertical CFL condition
+occur. In many practical applications these events may occur remote from the main area of
+interest or due to short-lived conditions such that the extra numerical diffusion or
+viscosity does not greatly affect the overall solution. With such applications, setting:
+\forcode{ln_zad_Aimp = .true.} should allow much longer model timesteps to be used whilst
+retaining the accuracy of the high order explicit schemes over most of the domain.
+
+\begin{table}[htbp]
+  \begin{center}
+    % \begin{tabular}{cp{70pt}cp{70pt}cp{70pt}cp{70pt}}
+    \begin{tabular}{r|ccc}
+      \hline
+      spatial discretization &   2nd order centered   & 3rd order upwind & 4th order compact  \\
+      advective CFL criterion     & 0.904 &   0.472  &   0.522    \\
+      \hline
+    \end{tabular}
+    \caption{
+      \protect\label{tab:zad_Aimp_CFLcrit}
+      The advective CFL criteria for a range of spatial discretizations for the Leap-Frog with Robert Asselin filter time-stepping
+      ($\nu=0.1$) as given in \citep{lemarie.debreu.ea_OM15}.
+    }
+  \end{center}
+\end{table}
+
+In particular, the advection scheme remains explicit everywhere except where and when
+local vertical velocities exceed a threshold set just below the explicit stability limit.
+Once the threshold is reached a tapered transition towards an implicit scheme is used by
+partitioning the vertical velocity into a part that can be treated explicitly and any
+excess that must be treated implicitly. The partitioning is achieved via a Courant-number
+dependent weighting algorithm as described in \citep{shchepetkin_OM15}.
+
+The local cell Courant number ($Cu$) used for this partitioning is:
+
+\begin{equation}
+  \label{eq:Eqn_zad_Aimp_Courant}
+  \begin{split}
+    Cu &= {2 \rdt \over e^n_{3t_{ijk}}} \bigg (\big [ \texttt{Max}(w^n_{ijk},0.0) - \texttt{Min}(w^n_{ijk+1},0.0) \big ]    \\
+       &\phantom{=} +\big [ \texttt{Max}(e_{{2_u}ij}e^n_{{3_{u}}ijk}u^n_{ijk},0.0) - \texttt{Min}(e_{{2_u}i-1j}e^n_{{3_{u}}i-1jk}u^n_{i-1jk},0.0) \big ]
+                     \big / e_{{1_t}ij}e_{{2_t}ij}            \\
+       &\phantom{=} +\big [ \texttt{Max}(e_{{1_v}ij}e^n_{{3_{v}}ijk}v^n_{ijk},0.0) - \texttt{Min}(e_{{1_v}ij-1}e^n_{{3_{v}}ij-1k}v^n_{ij-1k},0.0) \big ]
+                     \big / e_{{1_t}ij}e_{{2_t}ij} \bigg )    \\
+  \end{split}
+\end{equation}
+
+\noindent and the tapering algorithm follows \citep{shchepetkin_OM15} as:
+
+\begin{align}
+  \label{eq:Eqn_zad_Aimp_partition}
+Cu_{min} &= 0.15 \nonumber \\
+Cu_{max} &= 0.27 \nonumber \\
+Cu_{cut} &= 2Cu_{max} - Cu_{min} \nonumber \\
+Fcu    &= 4Cu_{max}*(Cu_{max}-Cu_{min}) \nonumber \\
+Cf &=
+     \begin{cases}
+        0.0                                                        &\text{if $Cu \leq Cu_{min}$} \\
+        (Cu - Cu_{min})^2 / (Fcu +  (Cu - Cu_{min})^2)             &\text{else if $Cu < Cu_{cut}$} \\
+        (Cu - Cu_{max}) / Cu                                       &\text{else}
+     \end{cases}
+\end{align}
+
+\noindent With these settings the coefficient ($Cf$) is shown in \autoref{fig:zad_Aimp_coeff}
+
+\begin{figure}[!t]
+  \begin{center}
+    \includegraphics[width=\textwidth]{Fig_ZDF_zad_Aimp_coeff}
+    \caption{
+      \protect\label{fig:zad_Aimp_coeff}
+      The value of the partitioning coefficient ($Cf$) used to partition vertical velocities into parts to
+      be treated implicitly and explicitly for a range of typical Courant numbers (\forcode{ln_zad_Aimp=.true.})
+    }
+  \end{center}
+\end{figure}
+
+\noindent The partitioning coefficient is used to determine the part of the vertical
+velocity that must be handled implicitly ($w_i$) and to subtract this from the total
+vertical velocity ($w_n$) to leave that which can continue to be handled explicitly:
+
+\begin{align}
+  \label{eq:Eqn_zad_Aimp_partition2}
+    w_{i_{ijk}} &= Cf_{ijk} w_{n_{ijk}}     \nonumber \\
+    w_{n_{ijk}} &= (1-Cf_{ijk}) w_{n_{ijk}}           
+\end{align}
+
+\noindent Note that the coefficient is such that the treatment is never fully implicit;
+the three cases from \autoref{eq:Eqn_zad_Aimp_partition} can be considered as:
+fully-explicit; mixed explicit/implicit and mostly-implicit.  The $w_i$ component is added
+to the implicit solvers for the vertical mixing in \mdl{dynzdf} and \mdl{trazdf} in a
+similar way to \citep{shchepetkin_OM15}.  This is sufficient for the flux-limited
+advection scheme (\forcode{ln_traadv_mus}) but further intervention is required when using
+the flux-corrected scheme (\forcode{ln_traadv_fct}). For these schemes the implicit
+upstream fluxes must be added to both the monotonic guess and to the higher order solution
+when calculating the antidiffusive fluxes. The implicit vertical fluxes are then removed
+since they are added by the implicit solver later on.
+
+The adaptive-implicit vertical advection option is new to NEMO at v4.0 and has yet to be 
+used in a wide range of simulations. The following test simulation, however, does illustrate
+the potential benefits and will hopefully encourage further testing and feedback from users:
+
+\subsection{Adaptive-implicit vertical advection in the OVERFLOW test-case}
+The \href{https://forge.ipsl.jussieu.fr/nemo/chrome/site/doc/NEMO/guide/html/test\_cases.html\#overflow}{OVERFLOW test case}
+provides a simple illustration of the adaptive-implicit advection in action. The example here differs from the basic test case
+by only a few extra physics choices namely:
+
+\begin{verbatim}
+     ln_dynldf_OFF = .false.
+     ln_dynldf_lap = .true.
+     ln_dynldf_hor = .true.
+     ln_zdfnpc     = .true.
+     ln_traadv_fct = .true.
+        nn_fct_h   =  4
+        nn_fct_v   =  4
+\end{verbatim}
+
+\noindent which were chosen to provide a slightly more stable and less noisy solution. The
+result when using the default value of \forcode{nn_rdt = 10.} without adaptive-implicit
+vertical velocity is illustrated in \autoref{fig:zad_Aimp_overflow_frames}. The mass of
+cold water, initially sitting on the shelf, moves down the slope and forms a
+bottom-trapped, dense plume. Even with these extra physics choices the model is close to
+stability limits and attempts with \forcode{nn_rdt = 30.} will fail after about 5.5 hours
+with excessively high horizontal velocities. This time-scale corresponds with the time the
+plume reaches the steepest part of the topography and, although detected as a horizontal
+CFL breach, the instability originates from a breach of the vertical CFL limit. This is a good
+candidate, therefore, for use of the adaptive-implicit vertical advection scheme.
+
+The results with \forcode{ln_zad_Aimp=.true.} and a variety of model timesteps are shown
+in \autoref{fig:zad_Aimp_overflow_all_rdt} (together with the equivalent frames from the
+base run).  In this simple example the use of the adaptive-implicit vertcal advection
+scheme has enabled a 12x increase in the model timestep without significantly altering the
+solution (although at this extreme the plume is more diffuse and has not travelled so far).
+Notably, the solution with and without the scheme is slightly different even with
+\forcode{nn_rdt = 10.}; suggesting that the base run was close enough to instability to
+trigger the scheme despite completing successfully.  To assist in diagnosing how active
+the scheme is, in both location and time, the 3D implicit and explicit components of the
+vertical velocity are available via XIOS as \texttt{wimp} and \texttt{wexp} respectively.
+Likewise, the partitioning coefficient ($Cf$) is also available as \texttt{wi\_cff}.
+
+\begin{figure}[!t]
+  \begin{center}
+    \includegraphics[width=\textwidth]{Fig_ZDF_zad_Aimp_overflow_frames}
+    \caption{
+      \protect\label{fig:zad_Aimp_overflow_frames}
+      A time-series of temperature vertical cross-sections for the OVERFLOW test case. These results are for the default
+      settings with \forcode{nn_rdt=10.0} and without adaptive implicit vertical advection (\forcode{ln_zad_Aimp=.false.}).
+    }
+  \end{center}
+\end{figure}
+
+\begin{figure}[!t]
+  \begin{center}
+    \includegraphics[width=\textwidth]{Fig_ZDF_zad_Aimp_overflow_all_rdt}
+    \caption{
+      \protect\label{fig:zad_Aimp_overflow_all_rdt}
+      Sample temperature vertical cross-sections from mid- and end-run using different values for \forcode{nn_rdt} 
+      and with or without adaptive implicit vertical advection. Without the adaptive implicit vertical advection only
+      the run with the shortest timestep is able to run to completion. Note also that the colour-scale has been
+      chosen to confirm that temperatures remain within the original range of 10$^o$ to 20$^o$.
+    }
+  \end{center}
+\end{figure}