The two unit vectors of the plane are denoted by
\[ \begin{align} \mathbf{i} &\coloneqq \left(\begin{array}{c} 1\\ 0 \end{array}\right),\\ \mathbf{j} &\coloneqq \left(\begin{array}{c} 0\\ 1 \end{array}\right)\\ \end{align} \]
Now one further defines a three-element generating system $\left(\mathbf{i}_1, \mathbf{i}_2 \mathbf{i}_3\right)$ by
\[ \begin{align} \mathbf{i}_1 &\coloneqq \left(\begin{array}{c} 1\\ 0 \end{array}\right) = \mathbf{i},\\ \mathbf{i}_2 &\coloneqq \left(\begin{array}{c} -\sin\left(30^\circ\right)\\ \cos\left(30^\circ\right) \end{array}\right) = \left(\begin{array}{c} -\frac{1}{2}\\ \frac{\sqrt{3}}{2} \end{array}\right) = -\frac{1}{2}\mathbf{i} + \frac{\sqrt{3}}{2}\mathbf{j},\\ \mathbf{i}_3 &\coloneqq \left(\begin{array}{c} -\frac{1}{2}\\ -\frac{\sqrt{3}}{2} \end{array}\right) = -\frac{1}{2}\mathbf{i} - \frac{\sqrt{3}}{2}\mathbf{j}, \end{align} \]
whose elements are each rotated by 120$^\circ$ relative to one another. We further define a generating system $\left(\mathbf{j}_1, \mathbf{j}_2 \mathbf{j}_3\right)$ rotated by 90$^\circ$ relative to it:
\[ \begin{align} \mathbf{j}_1 &\coloneqq \mathbf{k}\times\mathbf{i}_1 = \left(\begin{array}{c} 0\\ 0\\ 1 \end{array}\right)\times\left(\begin{array}{c} 1\\ 0\\ 0 \end{array}\right) = \left(\begin{array}{c} 0\\ 1 \end{array}\right),\\ \mathbf{j}_2 &\coloneqq \mathbf{k}\times\mathbf{i}_2 = \left(\begin{array}{c} 0\\ 0\\ 1 \end{array}\right)\times\left(\begin{array}{c} -\frac{1}{2}\\ \frac{\sqrt{3}}{2}\\ 0 \end{array}\right) = \left(\begin{array}{c} -\frac{\sqrt{3}}{2}\\ -\frac{1}{2} \end{array}\right) = -\frac{\sqrt{3}}{2}\mathbf{i} - \frac{1}{2}\mathbf{j},\\ \mathbf{j}_3 &\coloneqq \mathbf{k}\times\mathbf{i}_3 = \left(\begin{array}{c} 0\\ 0\\ 1 \end{array}\right)\times\left(\begin{array}{c} -\frac{1}{2}\\ -\frac{\sqrt{3}}{2}\\ 0 \end{array}\right) = \left(\begin{array}{c} \frac{\sqrt{3}}{2}\\ -\frac{1}{2} \end{array}\right) = \frac{\sqrt{3}}{2}\mathbf{i} - \frac{1}{2}\mathbf{j}. \end{align} \]
One observes
\[ \begin{align} \mathbf{i}_1 + \mathbf{i}_2 + \mathbf{i}_3 = \mathbf{0}, & {} & \mathbf{j}_1 + \mathbf{j}_2 + \mathbf{j}_3 = \mathbf{0}. \end{align} \]
For a two-dimensional vector $\mathbf{v}$ one can write
\[ \begin{align} \mathbf{v} = u\mathbf{i} + v\mathbf{j} \end{align} \]
where
\[ \begin{align} u = \mathbf{i}\cdot\mathbf{v}, & {} & v = \mathbf{j}\cdot\mathbf{v}. \end{align} \]
Since the $\mathbf{i}_k$, $\mathbf{j}_l$ are each pairwise linearly independent, one can write
\[ \begin{align} \mathbf{v} = u_k'\mathbf{i}_k + u_l'\mathbf{i}_l = v_k'\mathbf{j}_k + v_l'\mathbf{j}_l. \end{align} \]
Since the choice of $k$, $l$ is not unique and one could also use all three unit vectors, one can introduce another linear condition. One uses this freedom to require the notation
\[ \begin{align} \mathbf{v} &= \frac{2}{3}\left(u_1\mathbf{i}_1 + u_2\mathbf{i}_2 + u_3\mathbf{i}_3\right) = \frac{2}{3}\left(v_1\mathbf{j}_1 + v_2\mathbf{j}_2 + v_3\mathbf{j}_3\right) \end{align} \]
This is satisfied by
\[ \begin{align} u_k = \mathbf{i}_k\cdot\mathbf{v}, & {} & v_k = \mathbf{j}_k\cdot\mathbf{v} \end{align} \]
since it follows from this
\[ \begin{align} \mathbf{v} &= u\mathbf{i} + v\mathbf{j} = \frac{2}{3}\left(\frac{6}{4}u\mathbf{i} + \frac{6}{4}v\mathbf{j}\right)\nonumber\\ &= \frac{2}{3}\left(u\mathbf{i} + \frac{1}{4}u\mathbf{i} - \frac{\sqrt{3}}{4}u\mathbf{j} - \frac{\sqrt{3}}{4}v\mathbf{i} + \frac{3}{4}v\mathbf{j} + \frac{1}{4}u\mathbf{i} + \frac{\sqrt{3}}{4}u\mathbf{j} + \frac{\sqrt{3}}{4}v\mathbf{i} + \frac{3}{4}v\mathbf{j}\right)\nonumber\\ &= \frac{2}{3}\left(u\mathbf{i} + \left(-\frac{1}{2}u + \frac{\sqrt{3}}{2}v\right)\left(-\frac{1}{2}\mathbf{i} + \frac{\sqrt{3}}{2}\mathbf{j}\right) + \left(-\frac{1}{2}u - \frac{\sqrt{3}}{2}v\right)\left(-\frac{1}{2}\mathbf{i} - \frac{\sqrt{3}}{2}\mathbf{j}\right)\right)\nonumber\\ &= \frac{2}{3}\left(u_1\mathbf{i}_1 + u_2\mathbf{i}_2 + u_3\mathbf{i}_3\right), \end{align} \]
where
\[ \begin{align} u_1 &= u, & {} & u_2 = -\frac{1}{2}u + \frac{\sqrt{3}}{2}v, & {} & u_3 = -\frac{1}{2}u - \frac{\sqrt{3}}{2}v. \end{align} \]
was used. This applies analogously to the $v_k$. From this it follows
and analogously for the $v_k$. For the gradient $\nabla\alpha$ of a scalar field $\alpha$ one now obtains
\[ \begin{align} \nabla\alpha = \frac{2}{3}\left[\left(\mathbf{i}_1\cdot\nabla\alpha\right)\mathbf{i}_1 + \left(\mathbf{i}_2\cdot\nabla\alpha\right)\mathbf{i}_2 + \left(\mathbf{i}_3\cdot\nabla\alpha\right)\mathbf{i}_3\right]. \end{align} \]
Because of
\[ \begin{align} \mathbf{i}_k\cdot\nabla\alpha = \frac{\partial\alpha}{\partial x_k} \end{align} \]
it follows
With this one furthermore establishes
It holds that
\[ \begin{align} \mathbf{i}_1 &= \frac{1}{\sqrt{3}}\left(\mathbf{j}_3 - \mathbf{j}_2\right) \end{align} \]
and cyclically:
\[ \begin{align} \mathbf{i}_2 = \frac{1}{\sqrt{3}}\left(\mathbf{j}_1 - \mathbf{j}_3\right), & {} & \mathbf{i}_3 = \frac{1}{\sqrt{3}}\left(\mathbf{j}_2 - \mathbf{j}_1\right) \end{align} \]
From this it follows
\[ \begin{align} \nabla\alpha = \frac{2}{3\sqrt{3}}\left(\frac{\partial\alpha}{\partial x_1}\left(\mathbf{j}_3 - \mathbf{j}_2\right) + \frac{\partial\alpha}{\partial x_2}\left(\mathbf{j}_1 - \mathbf{j}_3\right) + \frac{\partial\alpha}{\partial x_3}\left(\mathbf{j}_2 - \mathbf{j}_1\right)\right). \end{align} \]
Rearranging yields
\[ \begin{align} \nabla\alpha = \frac{2}{3\sqrt{3}}\left(\left(\frac{\partial\alpha}{\partial x_2} - \frac{\partial\alpha}{\partial x_3}\right)\mathbf{j}_1 + \left(\frac{\partial\alpha}{\partial x_3} - \frac{\partial\alpha}{\partial x_1}\right)\mathbf{j}_2 + \left(\frac{\partial\alpha}{\partial x_1} - \frac{\partial\alpha}{\partial x_2}\right)\mathbf{j}_3\right). \end{align} \]
One has
\[ \begin{align} \mathbf{v} &= \frac{2}{3}\left(u_1\mathbf{i}_1 + u_2\mathbf{i}_2 + u_3\mathbf{i}_3\right) = \frac{2}{3}\left[u_1\left(\begin{array}{c}1\\0\end{array}\right) + u_2\left(\begin{array}{c}-\frac{1}{2}\\\frac{\sqrt{3}}{2}\end{array}\right) + u_3\left(\begin{array}{c}-\frac{1}{2}\\-\frac{\sqrt{3}}{2}\end{array}\right)\right] = \frac{2}{3}\left(\begin{array}{c}u_1 - \frac{u_2}{2} - \frac{u_3}{2}\\\frac{\sqrt{3}}{2}u_2 - \frac{\sqrt{3}}{2}u_3\end{array}\right)\nonumber\\ &= \left(\begin{array}{c}\frac{2}{3}u_1 - \frac{u_2}{3} - \frac{u_3}{3}\\\frac{1}{\sqrt{3}}u_2 - \frac{1}{\sqrt{3}}u_3\end{array}\right). \end{align} \]
From this one obtains for the divergence
\[ \begin{align} D \coloneqq \nabla\cdot\mathbf{v} = \frac{2}{3}\frac{u_1}{\partial x} - \frac{1}{3}\frac{u_2}{\partial x} - \frac{1}{3}\frac{u_3}{\partial x} + \frac{1}{\sqrt{3}}\frac{u_2}{\partial y} - \frac{1}{\sqrt{3}}\frac{u_3}{\partial y}. \end{align} \]
One has
\[ \begin{align} \frac{\partial}{\partial x} = \frac{\partial}{\partial x_1}, \end{align} \]
because of
\[ \begin{align} \mathbf{j} = \frac{1}{\sqrt{3}}\left(\mathbf{i}_2 - \mathbf{i}_3\right) \end{align} \]
one furthermore has
\[ \begin{align} \frac{\partial}{\partial y} = \frac{1}{\sqrt{3}}\left(\frac{\partial}{\partial x_2} - \frac{\partial}{\partial x_3}\right).\tag{27.31}\label{eq:ddy_hex} \end{align} \]
From this it follows
\[ \begin{align} D &= \frac{2}{3}\frac{\partial u_1}{\partial x} - \frac{1}{3}\frac{\partial u_2}{\partial x} - \frac{1}{3}\frac{\partial u_3}{\partial x} + \frac{1}{3}\left(\frac{\partial u_2}{\partial x_2} - \frac{\partial u_2}{\partial x_3}\right) - \frac{1}{3}\left(\frac{\partial u_3}{\partial x_2} - \frac{\partial u_3}{\partial x_3}\right)\nonumber\\ &= \frac{2}{3}\frac{\partial u_1}{\partial x_1} + \frac{1}{3}\frac{u_2}{\partial x_2} + \frac{1}{3}\frac{u_3}{\partial x_3} - \frac{1}{3}\frac{u_2}{\partial x_1} - \frac{1}{3}\frac{\partial u_3}{\partial x_1} - \frac{1}{3}\frac{\partial u_2}{\partial x_3} - \frac{1}{3}\frac{\partial u_3}{\partial x_2}\nonumber\\ &= \frac{2}{3}\left(\frac{\partial u_1}{\partial x_1} + \frac{\partial u_2}{\partial x_2} + \frac{\partial u_3}{\partial x_3}\right) - \frac{1}{3}\left(\frac{\partial u_2}{\partial x_1} + \frac{\partial u_2}{\partial x_2} + \frac{\partial u_2}{\partial x_3} + \frac{\partial u_3}{\partial x_1} + \frac{\partial u_3}{\partial x_2} + \frac{\partial u_3}{\partial x_3}\right). \end{align} \]
With Eq. (27.22) follows
Combining with Eq. (27.21), one further obtains
\[ \begin{align} \Delta\alpha = \frac{2}{3}\left(\frac{\partial^2\alpha}{\partial x_1^2} + \frac{\partial^2\alpha}{\partial x_2^2} + \frac{\partial^2\alpha}{\partial x_3^2}\right). \end{align} \]
One has
\[ \begin{align} u_1 &= \mathbf{v}\cdot\mathbf{i}_1 = \frac{2}{3}\left(v_1\mathbf{j}_1\cdot\mathbf{i}_1 + v_2\mathbf{j}_2\cdot\mathbf{i}_1 + v_3\mathbf{j}_3\cdot\mathbf{i}_1\right) = \frac{2}{3}\left(0\cdot v_1 - \frac{\sqrt{3}}{2}v_2 + \frac{\sqrt{3}}{2}v_3\right)\nonumber\\ &= \frac{1}{\sqrt{3}}\left(v_3 - v_2\right). \end{align} \]
By cyclically shifting the indices one obtains
\[ \begin{align} u_2 = \frac{1}{\sqrt{3}}\left(v_1 - v_3\right), & {} & u_3 = \frac{1}{\sqrt{3}}\left(v_2 - v_1\right). \end{align} \]
This gives another representation of the divergence:
\[ \begin{align} D = \frac{2}{3\sqrt{3}}\left[\left(\frac{\partial v_2}{\partial x_3} - \frac{\partial v_3}{\partial x_2}\right) + \left(\frac{\partial v_3}{\partial x_1} - \frac{\partial v_1}{\partial x_3}\right) + \left(\frac{\partial v_1}{\partial x_2} - \frac{\partial v_2}{\partial x_1}\right)\right]\tag{27.37}\label{eq:div_3elements2d_mod} \end{align} \]
One has
\[ \begin{align} \frac{\partial v_2}{\partial x_3} - \frac{\partial v_3}{\partial x_2} &= \frac{\partial v_2}{\partial x_3} + \frac{\partial v_2}{\partial x_2} - \frac{\partial v_2}{\partial x_2} - \frac{\partial v_3}{\partial x_2} = \left(\frac{\partial}{\partial x_3} + \frac{\partial}{\partial x_2}\right)v_2 - \frac{\partial}{\partial x_2}\left(v_2 + v_3\right)\nonumber\\ &= -\frac{\partial v_2}{\partial x_1} + \frac{\partial v_1}{\partial x_2} = \frac{\partial v_1}{\partial x_2} - \frac{\partial v_2}{\partial x_1}. \end{align} \]
By cyclically permuting this equation, one sees that all three bracketed summands in Eq. (27.37) are equal. From this follow
\[ \begin{align} D &= \frac{2}{\sqrt{3}}\left(\frac{\partial v_2}{\partial x_3} - \frac{\partial v_3}{\partial x_2}\right),\\ D &= \frac{2}{\sqrt{3}}\left(\frac{\partial v_3}{\partial x_1} - \frac{\partial v_1}{\partial x_3}\right),\\ D &= \frac{2}{\sqrt{3}}\left(\frac{\partial v_1}{\partial x_2} - \frac{\partial v_2}{\partial x_1}\right). \end{align} \]
One has
\[ \begin{align} \mathbf{v} &= \frac{2}{3}\left(v_1\mathbf{j}_1 + v_2\mathbf{j}_2 + v_3\mathbf{j}_3\right) = \frac{2}{3}\left[v_1\left(\begin{array}{c}0\\1\end{array}\right) + v_2\left(\begin{array}{c}-\frac{\sqrt{3}}{2}\\-\frac{1}{2}\end{array}\right) + v_3\left(\begin{array}{c}\frac{\sqrt{3}}{2}\\-\frac{1}{2}\end{array}\right)\right]\nonumber\\ &= \left(\begin{array}{c}\frac{1}{\sqrt{3}}v_3 - \frac{1}{\sqrt{3}}v_2\\\frac{2}{3}v_1 - \frac{v_2}{3} - \frac{v_3}{3}\end{array}\right). \end{align} \]
For the vorticity one obtains
\[ \begin{align} \zeta &= \frac{\partial v}{\partial x} - \frac{\partial u}{\partial y} = \frac{\partial}{\partial x}\left(\frac{2}{3}v_1 - \frac{v_2}{3} - \frac{v_3}{3}\right) - \frac{\partial}{\partial y}\left(\frac{1}{\sqrt{3}}v_3 - \frac{1}{\sqrt{3}}v_2\right)\nonumber\\ &= \frac{\partial}{\partial x_1}\left(\frac{2}{3}v_1 - \frac{v_2}{3} - \frac{v_3}{3}\right) - \frac{1}{\sqrt{3}}\left(\frac{\partial}{\partial x_2} - \frac{\partial}{\partial x_3}\right)\left(\frac{1}{\sqrt{3}}v_3 - \frac{1}{\sqrt{3}}v_2\right)\nonumber\\ &= \frac{2}{3}\frac{\partial v_1}{\partial x_1} + \frac{1}{3}\frac{\partial v_2}{\partial x_2} + \frac{1}{3}\frac{\partial v_3}{\partial x_3} - \frac{1}{3}\frac{\partial v_2}{\partial x_1} - \frac{1}{3}\frac{\partial v_3}{\partial x_1} - \frac{1}{3}\frac{\partial v_3}{\partial x_2} - \frac{1}{3}\frac{\partial v_2}{\partial x_3}\nonumber\\ &= \frac{2}{3}\left(\frac{\partial v_1}{\partial x_1} + \frac{\partial v_2}{\partial x_2} + \frac{\partial v_3}{\partial x_3}\right) - \frac{1}{3}\left(\frac{\partial v_2}{\partial x_1} + \frac{\partial v_2}{\partial x_2} + \frac{\partial v_2}{\partial x_3} + \frac{\partial v_3}{\partial x_1} + \frac{\partial v_3}{\partial x_2} + \frac{\partial v_3}{\partial x_3}\right). \end{align} \]
With Eq. (27.22) follows
\[ \begin{align} \zeta = \frac{2}{3}\left(\frac{\partial v_1}{\partial x_1} + \frac{\partial v_2}{\partial x_2} + \frac{\partial v_3}{\partial x_3}\right). \end{align} \]
One has
\[ \begin{align} v_1 &= \mathbf{v}\cdot\mathbf{j}_1 = \frac{2}{3}\left(u_1\mathbf{i}_1 + u_2\mathbf{i}_2 + u_3\mathbf{i}_3\right)\cdot\mathbf{j}_1 = \frac{2}{3}\left(u_1\mathbf{i}_1 + u_2\mathbf{i}_2 + u_3\mathbf{i}_3\right)\cdot\mathbf{j}\nonumber\\ &= \frac{2}{3}\left(\frac{\sqrt{3}}{2}u_2 - \frac{\sqrt{3}}{2}u_3\right) = \frac{1}{\sqrt{3}}\left(u_2 - u_3\right).\tag{27.45}\label{eq:turn_3elemets2d} \end{align} \]
By cyclically shifting the indices one obtains
\[ \begin{align} v_2 = \frac{1}{\sqrt{3}}\left(u_3 - u_1\right), & {} & v_3 = \frac{1}{\sqrt{3}}\left(u_1 - u_2\right). \end{align} \]
From this follows another representation of the vorticity:
\[ \begin{align} \zeta = \frac{2}{3\sqrt{3}}\left(\frac{\partial u_2}{\partial x_1} - \frac{\partial u_3}{\partial x_1} + \frac{\partial u_3}{\partial x_2} - \frac{\partial u_1}{\partial x_2} + \frac{\partial u_1}{\partial x_3} - \frac{\partial u_2}{\partial x_3}\right) \end{align} \]
By rearranging one obtains
\[ \begin{align} \zeta = \frac{2}{3\sqrt{3}}\left[\left(\frac{\partial u_3}{\partial x_2} - \frac{\partial u_2}{\partial x_3}\right) + \left(\frac{\partial u_1}{\partial x_3} - \frac{\partial u_3}{\partial x_1}\right) + \left(\frac{\partial u_2}{\partial x_1} - \frac{\partial u_1}{\partial x_2}\right)\right]. \end{align} \]
Analogously to the case of Eq. (27.37), the following hold here as well
\[ \begin{align} \zeta &= \frac{2}{\sqrt{3}}\left(\frac{\partial u_3}{\partial x_2} - \frac{\partial u_2}{\partial x_3}\right),\tag{27.49}\label{eq:hex_curl_symm_0}\\ \zeta &= \frac{2}{\sqrt{3}}\left(\frac{\partial u_1}{\partial x_3} - \frac{\partial u_3}{\partial x_1}\right),\tag{27.50}\label{eq:hex_curl_symm_1}\\ \zeta &= \frac{2}{\sqrt{3}}\left(\frac{\partial u_2}{\partial x_1} - \frac{\partial u_1}{\partial x_2}\right).\tag{27.51}\label{eq:hex_curl_symm_2} \end{align} \]
Using the fundamental theorem of vector analysis one can write the vector field $\mathbf{v}$ in the form
\[ \begin{align} \mathbf{v} = \mathbf{k}\times\nabla\psi + \nabla\chi \end{align} \]
where
\[ \begin{align} \mathbf{k}\times\nabla\psi &= \frac{2}{3}\left(\frac{\partial\psi}{\partial x_1}\mathbf{j}_1 + \frac{\partial\psi}{\partial x_2}\mathbf{j}_2 + \frac{\partial\psi}{\partial x_3}\mathbf{j}_3\right),\\ \nabla\chi &= \frac{2}{3}\left(\frac{\partial\chi}{\partial x_1}\mathbf{i}_1 + \frac{\partial\chi}{\partial x_2}\mathbf{i}_2 + \frac{\partial\chi}{\partial x_3}\mathbf{i}_3\right). \end{align} \]
From this follow
\[ \begin{align} u_1 &= \frac{1}{\sqrt{3}}\left(\frac{\partial\psi}{\partial x_3} - \frac{\partial\psi}{\partial x_2}\right) + \frac{\partial\chi}{\partial x_1} = -\frac{1}{\sqrt{3}}\left(\frac{\partial\psi}{\partial x_2} - \frac{\partial\psi}{\partial x_3}\right) + \frac{\partial\chi}{\partial x_1}\tag{27.55}\label{eq:trivariate_helmholtz_eq_0},\\ u_2 &= -\frac{1}{\sqrt{3}}\left(\frac{\partial\psi}{\partial x_3} - \frac{\partial\psi}{\partial x_1}\right) + \frac{\partial\chi}{\partial x_2}\tag{27.56}\label{eq:trivariate_helmholtz_eq_1},\\ u_3 &= -\frac{1}{\sqrt{3}}\left(\frac{\partial\psi}{\partial x_1} - \frac{\partial\psi}{\partial x_2}\right) + \frac{\partial\chi}{\partial x_3}\tag{27.57}\label{eq:trivariate_helmholtz_eq_2}. \end{align} \]
and
\[ \begin{align} v_1 &= \frac{\partial\psi}{\partial x_1} + \frac{1}{\sqrt{3}}\left(\frac{\partial\chi}{\partial x_2} - \frac{\partial\chi}{\partial x_3}\right)\tag{27.58}\label{eq:trivariate_helmholtz_eq_3},\\ v_2 &= \frac{\partial\psi}{\partial x_2} + \frac{1}{\sqrt{3}}\left(\frac{\partial\chi}{\partial x_3} - \frac{\partial\chi}{\partial x_1}\right)\tag{27.59}\label{eq:trivariate_helmholtz_eq_4},\\ v_3 &= \frac{\partial\psi}{\partial x_3} + \frac{1}{\sqrt{3}}\left(\frac{\partial\chi}{\partial x_1} - \frac{\partial\chi}{\partial x_2}\right)\tag{27.60}\label{eq:trivariate_helmholtz_eq_5}. \end{align} \]
Because of
\[ \begin{align} D = \Delta\chi, & {} & \zeta = \Delta\psi \end{align} \]
the following further hold
\[ \begin{align} \Delta u_1 &= -\frac{1}{\sqrt{3}}\left(\frac{\partial\zeta}{\partial x_2} - \frac{\partial\zeta}{\partial x_3}\right) + \frac{\partial D}{\partial x_1},\tag{27.62}\label{eq:hex_laplace_u_1}\\ \Delta u_2 &= -\frac{1}{\sqrt{3}}\left(\frac{\partial\zeta}{\partial x_3} - \frac{\partial\zeta}{\partial x_1}\right) + \frac{\partial D}{\partial x_2},\tag{27.63}\label{eq:hex_laplace_u_2}\\ \Delta u_3 &= -\frac{1}{\sqrt{3}}\left(\frac{\partial\zeta}{\partial x_1} - \frac{\partial\zeta}{\partial x_2}\right) + \frac{\partial D}{\partial x_3}\tag{27.64}\label{eq:hex_laplace_u_3}, \end{align} \]
\[ \begin{align} \Delta v_1 &= \frac{\partial\zeta}{\partial x_1} + \frac{1}{\sqrt{3}}\left(\frac{\partial D}{\partial x_2} - \frac{\partial D}{\partial x_3}\right),\\ \Delta v_2 &= \frac{\partial\zeta}{\partial x_2} + \frac{1}{\sqrt{3}}\left(\frac{\partial D}{\partial x_3} - \frac{\partial D}{\partial x_1}\right),\\ \Delta v_3 &= \frac{\partial\zeta}{\partial x_3} + \frac{1}{\sqrt{3}}\left(\frac{\partial D}{\partial x_1} - \frac{\partial D}{\partial x_2}\right). \end{align} \]
Let $\alpha$ be a field of the form
\[ \begin{align} \alpha\left(\mathbf{r}\right) = \exp\left(i\mathbf{k}\cdot\mathbf{r}\right). \end{align} \]
The distance between the centers of two hexagons is denoted by $d$ (from now on this quantity is referred to as the lattice constant). One now defines, for $1 \leq j \leq 3$, the central difference quotient in the $j$-direction by
\[ \begin{align} \delta_j\alpha &= \frac{\alpha\left(\mathbf{r} + \mathbf{i}_j\frac{d}{2}\right) - \alpha\left(\mathbf{r} - \mathbf{i}_j\frac{d}{2}\right)}{d}. \end{align} \]
For this one has
\[ \begin{align} \delta_j\alpha &= \frac{\exp\left[i\mathbf{k}\cdot\left(\mathbf{r} + \mathbf{i}_j\frac{d}{2}\right)\right] - \exp\left[i\mathbf{k}\cdot\left(\mathbf{r} - \mathbf{i}_j\frac{d}{2}\right)\right]}{d}\nonumber\\ &= \exp\left(i\mathbf{k}\cdot\mathbf{r}\right)\frac{\exp\left(i\mathbf{k}\cdot\mathbf{i}_j\frac{d}{2}\right) - \exp\left(-i\mathbf{k}\cdot\mathbf{i}_j\frac{d}{2}\right)}{d} = \alpha\left(\mathbf{r}\right)\frac{\exp\left(i\mathbf{k}\cdot\mathbf{i}_j\frac{d}{2}\right) - \exp\left(-i\mathbf{k}\cdot\mathbf{i}_j\frac{d}{2}\right)}{d}\nonumber\\ &= \alpha\left(\mathbf{r}\right)\frac{\exp\left(ik_j\frac{d}{2}\right) - \exp\left(-ik_j\frac{d}{2}\right)}{d} = \alpha\left(\mathbf{r}\right)\frac{2i}{d}\sin\left(\frac{k_jd}{2}\right). \end{align} \]
If one defines
\[ \begin{align} s_j \coloneqq \sin\left(\frac{k_jd}{2}\right), \end{align} \]
one can write this in the form
\[ \begin{align} \delta_j\alpha = \frac{2i}{d}s_j\alpha \end{align} \]
One now further defines a simple averaging operator in the $j$-direction by
\[ \begin{align} \newoverline{\alpha}^{(j)} & \coloneqq \frac{\alpha\left(\mathbf{r} + \mathbf{i}_j\frac{d}{2}\right) + \alpha\left(\mathbf{r} - \mathbf{i}_j\frac{d}{2}\right)}{2} = \frac{\exp\left[i\mathbf{k}\cdot\left(\mathbf{r} + \mathbf{i}_j\frac{d}{2}\right)\right] + \exp\left[i\mathbf{k}\cdot\left(\mathbf{r} - \mathbf{i}_j\frac{d}{2}\right)\right]}{2}\nonumber\\ &= \exp\left(i\mathbf{k}\cdot\mathbf{r}\right)\frac{\exp\left(i\mathbf{k}\cdot\mathbf{i}_j\frac{d}{2}\right) + \exp\left(-i\mathbf{k}\cdot\mathbf{i}_j\frac{d}{2}\right)}{2}\nonumber\\ &= \alpha\left(\mathbf{r}\right)\frac{\exp\left(ik_j\frac{d}{2}\right) + \exp\left(-ik_j\frac{d}{2}\right)}{2} = \alpha\left(\mathbf{r}\right)\cos\left(\frac{k_jd}{2}\right).\tag{27.73}\label{eq:hex_simple_average} \end{align} \]
If one defines
\[ \begin{align} c_j \coloneqq \cos\left(\frac{k_jd}{2}\right), \end{align} \]
one can write this in the form
\[ \begin{align} \newoverline{\alpha}^{(j)} = c_j\alpha \end{align} \]
The linearized shallow water equations (13.173) - (13.174) on the f-plane can be written with the geopotential $\Phi \coloneqq gh$ in the form
\[ \begin{align} \frac{\partial\Phi}{\partial t} + \Phi_0\nabla\cdot\mathbf{v} = 0, & {} & \frac{\partial\mathbf{v}}{\partial t} + f_0\mathbf{k}\times\mathbf{v} - \nabla\Phi = 0 \end{align} \]
Here $\Phi_0 \coloneqq gD$. With Eqs. (27.33) and (27.45), this leads to the evolution equations for the components with respect to the three-element generating system:
\[ \begin{align} \frac{\partial\Phi}{\partial t} + \frac{2}{3}\Phi_0\left(\frac{\partial u_1}{\partial x_1} + \frac{\partial u_2}{\partial x_2} + \frac{\partial u_3}{\partial x_3}\right) &= 0,\\ \frac{\partial u_1}{\partial t} - \frac{f_0}{\sqrt{3}}\left(u_2 - u_3\right) + \frac{\partial\Phi}{\partial x_1} &= 0,\\ \frac{\partial u_2}{\partial t} - \frac{f_0}{\sqrt{3}}\left(u_3 - u_1\right) + \frac{\partial\Phi}{\partial x_2} &= 0,\\ \frac{\partial u_3}{\partial t} - \frac{f_0}{\sqrt{3}}\left(u_1 - u_2\right) + \frac{\partial\Phi}{\partial x_3} &= 0. \end{align} \]
Spatial discretization yields
\[ \begin{align} \frac{\partial\Phi}{\partial t} + \frac{2}{3}\Phi_0\left(\delta_1u_1 + \delta_2u_2 + \delta_3u_3\right) &= 0,\tag{27.81}\label{eq:swe_lin_hex_0}\\ \frac{\partial u_1}{\partial t} - \frac{f_0}{\sqrt{3}}\left(\newoverline{u_2}^{3} - \newoverline{u_3}^{(2)}\right) + \delta_1\Phi &= 0,\tag{27.82}\label{eq:swe_lin_hex_1}\\ \frac{\partial u_2}{\partial t} - \frac{f_0}{\sqrt{3}}\left(\newoverline{u_3}^{(1)} - \newoverline{u_1}^{3}\right) + \delta_2\Phi &= 0,\tag{27.83}\label{eq:swe_lin_hex_2}\\ \frac{\partial u_3}{\partial t} - \frac{f_0}{\sqrt{3}}\left(\newoverline{u_1}^{(2)} - \newoverline{u_2}^{(1)}\right) + \delta_3\Phi &= 0.\tag{27.84}\label{eq:swe_lin_hex_3} \end{align} \]
The arithmetic mean of the two nearest velocity components in the desired direction was used to reconstruct the wind components for computing the Coriolis acceleration. If one assumes a monochromatic plane wave for all fields, i.e.,
\[ \begin{align} \psi\left(\mathbf{r}, t\right) = \newhat{\psi}\exp\left(i\mathbf{k}\cdot\mathbf{r} - i\omega t\right) \end{align} \]
for any field $\psi$ with a complex amplitude $\newhat{\psi}$, one obtains
\[ \begin{align} -i\omega\newhat{\Phi} + \frac{4i}{3d}\Phi_0\left(s_1\newhat{u}_1 + s_2\newhat{u}_2 + s_3\newhat{u}_3\right) &= 0,\\ -i\omega\newhat{u}_1 - \frac{f_0}{\sqrt{3}}\left(c_3\newhat{u}_2 - c_2\newhat{u}_3\right) + \frac{2is_1}{d}\newhat{\Phi} &= 0,\\ -i\omega\newhat{u}_2 - \frac{f_0}{\sqrt{3}}\left(c_1\newhat{u}_3 - c_3\newhat{u}_1\right) + \frac{2is_1}{d}\newhat{\Phi} &= 0,\\ -i\omega\newhat{u}_3 - \frac{f_0}{\sqrt{3}}\left(c_2\newhat{u}_1 - c_1\newhat{u}_2\right) + \frac{2is_1}{d}\newhat{\Phi} &= 0. \end{align} \]
It should be noted that there is a phase offset between the equations, arising from the fact that the individual equations hold at different locations. This offset has been canceled out. In matrix form one obtains
\[ \begin{align} \left(\begin{array}{cccc} -i\omega & \frac{4i}{3d}\Phi_0s_1 & \frac{4i}{3d}\Phi_0s_2 & \frac{4i}{3d}\Phi_0s_3 \\ \frac{2is_1}{d} & -i\omega & -\frac{f_0}{\sqrt{3}}c_3 & \frac{f_0}{\sqrt{3}}c_2\\ \frac{2is_2}{d} & \frac{f_0}{\sqrt{3}}c_3 & -i\omega & -\frac{f_0}{\sqrt{3}}c_1\\ \frac{2is_3}{d} & -\frac{f_0}{\sqrt{3}}c_2 & \frac{f_0}{\sqrt{3}}c_1 & -i\omega \end{array}\right)\left(\begin{array}{c} \newhat{\Phi}\\ \newhat{u}_1\\ \newhat{u}_2\\ \newhat{u}_3 \end{array}\right) = \mathbf{0}. \end{align} \]
Multiplying this by the imaginary unit $i$ gives
\[ \begin{align} \left(\begin{array}{cccc} \omega & -\frac{4}{3d}\Phi_0s_1 & -\frac{4}{3d}\Phi_0s_2 & -\frac{4}{3d}\Phi_0s_3 \\ -\frac{2s_1}{d} & \omega & -\frac{if_0}{\sqrt{3}}c_3 & \frac{if_0}{\sqrt{3}}c_2\\ -\frac{2s_2}{d} & \frac{if_0}{\sqrt{3}}c_3 & \omega & -\frac{if_0}{\sqrt{3}}c_1\\ -\frac{2s_3}{d} & -\frac{if_0}{\sqrt{3}}c_2 & \frac{if_0}{\sqrt{3}}c_1 & \omega \end{array}\right)\left(\begin{array}{c} \newhat{\Phi}\\ \newhat{u}_1\\ \newhat{u}_2\\ \newhat{u}_3 \end{array}\right) = \mathbf{0}. \end{align} \]
Nontrivial solutions exist if and only if the determinant of the coefficient matrix of this system of linear equations vanishes:
\[ \begin{align} \left|\begin{array}{cccc} \omega & -\frac{4}{3d}\Phi_0s_1 & -\frac{4}{3d}\Phi_0s_2 & -\frac{4}{3d}\Phi_0s_3 \\ -\frac{2s_1}{d} & \omega & -\frac{if_0}{\sqrt{3}}c_3 & \frac{if_0}{\sqrt{3}}c_2\\ -\frac{2s_2}{d} & \frac{if_0}{\sqrt{3}}c_3 & \omega & -\frac{if_0}{\sqrt{3}}c_1\\ -\frac{2s_3}{d} & -\frac{if_0}{\sqrt{3}}c_2 & \frac{if_0}{\sqrt{3}}c_1 & \omega \end{array}\right| &\stackrel{!}{=} 0\nonumber\\ \Leftrightarrow\omega\left|\begin{array}{ccc} \omega & -\frac{if_0}{\sqrt{3}}c_3 & \frac{if_0}{\sqrt{3}}c_2\\ \frac{if_0}{\sqrt{3}}c_3 & \omega & -\frac{if_0}{\sqrt{3}}c_1\\ -\frac{if_0}{\sqrt{3}}c_2 & \frac{if_0}{\sqrt{3}}c_1 & \omega \end{array} \right| + \frac{4}{3d}\Phi_0s_1\left|\begin{array}{ccc} -\frac{2s_1}{d} & -\frac{if_0}{\sqrt{3}}c_3 & \frac{if_0}{\sqrt{3}}c_2\\ -\frac{2s_2}{d} & \omega & -\frac{if_0}{\sqrt{3}}c_1\\ -\frac{2s_3}{d} & \frac{if_0}{\sqrt{3}}c_1 & \omega \end{array}\right| & \nonumber\\ - \frac{4}{3d}\Phi_0s_2\left|\begin{array}{ccc} -\frac{2s_1}{d} & \omega & \frac{if_0}{\sqrt{3}}c_2\\ -\frac{2s_2}{d} & \frac{if_0}{\sqrt{3}}c_3 & -\frac{if_0}{\sqrt{3}}c_1\\ -\frac{2s_3}{d} & -\frac{if_0}{\sqrt{3}}c_2 & \omega \end{array} \right| + \frac{4}{3d}\Phi_0s_3\left|\begin{array}{ccc} -\frac{2s_1}{d} & \omega & -\frac{if_0}{\sqrt{3}}c_3\\ -\frac{2s_2}{d} & \frac{if_0}{\sqrt{3}}c_3 & \omega\\ -\frac{2s_3}{d} & -\frac{if_0}{\sqrt{3}}c_2 & \frac{if_0}{\sqrt{3}}c_1 \end{array} \right| &= 0 \end{align} \] \[ \begin{align} \Leftrightarrow \omega^4 - \frac{\omega^2f_0^2}{3}\left(c_1^2 + c_2^2 + c_3^2\right) + \frac{4}{3d}\Phi_0s_1\left|\begin{array}{ccc} -\frac{2s_1}{d} & -\frac{if_0}{\sqrt{3}}c_3 & \frac{if_0}{\sqrt{3}}c_2\\ -\frac{2s_2}{d} & \omega & -\frac{if_0}{\sqrt{3}}c_1\\ -\frac{2s_3}{d} & \frac{if_0}{\sqrt{3}}c_1 & \omega \end{array}\right| & \nonumber\\ -\frac{4}{3d}\Phi_0s_2\left|\begin{array}{ccc} -\frac{2s_1}{d} & \omega & \frac{if_0}{\sqrt{3}}c_2\\ -\frac{2s_2}{d} & \frac{if_0}{\sqrt{3}}c_3 & -\frac{if_0}{\sqrt{3}}c_1\\ -\frac{2s_3}{d} & -\frac{if_0}{\sqrt{3}}c_2 & \omega \end{array}\right| + \frac{4}{3d}\Phi_0s_3\left|\begin{array}{ccc} -\frac{2s_1}{d} & \omega & -\frac{if_0}{\sqrt{3}}c_3\\ -\frac{2s_2}{d} & \frac{if_0}{\sqrt{3}}c_3 & \omega\\ -\frac{2s_3}{d} & -\frac{if_0}{\sqrt{3}}c_2 & \frac{if_0}{\sqrt{3}}c_1 \end{array} \right| &= 0 \end{align} \]
One has
\[ \begin{align} \frac{4}{3d}\Phi_0s_1\left|\begin{array}{ccc} -\frac{2s_1}{d} & -\frac{if_0}{\sqrt{3}}c_3 & \frac{if_0}{\sqrt{3}}c_2\\ -\frac{2s_2}{d} & \omega & -\frac{if_0}{\sqrt{3}}c_1\\ -\frac{2s_3}{d} & \frac{if_0}{\sqrt{3}}c_1 & \omega \end{array}\right| &= -\frac{8\omega^2}{3d^2}\Phi_0s_1^2 + \frac{8\Phi_0f_0^2}{9d^2}\left(s_1c_1s_3c_3 + s_1c_1s_2c_2 + s_1c_1s_1c_1\right) + \frac{8if_0}{3\sqrt{3}d}s_1\left(s_3c_2 - c_3s_2\right)\nonumber\\ &= -\frac{8\omega^2}{3d^2}\Phi_0s_1^2 + \frac{8\Phi_0f_0^2}{9d^2}s_1c_1\left(s_1c_1 + s_2c_2 + s_3c_3\right) - \frac{8if_0}{3\sqrt{3}d}s_1\left(s_2c_3 - s_3c_2\right). \end{align} \]
The remaining two partial determinants arise analogously by cyclically exchanging the indices:
\[ \begin{align} -\frac{4}{3d}\Phi_0s_2\left|\begin{array}{ccc} -\frac{2s_1}{d} & \omega & \frac{if_0}{\sqrt{3}}c_2\\ -\frac{2s_2}{d} & \frac{if_0}{\sqrt{3}}c_3 & -\frac{if_0}{\sqrt{3}}c_1\\ -\frac{2s_3}{d} & -\frac{if_0}{\sqrt{3}}c_2 & \omega \end{array}\right| &= -\frac{8\omega^2}{3d^2}\Phi_0s_2^2 + \frac{8\Phi_0f_0^2}{9d^2}s_2c_2\left(s_1c_1 + s_2c_2 + s_3c_3\right) - \frac{8if_0}{3\sqrt{3}d}s_2\left(s_3c_1 - s_1c_3\right),\nonumber\\ \frac{4}{3d}\Phi_0s_3\left|\begin{array}{ccc} -\frac{2s_1}{d} & \omega & -\frac{if_0}{\sqrt{3}}c_3\\ -\frac{2s_2}{d} & \frac{if_0}{\sqrt{3}}c_3 & \omega\\ -\frac{2s_3}{d} & -\frac{if_0}{\sqrt{3}}c_2 & \frac{if_0}{\sqrt{3}}c_1 \end{array} \right| &= -\frac{8\omega^2}{3d^2}\Phi_0s_3^2 + \frac{8\Phi_0f_0^2}{9d^2}s_3c_3\left(s_1c_1 + s_2c_2 + s_3c_3\right) - \frac{8if_0}{3\sqrt{3}d}s_3\left(s_1c_2 - s_2c_1\right). \end{align} \]
Thus, the characteristic polynomial finally reads:
For the case of well-resolved waves $\left|\mathbf{k}\right| \ll \frac{1}{d}$, one has
\[ \begin{align} c_j \approx 1, & {} & s_j \approx \frac{k_jd}{2}. \end{align} \]
Substituting this into Eq. (27.96), one obtains
\[ \begin{align} \omega^4 - \omega^2\left[f_0^2 + \frac{8\Phi_0}{3d^2}\left(\frac{k_1^2d^2}{4} + \frac{k_2^2d^2}{4} + \frac{k_3^2d^2}{4}\right)\right] + \frac{8\Phi_0f_0^2}{9d^2}\underbrace{\left(\frac{k_1d}{2} + \frac{k_1d}{2} + \frac{k_1d}{2}\right)^2}_{\approx 0} &\approx 0\nonumber\\ \Leftrightarrow\omega^4 - \omega^2\left[f_0^2 + \frac{2\Phi_0}{3}\left(k_1^2 + k_2^2 + k_3^2\right)\right] &\approx 0.\nonumber\\ \Leftrightarrow\omega^4 - \omega^2\left[f_0^2 + \Phi\mathbf{k}^2\right] &\approx 0. \end{align} \]
The last step took advantage of the fact that $k_1^2 + k_2^2 + k_3^2 = \frac{3}{2}\mathbf{k}^2$. This can be seen by assuming, without loss of generality, that $\mathbf{k} \parallel \mathbf{i}_1$, since then $k_1 = k$, $k_2 = -k\sin\left(30^\circ\right) = -\frac{k}{2}$ and $k_3 = -\frac{k}{2}$. This leads to the two inertia-gravity modes
\[ \begin{align} \omega^2 \approx f_0^2 + \Phi_0\left|\mathbf{k}\right|^2 \end{align} \]
as well as the geostrophic mode
\[ \begin{align} \omega \approx 0. \end{align} \]
However, for less well-resolved waves, the term $\propto\left(s_1c_1 + s_2c_2 + s_3c_3\right)^2$ causes the geostrophic mode to split into two submodes with $\omega \not= 0$. This is problematic and corresponds to the findings of Sect. 26.9.
In a planar triangular or hexagonal C-grid discretization, the velocity components are rotated by 120 or 60$^\circ$ with respect to each other, but the components are not defined at the same location. Therefore, the condition Eq. (27.18) does not hold straightforwardly for the vector components, but only after an interpolation to reference locations, for which the centers of the cells are used here. This interpolation is done using operators
\[ \begin{align} \newtilde{u}^{(i)} \end{align} \]
with $1 \leq i \leq 3$, where $u$ is any vector component. A separate averaging operator is permitted here for each spatial direction $i$. The conditions Eq. (27.18) can be expressed here in the form
\[ \begin{align} \newtilde{u}_1^{(1)} + \newtilde{u}_2^{(2)} + \newtilde{u}_3^{(3)} \hastobe 0, & {} & \newtilde{v}_1^{(1)} + \newtilde{v}_2^{(2)} + \newtilde{v}_3^{3} \hastobe 0\tag{27.102}\label{eq:linear_condition_trivariate_discrete} \end{align} \]
These equations are to be understood as requirements for the vector fields $\mathbf{v}$. If gradient fields $\nabla\alpha$ do not satisfy the condition Eq. (27.102), one has a problem, since the computation of gradient fields on the C-grid is fixed. Therefore, the averaging of Eq. (27.22) is used as a requirement for the averaging operator to be found:
For the Helmholtz decomposition, Eqs. (27.55) - (27.60), averaging in the respectively relevant spatial directions yields
\[ \begin{align} \newtilde{u}_1^{(1)} &= \frac{1}{\sqrt{3}}\left(\newtilde{\delta_3\psi}^{(1)} - \newtilde{\delta_2\psi}^{(1)}\right) + \newtilde{\delta_1\chi}^{(1)},\\ \newtilde{u}_2^{(2)} &= \frac{1}{\sqrt{3}}\left(\newtilde{\delta_1\psi}^{(2)} - \newtilde{\delta_3\psi}^{(2)}\right) + \newtilde{\delta_2\chi}^{(2)},\\ \newtilde{u}_3^{3} &= \frac{1}{\sqrt{3}}\left(\newtilde{\delta_2\psi}^{3} - \newtilde{\delta_1\psi}^{3}\right) + \newtilde{\delta_3\chi}^{3},\\ \newtilde{v}_1^{(1)} &= \frac{1}{\sqrt{3}}\left(\newtilde{\delta_2\chi}^{(1)} - \newtilde{\delta_3\chi}^{(1)}\right) + \newtilde{\delta_1\psi}^{(1)},\\ \newtilde{v}_2^{(2)} &= \frac{1}{\sqrt{3}}\left(\newtilde{\delta_3\chi}^{(2)} - \newtilde{\delta_1\chi}^{(2)}\right) + \newtilde{\delta_2\psi}^{(2)},\\ \newtilde{v}_3^{3} &= \frac{1}{\sqrt{3}}\left(\newtilde{\delta_1\chi}^{3} - \newtilde{\delta_2\chi}^{3}\right) + \newtilde{\delta_3\psi}^{3}. \end{align} \]
This does not yet satisfy Eqs. (27.102) - (27.102); rather, a preparatory averaging of the scalar fields is required:
\[ \begin{align} \newtilde{u}_1^{(1)} &= \frac{1}{\sqrt{3}}\left(\newtilde{\delta_3\newtilde{\psi}^{(2)}}^{(1)} - \newtilde{\delta_2\newtilde{\psi}^{3}}^{(1)}\right) + \newtilde{\delta_1\chi}^{(1)}\\ \newtilde{u}_2^{(2)} &= \frac{1}{\sqrt{3}}\left(\newtilde{\delta_1\newtilde{\psi}^{3}}^{(2)} - \newtilde{\delta_3\newtilde{\psi}^{(1)}}^{(2)}\right) + \newtilde{\delta_2\chi}^{(2)}\\ \newtilde{u}_3^{3} &= \frac{1}{\sqrt{3}}\left(\newtilde{\delta_2\newtilde{\psi}^{(1)}}^{3} - \newtilde{\delta_1\newtilde{\psi}^{(2)}}^{3}\right) + \newtilde{\delta_3\chi}^{3}\\ \newtilde{v}_1^{(1)} &= \frac{1}{\sqrt{3}}\left(\newtilde{\delta_2\newtilde{\chi}^{3}}^{(1)} - \newtilde{\delta_3\newtilde{\chi}^{(2)}}^{(1)}\right) + \newtilde{\delta_1\psi}^{(1)}\\ \newtilde{v}_2^{(2)} &= \frac{1}{\sqrt{3}}\left(\newtilde{\delta_3\newtilde{\chi}^{(1)}}^{(2)} - \newtilde{\delta_1\newtilde{\chi}^{3}}^{(2)}\right) + \newtilde{\delta_2\psi}^{(2)}\\ \newtilde{v}_3^{3} &= \frac{1}{\sqrt{3}}\left(\newtilde{\delta_1\newtilde{\chi}^{(2)}}^{3} - \newtilde{\delta_2\newtilde{\chi}^{(1)}}^{3}\right) + \newtilde{\delta_3\psi}^{3} \end{align} \]
This implies
\[ \begin{align} u_1 &= \frac{1}{\sqrt{3}}\left(\delta_3\newtilde{\psi}^{(2)} - \delta_2\newtilde{\psi}^{3}\right) + \delta_1\chi\tag{27.116}\label{eq:hex:helmholtz_0},\\ u_2 &= \frac{1}{\sqrt{3}}\left(\delta_1\newtilde{\psi}^{3} - \delta_3\newtilde{\psi}^{(1)}\right) + \delta_2\chi\tag{27.117}\label{eq:hex:helmholtz_1},\\ u_3 &= \frac{1}{\sqrt{3}}\left(\delta_2\newtilde{\psi}^{(1)} - \delta_1\newtilde{\psi}^{(2)}\right) + \delta_3\chi\tag{27.118}\label{eq:hex:helmholtz_2},\\ v_1 &= \frac{1}{\sqrt{3}}\left(\delta_2\newtilde{\chi}^{3} - \delta_3\newtilde{\chi}^{(2)}\right) + \delta_1\psi,\\ v_2 &= \frac{1}{\sqrt{3}}\left(\delta_3\newtilde{\chi}^{(1)} - \delta_1\newtilde{\chi}^{3}\right) + \delta_2\psi,\\ v_3 &= \frac{1}{\sqrt{3}}\left(\delta_1\newtilde{\chi}^{(2)} - \delta_2\newtilde{\chi}^{(1)}\right) + \delta_3\psi. \end{align} \]
Now the question arises what the averaging operator $\newtilde{\psi}^{(j)}$ must concretely look like. For the operator introduced in Eq. (27.73), one has
\[ \begin{align} &\newoverline{\delta_1\psi}^{(1)} + \newoverline{\delta_1\psi}^{(2)} + \newoverline{\delta_1\psi}^{3}\nonumber\\ &= \psi\left(\mathbf{r}\right)\frac{2i}{d}\left[\sin\left(\frac{k_1d}{2}\right)\cos\left(\frac{k_1d}{2}\right) + \sin\left(\frac{k_2d}{2}\right)\cos\left(\frac{k_2d}{2}\right) + \sin\left(\frac{k_3d}{2}\right)\cos\left(\frac{k_3d}{2}\right)\right]\nonumber\\ &= \psi\left(\mathbf{r}\right)\frac{i}{d}\left[2\sin\left(\frac{k_1d}{2}\right)\cos\left(\frac{k_1d}{2}\right) + 2\sin\left(\frac{k_2d}{2}\right)\cos\left(\frac{k_2d}{2}\right) + 2\sin\left(\frac{k_3d}{2}\right)\cos\left(\frac{k_3d}{2}\right)\right]\nonumber\\ &= \psi\left(\mathbf{r}\right)\frac{i}{d}\left[\sin\left(k_1d\right) + \sin\left(k_2d\right) + \sin\left(k_3d\right)\right] \stackrel{\text{in general}}{\not=} 0. \end{align} \]
This operator therefore does not satisfy the condition Eq. (27.103). Therefore, a second averaging operator is defined
\[ \begin{align} \newoverline{\psi}^{((1))} & \coloneqq \newoverline{\psi}^{(2, 3)} = \newoverline{\newoverline{\psi}^{3}}^{(2)} = \newoverline{\newoverline{\psi}^{(2)}}^{3} = \psi\left(\mathbf{r}\right)\cos\left(\frac{k_2d}{2}\right)\cos\left(\frac{k_3d}{2}\right).\text{ (and cyclic)} \end{align} \]
For this one has
\[ \begin{align} \newoverline{\delta_1\alpha}^{((1))} + \newoverline{\delta_2\alpha}^{((2))} + \newoverline{\delta_3\alpha}^{(3)} &= \alpha\left(\mathbf{r}\right)\frac{2i}{d}\Big[\sin\left(\frac{k_1d}{2}\right)\cos\left(\frac{k_2d}{2}\right)\cos\left(\frac{k_3d}{2}\right) + \cos\left(\frac{k_1d}{2}\right)\sin\left(\frac{k_2d}{2}\right)\cos\left(\frac{k_3d}{2}\right)\nonumber\\ & + \cos\left(\frac{k_1d}{2}\right)\cos\left(\frac{k_2d}{2}\right)\sin\left(\frac{k_3d}{2}\right)\Big]\nonumber\\ &= \alpha\left(\mathbf{r}\right)\frac{2i}{d}\Big[\sin\left(\frac{k_1d + k_2d}{2}\right)\cos\left(\frac{k_3d}{2}\right) + \cos\left(\frac{k_1d}{2}\right)\cos\left(\frac{k_2d}{2}\right)\sin\left(\frac{k_3d}{2}\right)\Big]\nonumber\\ &= \alpha\left(\mathbf{r}\right)\frac{2i}{d}\Big[\sin\left(\frac{k_1 + k_2 + k_3}{2}d\right) - \cos\left(\frac{k_1 + k_2}{2}d\right)\sin\left(\frac{k_3}{2}d\right)\nonumber\\ &+ \cos\left(\frac{k_1d}{2}\right)\cos\left(\frac{k_2d}{2}\right)\sin\left(\frac{k_3d}{2}\right)\Big]. \end{align} \]
With $k_1 + k_2 + k_3 = 0$ follows
\[ \begin{align} \newoverline{\delta_1\alpha}^{((1))} + \newoverline{\delta_2\alpha}^{((2))} + \newoverline{\delta_3\alpha}^{(3)} &= \alpha\left(\mathbf{r}\right)\frac{2i}{d}\Big[-\cos\left(\frac{k_1 + k_2}{2}d\right)\sin\left(\frac{k_3}{2}d\right) + \cos\left(\frac{k_1d}{2}\right)\cos\left(\frac{k_2d}{2}\right)\sin\left(\frac{k_3d}{2}\right)\Big]\nonumber\\ &= \alpha\left(\mathbf{r}\right)\frac{2i}{d}\sin\left(\frac{k_3}{2}d\right)\Big[\cos\left(\frac{k_1d}{2}\right)\cos\left(\frac{k_2d}{2}\right) - \cos\left(\frac{k_1 + k_2}{2}d\right)\Big]\nonumber\\ &= \alpha\left(\mathbf{r}\right)\frac{2i}{d}\sin\left(\frac{k_3}{2}d\right)\Big[\cos\left(\frac{k_1d}{2}\right)\cos\left(\frac{k_2d}{2}\right) - \cos\left(\frac{k_1d}{2}\right)\cos\left(\frac{k_2d}{2}\right)\nonumber\\ &+ \sin\left(\frac{k_1d}{2}\right)\sin\left(\frac{k_2d}{2}\right)\Big]\nonumber\\ &= \alpha\left(\mathbf{r}\right)\frac{2i}{d}\sin\left(\frac{k_1}{2}d\right)\sin\left(\frac{k_2d}{2}\right)\sin\left(\frac{k_3d}{2}\right) \stackrel{\text{in general}}{\not=} 0.\tag{27.125}\label{eq:thuburn_op_deriv_0} \end{align} \]
This operator also does not satisfy the condition Eq. (27.103). However, one continues the calculation, again using $k_1 + k_2 + k_3 = 0$:
\[ \begin{align} & \newoverline{\delta_1\alpha}^{((1))} + \newoverline{\delta_2\alpha}^{((2))} + \newoverline{\delta_3\alpha}^{(3)} = \alpha\left(\mathbf{r}\right)\frac{2i}{d}\sin\left(\frac{k_1}{2}d\right)\sin\left(\frac{k_2d}{2}\right)\sin\left(\frac{k_3d}{2}\right)\nonumber\\ &= \alpha\left(\mathbf{r}\right)\frac{2i}{d}\sin\left(\frac{k_1}{2}d\right)\sin\left(\frac{k_2d}{2}\right)\sin\left(-\frac{k_1 + k_2}{2}d\right) = -\alpha\left(\mathbf{r}\right)\frac{2i}{d}\sin\left(\frac{k_1}{2}d\right)\sin\left(\frac{k_2d}{2}\right)\sin\left(\frac{k_1 + k_2}{2}d\right)\nonumber\\ &= -\alpha\left(\mathbf{r}\right)\frac{2i}{d}\sin\left(\frac{k_1}{2}d\right)\sin\left(\frac{k_2d}{2}\right)\left[\sin\left(\frac{k_1}{2}d\right)\cos\left(\frac{k_2}{2}d\right) + \cos\left(\frac{k_1}{2}d\right)\sin\left(\frac{k_2}{2}d\right)\right]\nonumber\\ &= -\alpha\left(\mathbf{r}\right)\frac{2i}{d}\left[\sin\left(\frac{k_2}{2}d\right)\sin\left(\frac{k_1}{2}d\right)^2\cos\left(\frac{k_2}{2}d\right) + \sin\left(\frac{k_1}{2}d\right)\cos\left(\frac{k_1}{2}d\right)\sin\left(\frac{k_2}{2}d\right)^2\right]\nonumber\\ &= -\alpha\left(\mathbf{r}\right)\frac{i}{d}\left[\sin\left(k_2d\right)\sin\left(\frac{k_1}{2}d\right)^2 + \sin\left(k_1d\right)\sin\left(\frac{k_2}{2}d\right)^2\right]\nonumber\\ &= -\alpha\left(\mathbf{r}\right)\frac{i}{d}\left[\sin\left(k_2d\right) + \sin\left(k_1d\right) - \sin\left(k_2d\right)\cos\left(\frac{k_1}{2}d\right)^2 - \sin\left(k_1d\right)\cos\left(\frac{k_2}{2}d\right)^2\right]\nonumber \end{align} \] \[ \begin{align} &= -\alpha\left(\mathbf{r}\right)\frac{i}{d}\Big[\sin\left(k_2d\right) + \sin\left(k_1d\right) - \sin\left(k_2d\right)\cos\left(k_1\right) - \sin\left(k_1d\right)\cos\left(k_2d\right) + \sin\left(k_2d\right)\sin\left(\frac{k_1}{2}d\right)^2\nonumber\\ & + \sin\left(k_1d\right)\sin\left(\frac{k_2}{2}d\right)^2\Big]\nonumber\\ &= -\alpha\left(\mathbf{r}\right)\frac{i}{d}\Big[\sin\left(k_2d\right) + \sin\left(k_1d\right) - \sin\left[\left(k_1 + k_2\right)d\right] + \sin\left(k_2d\right)\sin\left(\frac{k_1}{2}d\right)^2 + \sin\left(k_1d\right)\sin\left(\frac{k_2}{2}d\right)^2\Big]\nonumber\\ &= -\alpha\left(\mathbf{r}\right)\frac{i}{d}\Big[\sin\left(k_2d\right) + \sin\left(k_1d\right) - \sin\left(-k_3d\right) + \sin\left(k_2d\right)\sin\left(\frac{k_1}{2}d\right)^2 + \sin\left(k_1d\right)\sin\left(\frac{k_2}{2}d\right)^2\Big]\nonumber\\ &= -\alpha\left(\mathbf{r}\right)\frac{i}{d}\Big[\sin\left(k_2d\right) + \sin\left(k_1d\right) + \sin\left(k_3d\right) + \sin\left(k_2d\right)\sin\left(\frac{k_1}{2}d\right)^2 + \sin\left(k_1d\right)\sin\left(\frac{k_2}{2}d\right)^2\Big]\nonumber\\ &= -\left(\newoverline{\delta_1\alpha}^{(1)} + \newoverline{\delta_1\alpha}^{(2)} + \newoverline{\delta_1\alpha}^{3}\right) - \alpha\left(\mathbf{r}\right)\frac{i}{d}\Big[\sin\left(k_2d\right)\sin\left(\frac{k_1}{2}d\right)^2 + \sin\left(k_1d\right)\sin\left(\frac{k_2}{2}d\right)^2\Big] \end{align} \]
One now computes
\[ \begin{align} &\sin\left(k_2d\right)\sin\left(\frac{k_1}{2}d\right)^2 + \sin\left(k_1d\right)\sin\left(\frac{k_2}{2}d\right)^2 = \sin\left(\frac{k_2}{2}d + \frac{k_2}{2}d\right)\sin\left(\frac{k_1}{2}d\right)^2 + \sin\left(\frac{k_1}{2}d + \frac{k_1}{2}d\right)\sin\left(\frac{k_2}{2}d\right)^2\nonumber\\ &= 2\sin\left(\frac{k_1}{2}d\right)^2\sin\left(\frac{k_2}{2}d\right)\cos\left(\frac{k_2}{2}d\right) + 2\sin\left(\frac{k_2}{2}d\right)^2\sin\left(\frac{k_1}{2}d\right)\cos\left(\frac{k_1}{2}d\right)\nonumber\\ &= \sin\left(\frac{k_1}{2}d\right)\sin\left(\frac{k_2}{2}d\right)\left[2\sin\left(\frac{k_1}{2}d\right)\cos\left(\frac{k_2}{2}d\right) + 2\cos\left(\frac{k_1}{2}d\right)\sin\left(\frac{k_2}{2}d\right)\right]\nonumber\\ &= \sin\left(\frac{k_1}{2}d\right)\sin\left(\frac{k_2}{2}d\right)\sin\left(\frac{k_1 + k_2}{2}d\right) = \sin\left(\frac{k_1}{2}d\right)\sin\left(\frac{k_2}{2}d\right)\sin\left(-\frac{k_3}{2}d\right)\nonumber\\ &= -\sin\left(\frac{k_1}{2}d\right)\sin\left(\frac{k_2}{2}d\right)\sin\left(\frac{k_3}{2}d\right). \end{align} \]
Thus follows
\[ \begin{align} & \newoverline{\delta_1\alpha}^{((1))} + \newoverline{\delta_2\alpha}^{((2))} + \newoverline{\delta_3\alpha}^{(3)} = -\left(\newoverline{\delta_1\alpha}^{(1)} + \newoverline{\delta_1\alpha}^{(2)} + \newoverline{\delta_1\alpha}^{3}\right) + \alpha\left(\mathbf{r}\right)\frac{i}{d}\sin\left(\frac{k_1}{2}d\right)\sin\left(\frac{k_2}{2}d\right)\sin\left(\frac{k_3}{2}d\right)\nonumber\\ &\stackrel{\href{#eq:thuburn_op_deriv_0}{\text{Eq. (27.125)}}}{=} -\left(\newoverline{\delta_1\alpha}^{(1)} + \newoverline{\delta_1\alpha}^{(2)} + \newoverline{\delta_1\alpha}^{3}\right) + \frac{1}{2}\left(\newoverline{\delta_1\alpha}^{((1))} + \newoverline{\delta_2\alpha}^{((2))} + \newoverline{\delta_3\alpha}^{(3)}\right)\nonumber\\ &\Rightarrow\newoverline{\delta_1\alpha}^{((1))} + \newoverline{\delta_2\alpha}^{((2))} + \newoverline{\delta_3\alpha}^{(3)} = -\frac{1}{2}\left(\newoverline{\delta_1\alpha}^{(1)} + \newoverline{\delta_1\alpha}^{(2)} + \newoverline{\delta_1\alpha}^{3}\right).\tag{27.128}\label{eq:thuburn_op_deriv_1} \end{align} \]
One now makes the ansatz
\[ \begin{align} \newtilde{\alpha}^{(j)} = \beta\newoverline{\alpha}^{(j)} + \left(1 - \beta\right)\newoverline{\alpha}^{((j))} \end{align} \]
with $0 \leq \beta \leq 0$. Substituting Eq. (27.128) here, one obtains
\[ \begin{align} \newtilde{\delta_1\alpha}^{(1)} + \newtilde{\delta_2\alpha}^{(2)} + \newtilde{\delta_3\alpha}^{3} = \left(\beta - \frac{1}{2}\left(1 - \beta\right)\right)\left(\newoverline{\delta_1\alpha}^{(1)} + \newoverline{\delta_1\alpha}^{(2)} + \newoverline{\delta_1\alpha}^{3}\right) \stackrel{!}{=} 0. \end{align} \]
This is fulfilled for
\[ \begin{align} \beta - \frac{1}{2}\left(1 - \beta\right) &= 0 \Leftrightarrow \beta - \frac{1}{2} + \frac{\beta}{2} = \frac{3\beta}{2} - \frac{1}{2} = 0\nonumber\\ \Leftrightarrow \beta = \frac{1}{3}. \end{align} \]
Thus, one has
This operator is called the Thuburn operator because it was first used by Thuburn in [38].
Applying the Thuburn operator in the linearized shallow water equations on the hexagonal grid, Eqs. (27.81) - (27.84), leads to
\[ \begin{align} \frac{\partial\Phi}{\partial t} + \frac{2}{3}\Phi_0\left(\delta_1u_1 + \delta_2u_2 + \delta_3u_3\right) &= 0,\\ \frac{\partial u_1}{\partial t} - \frac{f_0}{\sqrt{3}}\left(\newtilde{u}_2^{(3)} - \newtilde{u}_3^{(2)}\right) + \delta_1\Phi &= 0,\tag{27.134}\label{eq:x_1_momentum_thuburn_c-grid}\\ \frac{\partial u_2}{\partial t} - \frac{f_0}{\sqrt{3}}\left(\newtilde{u}_3^{(1)} - \newtilde{u}_1^{(3)}\right) + \delta_2\Phi &= 0,\\ \frac{\partial u_3}{\partial t} - \frac{f_0}{\sqrt{3}}\left(\newtilde{u}_1^{(2)} - \newtilde{u}_2^{(1)}\right) + \delta_3\Phi &= 0. \end{align} \]
From this it follows immediately
\[ \begin{align} \frac{\partial}{\partial t}\left(\newtilde{u}_1^{(1)} + \newtilde{u}_2^{(2)} + \newtilde{u}_3^{(3)}\right) = 0.\tag{27.137}\label{eq:thuburn_con_swe_c-grid} \end{align} \]
The dispersion relation of this system of equations can be calculated analogously to the derivation of Eq. (27.96) by performing the substitution
\[ \begin{align} c_1 \to a_1 \coloneqq\frac{c_1 + 2c_2c_3}{3}\text{ (and cyclic)} \end{align} \]
This leads to
\[ \begin{align} \omega^4 - \omega^2\left[\frac{f_0^2}{3}\left(a_1^2 + a_2^2 + a_3^2\right) + \frac{8\Phi_0}{3d^2}\left(s_1^2 + s_2^2 + s_3^2\right)\right] + \frac{8\Phi_0f_0^2}{9d^2}\left(s_1a_1 + s_2a_2 + s_3a_3\right)^2 = 0.\tag{27.139}\label{eq:disp_rel_hex_c_mod_0} \end{align} \]
If one sets $\alpha = \exp\left(i\mathbf{k}\cdot\mathbf{r}\right)$ in Eq. (27.103), one obtains
\[ \begin{align} \frac{2i}{d}\left(a_1s_1 + a_2s_2 + a_3s_3\right)\alpha = 0 \Leftrightarrow a_1s_1 + a_2s_2 + a_3s_3 = 0. \end{align} \]
Thus, Eq. (27.139) becomes
\[ \begin{align} \omega^4 - \omega^2\left[\frac{f_0^2}{3}\left(a_1^2 + a_2^2 + a_3^2\right) + \frac{8\Phi_0}{3d^2}\left(s_1^2 + s_2^2 + s_3^2\right)\right] = 0.\tag{27.141}\label{eq:thuburn_leads_to_geostrophic_mode} \end{align} \]
If one uses Eq. (27.132) to reconstruct the tangential velocity components, a stationary geostrophic mode thus arises.
On a regular hexagonal lattice with lattice constant $d$, a triangle has area
\[ \begin{align} A_v = \frac{d}{2}l. \end{align} \]
Here one has
\[ \begin{align} l = \cos\left(30^\circ\right)d = \frac{\sqrt{3}}{2}d, \end{align} \]
from which
\[ \begin{align} A_v = \frac{d^2\sqrt{3}}{4}\tag{27.144}\label{eq:hex_face_triangle} \end{align} \]
follows. Thus, on a triangle (superscript $t$) pointing upward (index $u$), one can compute the divergence according to Gauss's theorem in the form
\[ \begin{align} D^{(t)}_u = -\frac{dv_1 + dv_2 + dv_3}{A_v} = -\frac{4}{\sqrt{3}d}\left(v_1 + v_2 + v_3\right) \end{align} \]
The minus sign arises from the fact that in such a triangle the coordinate axes of the $j$-system point inward. For a triangle pointing downward (index $l$), it is the other way around; one has
\[ \begin{align} D^{(t)}_l = \frac{dv_1 + dv_2 + dv_3}{A_v} = \frac{4}{\sqrt{3}d}\left(v_1 + v_2 + v_3\right). \end{align} \]
This can be summarized in the form
\[ \begin{align} D^{(t)}_{u, l} = \mp\frac{dv_1 + dv_2 + dv_3}{A_v} = \mp\frac{4}{\sqrt{3}d}\left(v_1 + v_2 + v_3\right) = \mp\frac{4}{\sqrt{3}d}\sum_{e \in t}v_e \end{align} \]
where $e$ denotes the edges of the triangle. Analogously, for the relative vorticity, Stokes' theorem gives
\[ \begin{align} \zeta^{(t)}_{u, l} = \pm\frac{du_1 + du_2 + du_3}{A_v} = \pm\frac{4}{\sqrt{3}d}\sum_{e \in t}u_e. \end{align} \]
Let $h$ now denote a hexagon. Then $u \in h$ denote the three upward-pointing triangles that overlap with the hexagon $h$ under consideration and similarly $l \in h$ denote the three downward-pointing triangles that overlap with $h$. With these notations one can write
\[ \begin{align} \frac{\sqrt{3}d}{4}\left(\sum_{l \in h}D_l^{(t)} - \sum_{u \in h}D_u^{(t)}\right) &= \sum_{l \in h}\sum_{e \in l}v_e + \sum_{u \in h}\sum_{e \in u}v_e = \sum_{l \in h}\left(\sum_{e \in l, e \in h}v_e + \sum_{e \in l, e \not\in h}v_e\right) + \sum_{u \in h}\left(\sum_{e \in u, e \in h}v_e + \sum_{e \in u, e \not\in h}v_e\right) \end{align} \]
In the last transformation step, the edges of the triangles were divided into those that are also edges of the hexagon under consideration ($e \in h$) and those that are not ($e \not\in h$). Factoring out yields
\[ \begin{align} \frac{\sqrt{3}d}{4}\left(\sum_{l \in h}D_l^{(t)} - \sum_{u \in h}D_u^{(t)}\right) &= \sum_{l \in h}\sum_{e \in l, e \in h}v_e + \sum_{l \in h}\sum_{e \in l, e \not\in h}v_e + \sum_{u \in h}\sum_{e \in u, e \in h}v_e + \sum_{u \in h}\sum_{e \in u, e \not\in h}v_e. \end{align} \]
As one can see, all velocity components contribute with a positive sign to $\sum_{l \in h}D_l^{(t)} - \sum_{u \in h}D_u^{(t)}$. In the sums $\sum_{l \in h}\sum_{e \in l, e \in h}$ and $\sum_{u \in h}\sum_{e \in u, e \in h}$, each edge of the hexagon $h$ occurs exactly once. Thus one has
\[ \begin{align} \frac{\sqrt{3}d}{4}\left(\sum_{l \in h}D_l^{(t)} - \sum_{u \in h}D_u^{(t)}\right) &= 2\sum_{e \in h}v_e + \sum_{l \in h}\sum_{e \in l, e \not\in h}v_e + \sum_{u \in h}\sum_{e \in u, e \not\in h}v_e. \end{align} \]
Those edges that do not border the hexagon under consideration each contribute exactly once to $\sum_{l \in h}D_l^{(t)} - \sum_{u \in h}D_u^{(t)}$:
\[ \begin{align} \frac{\sqrt{3}d}{4}\left(\sum_{l \in h}D_l^{(t)} - \sum_{u \in h}D_u^{(t)}\right) &= 2\sum_{e \in h}v_e + \sum_{t \in h}\sum_{e \in t, e \not\in h}v_e \end{align} \]
At this point, consider the quantity $\newtilde{v}_1^{(1)}$, which is defined at the center of the hexagon under consideration. It is a superposition of four velocity components, two of which lie on the $\mathbf{i}_1$-axis and two perpendicular to it. The former contribute to $\newtilde{v}_1^{(1)}$ with weight $1/3$, the latter with weight $1/6.$ This can also be transferred cyclically to the other spatial directions. Thus one has
\[ \begin{align} \frac{\sqrt{3}d}{4}\left(\sum_{l \in h}D_l^{(t)} - \sum_{u \in h}D_u^{(t)}\right) &= 6\left(\newtilde{v}_1^{(1)} + \newtilde{v}_2^{(2)} + \newtilde{v}_3^{3}\right).\tag{27.153}\label{eq:checkerboard_triangular} \end{align} \]
If the condition $\newtilde{v}_1^{(1)} + \newtilde{v}_2^{(2)} + \newtilde{v}_3^{3} = 0$ is not met, then the mean of the divergence of the three upward-pointing triangles is in general not equal to that of the downward-pointing triangles. This leads to the creation of a so-called checkerboard pattern in the divergence. This noise also propagates from there into other fields. It is therefore important to always adhere to the condition $\newtilde{v}_1^{(1)} + \newtilde{v}_2^{(2)} + \newtilde{v}_3^{3} = 0$.
On the hexagonal grid, the reasoning can be transferred completely analogously to the relative vorticity, which leads to
\[ \begin{align} \frac{\sqrt{3}d}{4}\left(\sum_{u \in h}\zeta_u^{(t)} - \sum_{l \in h}\zeta_l^{(t)}\right) &= 6\left(\newtilde{u}_1^{(1)} + \newtilde{u}_2^{(2)} + \newtilde{u}_3^{3}\right) \end{align} \]
Failure to comply with the condition $\newtilde{u}_1^{(1)} + \newtilde{u}_2^{(2)} + \newtilde{u}_3^{3} = 0$ also leads to a checkerboard pattern on the hexagonal grid (but here in the vorticity) and must therefore also be taken into account here.
In the previous section it was shown that a checkerboard pattern exists on the hexagonal grid in the vorticity calculated on triangles if the condition $\newtilde{u}_1^{(1)} + \newtilde{u}_2^{(2)} + \newtilde{u}_3^{3} = 0$ for the wind field $\mathbf{v}$ is not met. However, in order for this condition for the wind field to be fulfilled for all time steps, this must apply to the wind field in the initial state as well as to all forcings $\mathbf{w}$. However, simulations are usually carried out with momentum diffusion to avoid a build-up of kinetic energy on the grid scale. These diffusive terms break down gradients on the grid scale and can therefore also break down an existing checkerboard pattern as long as it is not too strong. Accordingly, the requirement $\newtilde{u}_1^{(1)} + \newtilde{u}_2^{(2)} + \newtilde{u}_3^{3} = 0$ for the initial wind field $\mathbf{v}$ can be dispensed with. However, the following is important:
For all velocity tendencies $\mathbf{w}$, $\newtilde{w}_1^{(1)} + \newtilde{w}_2^{(2)} + \newtilde{w}_3^{3} = 0$ must hold.
The analogous facts apply on the triangular grid.
If one applies the Laplace operator to Eqs. (27.116) - (27.118), one obtains
\[ \begin{align} \Delta u_1 &= \frac{1}{\sqrt{3}}\left(\delta_3\Delta\newtilde{\psi}^{(2)} - \delta_2\Delta \newtilde{\psi}^{3}\right) + \delta_1\Delta\chi,\tag{27.155}\label{eq:rhombus_curl_deriv_0}\\ \Delta u_2 &= \frac{1}{\sqrt{3}}\left(\delta_1\Delta\newtilde{\psi}^{3} - \delta_3\Delta \newtilde{\psi}^{(1)}\right) + \delta_2\Delta\chi,\tag{27.156}\label{eq:rhombus_curl_deriv_1}\\ \Delta u_3 &= \frac{1}{\sqrt{3}}\left(\delta_2\Delta\newtilde{\psi}^{(1)} - \delta_1\Delta \newtilde{\psi}^{(2)}\right) + \delta_3\Delta\chi.\tag{27.157}\label{eq:rhombus_curl_deriv_2} \end{align} \]
Since discretized partial derivatives and averaging operators commute, it follows immediately
\[ \begin{align} \newtilde{\Delta u_1}^{(1)} + \newtilde{\Delta u_2}^{(2)} + \newtilde{\Delta u_3}^{(3)} = 0. \end{align} \]
If one discretizes Eqs. (27.62) - (27.64), one obtains
\[ \begin{align} \Delta u_1 &= \frac{1}{\sqrt{3}}\left(\delta_3\zeta - \delta_2\zeta\right) + \delta_1 D,\tag{27.159}\label{eq:rhombus_curl_deriv_3}\\ \Delta u_2 &= \frac{1}{\sqrt{3}}\left(\delta_1\zeta - \delta_3\zeta\right) + \delta_2 D,\tag{27.160}\label{eq:rhombus_curl_deriv_4}\\ \Delta u_3 &= \frac{1}{\sqrt{3}}\left(\delta_2\zeta - \delta_1\zeta\right) + \delta_3 D.\tag{27.161}\label{eq:rhombus_curl_deriv_5} \end{align} \]
By comparison, one finds that there are three types of vorticities in the $\mathbf{i}$ system:
\[ \begin{align} \zeta_1 \coloneqq \Delta\newtilde{\psi}^{(1)}, & {} & \zeta_2 \coloneqq \Delta\newtilde{\psi}^{(2)}, & {} & \zeta_3 \coloneqq \Delta\newtilde{\psi}^{3} \end{align} \]
It therefore holds that
\[ \begin{align} \Delta u_1 &= \frac{1}{\sqrt{3}}\left(\delta_3\zeta_2 - \delta_2\zeta_3\right) + \delta_1 D,\tag{27.163}\label{eq:rhombus_curl_deriv_6}\\ \Delta u_2 &= \frac{1}{\sqrt{3}}\left(\delta_1\zeta_3 - \delta_3\zeta_1\right) + \delta_2 D,\tag{27.164}\label{eq:rhombus_curl_deriv_7}\\ \Delta u_3 &= \frac{1}{\sqrt{3}}\left(\delta_2\zeta_1 - \delta_1\zeta_2\right) + \delta_3 D.\tag{27.165}\label{eq:rhombus_curl_deriv_8} \end{align} \]
For the divergence, according to Eq. (27.33), one has
\[ \begin{align} D &= \frac{2}{3}\left(\delta_1u_1 + \delta_2u_2 + \delta_3u_3\right). \end{align} \]
Substituting this into Eqs. (27.163) - (27.165), one obtains
\[ \begin{align} \Delta u_1 &= \frac{1}{\sqrt{3}}\left(\delta_3\zeta_2 - \delta_2\zeta_3\right) + \frac{2}{3}\left(\delta_1\delta_1u_1 + \delta_1\delta_2u_2 + \delta_1\delta_3u_3\right),\tag{27.167}\label{eq:hex_laplace_vec_0}\\ \Delta u_2 &= \frac{1}{\sqrt{3}}\left(\delta_1\zeta_3 - \delta_3\zeta_1\right) + \frac{2}{3}\left(\delta_2\delta_1u_1 + \delta_2\delta_2u_2 + \delta_2\delta_3u_3\right),\tag{27.168}\label{eq:hex_laplace_vec_1}\\ \Delta u_3 &= \frac{1}{\sqrt{3}}\left(\delta_2\zeta_1 - \delta_1\zeta_2\right) + \frac{2}{3}\left(\delta_3\delta_1u_1 + \delta_3\delta_2u_2 + \delta_3\delta_3u_3\right).\tag{27.169}\label{eq:hex_laplace_vec_2} \end{align} \]
For the vorticities $\zeta_i$, motivated by Eqs. (27.49) - (27.51), one makes the ansätze
\[ \begin{align} \zeta_1 &= \frac{2}{\sqrt{3}}\left(\delta_2u_3 - \delta_3u_2\right),\tag{27.170}\label{eq:hex_curl_rhombus_0}\\ \zeta_2 &= \frac{2}{\sqrt{3}}\left(\delta_3u_1 - \delta_1u_3\right),\tag{27.171}\label{eq:hex_curl_rhombus_1}\\ \zeta_3 &= \frac{2}{\sqrt{3}}\left(\delta_1u_2 - \delta_2u_1\right).\tag{27.172}\label{eq:hex_curl_rhombus_2} \end{align} \]
Substituting Eq. (27.170) into Eq. (27.167), one obtains
\[ \begin{align} \Delta u_1 &= \frac{2}{3}\left(\delta_3\delta_3u_1 - \delta_3\delta_1u_3 - \delta_2\delta_1u_2 + \delta_2\delta_2u_1\right) + \frac{2}{3}\left(\delta_1\delta_1u_1 + \delta_1\delta_2u_2 + \delta_1\delta_3u_3\right) = \frac{2}{3}\left(\delta_1\delta_1u_1 + \delta_2\delta_2u_1 + \delta_3\delta_3u_1\right). \end{align} \]
By cyclic permutations, the analogous statement is obtained for the other components as well. The ansätze of Eqs. (27.170) - (27.172) are therefore justified. Two neighboring triangles together form a rhombus (a parallelogram); by Eq. (27.144), such a rhombus has the area
\[ \begin{align} A_r = 2A_v = 2\frac{d^2\sqrt{3}}{4} = \frac{d^2\sqrt{3}}{2}. \end{align} \]
Thus, $\zeta_i$ corresponds to the vorticity calculated using Stokes' theorem over a rhombus whose short axis of symmetry is parallel to the $i$-axis. On the hexagonal grid, the rhombuses form the natural sets on which Stokes' theorem must be evaluated. This is a deviation from the otherwise intuitive assumption that the vorticity must be calculated via the dual grid (the triangular grid).
As usual, divergence and vorticity appear in reversed roles on the triangular grid. There, the vorticity has the unique value
\[ \begin{align} \zeta = \frac{2}{3}\left(\delta_1v_1 + \delta_2v_2 + \delta_3v_3\right) \end{align} \]
However, three types of divergences occur there, each of which is calculated using Gauss's theorem over the respective rhombuses:
\[ \begin{align} D_1 &= \frac{2}{\sqrt{3}}\left(\delta_3v_2 - \delta_2v_3\right)\\ D_2 &= \frac{2}{\sqrt{3}}\left(\delta_1v_3 - \delta_3v_1\right)\\ D_3 &= \frac{2}{\sqrt{3}}\left(\delta_2v_1 - \delta_1v_2\right) \end{align} \]
This is also underscored by the discretization of the vector Laplace operator according to Eq. (B.54), which reads
\[ \begin{align} \Delta\mathbf{v} = \nabla\left(\nabla\cdot\mathbf{v}\right) - \nabla\times\left(\nabla\times\mathbf{v}\right) \end{align} \]
Applying this to the vector component $u_1$ gives
\[ \begin{align} \Delta u_1 = \delta_1D_h - \delta_{i, \perp}\zeta_t. \end{align} \]
Here $D_h$ is the divergence on hexagons and $\zeta_t$ is the vorticity on triangles. $\delta_{i, \perp}$ denotes the central difference quotient in the $\mathbf{j}_1-$direction. For this one has
\[ \begin{align} -\delta_{i, \perp}\zeta_t &= \frac{\sqrt{3}}{d}\left(\zeta_{t, l} - \zeta_{t, u}\right). \end{align} \]
From Eq. (27.163) follows
\[ \begin{align} \Delta u_1 &= \delta_1D_h + \frac{\sqrt{3}}{d}\left(\zeta_{t, l} - \zeta_{t, u}\right) = \delta_1 D_h + \frac{1}{\sqrt{3}}\left(\delta_3\zeta_2 - \delta_2\zeta_3\right)\nonumber\\ & \Leftrightarrow \sqrt{3}\frac{\zeta_{t, l} - \zeta_{t, u}}{d} = \frac{\sqrt{3}}{3}\left(\delta_3\zeta_2 - \delta_2\zeta_3\right) \Leftrightarrow \sqrt{3}\frac{\zeta_{t, l} - \zeta_{t, u}}{d} = \frac{\sqrt{3}}{3}\frac{\zeta_{2, l} - \zeta_{2, u} - \zeta_{3, u} + \zeta_{3, l}}{d}\nonumber\\ & \Leftrightarrow \sqrt{3}\frac{\zeta_{t, l} - \zeta_{t, u}}{d} = \frac{\sqrt{3}}{3}\frac{\zeta_{2, l} + \zeta_{3, l} - \left(\zeta_{2, u} - \zeta_{3, u}\right)}{d}\nonumber\\ & \Leftrightarrow \sqrt{3}\frac{\zeta_{t, l} - \zeta_{t, u}}{d} = \frac{\sqrt{3}}{3}\frac{\zeta_{1, l} + \zeta_{2, l} + \zeta_{3, l} - \left(\zeta_{1, u} + \zeta_{2, u} - \zeta_{3, u}\right)}{d}. \end{align} \]
The last equivalence transformation follows from the fact $\zeta_{1, l} = \zeta_{1, u}$. This is satisfied by the definitions
\[ \begin{align} \zeta_{t, u} &\coloneqq \frac{\zeta_{1, u} + \zeta_{2, u} + \zeta_{3, u}}{3},\\ \zeta_{t, l} &\coloneqq \frac{\zeta_{1, l} + \zeta_{2, l} + \zeta_{3, l}}{3} \end{align} \]
The rotation on a triangle can therefore be calculated as the average of the three rhombuses that overlap with this triangle.
On the hexagonal C-grid, gradient fields $\mathbf{u} = \nabla\psi$ satisfy the condition
\[ \begin{align} \newtilde{u}_1^{(1)} + \newtilde{u}_2^{(2)} + \newtilde{u}_3^{3} = 0. \end{align} \]
On the triangular C-grid, a gradient field $\mathbf{v} = \nabla\psi$ naturally has components parallel to the basis elements $\mathbf{j}_i$. This does not satisfy
\[ \begin{align} \newtilde{v}_1^{(1)} + \newtilde{v}_2^{(2)} + \newtilde{v}_3^{3} = 0 \end{align} \]
as one can easily verify. Therefore, according to Eq. (27.153), gradient fields on the triangular C-grid produce a checkerboard pattern in the divergence, which would require unphysical numerical stabilizers such as divergence damping. This grid was therefore rightly excluded in Sect. 26.9.1. This finding was first put forward by Gassmann in [24].