\documentclass[../tex_main/NEMO_manual]{subfiles}
\begin{document}
% ================================================================
% Chapter Ñ Appendix B : Diffusive Operators
% ================================================================
\chapter{Appendix B : Diffusive Operators}
\label{Apdx_B}
\minitoc


\newpage
$\ $\newline    % force a new ligne

% ================================================================
% Horizontal/Vertical 2nd Order Tracer Diffusive Operators
% ================================================================
\section{Horizontal/Vertical $2^{nd}$ order tracer diffusive operators}
\label{Apdx_B_1}

\subsubsection*{In z-coordinates}
In $z$-coordinates, the horizontal/vertical second order tracer diffusion operator
is given by:
\begin{eqnarray} \label{Apdx_B1}
 &D^T = \frac{1}{e_1 \, e_2}      \left[
  \left. \frac{\partial}{\partial i} \left( 	\frac{e_2}{e_1}A^{lT} \;\left. \frac{\partial T}{\partial i} \right|_z   \right)   \right|_z      \right.
  						     \left.
+ \left. \frac{\partial}{\partial j} \left(  \frac{e_1}{e_2}A^{lT} \;\left. \frac{\partial T}{\partial j} \right|_z   \right)   \right|_z      \right]
+ \frac{\partial }{\partial z}\left( {A^{vT} \;\frac{\partial T}{\partial z}} \right)
\end{eqnarray}

\subsubsection*{In generalized vertical coordinates}
In $s$-coordinates, we defined the slopes of $s$-surfaces, $\sigma_1$ and
$\sigma_2$ by \eqref{Apdx_A_s_slope} and the vertical/horizontal ratio of diffusion
coefficient by $\epsilon = A^{vT} / A^{lT}$. The diffusion operator is given by:

\begin{equation} \label{Apdx_B2}
D^T = \left. \nabla \right|_s \cdot
           \left[ A^{lT} \;\Re \cdot \left. \nabla \right|_s T  \right] \\
\;\;\text{where} \;\Re =\left( {{\begin{array}{*{20}c}
 1 \hfill & 0 \hfill & {-\sigma _1 } \hfill \\
 0 \hfill & 1 \hfill & {-\sigma _2 } \hfill \\
 {-\sigma _1 } \hfill & {-\sigma _2 } \hfill & {\varepsilon +\sigma _1
^2+\sigma _2 ^2} \hfill \\
\end{array} }} \right)
\end{equation}
or in expanded form:
\begin{subequations}
\begin{align*} {\begin{array}{*{20}l}
D^T=& \frac{1}{e_1\,e_2\,e_3 }\;\left[ {\ \ \ \ e_2\,e_3\,A^{lT} \;\left.
{\frac{\partial }{\partial i}\left( {\frac{1}{e_1}\;\left. {\frac{\partial T}{\partial i}} \right|_s -\frac{\sigma _1 }{e_3 }\;\frac{\partial T}{\partial s}} \right)} \right|_s } \right.  \\
&\qquad \quad \ \ \ +e_1\,e_3\,A^{lT} \;\left. {\frac{\partial }{\partial j}\left( {\frac{1}{e_2 }\;\left. {\frac{\partial T}{\partial j}} \right|_s -\frac{\sigma _2 }{e_3 }\;\frac{\partial T}{\partial s}} \right)} \right|_s \\
&\qquad \quad \ \ \ + e_1\,e_2\,A^{lT} \;\frac{\partial }{\partial s}\left( {-\frac{\sigma _1 }{e_1 }\;\left. {\frac{\partial T}{\partial i}} \right|_s -\frac{\sigma _2 }{e_2 }\;\left. {\frac{\partial T}{\partial j}} \right|_s } \right.
 \left. {\left. {+\left( {\varepsilon +\sigma _1^2+\sigma _2 ^2} \right)\;\frac{1}{e_3 }\;\frac{\partial T}{\partial s}} \right)\;\;} \right]
\end{array} }
\end{align*}
\end{subequations}

Equation \eqref{Apdx_B2} is obtained from \eqref{Apdx_B1} without any
additional assumption. Indeed, for the special case $k=z$ and thus $e_3 =1$,
we introduce an arbitrary vertical coordinate $s = s (i,j,z)$ as in Appendix~\ref{Apdx_A}
and use \eqref{Apdx_A_s_slope} and \eqref{Apdx_A_s_chain_rule}.
Since no cross horizontal derivative $\partial _i \partial _j $ appears in
\eqref{Apdx_B1}, the ($i$,$z$) and ($j$,$z$) planes are independent.
The derivation can then be demonstrated for the ($i$,$z$)~$\to$~($j$,$s$)
transformation without any loss of generality:

\begin{subequations}
\begin{align*} {\begin{array}{*{20}l}
D^T&=\frac{1}{e_1\,e_2} \left. {\frac{\partial }{\partial i}\left( {\frac{e_2}{e_1}A^{lT}\;\left. {\frac{\partial T}{\partial i}} \right|_z } \right)} \right|_z
                     +\frac{\partial }{\partial z}\left( {A^{vT}\;\frac{\partial T}{\partial z}} \right)     \\
 \\
%
&=\frac{1}{e_1\,e_2 }\left[ {\left. {\;\frac{\partial }{\partial i}\left( {\frac{e_2}{e_1}A^{lT}\;\left( {\left. {\frac{\partial T}{\partial i}} \right|_s
                                                    -\frac{e_1\,\sigma _1 }{e_3 }\frac{\partial T}{\partial s}} \right)} \right)} \right|_s } \right. \\
& \qquad \qquad \left. { -\frac{e_1\,\sigma _1 }{e_3 }\frac{\partial }{\partial s}\left( {\frac{e_2 }{e_1 }A^{lT}\;\left. {\left( {\left. {\frac{\partial T}{\partial i}} \right|_s -\frac{e_1 \,\sigma _1 }{e_3 }\frac{\partial T}{\partial s}} \right)} \right|_s } \right)\;} \right]
\shoveright{ +\frac{1}{e_3 }\frac{\partial }{\partial s}\left[ {\frac{A^{vT}}{e_3 }\;\frac{\partial T}{\partial s}} \right]}  \qquad \qquad \qquad \\
 \\
%
&=\frac{1}{e_1 \,e_2 \,e_3 }\left[ {\left. {\;\;\frac{\partial }{\partial i}\left( {\frac{e_2 \,e_3 }{e_1 }A^{lT}\;\left. {\frac{\partial T}{\partial i}} \right|_s } \right)} \right|_s -\left. {\frac{e_2 }{e_1}A^{lT}\;\frac{\partial e_3 }{\partial i}} \right|_s \left. {\frac{\partial T}{\partial i}} \right|_s } \right. \\
&  \qquad \qquad \quad \left. {-e_3 \frac{\partial }{\partial i}\left( {\frac{e_2 \,\sigma _1 }{e_3 }A^{lT}\;\frac{\partial T}{\partial s}} \right)} \right|_s -e_1 \,\sigma _1 \frac{\partial }{\partial s}\left( {\frac{e_2 }{e_1 }A^{lT}\;\left. {\frac{\partial T}{\partial i}} \right|_s } \right) \\
&  \qquad \qquad \quad \shoveright{ -e_1 \,\sigma _1 \frac{\partial }{\partial s}\left( {-\frac{e_2 \,\sigma _1 }{e_3 }A^{lT}\;\frac{\partial T}{\partial s}} \right)\;\,\left. {+\frac{\partial }{\partial s}\left( {\frac{e_1 \,e_2 }{e_3 }A^{vT}\;\frac{\partial T}{\partial s}} \right)\quad} \right] }\\
\end{array} } 		\\
%
 {\begin{array}{*{20}l}
\intertext{Noting that $\frac{1}{e_1} \left. \frac{\partial e_3 }{\partial i} \right|_s = \frac{\partial \sigma _1 }{\partial s}$, it becomes:}
%
& =\frac{1}{e_1\,e_2\,e_3 }\left[ {\left. {\;\;\;\frac{\partial }{\partial i}\left( {\frac{e_2\,e_3 }{e_1}\,A^{lT}\;\left. {\frac{\partial T}{\partial i}} \right|_s } \right)} \right|_s \left. -\, {e_3 \frac{\partial }{\partial i}\left( {\frac{e_2 \,\sigma _1 }{e_3 }A^{lT}\;\frac{\partial T}{\partial s}} \right)} \right|_s } \right. \\
& \qquad \qquad \quad -e_2 A^{lT}\;\frac{\partial \sigma _1 }{\partial s}\left. {\frac{\partial T}{\partial i}} \right|_s -e_1 \,\sigma_1 \frac{\partial }{\partial s}\left( {\frac{e_2 }{e_1 }A^{lT}\;\left. {\frac{\partial T}{\partial i}} \right|_s } \right) \\
& \qquad \qquad \quad\shoveright{ \left. { +e_1 \,\sigma _1 \frac{\partial }{\partial s}\left( {\frac{e_2 \,\sigma _1 }{e_3 }A^{lT}\;\frac{\partial T}{\partial s}} \right)+\frac{\partial }{\partial s}\left( {\frac{e_1 \,e_2 }{e_3 }A^{vT}\;\frac{\partial T}{\partial z}} \right)\;\;\;} \right] }\\
\\
&=\frac{1}{e_1 \,e_2 \,e_3 } \left[ {\left. {\;\;\;\frac{\partial }{\partial i} \left( {\frac{e_2 \,e_3 }{e_1 }A^{lT}\;\left. {\frac{\partial T}{\partial i}} \right|_s } \right)} \right|_s \left. {-\frac{\partial }{\partial i}\left( {e_2 \,\sigma _1 A^{lT}\;\frac{\partial T}{\partial s}} \right)} \right|_s } \right. \\
& \qquad \qquad \quad \left. {+\frac{e_2 \,\sigma _1 }{e_3}A^{lT}\;\frac{\partial T}{\partial s} \;\frac{\partial e_3 }{\partial i}}  \right|_s -e_2 A^{lT}\;\frac{\partial \sigma _1 }{\partial s}\left. {\frac{\partial T}{\partial i}} \right|_s \\
& \qquad \qquad \quad-e_2 \,\sigma _1 \frac{\partial}{\partial s}\left( {A^{lT}\;\left. {\frac{\partial T}{\partial i}} \right|_s } \right)+\frac{\partial }{\partial s}\left( {\frac{e_1 \,e_2 \,\sigma _1 ^2}{e_3 }A^{lT}\;\frac{\partial T}{\partial s}} \right) \\
& \qquad \qquad \quad\shoveright{ \left. {-\frac{\partial \left( {e_1 \,e_2 \,\sigma _1 } \right)}{\partial s} \left( {\frac{\sigma _1 }{e_3}A^{lT}\;\frac{\partial T}{\partial s}} \right) + \frac{\partial }{\partial s}\left( {\frac{e_1 \,e_2 }{e_3 }A^{vT}\;\frac{\partial T}{\partial s}} \right)\;\;\;} \right]}
\end{array} } \\
{\begin{array}{*{20}l}
%
\intertext{using the same remark as just above, it becomes:}
%
&= \frac{1}{e_1 \,e_2 \,e_3 } \left[ {\left. {\;\;\;\frac{\partial }{\partial i} \left( {\frac{e_2 \,e_3 }{e_1 }A^{lT}\;\left. {\frac{\partial T}{\partial i}} \right|_s -e_2 \,\sigma _1 A^{lT}\;\frac{\partial T}{\partial s}} \right)} \right|_s } \right.\;\;\; \\
& \qquad \qquad \quad+\frac{e_1 \,e_2 \,\sigma _1 }{e_3 }A^{lT}\;\frac{\partial T}{\partial s}\;\frac{\partial \sigma _1 }{\partial s} - \frac {\sigma _1 }{e_3} A^{lT} \;\frac{\partial \left( {e_1 \,e_2 \,\sigma _1 } \right)}{\partial s}\;\frac{\partial T}{\partial s} \\
& \qquad \qquad \quad-e_2 \left( {A^{lT}\;\frac{\partial \sigma _1 }{\partial s}\left. {\frac{\partial T}{\partial i}} \right|_s +\frac{\partial }{\partial s}\left( {\sigma _1 A^{lT}\;\left. {\frac{\partial T}{\partial i}} \right|_s } \right)-\frac{\partial \sigma _1 }{\partial s}\;A^{lT}\;\left. {\frac{\partial T}{\partial i}} \right|_s } \right) \\
& \qquad \qquad \quad\shoveright{\left. {+\frac{\partial }{\partial s}\left( {\frac{e_1 \,e_2 \,\sigma _1 ^2}{e_3 }A^{lT}\;\frac{\partial T}{\partial s}+\frac{e_1 \,e_2}{e_3 }A^{vT}\;\frac{\partial T}{\partial s}} \right)\;\;\;} \right] }
 \end{array} } \\
{\begin{array}{*{20}l}
%
\intertext{Since the horizontal scale factors do not depend on the vertical coordinate,
the last term of the first line and the first term of the last line cancel, while
the second line reduces to a single vertical derivative, so it becomes:}
%
& =\frac{1}{e_1 \,e_2 \,e_3 }\left[ {\left. {\;\;\;\frac{\partial }{\partial i}\left( {\frac{e_2 \,e_3 }{e_1 }A^{lT}\;\left. {\frac{\partial T}{\partial i}} \right|_s -e_2 \,\sigma _1 \,A^{lT}\;\frac{\partial T}{\partial s}} \right)} \right|_s } \right. \\
& \qquad \qquad \quad \shoveright{ \left. {+\frac{\partial }{\partial s}\left( {-e_2 \,\sigma _1 \,A^{lT}\;\left. {\frac{\partial T}{\partial i}} \right|_s +A^{lT}\frac{e_1 \,e_2 }{e_3 }\;\left( {\varepsilon +\sigma _1 ^2} \right)\frac{\partial T}{\partial s}} \right)\;\;\;} \right]}
 \\
%
\intertext{in other words, the horizontal/vertical Laplacian operator in the ($i$,$s$) plane takes the following form:}
\end{array} } \\
%
{\frac{1}{e_1\,e_2\,e_3}}
\left( {{\begin{array}{*{30}c}
{\left. {\frac{\partial \left( {e_2 e_3 \bullet } \right)}{\partial i}} \right|_s } \hfill \\
{\frac{\partial \left( {e_1 e_2 \bullet } \right)}{\partial s}} \hfill \\
\end{array}}}\right)
\cdot \left[ {A^{lT}
\left( {{\begin{array}{*{30}c}
 {1} \hfill & {-\sigma_1 } \hfill \\
 {-\sigma_1} \hfill & {\varepsilon + \sigma_1^2} \hfill \\
\end{array} }} \right)
\cdot
\left( {{\begin{array}{*{30}c}
{\frac{1}{e_1 }\;\left. {\frac{\partial \bullet }{\partial i}} \right|_s } \hfill \\
{\frac{1}{e_3 }\;\frac{\partial \bullet }{\partial s}} \hfill \\
\end{array}}}       \right) \left( T \right)} \right]
\end{align*}
\end{subequations}
\addtocounter{equation}{-2}

% ================================================================
% Isopycnal/Vertical 2nd Order Tracer Diffusive Operators
% ================================================================
\section{Iso/Diapycnal $2^{nd}$ order tracer diffusive operators}
\label{Apdx_B_2}

\subsubsection*{In z-coordinates}

The iso/diapycnal diffusive tensor $\textbf {A}_{\textbf I}$ expressed in the ($i$,$j$,$k$)
curvilinear coordinate system in which the equations of the ocean circulation model are
formulated, takes the following form \citep{Redi_JPO82}:

\begin{equation} \label{Apdx_B3}
\textbf {A}_{\textbf I} = \frac{A^{lT}}{\left( {1+a_1 ^2+a_2 ^2} \right)}
\left[ {{\begin{array}{*{20}c}
 {1+a_1 ^2} \hfill & {-a_1 a_2 } \hfill & {-a_1 } \hfill \\
 {-a_1 a_2 } \hfill & {1+a_2 ^2} \hfill & {-a_2 } \hfill \\
 {-a_1 } \hfill & {-a_2 } \hfill & {\varepsilon +a_1 ^2+a_2 ^2} \hfill \\
\end{array} }} \right]
\end{equation}
where ($a_1$, $a_2$) are the isopycnal slopes in ($\textbf{i}$,
$\textbf{j}$) directions, relative to geopotentials:
\begin{equation*}
a_1 =\frac{e_3 }{e_1 }\left( {\frac{\partial \rho }{\partial i}} \right)\left( {\frac{\partial \rho }{\partial k}} \right)^{-1}
\qquad , \qquad
a_2 =\frac{e_3 }{e_2 }\left( {\frac{\partial \rho }{\partial j}}
\right)\left( {\frac{\partial \rho }{\partial k}} \right)^{-1}
\end{equation*}

In practice, isopycnal slopes are generally less than $10^{-2}$ in the ocean, so
$\textbf {A}_{\textbf I}$ can be simplified appreciably \citep{Cox1987}:
\begin{subequations} \label{Apdx_B4}
\begin{equation} \label{Apdx_B4a}
{\textbf{A}_{\textbf{I}}} \approx A^{lT}\;\Re\;\text{where} \;\Re =
\left[ {{\begin{array}{*{20}c}
 1 \hfill & 0 \hfill & {-a_1 } \hfill \\
 0 \hfill & 1 \hfill & {-a_2 } \hfill \\
 {-a_1 } \hfill & {-a_2 } \hfill & {\varepsilon +a_1 ^2+a_2 ^2} \hfill \\
\end{array} }} \right],
\end{equation}
and the iso/dianeutral diffusive operator in $z$-coordinates is then
\begin{equation}\label{Apdx_B4b}
 D^T = \left. \nabla \right|_z \cdot
           \left[ A^{lT} \;\Re \cdot \left. \nabla \right|_z T  \right]. \\
\end{equation}
\end{subequations}


Physically, the full tensor \eqref{Apdx_B3}
represents strong isoneutral diffusion on a plane parallel to the isoneutral
surface and weak dianeutral diffusion perpendicular to this plane.
However, the approximate `weak-slope' tensor \eqref{Apdx_B4a} represents strong
diffusion along the isoneutral surface, with weak
\emph{vertical}  diffusion -- the principal axes of the tensor are no
longer orthogonal. This simplification also decouples
the ($i$,$z$) and ($j$,$z$) planes of the tensor. The weak-slope operator therefore takes the same
form, \eqref{Apdx_B4}, as \eqref{Apdx_B2}, the diffusion operator for geopotential
diffusion written in non-orthogonal $i,j,s$-coordinates. Written out
explicitly,

\begin{multline} \label{Apdx_B_ldfiso}
 D^T=\frac{1}{e_1 e_2 }\left\{
 {\;\frac{\partial }{\partial i}\left[ {A_h \left( {\frac{e_2}{e_1}\frac{\partial T}{\partial i}-a_1 \frac{e_2}{e_3}\frac{\partial T}{\partial k}} \right)} \right]}
 {+\frac{\partial}{\partial j}\left[ {A_h \left( {\frac{e_1}{e_2}\frac{\partial T}{\partial j}-a_2 \frac{e_1}{e_3}\frac{\partial T}{\partial k}} \right)} \right]\;} \right\} \\
\shoveright{+\frac{1}{e_3 }\frac{\partial }{\partial k}\left[ {A_h \left( {-\frac{a_1 }{e_1 }\frac{\partial T}{\partial i}-\frac{a_2 }{e_2 }\frac{\partial T}{\partial j}+\frac{\left( {a_1 ^2+a_2 ^2+\varepsilon} \right)}{e_3 }\frac{\partial T}{\partial k}} \right)} \right]}. \\
\end{multline}


The isopycnal diffusion operator \eqref{Apdx_B4},
\eqref{Apdx_B_ldfiso} conserves tracer quantity and dissipates its
square. The demonstration of the first property is trivial as \eqref{Apdx_B4} is the divergence
of fluxes. Let us demonstrate the second one:
\begin{equation*}
\iiint\limits_D T\;\nabla .\left( {\textbf{A}}_{\textbf{I}} \nabla T \right)\,dv
          = -\iiint\limits_D \nabla T\;.\left( {\textbf{A}}_{\textbf{I}} \nabla T \right)\,dv,
\end{equation*}
and since
\begin{subequations}
\begin{align*} {\begin{array}{*{20}l}
\nabla T\;.\left( {{\rm {\bf A}}_{\rm {\bf I}} \nabla T}
\right)&=A^{lT}\left[ {\left( {\frac{\partial T}{\partial i}} \right)^2-2a_1
\frac{\partial T}{\partial i}\frac{\partial T}{\partial k}+\left(
{\frac{\partial T}{\partial j}} \right)^2} \right. \\
&\qquad \qquad \qquad
{ \left. -\,{2a_2 \frac{\partial T}{\partial j}\frac{\partial T}{\partial k}+\left( {a_1 ^2+a_2 ^2+\varepsilon} \right)\left( {\frac{\partial T}{\partial k}} \right)^2} \right]} \\
&=A_h \left[ {\left( {\frac{\partial T}{\partial i}-a_1 \frac{\partial
          T}{\partial k}} \right)^2+\left( {\frac{\partial T}{\partial
          j}-a_2 \frac{\partial T}{\partial k}} \right)^2}
  +\varepsilon \left(\frac{\partial T}{\partial k}\right) ^2\right]      \\
& \geq 0
\end{array} }
\end{align*}
\end{subequations}
\addtocounter{equation}{-1}
 the property becomes obvious.

\subsubsection*{In generalized vertical coordinates}

Because the weak-slope operator \eqref{Apdx_B4}, \eqref{Apdx_B_ldfiso} is decoupled
in the ($i$,$z$) and ($j$,$z$) planes, it may be transformed into
generalized $s$-coordinates in the same way as \eqref{Apdx_B_1} was transformed into
\eqref{Apdx_B_2}. The resulting operator then takes the simple form

\begin{equation} \label{Apdx_B_ldfiso_s}
D^T = \left. \nabla \right|_s \cdot
           \left[ A^{lT} \;\Re \cdot \left. \nabla \right|_s T  \right] \\
\;\;\text{where} \;\Re =\left( {{\begin{array}{*{20}c}
 1 \hfill & 0 \hfill & {-r _1 } \hfill \\
 0 \hfill & 1 \hfill & {-r _2 } \hfill \\
 {-r _1 } \hfill & {-r _2 } \hfill & {\varepsilon +r _1
^2+r _2 ^2} \hfill \\
\end{array} }} \right),
\end{equation}

where ($r_1$, $r_2$) are the isopycnal slopes in ($\textbf{i}$,
$\textbf{j}$) directions, relative to $s$-coordinate surfaces:
\begin{equation*}
r_1 =\frac{e_3 }{e_1 }\left( {\frac{\partial \rho }{\partial i}} \right)\left( {\frac{\partial \rho }{\partial s}} \right)^{-1}
\qquad , \qquad
r_2 =\frac{e_3 }{e_2 }\left( {\frac{\partial \rho }{\partial j}}
\right)\left( {\frac{\partial \rho }{\partial s}} \right)^{-1}.
\end{equation*}

To prove  \eqref{Apdx_B5}  by direct re-expression of \eqref{Apdx_B_ldfiso} is
straightforward, but laborious. An easier way is first to note (by reversing the
derivation of \eqref{Apdx_B_2} from \eqref{Apdx_B_1} ) that the
weak-slope operator may be \emph{exactly} reexpressed in 
non-orthogonal $i,j,\rho$-coordinates as

\begin{equation} \label{Apdx_B5}
D^T = \left. \nabla \right|_\rho \cdot
           \left[ A^{lT} \;\Re \cdot \left. \nabla \right|_\rho T  \right] \\
\;\;\text{where} \;\Re =\left( {{\begin{array}{*{20}c}
 1 \hfill & 0 \hfill &0 \hfill \\
 0 \hfill & 1 \hfill & 0 \hfill \\
0 \hfill & 0 \hfill & \varepsilon \hfill \\
\end{array} }} \right).
\end{equation}
Then direct transformation from $i,j,\rho$-coordinates to
$i,j,s$-coordinates gives \eqref{Apdx_B_ldfiso_s} immediately.

Note that the weak-slope approximation is only made in
transforming from the (rotated,orthogonal) isoneutral axes to the
non-orthogonal $i,j,\rho$-coordinates. The further transformation
into $i,j,s$-coordinates is exact, whatever the steepness of
the  $s$-surfaces, in the same way as the transformation of
horizontal/vertical Laplacian diffusion in $z$-coordinates,
\eqref{Apdx_B_1} onto $s$-coordinates is exact, however steep the $s$-surfaces.


% ================================================================
% Lateral/Vertical Momentum Diffusive Operators
% ================================================================
\section{Lateral/Vertical momentum diffusive operators}
\label{Apdx_B_3}

The second order momentum diffusion operator (Laplacian) in the $z$-coordinate
is found by applying \eqref{Eq_PE_lap_vector}, the expression for the Laplacian
of a vector,  to the horizontal velocity vector :
\begin{align*}
\Delta {\textbf{U}}_h
&=\nabla \left( {\nabla \cdot {\textbf{U}}_h } \right)-
\nabla \times \left( {\nabla \times {\textbf{U}}_h } \right)    \\
\\
&=\left( {{\begin{array}{*{20}c}
 {\frac{1}{e_1 }\frac{\partial \chi }{\partial i}} \hfill \\
 {\frac{1}{e_2 }\frac{\partial \chi }{\partial j}} \hfill \\
 {\frac{1}{e_3 }\frac{\partial \chi }{\partial k}} \hfill \\
\end{array} }} \right)-\left( {{\begin{array}{*{20}c}
 {\frac{1}{e_2 }\frac{\partial \zeta }{\partial j}-\frac{1}{e_3
}\frac{\partial }{\partial k}\left( {\frac{1}{e_3 }\frac{\partial
u}{\partial k}} \right)} \hfill \\
 {\frac{1}{e_3 }\frac{\partial }{\partial k}\left( {-\frac{1}{e_3
}\frac{\partial v}{\partial k}} \right)-\frac{1}{e_1 }\frac{\partial \zeta
}{\partial i}} \hfill \\
 {\frac{1}{e_1 e_2 }\left[ {\frac{\partial }{\partial i}\left( {\frac{e_2
}{e_3 }\frac{\partial u}{\partial k}} \right)-\frac{\partial }{\partial
j}\left( {-\frac{e_1 }{e_3 }\frac{\partial v}{\partial k}} \right)} \right]}
\hfill \\
\end{array} }} \right)
\\
\\
&=\left( {{\begin{array}{*{20}c}
{\frac{1}{e_1 }\frac{\partial \chi }{\partial i}-\frac{1}{e_2 }\frac{\partial \zeta }{\partial j}} \\
{\frac{1}{e_2 }\frac{\partial \chi }{\partial j}+\frac{1}{e_1 }\frac{\partial \zeta }{\partial i}} \\
0 \\
\end{array} }} \right)
+\frac{1}{e_3 }
\left( {{\begin{array}{*{20}c}
{\frac{\partial }{\partial k}\left( {\frac{1}{e_3 }\frac{\partial u}{\partial k}} \right)} \\
{\frac{\partial }{\partial k}\left( {\frac{1}{e_3 }\frac{\partial v}{\partial k}} \right)} \\
{\frac{\partial \chi }{\partial k}-\frac{1}{e_1 e_2 }\left( {\frac{\partial ^2\left( {e_2 \,u} \right)}{\partial i\partial k}+\frac{\partial ^2\left( {e_1 \,v} \right)}{\partial j\partial k}} \right)} \\
\end{array} }} \right)
\end{align*}
Using \eqref{Eq_PE_div}, the definition of the horizontal divergence, the third
componant of the second vector is obviously zero and thus :
\begin{equation*}
\Delta {\textbf{U}}_h = \nabla _h \left( \chi \right) - \nabla _h \times \left( \zeta \right) + \frac {1}{e_3 } \frac {\partial }{\partial k} \left( {\frac {1}{e_3 } \frac{\partial {\textbf{ U}}_h }{\partial k}} \right)
\end{equation*}

Note that this operator ensures a full separation between the vorticity and horizontal
divergence fields (see Appendix~\ref{Apdx_C}). It is only equal to a Laplacian
applied to each component in Cartesian coordinates, not on the sphere.

The horizontal/vertical second order (Laplacian type) operator used to diffuse
horizontal momentum in the $z$-coordinate therefore takes the following form :
\begin{equation} \label{Apdx_B_Lap_U}
 {\textbf{D}}^{\textbf{U}} =
 	  \nabla _h \left( {A^{lm}\;\chi } \right)
	- \nabla _h \times \left( {A^{lm}\;\zeta \;{\textbf{k}}} \right)
	+ \frac{1}{e_3 }\frac{\partial }{\partial k}\left( {\frac{A^{vm}\;}{e_3 }
				\frac{\partial {\rm {\bf U}}_h }{\partial k}} \right) \\
\end{equation}
that is, in expanded form:
\begin{align*}
D^{\textbf{U}}_u
& = \frac{1}{e_1} \frac{\partial \left( {A^{lm}\chi   } \right)}{\partial i}
     -\frac{1}{e_2} \frac{\partial \left( {A^{lm}\zeta } \right)}{\partial j}
     +\frac{1}{e_3} \frac{\partial u}{\partial k}      \\
D^{\textbf{U}}_v
& = \frac{1}{e_2 }\frac{\partial \left( {A^{lm}\chi   } \right)}{\partial j}
     +\frac{1}{e_1 }\frac{\partial \left( {A^{lm}\zeta } \right)}{\partial i}
     +\frac{1}{e_3} \frac{\partial v}{\partial k}
\end{align*}

Note Bene: introducing a rotation in \eqref{Apdx_B_Lap_U} does not lead to a
useful expression for the iso/diapycnal Laplacian operator in the $z$-coordinate.
Similarly, we did not found an expression of practical use for the geopotential
horizontal/vertical Laplacian operator in the $s$-coordinate. Generally,
\eqref{Apdx_B_Lap_U} is used in both $z$- and $s$-coordinate systems, that is
a Laplacian diffusion is applied on momentum along the coordinate directions.
\end{document}
