[2282] | 1 | % ================================================================ |
---|
| 2 | % Iso-neutral diffusion : |
---|
| 3 | % ================================================================ |
---|
[2285] | 4 | \chapter{Griffies's iso-neutral diffusion} |
---|
[2282] | 5 | \label{Apdx_C} |
---|
| 6 | \minitoc |
---|
| 7 | |
---|
[2285] | 8 | \section{Griffies's formulation of iso-neutral diffusion} |
---|
[2282] | 9 | |
---|
| 10 | \subsection{Introduction} |
---|
| 11 | We define a scheme that get its inspiration from the scheme of |
---|
| 12 | \citet{Griffies_al_JPO98}, but is formulated within the \NEMO |
---|
| 13 | framework, using scale factors rather than grid-sizes. |
---|
| 14 | |
---|
| 15 | The off-diagonal terms of the small angle diffusion tensor |
---|
| 16 | \eqref{Eq_PE_iso_tensor} produce skew-fluxes along the |
---|
| 17 | i- and j-directions resulting from the vertical tracer gradient: |
---|
| 18 | \begin{align} |
---|
| 19 | \label{eq:i13c} |
---|
| 20 | &+\kappa r_1\frac{1}{e_3}\frac{\partial T}{\partial k},\qquad+\kappa r_2\frac{1}{e_3}\frac{\partial T}{\partial k}\\ |
---|
| 21 | \intertext{and in the k-direction resulting from the lateral tracer gradients} |
---|
| 22 | \label{eq:i31c} |
---|
| 23 | & \kappa r_1\frac{1}{e_1}\frac{\partial T}{\partial i}+\kappa r_2\frac{1}{e_1}\frac{\partial T}{\partial i} |
---|
| 24 | \end{align} |
---|
| 25 | where \eqref{Eq_PE_iso_slopes} |
---|
| 26 | \begin{align*} |
---|
| 27 | r_1 &=-\frac{e_3 }{e_1 } \left( \frac{\partial \rho }{\partial i} |
---|
| 28 | \right) |
---|
| 29 | \left( {\frac{\partial \rho }{\partial k}} \right)^{-1} \\ |
---|
| 30 | &=-\frac{e_3 }{e_1 } \left( -\alpha\frac{\partial T }{\partial i} + |
---|
| 31 | \beta\frac{\partial S }{\partial i} \right) \left( |
---|
| 32 | -\alpha\frac{\partial T }{\partial k} + \beta\frac{\partial S |
---|
| 33 | }{\partial k} \right)^{-1} |
---|
| 34 | \end{align*} |
---|
| 35 | is the i-component of the slope of the isoneutral surface relative to the computational |
---|
| 36 | surface, and $r_2$ is the j-component. |
---|
| 37 | |
---|
| 38 | The extra vertical diffusive flux associated with the $_{33}$ |
---|
| 39 | component of the small angle diffusion tensor is |
---|
| 40 | \begin{equation} |
---|
| 41 | \label{eq:i33c} |
---|
| 42 | -\kappa(r_1^2 +r_2^2) \frac{1}{e_3}\frac{\partial T}{\partial k}. |
---|
| 43 | \end{equation} |
---|
| 44 | |
---|
| 45 | Since there are no cross terms involving $r_1$ and $r_2$ in the above, we can |
---|
| 46 | consider the isoneutral diffusive fluxes separately in the i-k and j-k |
---|
| 47 | planes, just adding together the vertical components from each |
---|
| 48 | plane. The following description will describe the fluxes on the i-k |
---|
| 49 | plane. |
---|
| 50 | |
---|
| 51 | There is no natural discretization for the i-component of the |
---|
| 52 | skew-flux, \eqref{eq:i13c}, as |
---|
| 53 | although it must be evaluated at u-points, it involves vertical |
---|
| 54 | gradients (both for the tracer and the slope $r_1$), defined at |
---|
| 55 | w-points. Similarly, the vertical skew flux, \eqref{eq:i31c}, is evaluated at |
---|
| 56 | w-points but involves horizontal gradients defined at u-points. |
---|
| 57 | |
---|
| 58 | \subsection{The standard discretization} |
---|
| 59 | The straightforward approach to discretize the lateral skew flux |
---|
| 60 | \eqref{eq:i13c} from tracer cell $i,k$ to $i+1,k$, introduced in 1995 |
---|
| 61 | into OPA, \eqref{Eq_tra_ldf_iso}, is to calculate a mean vertical |
---|
| 62 | gradient at the u-point from the average of the four surrounding |
---|
| 63 | vertical tracer gradients, and multiply this by a mean slope at the |
---|
| 64 | u-point, calculated from the averaged surrounding vertical density |
---|
| 65 | gradients. The total area-integrated skew-flux from tracer cell $i,k$ |
---|
| 66 | to $i+1,k$, noting that the $e_{{3u}_{i+1/2}^k}$ in the area |
---|
| 67 | $e_{{3u}_{i+1/2}^k}e_{{2u}_{i+1/2}^k}$ at the u-point cancels out with |
---|
| 68 | the $1/e_{{3u}_{i+1/2}^k}$ associated with the vertical tracer |
---|
| 69 | gradient, is then \eqref{Eq_tra_ldf_iso} |
---|
| 70 | \begin{equation*} |
---|
| 71 | \left(F_u^{\mathrm{skew}} \right)_{i+\hhalf}^k = \kappa _{i+\hhalf}^k |
---|
| 72 | e_{{2u}_{i+1/2}^k} \overline{\overline |
---|
| 73 | r_1} ^{\,i,k}\,\overline{\overline{\delta_k T}}^{\,i,k}, |
---|
| 74 | \end{equation*} |
---|
| 75 | where |
---|
| 76 | \begin{equation*} |
---|
| 77 | \overline{\overline |
---|
| 78 | r_1} ^{\,i,k} = -\frac{e_{{3u}_{i+1/2}^k}}{e_{{1u}_{i+1/2}^k}} |
---|
| 79 | \frac{\delta_{i+1/2} [\rho]}{\overline{\overline{\delta_k \rho}}^{\,i,k}}. |
---|
| 80 | \end{equation*} |
---|
| 81 | |
---|
| 82 | Unfortunately the resulting combination $\overline{\overline{\delta_k |
---|
| 83 | \bullet}}^{\,i,k}$ of a $k$ average and a $k$ difference %of the tracer |
---|
| 84 | reduces to $\bullet_{k+1}-\bullet_{k-1}$, so two-grid-point oscillations are |
---|
| 85 | invisible to this discretization of the iso-neutral operator. These |
---|
| 86 | \emph{computational modes} will not be damped by this operator, and |
---|
| 87 | may even possibly be amplified by it. Consequently, applying this |
---|
| 88 | operator to a tracer does not guarantee the decrease of its |
---|
| 89 | global-average variance. To correct this, we introduced a smoothing of |
---|
| 90 | the slopes of the iso-neutral surfaces (see \S\ref{LDF}). This |
---|
| 91 | technique works fine for $T$ and $S$ as they are active tracers |
---|
| 92 | ($i.e.$ they enter the computation of density), but it does not work |
---|
| 93 | for a passive tracer. |
---|
| 94 | \subsection{Expression of the skew-flux in terms of triad slopes} |
---|
| 95 | \citep{Griffies_al_JPO98} introduce a different discretization of the |
---|
| 96 | off-diagonal terms that nicely solves the problem. |
---|
| 97 | % Instead of multiplying the mean slope calculated at the u-point by |
---|
| 98 | % the mean vertical gradient at the u-point, |
---|
| 99 | % >>>>>>>>>>>>>>>>>>>>>>>>>>>> |
---|
| 100 | \begin{figure}[h] \begin{center} |
---|
| 101 | \includegraphics[width=0.90\textwidth]{./TexFiles/Figures/Fig_triad_fluxes} |
---|
[2376] | 102 | \caption{ \label{Fig_ISO_triad} |
---|
| 103 | (a) Arrangement of triads $S_i$ and tracer gradients to |
---|
| 104 | give lateral tracer flux from box $i,k$ to $i+1,k$ |
---|
| 105 | (b) Triads $S'_i$ and tracer gradients to give vertical tracer flux from |
---|
| 106 | box $i,k$ to $i,k+1$.} |
---|
[2282] | 107 | \label{Fig_triad} |
---|
| 108 | \end{center} \end{figure} |
---|
| 109 | % >>>>>>>>>>>>>>>>>>>>>>>>>>>> |
---|
| 110 | They get the skew flux from the products of the vertical gradients at |
---|
| 111 | each w-point surrounding the u-point with the corresponding `triad' |
---|
| 112 | slope calculated from the lateral density gradient across the u-point |
---|
| 113 | divided by the vertical density gradient at the same w-point as the |
---|
| 114 | tracer gradient. See Fig.~\ref{Fig_triad}a, where the thick lines |
---|
| 115 | denote the tracer gradients, and the thin lines the corresponding |
---|
| 116 | triads, with slopes $s_1, \dotsc s_4$. The total area-integrated |
---|
| 117 | skew-flux from tracer cell $i,k$ to $i+1,k$ |
---|
| 118 | \begin{multline} |
---|
| 119 | \label{eq:i13} |
---|
| 120 | \left( F_u^{13} \right)_{i+\frac{1}{2}}^k = \kappa _{i+1}^k a_1 s_1 |
---|
| 121 | \delta _{k+\frac{1}{2}} \left[ T^{i+1} |
---|
| 122 | \right]/e_{{3w}_{i+1}}^{k+\frac{1}{2}} + \kappa _i^k a_2 s_2 \delta |
---|
| 123 | _{k+\frac{1}{2}} \left[ T^i |
---|
| 124 | \right]/e_{{3w}_{i+1}}^{k+\frac{1}{2}} \\ |
---|
| 125 | +\kappa _{i+1}^k a_3 s_3 \delta _{k-\frac{1}{2}} \left[ T^{i+1} |
---|
| 126 | \right]/e_{{3w}_{i+1}}^{k+\frac{1}{2}} +\kappa _i^k a_4 s_4 \delta |
---|
| 127 | _{k-\frac{1}{2}} \left[ T^i \right]/e_{{3w}_{i+1}}^{k+\frac{1}{2}}, |
---|
| 128 | \end{multline} |
---|
| 129 | where the contributions of the triad fluxes are weighted by areas |
---|
| 130 | $a_1, \dotsc a_4$, and $\kappa$ is now defined at the tracer points |
---|
| 131 | rather than the u-points. This discretization gives a much closer |
---|
| 132 | stencil, and disallows the two-point computational modes. |
---|
| 133 | |
---|
| 134 | The vertical skew flux \eqref{eq:i31c} from tracer cell $i,k$ to $i,k+1$ at the |
---|
| 135 | w-point $i,k+\hhalf$ is constructed similarly (Fig.~\ref{Fig_triad}b) |
---|
| 136 | by multiplying lateral tracer gradients from each of the four |
---|
| 137 | surrounding u-points by the appropriate triad slope: |
---|
| 138 | \begin{multline} |
---|
| 139 | \label{eq:i31} |
---|
| 140 | \left( F_w^{31} \right) _i ^{k+\frac{1}{2}} = \kappa_i^{k+1} a_{1}' |
---|
| 141 | s_{1}' \delta _{i-\frac{1}{2}} \left[ T^{k+1} \right]/{e_{3u}}_{i-\frac{1}{2}}^{k+1} |
---|
| 142 | +\kappa_i^{k+1} a_{2}' s_{2}' \delta _{i+\frac{1}{2}} \left[ T^{k+1} \right]/{e_{3u}}_{i+\frac{1}{2}}^{k+1}\\ |
---|
| 143 | + \kappa_i^k a_{3}' s_{3}' \delta _{i-\frac{1}{2}} \left[ T^k\right]/{e_{3u}}_{i-\frac{1}{2}}^k |
---|
| 144 | +\kappa_i^k a_{4}' s_{4}' \delta _{i+\frac{1}{2}} \left[ T^k \right]/{e_{3u}}_{i+\frac{1}{2}}^k. |
---|
| 145 | \end{multline} |
---|
| 146 | |
---|
| 147 | We notate the triad slopes in terms of the `anchor point' $i,k$ |
---|
| 148 | (appearing in both the vertical and lateral gradient), and the u- and |
---|
| 149 | w-points $(i+i_p,k)$, $(i,k+k_p)$ at the centres of the `arms' of the |
---|
| 150 | triad as follows (see also Fig.~\ref{Fig_triad}): |
---|
| 151 | \begin{equation} |
---|
| 152 | \label{Gf_slopes} |
---|
| 153 | _i^k \mathbb{R}_{i_p}^{k_p} |
---|
| 154 | =\frac{ {e_{3w}}_{\,i}^{\,k+k_p}} { {e_{1u}}_{\,i+i_p}^{\,k}} |
---|
| 155 | \ |
---|
| 156 | \frac |
---|
| 157 | {\left(\alpha / \beta \right)_i^k \ \delta_{i + i_p}[T^k] - \delta_{i + i_p}[S^k] } |
---|
| 158 | {\left(\alpha / \beta \right)_i^k \ \delta_{k+k_p}[T^i ] - \delta_{k+k_p}[S^i ] }. |
---|
| 159 | \end{equation} |
---|
| 160 | In calculating the slopes of the local neutral |
---|
| 161 | surfaces, the expansion coefficients $\alpha$ and $\beta$ are |
---|
| 162 | evaluated at the anchor points of the triad \footnote{Note that in \eqref{Gf_slopes} we use the ratio $\alpha / \beta$ |
---|
| 163 | instead of multiplying the temperature derivative by $\alpha$ and the |
---|
| 164 | salinity derivative by $\beta$. This is more efficient as the ratio |
---|
| 165 | $\alpha / \beta$ can to be evaluated directly}, while the metrics are |
---|
| 166 | calculated at the u- and w-points on the arms. |
---|
| 167 | |
---|
| 168 | % >>>>>>>>>>>>>>>>>>>>>>>>>>>> |
---|
| 169 | \begin{figure}[h] \begin{center} |
---|
| 170 | \includegraphics[width=0.60\textwidth]{./TexFiles/Figures/Fig_qcells} |
---|
[2376] | 171 | \caption{ \label{Fig_ISO_triad_notation} |
---|
| 172 | Triad notation for quarter cells.T-cells are inside |
---|
[2282] | 173 | boxes, while the $i+\half,k$ u-cell is shaded in green and the |
---|
| 174 | $i,k+\half$ w-cell is shaded in pink.} |
---|
| 175 | \label{qcells} |
---|
| 176 | \end{center} \end{figure} |
---|
| 177 | % >>>>>>>>>>>>>>>>>>>>>>>>>>>> |
---|
| 178 | |
---|
| 179 | Each triad $\{_i^k\:_{i_p}^{k_p}\}$ is associated (Fig.~\ref{qcells}) with the quarter |
---|
| 180 | cell that is the intersection of the $i,k$ T-cell, the $i+i_p,k$ |
---|
| 181 | u-cell and the $i,k+k_p$ w-cell. Expressing the slopes $s_i$ and |
---|
| 182 | $s'_i$ in \eqref{eq:i13} and \eqref{eq:i31} in this notation, we have |
---|
| 183 | e.g.\ $s_1=s'_1={\:}_i^k \mathbb{R}_{1/2}^{1/2}$. Each triad slope $_i^k |
---|
| 184 | \mathbb{R}_{i_p}^{k_p}$ is used once (as an $s$) to calculate the |
---|
| 185 | lateral flux along its u-arm, at $(i+i_p,k)$, and then again as an |
---|
| 186 | $s'$ to calculate the vertical flux along its w-arm at |
---|
| 187 | $(i,k+k_p)$. Each vertical area $a_i$ used to calculate the lateral |
---|
| 188 | flux and horizontal area $a'_i$ used to calculate the vertical flux |
---|
| 189 | can also be identified as the area across the u- and w-arms of a |
---|
| 190 | unique triad, and we can notate these areas, similarly to the triad |
---|
| 191 | slopes, as $_i^k{\mathbb{A}_u}_{i_p}^{k_p}$, |
---|
| 192 | $_i^k{\mathbb{A}_w}_{i_p}^{k_p}$, where e.g. in \eqref{eq:i13} |
---|
| 193 | $a_{1}={\:}_i^k{\mathbb{A}_u}_{1/2}^{1/2}$, and in \eqref{eq:i31} |
---|
| 194 | $a'_{1}={\:}_i^k{\mathbb{A}_w}_{1/2}^{1/2}$. |
---|
| 195 | |
---|
| 196 | \subsection{The full triad fluxes} |
---|
| 197 | A key property of isoneutral diffusion is that it should not affect |
---|
| 198 | the (locally referenced) density. In particular there should be no |
---|
| 199 | lateral or vertical density flux. The lateral density flux disappears so long as the |
---|
| 200 | area-integrated lateral diffusive flux from tracer cell $i,k$ to |
---|
| 201 | $i+1,k$ coming from the $_{11}$ term of the diffusion tensor takes the |
---|
| 202 | form |
---|
| 203 | \begin{equation} |
---|
| 204 | \label{eq:i11} |
---|
| 205 | \left( F_u^{11} \right) _{i+\frac{1}{2}} ^{k} = |
---|
| 206 | - \left( \kappa_i^{k+1} a_{1} + \kappa_i^{k+1} a_{2} + \kappa_i^k |
---|
| 207 | a_{3} + \kappa_i^k a_{4} \right) |
---|
| 208 | \frac{\delta _{i+1/2} \left[ T^k\right]}{{e_{1u}}_{\,i+1/2}^{\,k}}, |
---|
| 209 | \end{equation} |
---|
| 210 | where the areas $a_i$ are as in \eqref{eq:i13}. In this case, |
---|
| 211 | separating the total lateral flux, the sum of \eqref{eq:i13} and |
---|
| 212 | \eqref{eq:i11}, into triad components, a lateral tracer |
---|
| 213 | flux |
---|
| 214 | \begin{equation} |
---|
| 215 | \label{eq:latflux-triad} |
---|
| 216 | _i^k {\mathbb{F}_u}_{i_p}^{k_p} (T) = -A_i^k{ \:}_i^k{\mathbb{A}_u}_{i_p}^{k_p} |
---|
| 217 | \left( |
---|
| 218 | \frac{ \delta_{i+ i_p}[T^k] }{ {e_{1u}}_{\,i + i_p}^{\,k} } |
---|
| 219 | -\ {_i^k\mathbb{R}_{i_p}^{k_p}} \ |
---|
| 220 | \frac{ \delta_{k+k_p} [T^i] }{{e_{3w}}_{\,i}^{\,k+k_p} } |
---|
| 221 | \right) |
---|
| 222 | \end{equation} |
---|
| 223 | can be identified with each triad. Then, because the |
---|
| 224 | same metric factors ${e_{3w}}_{\,i}^{\,k+k_p}$ and |
---|
| 225 | ${e_{1u}}_{\,i+i_p}^{\,k}$ are employed for both the density gradients |
---|
| 226 | in $ _i^k \mathbb{R}_{i_p}^{k_p}$ and the tracer gradients, the lateral |
---|
| 227 | density flux associated with each triad separately disappears. |
---|
| 228 | \begin{equation} |
---|
| 229 | \label{eq:latflux-rho} |
---|
| 230 | {\mathbb{F}_u}_{i_p}^{k_p} (\rho)=-\alpha _i^k {\:}_i^k {\mathbb{F}_u}_{i_p}^{k_p} (T) + \beta_i^k {\:}_i^k {\mathbb{F}_u}_{i_p}^{k_p} (S)=0 |
---|
| 231 | \end{equation} |
---|
| 232 | Thus the total flux $\left( F_u^{31} \right) ^i _{i,k+\frac{1}{2}} + |
---|
| 233 | \left( F_u^{11} \right) ^i _{i,k+\frac{1}{2}}$ from tracer cell $i,k$ |
---|
| 234 | to $i+1,k$ must also vanish since it is a sum of four such triad fluxes. |
---|
| 235 | |
---|
| 236 | The squared slope $r_1^2$ in the expression \eqref{eq:i33c} for the |
---|
| 237 | $_{33}$ component is also expressed in terms of area-weighted |
---|
| 238 | squared triad slopes, so the area-integrated vertical flux from tracer |
---|
| 239 | cell $i,k$ to $i,k+1$ resulting from the $r_1^2$ term is |
---|
| 240 | \begin{equation} |
---|
| 241 | \label{eq:i33} |
---|
| 242 | \left( F_w^{33} \right) _i^{k+\frac{1}{2}} = |
---|
| 243 | - \left( \kappa_i^{k+1} a_{1}' s_{1}'^2 |
---|
| 244 | + \kappa_i^{k+1} a_{2}' s_{2}'^2 |
---|
| 245 | + \kappa_i^k a_{3}' s_{3}'^2 |
---|
| 246 | + \kappa_i^k a_{4}' s_{4}'^2 \right)\delta_{k+\frac{1}{2}} \left[ T^{i+1} \right], |
---|
| 247 | \end{equation} |
---|
| 248 | where the areas $a'$ and slopes $s'$ are the same as in |
---|
| 249 | \eqref{eq:i31}. |
---|
| 250 | Then, separating the total vertical flux, the sum of \eqref{eq:i31} and |
---|
| 251 | \eqref{eq:i33}, into triad components, a vertical flux |
---|
| 252 | \begin{align} |
---|
| 253 | \label{eq:vertflux-triad} |
---|
| 254 | _i^k {\mathbb{F}_w}_{i_p}^{k_p} (T) |
---|
| 255 | &= A_i^k{\: }_i^k{\mathbb{A}_w}_{i_p}^{k_p} |
---|
| 256 | \left( |
---|
| 257 | {_i^k\mathbb{R}_{i_p}^{k_p}}\frac{ \delta_{i+ i_p}[T^k] }{ {e_{1u}}_{\,i + i_p}^{\,k} } |
---|
| 258 | -\ \left({_i^k\mathbb{R}_{i_p}^{k_p}}\right)^2 \ |
---|
| 259 | \frac{ \delta_{k+k_p} [T^i] }{{e_{3w}}_{\,i}^{\,k+k_p} } |
---|
| 260 | \right) \\ |
---|
| 261 | &= - \left(\left.{ }_i^k{\mathbb{A}_w}_{i_p}^{k_p}\right/{ }_i^k{\mathbb{A}_u}_{i_p}^{k_p}\right) |
---|
| 262 | {_i^k\mathbb{R}_{i_p}^{k_p}}{\: }_i^k{\mathbb{F}_u}_{i_p}^{k_p} (T) \label{eq:vertflux-triad2} |
---|
| 263 | \end{align} |
---|
| 264 | may be associated with each triad. Each vertical density flux $_i^k {\mathbb{F}_w}_{i_p}^{k_p} (\rho)$ |
---|
| 265 | associated with a triad then separately disappears (because the |
---|
| 266 | lateral flux $_i^k{\mathbb{F}_u}_{i_p}^{k_p} (\rho)$ |
---|
| 267 | disappears). Consequently the total vertical density flux $\left( F_w^{31} \right)_i ^{k+\frac{1}{2}} + |
---|
| 268 | \left( F_w^{33} \right)_i^{k+\frac{1}{2}}$ from tracer cell $i,k$ |
---|
| 269 | to $i,k+1$ must also vanish since it is a sum of four such triad |
---|
| 270 | fluxes. |
---|
| 271 | |
---|
| 272 | We can explicitly identify (Fig.~\ref{qcells}) the triads associated with the $s_i$, $a_i$, and $s'_i$, $a'_i$ used in the definition of |
---|
| 273 | the u-fluxes and w-fluxes in |
---|
| 274 | \eqref{eq:i31}, \eqref{eq:i13}, \eqref{eq:i11} \eqref{eq:i33} and |
---|
| 275 | Fig.~\ref{Fig_triad} to write out the iso-neutral fluxes at $u$- and |
---|
| 276 | $w$-points as sums of the triad fluxes that cross the $u$- and $w$-faces: |
---|
| 277 | %(Fig.~\ref{Fig_ISO_triad}): |
---|
| 278 | \begin{flalign} \label{Eq_iso_flux} \textbf{F}_{iso}(T) &\equiv |
---|
| 279 | \sum_{\substack{i_p,\,k_p}} |
---|
| 280 | \begin{pmatrix} |
---|
| 281 | {_{i+1/2-i_p}^k {\mathbb{F}_u}_{i_p}^{k_p} } (T) \\ |
---|
| 282 | \\ |
---|
| 283 | {_i^{k+1/2-k_p} {\mathbb{F}_w}_{i_p}^{k_p} } (T) \\ |
---|
| 284 | \end{pmatrix}. |
---|
| 285 | \end{flalign} |
---|
| 286 | \subsection{Ensuring the scheme cannot increase tracer variance} |
---|
| 287 | \label{sec:variance} |
---|
| 288 | |
---|
| 289 | We now require that this operator cannot increase the |
---|
| 290 | globally-integrated tracer variance. |
---|
| 291 | %This changes according to |
---|
| 292 | % \begin{align*} |
---|
| 293 | % &\int_D D_l^T \; T \;dv \equiv \sum_{i,k} \left\{ T \ D_l^T \ b_T \right\} \\ |
---|
| 294 | % &\equiv + \sum_{i,k} \sum_{\substack{i_p,\,k_p}} \left\{ |
---|
| 295 | % \delta_{i} \left[{_{i+1/2-i_p}^k {\mathbb{F}_u }_{i_p}^{k_p}} \right] |
---|
| 296 | % + \delta_{k} \left[ {_i^{k+1/2-k_p} {\mathbb{F}_w}_{i_p}^{k_p}} \right] \ T \right\} \\ |
---|
| 297 | % &\equiv - \sum_{i,k} \sum_{\substack{i_p,\,k_p}} \left\{ |
---|
| 298 | % {_{i+1/2-i_p}^k {\mathbb{F}_u }_{i_p}^{k_p}} \ \delta_{i+1/2} [T] |
---|
| 299 | % + {_i^{k+1/2-k_p} {\mathbb{F}_w}_{i_p}^{k_p}} \ \delta_{k+1/2} [T] \right\} \\ |
---|
| 300 | % \end{align*} |
---|
| 301 | Each triad slope $_i^k\mathbb{R}_{i_p}^{k_p}$ drives a lateral flux |
---|
| 302 | $_i^k{\mathbb{F}_u}_{i_p}^{k_p} (T)$ across the u-point $i+i_p,k$ and |
---|
| 303 | a vertical flux $_i^k{\mathbb{F}_w}_{i_p}^{k_p} (T)$ across the |
---|
| 304 | w-point $i,k+k_p$. The lateral flux drives a net rate of change of |
---|
| 305 | variance at points $i+i_p-\half,k$ and $i+i_p+\half,k$ of |
---|
| 306 | \begin{multline} |
---|
| 307 | {b_T}_{i+i_p-1/2}^k\left(\frac{\partial T}{\partial t}T\right)_{i+i_p-1/2}^k+ |
---|
| 308 | \quad {b_T}_{i+i_p+1/2}^k\left(\frac{\partial T}{\partial |
---|
| 309 | t}T\right)_{i+i_p+1/2}^k \\ |
---|
| 310 | \begin{split} |
---|
| 311 | &= -T_{i+i_p-1/2}^k{\;} _i^k{\mathbb{F}_u}_{i_p}^{k_p} (T) \quad + \quad T_{i+i_p+1/2}^k |
---|
| 312 | {\;}_i^k{\mathbb{F}_u}_{i_p}^{k_p} (T) \\ |
---|
| 313 | &={\;} _i^k{\mathbb{F}_u}_{i_p}^{k_p} (T)\,\delta_{i+ i_p}[T^k], \label{locFdtdx} |
---|
| 314 | \end{split} |
---|
| 315 | \end{multline} |
---|
| 316 | while the vertical flux similarly drives a net rate of change of variance at points $i,k+k_p-\half$ and |
---|
| 317 | $i,k+k_p+\half$ of |
---|
| 318 | \begin{equation} |
---|
| 319 | \label{locFdtdz} |
---|
| 320 | _i^k{\mathbb{F}_w}_{i_p}^{k_p} (T) \,\delta_{k+ k_p}[T^i]. |
---|
| 321 | \end{equation} |
---|
| 322 | The total variance tendency driven by the triad is the sum of these |
---|
| 323 | two. Expanding $_i^k{\mathbb{F}_u}_{i_p}^{k_p} (T)$ and |
---|
| 324 | $_i^k{\mathbb{F}_w}_{i_p}^{k_p} (T)$ with \eqref{eq:latflux-triad} and |
---|
| 325 | \eqref{eq:vertflux-triad}, it is |
---|
| 326 | \begin{multline*} |
---|
| 327 | -A_i^k\left \{ |
---|
| 328 | { } _i^k{\mathbb{A}_u}_{i_p}^{k_p} |
---|
| 329 | \left( |
---|
| 330 | \frac{ \delta_{i+ i_p}[T^k] }{ {e_{1u}}_{\,i + i_p}^{\,k} } |
---|
| 331 | - {_i^k\mathbb{R}_{i_p}^{k_p}} \ |
---|
| 332 | \frac{ \delta_{k+k_p} [T^i] }{{e_{3w}}_{\,i}^{\,k+k_p} }\right)\,\delta_{i+ i_p}[T^k] \right.\\ |
---|
| 333 | - \left. { } _i^k{\mathbb{A}_w}_{i_p}^{k_p} |
---|
| 334 | \left( |
---|
| 335 | \frac{ \delta_{i+ i_p}[T^k] }{ {e_{1u}}_{\,i + i_p}^{\,k} } |
---|
| 336 | -{\:}_i^k\mathbb{R}_{i_p}^{k_p} |
---|
| 337 | \frac{ \delta_{k+k_p} [T^i] }{{e_{3w}}_{\,i}^{\,k+k_p} } |
---|
| 338 | \right) {\,}_i^k\mathbb{R}_{i_p}^{k_p}\delta_{k+ k_p}[T^i] |
---|
| 339 | \right \}. |
---|
| 340 | \end{multline*} |
---|
| 341 | The key point is then that if we require |
---|
| 342 | $_i^k{\mathbb{A}_u}_{i_p}^{k_p}$ and $_i^k{\mathbb{A}_w}_{i_p}^{k_p}$ |
---|
| 343 | to be related to a triad volume $_i^k\mathbb{V}_{i_p}^{k_p}$ by |
---|
| 344 | \begin{equation} |
---|
| 345 | \label{eq:V-A} |
---|
| 346 | _i^k\mathbb{V}_{i_p}^{k_p} |
---|
| 347 | ={\;}_i^k{\mathbb{A}_u}_{i_p}^{k_p}\,{e_{1u}}_{\,i + i_p}^{\,k} |
---|
| 348 | ={\;}_i^k{\mathbb{A}_w}_{i_p}^{k_p}\,{e_{3w}}_{\,i}^{\,k + k_p}, |
---|
| 349 | \end{equation} |
---|
| 350 | the variance tendency reduces to the perfect square |
---|
| 351 | \begin{equation} |
---|
| 352 | \label{eq:perfect-square} |
---|
| 353 | -A_i^k{\:} _i^k\mathbb{V}_{i_p}^{k_p} |
---|
| 354 | \left( |
---|
| 355 | \frac{ \delta_{i+ i_p}[T^k] }{ {e_{1u}}_{\,i + i_p}^{\,k} } |
---|
| 356 | -{\:}_i^k\mathbb{R}_{i_p}^{k_p} |
---|
| 357 | \frac{ \delta_{k+k_p} [T^i] }{{e_{3w}}_{\,i}^{\,k+k_p} } |
---|
| 358 | \right)^2\leq 0. |
---|
| 359 | \end{equation} |
---|
| 360 | Thus, the constraint \eqref{eq:V-A} ensures that the fluxes associated |
---|
| 361 | with a given slope triad $_i^k\mathbb{R}_{i_p}^{k_p}$ do not increase |
---|
| 362 | the net variance. Since the total fluxes are sums of such fluxes from |
---|
| 363 | the various triads, this constraint, applied to all triads, is |
---|
| 364 | sufficient to ensure that the globally integrated variance does not |
---|
| 365 | increase. |
---|
| 366 | |
---|
| 367 | The expression \eqref{eq:V-A} can be interpreted as a discretization |
---|
| 368 | of the global integral |
---|
| 369 | \begin{equation} |
---|
| 370 | \label{eq:cts-var} |
---|
| 371 | \frac{\partial}{\partial t}\int\half T^2\, dV = |
---|
| 372 | \int\mathbf{F}\cdot\nabla T\, dV, |
---|
| 373 | \end{equation} |
---|
| 374 | where, within each triad volume $_i^k\mathbb{V}_{i_p}^{k_p}$, the |
---|
| 375 | lateral and vertical fluxes/unit area |
---|
| 376 | \[ |
---|
| 377 | \mathbf{F}=\left( |
---|
| 378 | _i^k{\mathbb{F}_u}_{i_p}^{k_p} (T)/_i^k{\mathbb{A}_u}_{i_p}^{k_p}, |
---|
| 379 | {\:}_i^k{\mathbb{F}_w}_{i_p}^{k_p} (T)/_i^k{\mathbb{A}_w}_{i_p}^{k_p} |
---|
| 380 | \right) |
---|
| 381 | \] |
---|
| 382 | and the gradient |
---|
| 383 | \[\nabla T = \left( |
---|
| 384 | \delta_{i+ i_p}[T^k] / {e_{1u}}_{\,i + i_p}^{\,k}, |
---|
| 385 | \delta_{k+ k_p}[T^i] / {e_{3w}}_{\,i}^{\,k + k_p} |
---|
| 386 | \right) |
---|
| 387 | \] |
---|
| 388 | \subsection{Triad volumes in Griffes's scheme and in \NEMO} |
---|
| 389 | To complete the discretization we now need only specify the triad |
---|
| 390 | volumes $_i^k\mathbb{V}_{i_p}^{k_p}$. \citet{Griffies_al_JPO98} identify |
---|
| 391 | these $_i^k\mathbb{V}_{i_p}^{k_p}$ as the volumes of the quarter |
---|
| 392 | cells, defined in terms of the distances between T, u,f and |
---|
| 393 | w-points. This is the natural discretization of |
---|
| 394 | \eqref{eq:cts-var}. The \NEMO model, however, operates with scale |
---|
| 395 | factors instead of grid sizes, and scale factors for the quarter |
---|
| 396 | cells are not defined. Instead, therefore we simply choose |
---|
| 397 | \begin{equation} |
---|
| 398 | \label{eq:V-NEMO} |
---|
| 399 | _i^k\mathbb{V}_{i_p}^{k_p}=\quarter {b_u}_{i+i_p}^k, |
---|
| 400 | \end{equation} |
---|
| 401 | as a quarter of the volume of the u-cell inside which the triad |
---|
| 402 | quarter-cell lies. This has the nice property that when the slopes |
---|
| 403 | $\mathbb{R}$ vanish, the lateral flux from tracer cell $i,k$ to |
---|
| 404 | $i+1,k$ reduces to the classical form |
---|
| 405 | \begin{equation} |
---|
| 406 | \label{eq:lat-normal} |
---|
| 407 | -\overline{A}_{\,i+1/2}^k\; |
---|
| 408 | \frac{{b_u}_{i+1/2}^k}{{e_{1u}}_{\,i + i_p}^{\,k}} |
---|
| 409 | \;\frac{\delta_{i+ 1/2}[T^k] }{{e_{1u}}_{\,i + i_p}^{\,k}} |
---|
| 410 | = -\overline{A}_{\,i+1/2}^k\;\frac{{e_{1w}}_{\,i + 1/2}^{\,k}\:{e_{1v}}_{\,i + 1/2}^{\,k}}{{e_{1u}}_{\,i + 1/2}^{\,k}}. |
---|
| 411 | \end{equation} |
---|
| 412 | In fact if the diffusive coefficient is defined at u-points, so that |
---|
| 413 | we employ $A_{i+i_p}^k$ instead of $A_i^k$ in the definitions of the |
---|
| 414 | triad fluxes \eqref{eq:latflux-triad} and \eqref{eq:vertflux-triad}, |
---|
| 415 | we can replace $\overline{A}_{\,i+1/2}^k$ by $A_{i+1/2}^k$ in the above. |
---|
| 416 | |
---|
| 417 | \subsection{Summary of the scheme} |
---|
| 418 | The divergence of the expression \eqref{Eq_iso_flux} for the fluxes gives the iso-neutral diffusion tendency at |
---|
| 419 | each tracer point: |
---|
| 420 | \begin{equation} \label{Gf_operator} D_l^T = \frac{1}{b_T} |
---|
| 421 | \sum_{\substack{i_p,\,k_p}} \left\{ \delta_{i} \left[{_{i+1/2-i_p}^k |
---|
| 422 | {\mathbb{F}_u }_{i_p}^{k_p}} \right] + \delta_{k} \left[ |
---|
| 423 | {_i^{k+1/2-k_p} {\mathbb{F}_w}_{i_p}^{k_p}} \right] \right\} |
---|
| 424 | \end{equation} |
---|
| 425 | where $b_T= e_{1T}\,e_{2T}\,e_{3T}$ is the volume of $T$-cells. |
---|
| 426 | The diffusion scheme satisfies the following six properties: |
---|
| 427 | \begin{description} |
---|
| 428 | \item[$\bullet$ horizontal diffusion] The discretization of the |
---|
| 429 | diffusion operator recovers \eqref{eq:lat-normal} the traditional five-point Laplacian in |
---|
| 430 | the limit of flat iso-neutral direction : |
---|
| 431 | \begin{equation} \label{Gf_property1a} D_l^T = \frac{1}{b_T} \ |
---|
| 432 | \delta_{i} \left[ \frac{e_{2u}\,e_{3u}}{e_{1u}} \; |
---|
| 433 | \overline{A}^{\,i} \; \delta_{i+1/2}[T] \right] \qquad |
---|
| 434 | \text{when} \quad { _i^k \mathbb{R}_{i_p}^{k_p} }=0 |
---|
| 435 | \end{equation} |
---|
| 436 | |
---|
| 437 | \item[$\bullet$ implicit treatment in the vertical] Only tracer values |
---|
| 438 | associated with a single water column appear in the expression |
---|
| 439 | \eqref{eq:i33} for the $_{33}$ fluxes, vertical fluxes driven by |
---|
| 440 | vertical gradients. This is of paramount importance since it means |
---|
| 441 | that an implicit in time algorithm can be used to solve the vertical |
---|
| 442 | diffusion equation. This is a necessity since the vertical eddy |
---|
| 443 | diffusivity associated with this term, |
---|
| 444 | \begin{equation} |
---|
| 445 | \frac{1}{b_w}\sum_{\substack{i_p, \,k_p}} \left\{ |
---|
| 446 | {\:}_i^k\mathbb{V}_{i_p}^{k_p} \:A_i^k \: \left(_i^k \mathbb{R}_{i_p}^{k_p}\right)^2 |
---|
| 447 | \right\} = |
---|
| 448 | \frac{1}{4b_w}\sum_{\substack{i_p, \,k_p}} \left\{ |
---|
| 449 | {b_u}_{i+i_p}^k\:A_i^k \: \left(_i^k \mathbb{R}_{i_p}^{k_p}\right)^2 |
---|
| 450 | \right\}, |
---|
| 451 | \end{equation} |
---|
| 452 | (where $b_w= e_{1w}\,e_{2w}\,e_{3w}$ is the volume of $w$-cells) can be quite large. |
---|
| 453 | |
---|
| 454 | \item[$\bullet$ pure iso-neutral operator] The iso-neutral flux of |
---|
| 455 | locally referenced potential density is zero. See |
---|
| 456 | \eqref{eq:latflux-rho} and \eqref{eq:vertflux-triad2}. |
---|
| 457 | |
---|
| 458 | \item[$\bullet$ conservation of tracer] The iso-neutral diffusion |
---|
| 459 | conserves tracer content, $i.e.$ |
---|
| 460 | \begin{equation} \label{Gf_property1} \sum_{i,j,k} \left\{ D_l^T \ |
---|
| 461 | b_T \right\} = 0 |
---|
| 462 | \end{equation} |
---|
| 463 | This property is trivially satisfied since the iso-neutral diffusive |
---|
| 464 | operator is written in flux form. |
---|
| 465 | |
---|
| 466 | \item[$\bullet$ no increase of tracer variance] The iso-neutral diffusion |
---|
| 467 | does not increase the tracer variance, $i.e.$ |
---|
| 468 | \begin{equation} \label{Gf_property1} \sum_{i,j,k} \left\{ T \ D_l^T |
---|
| 469 | \ b_T \right\} \leq 0 |
---|
| 470 | \end{equation} |
---|
| 471 | The property is demonstrated in |
---|
| 472 | \ref{sec:variance} above. It is a key property for a diffusion |
---|
| 473 | term. It means that it is also a dissipation term, $i.e.$ it is a |
---|
| 474 | diffusion of the square of the quantity on which it is applied. It |
---|
| 475 | therefore ensures that, when the diffusivity coefficient is large |
---|
| 476 | enough, the field on which it is applied become free of grid-point |
---|
| 477 | noise. |
---|
| 478 | |
---|
| 479 | \item[$\bullet$ self-adjoint operator] The iso-neutral diffusion |
---|
| 480 | operator is self-adjoint, $i.e.$ |
---|
| 481 | \begin{equation} \label{Gf_property1} \sum_{i,j,k} \left\{ S \ D_l^T |
---|
| 482 | \ b_T \right\} = \sum_{i,j,k} \left\{ D_l^S \ T \ b_T \right\} |
---|
| 483 | \end{equation} |
---|
| 484 | In other word, there is no need to develop a specific routine from |
---|
| 485 | the adjoint of this operator. We just have to apply the same |
---|
| 486 | routine. This property can be demonstrated similarly to the proof of |
---|
| 487 | the `no increase of tracer variance' property. The contribution by a |
---|
| 488 | single triad towards the left hand side of \eqref{Gf_property1}, can |
---|
| 489 | be found by replacing $\delta[T]$ by $\delta[S]$ in \eqref{locFdtdx} |
---|
| 490 | and \eqref{locFdtdx}. This results in a term similar to |
---|
| 491 | \eqref{eq:perfect-square}, |
---|
| 492 | \begin{equation} |
---|
| 493 | \label{eq:TScovar} |
---|
| 494 | -A_i^k{\:} _i^k\mathbb{V}_{i_p}^{k_p} |
---|
| 495 | \left( |
---|
| 496 | \frac{ \delta_{i+ i_p}[T^k] }{ {e_{1u}}_{\,i + i_p}^{\,k} } |
---|
| 497 | -{\:}_i^k\mathbb{R}_{i_p}^{k_p} |
---|
| 498 | \frac{ \delta_{k+k_p} [T^i] }{{e_{3w}}_{\,i}^{\,k+k_p} } |
---|
| 499 | \right) |
---|
| 500 | \left( |
---|
| 501 | \frac{ \delta_{i+ i_p}[S^k] }{ {e_{1u}}_{\,i + i_p}^{\,k} } |
---|
| 502 | -{\:}_i^k\mathbb{R}_{i_p}^{k_p} |
---|
| 503 | \frac{ \delta_{k+k_p} [S^i] }{{e_{3w}}_{\,i}^{\,k+k_p} } |
---|
| 504 | \right). |
---|
| 505 | \end{equation} |
---|
| 506 | This is symmetrical in $T $ and $S$, so exactly the same term arises |
---|
| 507 | from the discretization of this triad's contribution towards the |
---|
| 508 | RHS of \eqref{Gf_property1}. |
---|
| 509 | \end{description} |
---|
| 510 | |
---|
| 511 | |
---|
| 512 | $\ $\newline %force an empty line |
---|
| 513 | % ================================================================ |
---|
| 514 | % Skew flux formulation for Eddy Induced Velocity : |
---|
| 515 | % ================================================================ |
---|
| 516 | \section{ Eddy induced velocity and Skew flux formulation} |
---|
| 517 | |
---|
| 518 | When Gent and McWilliams's [1990] diffusion is used (\key{traldf\_eiv} defined), |
---|
| 519 | an additional advection term is added. The associated velocity is the so called |
---|
| 520 | eddy induced velocity, the formulation of which depends on the slopes of iso- |
---|
| 521 | neutral surfaces. Contrary to the case of iso-neutral mixing, the slopes used |
---|
| 522 | here are referenced to the geopotential surfaces, $i.e.$ \eqref{Eq_ldfslp_geo} |
---|
| 523 | is used in $z$-coordinate, and the sum \eqref{Eq_ldfslp_geo} |
---|
| 524 | + \eqref{Eq_ldfslp_iso} in $z^*$ or $s$-coordinates. |
---|
| 525 | |
---|
| 526 | The eddy induced velocity is given by: |
---|
| 527 | \begin{equation} \label{Eq_eiv_v} |
---|
| 528 | \begin{split} |
---|
| 529 | u^* & = - \frac{1}{e_{2}e_{3}}\; \partial_k \left( e_{2} \, A_{e} \; r_i \right) \\ |
---|
| 530 | v^* & = - \frac{1}{e_{1}e_{3}}\; \partial_k \left( e_{1} \, A_{e} \; r_j \right) \\ |
---|
| 531 | w^* & = \frac{1}{e_{1}e_{2}}\; \left\{ \partial_i \left( e_{2} \, A_{e} \; r_i \right) |
---|
| 532 | + \partial_j \left( e_{1} \, A_{e} \;r_j \right) \right\} \\ |
---|
| 533 | \end{split} |
---|
| 534 | \end{equation} |
---|
| 535 | where $A_{e}$ is the eddy induced velocity coefficient, and $r_i$ and $r_j$ the slopes between the iso-neutral and the geopotential surfaces. |
---|
| 536 | |
---|
| 537 | The traditional way to implement this additional advection is to add it to the Eulerian |
---|
| 538 | velocity prior to computing the tracer advection. This allows us to take advantage of |
---|
| 539 | all the advection schemes offered for the tracers (see \S\ref{TRA_adv}) and not just |
---|
| 540 | a $2^{nd}$ order advection scheme. This is particularly useful for passive tracers |
---|
| 541 | where \emph{positivity} of the advection scheme is of paramount importance. |
---|
| 542 | |
---|
| 543 | \citet{Griffies_JPO98} introduces another way to implement the eddy induced advection, |
---|
| 544 | the so-called skew form. It is based on a transformation of the advective fluxes |
---|
| 545 | using the non-divergent nature of the eddy induced velocity. |
---|
| 546 | For example in the (\textbf{i},\textbf{k}) plane, the tracer advective fluxes can be |
---|
| 547 | transformed as follows: |
---|
| 548 | \begin{flalign*} |
---|
| 549 | \begin{split} |
---|
| 550 | \textbf{F}_{eiv}^T = |
---|
| 551 | \begin{pmatrix} |
---|
| 552 | {e_{2}\,e_{3}\; u^*} \\ |
---|
| 553 | {e_{1}\,e_{2}\; w^*} \\ |
---|
| 554 | \end{pmatrix} \; T |
---|
| 555 | &= |
---|
| 556 | \begin{pmatrix} |
---|
| 557 | { - \partial_k \left( e_{2} \, A_{e} \; r_i \right) \; T \;} \\ |
---|
| 558 | {+ \partial_i \left( e_{2} \, A_{e} \; r_i \right) \; T \;} \\ |
---|
| 559 | \end{pmatrix} \\ |
---|
| 560 | &= |
---|
| 561 | \begin{pmatrix} |
---|
| 562 | { - \partial_k \left( e_{2} \, A_{e} \; r_i \; T \right) \;} \\ |
---|
| 563 | {+ \partial_i \left( e_{2} \, A_{e} \; r_i \; T \right) \;} \\ |
---|
| 564 | \end{pmatrix} |
---|
| 565 | + |
---|
| 566 | \begin{pmatrix} |
---|
| 567 | {+ e_{2} \, A_{e} \; r_i \; \partial_k T} \\ |
---|
| 568 | { - e_{2} \, A_{e} \; r_i \; \partial_i T} \\ |
---|
| 569 | \end{pmatrix} |
---|
| 570 | \end{split} |
---|
| 571 | \end{flalign*} |
---|
| 572 | and since the eddy induces velocity field is no-divergent, we end up with the skew |
---|
| 573 | form of the eddy induced advective fluxes: |
---|
| 574 | \begin{equation} \label{Eq_eiv_skew_continuous} |
---|
| 575 | \textbf{F}_{eiv}^T = \begin{pmatrix} |
---|
| 576 | {+ e_{2} \, A_{e} \; r_i \; \partial_k T} \\ |
---|
| 577 | { - e_{2} \, A_{e} \; r_i \; \partial_i T} \\ |
---|
| 578 | \end{pmatrix} |
---|
| 579 | \end{equation} |
---|
| 580 | The tendency associated with eddy induced velocity is then simply the divergence |
---|
| 581 | of the \eqref{Eq_eiv_skew_continuous} fluxes. It naturally conserves the tracer |
---|
| 582 | content, as it is expressed in flux form. It also preserve the tracer variance. |
---|
| 583 | |
---|
| 584 | The discrete form of \eqref{Eq_eiv_skew_continuous} using the slopes \eqref{Gf_slopes} and defining $A_e$ at $T$-point is then given by: |
---|
| 585 | \begin{flalign} \label{Eq_eiv_skew} \begin{split} |
---|
| 586 | \textbf{F}_{eiv}^T \equiv |
---|
| 587 | \sum_{\substack{i_p,\,k_p}} \begin{pmatrix} |
---|
| 588 | +{e_{2u}}_{i+1/2-i_p}^{k} \ {A_{e}}_{i+1/2-i_p}^{k} |
---|
| 589 | \ \ { _{i+1/2-i_p}^k \mathbb{R}_{i_p}^{k_p} } \ \ \delta_{k+k_p}[T_{i+1/2-i_p}] \\ |
---|
| 590 | \\ |
---|
| 591 | - {e_{2u}}_i^{k+1/2-k_p} \ {A_{e}}_i^{k+1/2-k_p} |
---|
| 592 | \ \ { _i^{k+1/2-k_p} \mathbb{R}_{i_p}^{k_p} } \ \delta_{i+i_p}[T^{k+1/2-k_p}] \\ |
---|
| 593 | \end{pmatrix} |
---|
| 594 | \end{split} \end{flalign} |
---|
| 595 | Note that \eqref{Eq_eiv_skew} is valid in $z$-coordinate with or without partial cells. In $z^*$ or $s$-coordinate, the the slope between the level and the geopotential surfaces must be added to $\mathbb{R}$ for the discret form to be exact. |
---|
| 596 | |
---|
| 597 | Such a choice of discretisation is consistent with the iso-neutral operator as it uses the same definition for the slopes. It also ensures the conservation of the tracer variance (see Appendix \ref{Apdx_eiv_skew}), $i.e.$ it does not include a diffusive component but is a `pure' advection term. |
---|
| 598 | |
---|
| 599 | |
---|
| 600 | |
---|
| 601 | |
---|
| 602 | \newpage %force an empty line |
---|
| 603 | % ================================================================ |
---|
| 604 | % Discrete Invariants of the skew flux formulation |
---|
| 605 | % ================================================================ |
---|
| 606 | \subsection{Discrete Invariants of the skew flux formulation} |
---|
| 607 | \label{Apdx_eiv_skew} |
---|
| 608 | |
---|
| 609 | |
---|
| 610 | Demonstration for the conservation of the tracer variance in the (\textbf{i},\textbf{j}) plane. |
---|
| 611 | |
---|
| 612 | This have to be moved in an Appendix. |
---|
| 613 | |
---|
| 614 | The continuous property to be demonstrated is : |
---|
| 615 | \begin{align*} |
---|
| 616 | \int_D \nabla \cdot \textbf{F}_{eiv}(T) \; T \;dv \equiv 0 |
---|
| 617 | \end{align*} |
---|
| 618 | The discrete form of its left hand side is obtained using \eqref{Eq_eiv_skew} |
---|
| 619 | \begin{align*} |
---|
| 620 | \sum\limits_{i,k} \sum_{\substack{i_p,\,k_p}} \Biggl\{ \;\; |
---|
| 621 | \delta_i &\left[ |
---|
| 622 | {e_{2u}}_{i+i_p+1/2}^{k} \;\ \ {A_{e}}_{i+i_p+1/2}^{k} |
---|
| 623 | \ \ \ { _{i+i_p+1/2}^k \mathbb{R}_{-i_p}^{k_p} } \quad \delta_{k+k_p}[T_{i+i_p+1/2}] |
---|
| 624 | \right] \; T_i^k \\ |
---|
| 625 | - \delta_k &\left[ |
---|
| 626 | {e_{2u}}_i^{k+k_p+1/2} \ \ {A_{e}}_i^{k+k_p+1/2} |
---|
| 627 | \ \ { _i^{k+k_p+1/2} \mathbb{R}_{i_p}^{-k_p} } \ \ \delta_{i+i_p}[T^{k+k_p+1/2}] |
---|
| 628 | \right] \; T_i^k \ \Biggr\} |
---|
| 629 | \end{align*} |
---|
| 630 | apply the adjoint of delta operator, it becomes |
---|
| 631 | \begin{align*} |
---|
| 632 | \sum\limits_{i,k} \sum_{\substack{i_p,\,k_p}} \Biggl\{ \;\; |
---|
| 633 | &\left( |
---|
| 634 | {e_{2u}}_{i+i_p+1/2}^{k} \;\ \ {A_{e}}_{i+i_p+1/2}^{k} |
---|
| 635 | \ \ \ { _{i+i_p+1/2}^k \mathbb{R}_{-i_p}^{k_p} } \quad \delta_{k+k_p}[T_{i+i_p+1/2}] |
---|
| 636 | \right) \; \delta_{i+1/2}[T^{k}] \\ |
---|
| 637 | - &\left( |
---|
| 638 | {e_{2u}}_i^{k+k_p+1/2} \ \ {A_{e}}_i^{k+k_p+1/2} |
---|
| 639 | \ \ { _i^{k+k_p+1/2} \mathbb{R}_{i_p}^{-k_p} } \ \ \delta_{i+i_p}[T^{k+k_p+1/2}] |
---|
| 640 | \right) \; \delta_{k+1/2}[T_{i}] \ \Biggr\} |
---|
| 641 | \end{align*} |
---|
| 642 | Expending the summation on $i_p$ and $k_p$, it becomes: |
---|
| 643 | \begin{align*} |
---|
| 644 | \begin{matrix} |
---|
| 645 | &\sum\limits_{i,k} \Bigl\{ |
---|
| 646 | &+{e_{2u}}_{i+1}^{k} &{A_{e}}_{i+1 }^{k} |
---|
| 647 | &\ {_{i+1}^k \mathbb{R}_{- 1/2}^{-1/2}} &\delta_{k-1/2}[T_{i+1}] &\delta_{i+1/2}[T^{k}] &\\ |
---|
| 648 | &&+{e_{2u}}_i^{k\ \ \ \:} &{A_{e}}_{i}^{k\ \ \ \:} |
---|
| 649 | &\ {\ \ \;_i^k \mathbb{R}_{+1/2}^{-1/2}} &\delta_{k-1/2}[T_{i\ \ \ \;}] &\delta_{i+1/2}[T^{k}] &\\ |
---|
| 650 | &&+{e_{2u}}_{i+1}^{k} &{A_{e}}_{i+1 }^{k} |
---|
| 651 | &\ {_{i+1}^k \mathbb{R}_{- 1/2}^{+1/2}} &\delta_{k+1/2}[T_{i+1}] &\delta_{i+1/2}[T^{k}] &\\ |
---|
| 652 | &&+{e_{2u}}_i^{k\ \ \ \:} &{A_{e}}_{i}^{k\ \ \ \:} |
---|
| 653 | &\ {\ \ \;_i^k \mathbb{R}_{+1/2}^{+1/2}} &\delta_{k+1/2}[T_{i\ \ \ \;}] &\delta_{i+1/2}[T^{k}] &\\ |
---|
| 654 | % |
---|
| 655 | &&-{e_{2u}}_i^{k+1} &{A_{e}}_i^{k+1} |
---|
| 656 | &{_i^{k+1} \mathbb{R}_{-1/2}^{- 1/2}} &\delta_{i-1/2}[T^{k+1}] &\delta_{k+1/2}[T_{i}] &\\ |
---|
| 657 | &&-{e_{2u}}_i^{k\ \ \ \:} &{A_{e}}_i^{k\ \ \ \:} |
---|
| 658 | &{\ \ \;_i^k \mathbb{R}_{-1/2}^{+1/2}} &\delta_{i-1/2}[T^{k\ \ \ \:}] &\delta_{k+1/2}[T_{i}] &\\ |
---|
| 659 | &&-{e_{2u}}_i^{k+1 } &{A_{e}}_i^{k+1} |
---|
| 660 | &{_i^{k+1} \mathbb{R}_{+1/2}^{- 1/2}} &\delta_{i+1/2}[T^{k+1}] &\delta_{k+1/2}[T_{i}] &\\ |
---|
| 661 | &&-{e_{2u}}_i^{k\ \ \ \:} &{A_{e}}_i^{k\ \ \ \:} |
---|
| 662 | &{\ \ \;_i^k \mathbb{R}_{+1/2}^{+1/2}} &\delta_{i+1/2}[T^{k\ \ \ \:}] &\delta_{k+1/2}[T_{i}] |
---|
| 663 | &\Bigr\} \\ |
---|
| 664 | \end{matrix} |
---|
| 665 | \end{align*} |
---|
| 666 | The two terms associated with the triad ${_i^k \mathbb{R}_{+1/2}^{+1/2}}$ are the |
---|
| 667 | same but of opposite signs, they cancel out. |
---|
| 668 | Exactly the same thing occurs for the triad ${_i^k \mathbb{R}_{-1/2}^{-1/2}}$. |
---|
| 669 | The two terms associated with the triad ${_i^k \mathbb{R}_{+1/2}^{-1/2}}$ are the |
---|
| 670 | same but both of opposite signs and shifted by 1 in $k$ direction. When summing over $k$ |
---|
| 671 | they cancel out with the neighbouring grid points. |
---|
| 672 | Exactly the same thing occurs for the triad ${_i^k \mathbb{R}_{-1/2}^{+1/2}}$ in the |
---|
| 673 | $i$ direction. Therefore the sum over the domain is zero, $i.e.$ the variance of the |
---|
| 674 | tracer is preserved by the discretisation of the skew fluxes. |
---|
| 675 | |
---|
| 676 | %%% Local Variables: |
---|
| 677 | %%% TeX-master: "../../NEMO_book.tex" |
---|
| 678 | %%% End: |
---|