User:Egm6936.f09/Vorticity equation

Vorticity equation for an incompressible flow is obtained when the curl of the Navier-Stokes equation, Eq.(15), is taken. The vorticity equation is given as


 * {| style="width:100%" border="0"

$$  \displaystyle \frac{D \boldsymbol \omega}{Dt} = \nabla \mathbf u \cdot \boldsymbol \omega + \nu \nabla^2 \boldsymbol \omega $$     (1a)
 * style="width:95%" |
 * style="width:95%" |
 * 
 * }

In component form,


 * {| style="width:100%" border="0"

$$  \displaystyle
 * style="width:95%" |
 * style="width:95%" |

\frac{D \omega_i}{Dt} = \frac{\partial u_i}{\partial x_j} \omega_j + \nu \nabla^2 \omega_i $$     (1b)
 * 
 * }

where the component form of the gradient of the velocity is given in Eq.(21), which might be new for especially fluid mechanics people. To get used to two different conventions, the reader is referred to the article Gradient of vector: Two tensor conventions. In Eq.(1b), $$ \displaystyle \nu $$ is the kinematic viscosity. Vorticity equation is used to have a more convenient explanation of the physical phenomena. The left hand side represents the change of the vorticity and the second term on the right hand side is the diffusion of the vorticity. The first term on the right hand side represents the change of vorticity due to the vortex-line stretching (see Meaning of vortex-line stretching). Vortex lines are the axis of rotation and they are tangent to the local vorticity vector. Vorticity vector reflects the rotation of the fluid particles or the circulation and is defined as


 * {| style="width:100%" border="0"

$$  \displaystyle \boldsymbol \omega := \nabla \times \mathbf u $$ (1c)
 * style="width:95%" |
 * style="width:95%" |
 * 
 * }

The vorticity video (Experiments were carried out by A. H. Shapiro) from National Committee for Fluid Mechanics Films (CFM), 1969, which can be found in the open-access collection of education videos and film notes by the NCFM is very instructive to understand the meaning of the terms vorticity, rotation, and circulation. See also Schwarzchild 2010 on vortex lines in turbulent superfluid.

Examples 2 and 3 in Tritton (1988, p.82) are highly instructive and explain the distinction of the rotation, which is changing orientation in space and closed path motion (this type of flow is called irrotational). Example 2 illustrates solid body rotation, an example of closed path motion. Here, the streamlines are curved. However, the flow is irrotational and thus does not have vorticity. Example 3 is laminar flow in a channel with straight vertical walls. The only component of the velocity is in the axial direction, $$\displaystyle v_x $$, which vanishes on the wall at $$\displaystyle  y=0 $$ due to no-slip condition and this axial velocity depends solely on the lateral position, i.e. $$\displaystyle  v_x = v_x(y) $$. Each fluid element on one lateral position stays on the same lateral position along the axis of the flow. Therefore, the streamlines are straight. However, as the magnitude of the velocity is larger at some lateral position compared to another position, e.g. $$\displaystyle v_x(y=y_1)>v_x(y=y_2) $$ where $$\displaystyle  y_1>y_2 $$, the fluid particle rotates while flowing down the channel and the flow is rotational. Therefore, the flow does not have to be curved to have vorticity (see Fig. 6.8, p.83, Tritton, 1988 ). Vorticity is used to better understand the characteristics of some flows such as turbulence. Vorticity or an irregular distribution of the vortices is one of the main characteristics of turbulence as explained in the Film Notes for Turbulence by R.W. Stewart 1969.

Now, the aim is to derive the vorticity equation and then explain the meaning of the vortex-line stretching.



Derivation of the vorticity equation
This derivation does not use any vector identities. It is built on the first principles. In other words, it is implicitly deriving the vector identities used in the derivation Using vector identities. Therefore, it is longer than that section.

We will first derive the Navier-Stokes (NS) equation from the balance of the linear momentum equation, which is


 * {| style="width:100%" border="0"

$$  \displaystyle \rho \frac {D \mathbf u}  {Dt} =  {\rm div} \boldsymbol \tau +  \mathbf b   \Longrightarrow \rho \frac {D u_j} {Dt} \mathbf e_j =  \frac {\partial \tau_{ij}} {\partial x_i} \mathbf e_j +  b_j \mathbf e_j $$     (1)
 * style="width:95%" |
 * style="width:95%" |
 * 
 * }

Here, $$\displaystyle \rho $$ is the density, $$\displaystyle \mathbf u $$ the velocity field, $$\displaystyle \boldsymbol \tau $$ the stress tensor, $$\displaystyle \mathbf b $$ body force, For constant property Newtonian fluid


 * {| style="width:100%" border="0"

$$  \displaystyle \boldsymbol \tau =  - p   \mathbf I   + 2 \mu \boldsymbol \gamma (\mathbf u) $$ ||  (2)
 * style="width:95%" |
 * style="width:95%" |
 * }

Here, $$\displaystyle \mathbf I $$ is the identity tensor, $$\displaystyle \mu $$ the dynamic viscosity, and $$\displaystyle \boldsymbol \gamma $$ the strain rate defined as

$$\displaystyle \boldsymbol \gamma := \frac{1}{2} \left[ \frac{\partial u_i} {\partial x_j} + \frac{\partial u_j} {\partial x_i} \right ] $$

Note about the minus sign in front of the pressure term. $$\displaystyle p>0 $$, since stress is compressive and thus negative. The stress tensor in component form can be written as


 * {| style="width:100%" border="0"

$$  \displaystyle \tau_{ij} =     -p \delta_{ij} +      \mu \left[ \frac{\partial u_i} {\partial x_j} + \frac{\partial u_j} {\partial x_i} \right ] $$     (3)
 * style="width:95%" |
 * style="width:95%" |
 * 
 * }

The stress tensor can be decomposed into the spherical (isotropic) part, denoted as $$\displaystyle \mathbf s $$, and the deviatoric part, $$\displaystyle \mathbf d $$ as


 * {| style="width:100%" border="0"

$$  \displaystyle \boldsymbol \tau =  \mathbf s +    \mathbf d $$ (4)
 * style="width:95%" |
 * style="width:95%" |
 * 
 * }

and


 * {| style="width:100%" border="0"

$$  \displaystyle \mathbf s = \bar \tau  \mathbf I $$ (5)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

where $$\displaystyle \bar \tau $$ is the mean stress defined in terms of the trace of the stress tensor, $$\displaystyle Tr (\boldsymbol \tau) $$, as

<span id="(6)">
 * {| style="width:100%" border="0"

$$ \displaystyle \bar \tau = \frac{1}{3}  \tau_{ii} = \frac{1}{3} Tr (\boldsymbol \tau) $$ (6)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

<span id="(7)">
 * {| style="width:100%" border="0"

$$ \displaystyle Tr (\boldsymbol \tau) = \tau_{ii} = -p \underbrace{\delta_{ii}}_{=3} + 2 \mu \underbrace{\frac{\partial u_i} {\partial x_i}}_{= {\rm div} \mathbf u} $$ (7)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

For an incompressible flow: $$\displaystyle {\rm div} \mathbf u = \mathbf 0 $$. Then,

<span id="(8)">
 * {| style="width:100%" border="0"

$$\displaystyle Tr (\boldsymbol \tau) = -3p $$ (8)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

<span id="(9)">
 * {| style="width:100%" border="0"

$$\displaystyle p = -\frac{1}{3} Tr (\boldsymbol \tau) $$ (9)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Thus $$\displaystyle p $$ is the mean stress $$\displaystyle \bar \tau $$ with a minus sign. Hence,

<span id="(10)">
 * {| style="width:100%" border="0"

$$\displaystyle \mathbf s = -p \mathbf I $$ (10)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

<span id="(11)">
 * {| style="width:100%" border="0"

$$\displaystyle \mathbf d = \boldsymbol \tau - \mathbf s = 2 \mu   \boldsymbol \gamma (\mathbf u) $$ (11)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Body force $$\displaystyle \mathbf b $$ : Assume $$\displaystyle \mathbf b $$ is the gradient of a potential $$\displaystyle \Psi $$ (like gravitational potential), then $$\displaystyle \mathbf b = - \rho \, {\rm grad} \Psi $$. Using Eq.(2), Eq.(1) becomes

<span id="(12)">
 * {| style="width:100%" border="0"

$$  \displaystyle \rho \frac {D \mathbf u}  {Dt} =  {\rm div} (-p \mathbf I)  + 2 \mu \, {\rm div} \boldsymbol \gamma (\mathbf u) - \rho \, {\rm grad} \Psi $$ (12)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

But,

<span id="(13)">
 * {| style="width:100%" border="0"

$$  \displaystyle {\rm div} (p \mathbf I) = \left ( \frac{\partial } {\partial x_i} p \, \delta_{ij} \right ) \mathbf e_j = \frac{\partial p} {\partial x_j} \mathbf e_j = {\rm grad}\, p  $$ (13)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Now, let’s focus on the $$\displaystyle {\rm div} \boldsymbol \gamma (\mathbf u) $$ term.

<span id="(14)">
 * {| style="width:100%" border="0"

$$  \displaystyle {\rm div} \boldsymbol \gamma (\mathbf u) = \frac{\partial } {\partial x_i} \gamma_{ij} \, \mathbf e_j = \frac{1}{2} \frac{\partial } {\partial x_i} \left ( \frac{\partial u_i} {\partial x_j} + \frac{\partial u_j} {\partial x_i} \right) \mathbf e_j = \frac{1}{2} \left [ \frac{\partial } {\partial x_j} \underbrace{ \left(\frac{\partial u_i} {\partial x_i} \right) }_{ {\rm div} \mathbf u = \mathbf 0} + \frac{\partial^2 u_j} {\partial x_i \partial x_i}
 * style="width:95%" |
 * style="width:95%" |

\right] \mathbf e_j = \frac{1}{2} \nabla^2 \mathbf u

$$ (14)
 * <p style="text-align:right">
 * }

Therefore, for an incompressible flow and assuming constant properties

<span id="(15)">
 * {| style="width:100%" border="0"

$$  \displaystyle \rho \frac {D \mathbf u}  {Dt} =  -{\rm grad} \,  P     + \mu \nabla^2 \mathbf u $$ (15)
 * style="width:92%; padding:10px; border:2px solid #8888aa" |
 * style="width:92%; padding:10px; border:2px solid #8888aa" |
 * <p style="text-align:right">
 * }

where the modified pressure is defined as $$\displaystyle P = p + \rho \Psi $$. Eq.(15) is known as Navier-Stokes (NS) equation. We will take the curl of the NS equation to arrive to the vorticity transport equation, given as Eq.(2.60) (Pope, 2000, p.126)

<span id="(16)">
 * {| style="width:100%" border="0"

$$  \displaystyle \frac{D \boldsymbol \omega}{Dt} = \boldsymbol \omega \cdot \overset{\rightarrow}{\nabla} \mathbf u + \nu \nabla^2 \boldsymbol \omega $$ (16)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Here, the gradient of the velocity vector is represented with a right arrow on the nabla operator, which is expressed in component form as

<span id="(16b)">
 * {| style="width:100%" border="0"

$$  \displaystyle \overset{\rightarrow}{\nabla} \mathbf u = \frac{\partial u_i}{\partial x_j} {\mathbf e_j} \otimes {\mathbf e_i} $$ (16b)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

This definition is the transpose of our definition given in Eq.(21). To have a clear understanding of the distinctions between the two component forms of the gradient of a vector, the reader is referred to the article Gradient of vector: Two tensor conventions. In Eq.(16), $$\displaystyle \boldsymbol \omega $$ is the vorticity, defined as $$\displaystyle \nabla \times \mathbf u = e_{ijk} \partial_i u_j \mathbf e_k $$. Note that our definition of the gradient is the transpose of Pope’s, the vorticity equation in terms of our gradient is given in Eq.(1).

The curl of the NS equation, i.e., Eq.(15), becomes

<span id="(18)">
 * {| style="width:100%" border="0"

$$  \displaystyle \nabla \times \frac {D \mathbf u}  {Dt} =  \nabla \times \left( -\frac {1}{\rho} {\rm grad} \, P \right) +  \nabla \times \left(  \mu \nabla^2 \mathbf u \right) $$ (18)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Let’s first focus on the left hand side, using the definition of the material time derivative

<span id="(19)">
 * {| style="width:100%" border="0"

$$  \displaystyle \nabla \times \frac {D \mathbf u}  {Dt} =  \nabla \times \left( \nabla \mathbf u \cdot \mathbf  u \right) +  \nabla \times \frac{\partial \mathbf u} {\partial t}
 * style="width:95%" |
 * style="width:95%" |

$$ (19)
 * <p style="text-align:right">
 * }

Here, the last terms is

<span id="(20)">
 * {| style="width:100%" border="0"

$$  \displaystyle \nabla \times \frac{\partial \mathbf u} {\partial t} = \frac{\partial \boldsymbol \omega} {\partial t}
 * style="width:92%; padding:10px; border:2px solid #8888aa" |
 * style="width:92%; padding:10px; border:2px solid #8888aa" |

$$ (20)
 * <p style="text-align:right">
 * }

Note that

<span id="(21)">
 * {| style="width:100%" border="0"

$$  \displaystyle \nabla \mathbf u = \frac{\partial u_i} {\partial x_j} \mathbf e_i \otimes \mathbf e_j $$ (21)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

and

<span id="(22)">
 * {| style="width:100%" border="0"

$$  \displaystyle \nabla \mathbf u \cdot \mathbf u = \frac{\partial u_i} {\partial x_j} u_j \mathbf e_i $$ (22)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

So, the second term coming from the material time derivative is

<span id="(23)">
 * {| style="width:100%" border="0"

$$  \displaystyle \nabla \times \left( \nabla \mathbf u \cdot \mathbf  u \right) = e_{ijk} \partial_i \left[ (\partial_p u_j) u_p \right] \mathbf e_k = e_{ijk} \left[ (\partial_i \partial_p u_j) u_p + (\partial_p u_j)( \partial_i u_p) \right] \mathbf e_k $$ (23)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Let

<span id="(24)">
 * {| style="width:100%" border="0"

$$  \displaystyle \boldsymbol \alpha = e_{ijk} (\partial_i \partial_p u_j) u_p  \mathbf e_k = u_p \partial_p \underbrace{ ( e_{ijk} \partial_i u_j) \mathbf e_k}_{\nabla \times \mathbf u = \boldsymbol \omega} $$ (24)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

which gives

<span id="(25)">
 * {| style="width:100%" border="0"

$$  \displaystyle \boldsymbol \alpha = \nabla \boldsymbol \omega \cdot \mathbf u
 * style="width:92%; padding:10px; border:2px solid #8888aa" |
 * style="width:92%; padding:10px; border:2px solid #8888aa" |

$$ (25)
 * <p style="text-align:right">
 * }

Let

<span id="(26)">
 * {| style="width:100%" border="0"

$$  \displaystyle \boldsymbol \beta = e_{ijk} (\partial_p u_j)( \partial_i u_p) \mathbf e_k $$ (26)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

To facilitate the writing, let’s define Jacobian as

<span id="(27)">
 * {| style="width:100%" border="0"

$$  \displaystyle \mathbf J = {\rm grad} \ \mathbf u = \frac{\partial u_j} {\partial x_i} \mathbf e_j\otimes \mathbf e_i $$ (27)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Or, in component form

<span id="(28)">
 * {| style="width:100%" border="0"

$$  \displaystyle J_i^j = \frac{\partial u_j} {\partial x_i} $$ (28)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Here, the superscript represents the row, and the subscript the column index. Also, let

<span id="(29)">
 * {| style="width:100%" border="0"

$$  \displaystyle \mathbf A = A_i^j \mathbf e_j\otimes \mathbf e_i = \mathbf J^2 $$ (29)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Then, $$\displaystyle \boldsymbol \beta $$ can be rewritten as

<span id="(30)">
 * {| style="width:100%" border="0"

$$  \displaystyle \boldsymbol \beta = e_{ijk} A_i^j \mathbf e_k = \mathbf e_1 (A_2^3 - A_3^2) + \mathbf e_2 (A_3^1 - A_1^3) + \mathbf e_3 (A_1^2 - A_2^1) $$ (30)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

We can write the first component of $$\displaystyle \boldsymbol \beta $$ in terms of the Jacobian as

<span id="(31)">
 * {| style="width:100%" border="0"

$$  \displaystyle A_2^3 = \frac{\partial u_3} {\partial x_p} \frac{\partial u_p} {\partial x_2} = J_p^3 J_2^p = J_1^3 J_2^1 + J_2^3 J_2^2 + J_3^3 J_2^3 $$ (31)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

<span id="(32)">
 * {| style="width:100%" border="0"

$$  \displaystyle A_3^2 = \frac{\partial u_3} {\partial x_p} \frac{\partial u_p} {\partial x_2} = J_p^2 J_3^p = J_1^2 J_3^1 + J_2^2 J_3^2 + J_3^2 J_3^3 $$ (32)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Then,

<span id="(33)">
 * {| style="width:100%" border="0"

$$  \displaystyle A_2^3 - A_3^2 = \frac{\partial u_3} {\partial x_p} \frac{\partial u_p} {\partial x_2} = J_1^3 J_2^1 - J_3^1 J_1^2 + (J_2^2 + J_3^3)( J_2^3 - J_3^2) $$ (33)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Similar expressions can be written for the other components of $$\displaystyle \boldsymbol \beta $$. To obtain the vorticity equation, we want to relate $$\displaystyle \boldsymbol \beta $$ to $$\displaystyle \boldsymbol \gamma $$, which is given as $$\displaystyle \nabla \mathbf u \cdot \boldsymbol \omega$$. Let’s first write the vorticity in ters of the Jacobian, then we will obtain $$\displaystyle \boldsymbol \gamma $$.

<span id="(34)">
 * {| style="width:100%" border="0"

$$  \displaystyle \boldsymbol \omega = \nabla \times \mathbf u = \underbrace{ (\partial_ 2 u_3 - \partial_3 u_2)}_{ J_2^3 - J_3^2} \mathbf e_1 + \underbrace{ (\partial_3 u_1 - \partial_1 u_3)}_{ J_3^1 J_1^3}\mathbf e_2 + \underbrace{ (\partial_1 u_2 - \partial_2 u_1) }_{ J_1^2 J_2^1} \mathbf e_3 $$ (34)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

<span id="(35)">
 * {| style="width:100%" border="0"

$$  \displaystyle \boldsymbol \gamma = \nabla \mathbf u \cdot \boldsymbol \omega = \left( \frac{\partial u_i} {\partial x_j} \mathbf e_i\otimes \mathbf e_j \right) \cdot (\omega_k \mathbf e_k) = \frac{\partial u_i} {\partial x_j} \omega_j \mathbf e_i = \underbrace{ \frac{\partial u_i} {\partial x_j}}_{J_j^i} \left(e_{klj} \partial_k u_l \right) \mathbf e_i = a_i \mathbf e_i $$ (35)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

The components of $$\displaystyle \boldsymbol \gamma $$ can be calculated using matrix algebra as

<span id="(36)">
 * {| style="width:100%" border="0"

$$  \displaystyle \begin{vmatrix} J_1^1 & J_2^1 & J_3^1 \\ J_1^2 & J_2^2 & J_3^2 \\ J_1^3 & J_2^3 & J_3^3 \\ \end{vmatrix} \begin{vmatrix} J_2^3 - J_3^2 \\ J_3^1 - J_1^3 \\ J_1^2 - J_2^1 \\ \end{vmatrix}
 * style="width:95%" |
 * style="width:95%" |

$$ (36)
 * <p style="text-align:right">
 * }

Let’s compute the first component of $$\displaystyle \boldsymbol \gamma $$, i.e., $$\displaystyle a_1 $$

<span id="(37)">
 * {| style="width:100%" border="0"

$$  \displaystyle a_1 = J_1^1 (J_2^3 - J_3^2) + J_2^1 (J_3^1 - J_1^3) + J_3^1 (J_1^2 - J_2^1) $$ (37)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

But, $$\displaystyle J_1^1 = -(J_2^2+J_3^3) $$ as $$\displaystyle J_i^i = 0 = {\rm div} \mathbf u $$. Therefore,

<span id="(38)">
 * {| style="width:100%" border="0"

$$  \displaystyle a_1 = -(J_2^2+J_3^3) (J_2^3 - J_3^2) + J_3^1 J_1^2 - J_2^1 J_1^3 $$ (38)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Now, we can relate the components of $$\displaystyle \boldsymbol \beta $$ to those of $$\displaystyle \boldsymbol \gamma $$ as

<span id="(39)"> :{| style="width:100%" border="0" $$  \displaystyle A_2^3 - A_3^2 = - a_1 $$ (39)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Similarly

<span id="(40)">
 * {| style="width:100%" border="0"

$$  \displaystyle A_3^1 - A_1^3 = - a_2 $$ (40)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

and

<span id="(41)">
 * {| style="width:100%" border="0"

$$  \displaystyle A_1^2 - A_2^1 = - a_3 $$ (41)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Hence,

<span id="(42)"> :{| style="width:100%" border="0" $$  \displaystyle \boldsymbol \beta = - \boldsymbol \gamma = -\nabla \mathbf u \cdot \boldsymbol \omega $$ (42)
 * style="width:92%; padding:10px; border:2px solid #8888aa" |
 * style="width:92%; padding:10px; border:2px solid #8888aa" |
 * <p style="text-align:right">
 * }

whose component form is given in Eq.(35). Now, we will work with the right hand side of the NS. The first terms is the pressure one

<span id="(43)"> :{| style="width:100%" border="0" $$  \displaystyle \nabla \times (\nabla P) = \mathbf 0 $$ (43)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

and as the density is constant

<span id="(44)"> :{| style="width:100%" border="0" $$  \displaystyle \nabla \times (\frac {1}{\rho} \nabla P) = \mathbf 0 $$ (44)
 * style="width:92%; padding:10px; border:2px solid #8888aa" |
 * style="width:92%; padding:10px; border:2px solid #8888aa" |
 * <p style="text-align:right">
 * }

Next,

<span id="(45)">
 * {| style="width:100%" border="0"

$$  \displaystyle \nabla \times (\nabla^2 \mathbf u) = {\rm curl} ({\rm div} ({\rm grad^T}\mathbf u)) = e_{ijk} \partial_i (\partial_p \partial_p u_j) \mathbf e_k = \partial_p \partial_p \underbrace{(e_{ijk} \partial_i u_j \mathbf e_k )}_{ \nabla \times \mathbf u = \boldsymbol \omega} = \nabla^2 \boldsymbol \omega $$ (45)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

As the viscosity and the density are constant

<span id="(46)">
 * {| style="width:100%" border="0"

$$  \displaystyle \nabla \times (\frac{\mu}{\rho} \nabla^2 \mathbf u) = \nu \nabla^2 \boldsymbol \omega $$ (46)
 * style="width:92%; padding:10px; border:2px solid #8888aa" |
 * style="width:92%; padding:10px; border:2px solid #8888aa" |
 * <p style="text-align:right">
 * }

Combining all terms, we obtain the evolution equation for the vorticity for incompressible flows as

<span id="(47)">
 * {| style="width:100%" border="0"

$$  \displaystyle \frac{D \boldsymbol \omega}{Dt} = \nabla \mathbf u \cdot \boldsymbol \omega + \nu \nabla^2 \boldsymbol \omega $$ (47)
 * style="width:92%; padding:10px; border:2px solid #8888aa" |
 * style="width:92%; padding:10px; border:2px solid #8888aa" |
 * <p style="text-align:right">
 * }

The vorticity equation in component form is

<span id="(48)">
 * {| style="width:100%" border="0"

$$  \displaystyle \frac{D \omega_i}{Dt} = \frac{\partial u_i}{\partial x_j} \omega_j + \nu \nabla^2 \omega_i $$ (48)
 * style="width:92%; padding:10px; border:2px solid #8888aa" |
 * style="width:92%; padding:10px; border:2px solid #8888aa" |
 * <p style="text-align:right">
 * }

Here, the first term is the time rate of change of vorticity and it represents the vorticity production when it is positive. $$\displaystyle \nabla \mathbf u \cdot \boldsymbol \omega $$ is the vorticity production due to vortex-line stretching and $$\displaystyle  \nu \nabla^2 \boldsymbol \omega $$ is the vorticity diffusion.

<span id="(An alternative method)">

An alternative method
The Jacobian $$\displaystyle \mathbf J $$ can be decomposed as

<span id="(49)"> :{| style="width:100%" border="0" $$  \displaystyle J_j^i = \underbrace{\frac{1}{2}(J_j^i + J_i^j)}_{S_j^i} - \underbrace{\frac{1}{2}(J_i^j - J_j^i)}_{\Omega_j^i} $$ (49)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

So,

<span id="(50)">
 * {| style="width:100%" border="0"

$$  \displaystyle \mathbf J = \mathbf S + \boldsymbol \Omega
 * style="width:92%; padding:10px; border:2px solid #8888aa" |
 * style="width:92%; padding:10px; border:2px solid #8888aa" |

$$ (50)
 * <p style="text-align:right">
 * }

where $$\displaystyle \mathbf S $$ is a small strain tensor and $$\displaystyle  \boldsymbol  \Omega$$ is a small rotation tensor. The vorticity can be expressed as

<span id="(51)"> :{| style="width:100%" border="0" $$  \displaystyle \boldsymbol \omega = \nabla \times \mathbf u = e_{ijk} \underbrace{\partial_i u_j}_{J_i^j} \mathbf e_k = \underbrace{(J_2^3 - J_3^2)}_{2 \Omega_2^3} \mathbf e_1 + \underbrace{(J_3^1 - J_1^3)}_{2 \Omega_3^1} \mathbf e_2 + \underbrace{(J_1^2 - J_2^1)}_{2 \Omega_1^2} \mathbf e_3 = 2 \Omega_i^j \mathbf e_k $$ (51)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Let’s define an extraction operator, $$\displaystyle \zeta$$. When it operates on the Jacobian it extracts the pseudo vector $$\displaystyle \boldsymbol \omega $$ as

<span id="(52)"> :{| style="width:100%" border="0" $$  \displaystyle \zeta \mathbf J = e_{ijk} J_i^j \mathbf e_k = \boldsymbol \omega $$ (52)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

<span id="(53)"> :{| style="width:100%" border="0" $$  \displaystyle \zeta \mathbf J = \zeta \mathbf S + \zeta \boldsymbol \Omega = \zeta \boldsymbol \Omega = \boldsymbol \omega $$ (53)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Note that as $$\displaystyle \mathbf  S$$  is symmetric $$\displaystyle  \zeta \mathbf S = \mathbf  0$$. Our goal is to find an expression for $$\displaystyle \boldsymbol \beta  $$, which is expressed as

<span id="(54)"> :{| style="width:100%" border="0" $$  \displaystyle \boldsymbol \beta = e_{ijk} \underbrace{(\partial_p u_j)}_{J_p^j} \underbrace{(\partial_i u_p)}_{J_i^p} \mathbf e_k = \zeta \mathbf J^2 $$ (54) and
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

<span id="(55)"> :{| style="width:100%" border="0" $$  \displaystyle \mathbf J^2 = (\mathbf S + \boldsymbol \Omega) (\mathbf S + \boldsymbol \Omega) = \mathbf S^2 + \boldsymbol \Omega^2 + \boldsymbol \Omega \mathbf S + \mathbf S \boldsymbol \Omega $$ (55)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Here,

<span id="(56)"> :{| style="width:100%" border="0" $$  \displaystyle (\mathbf S^2)^T = \mathbf S^T \mathbf S^T = \mathbf S^2 $$ (56)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

<span id="(57)"> :{| style="width:100%" border="0" $$  \displaystyle (\boldsymbol \Omega^2)^T = \boldsymbol \Omega^T \boldsymbol \Omega^T = (-\boldsymbol \Omega) (-\boldsymbol \Omega) = \boldsymbol \Omega^2 $$ (57)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

As both $$\displaystyle \mathbf S^2 $$ and $$\displaystyle \boldsymbol \Omega^2 $$ are symmetric

<span id="(58)"> :{| style="width:100%" border="0" $$  \displaystyle \zeta \mathbf S^2 = \zeta \boldsymbol \Omega^2 = \mathbf 0 $$ (58)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

<span id="(59)"> :{| style="width:100%" border="0" $$  \displaystyle (\boldsymbol \Omega \mathbf S)^T = \mathbf S^T \boldsymbol \Omega^T = \mathbf S (-\boldsymbol \Omega) = - \mathbf S \boldsymbol \Omega \Longrightarrow \mathbf S \boldsymbol \Omega = - (\boldsymbol \Omega \mathbf S)^T $$ (59)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Eq.(59) implies

<span id="(60)"> :{| style="width:100%" border="0" $$  \displaystyle \boldsymbol \Omega \mathbf S + \mathbf S \boldsymbol \Omega = \boldsymbol \Omega \mathbf S - (\boldsymbol \Omega \mathbf S)^T = \mathbf S \boldsymbol \Omega - (\mathbf S \boldsymbol \Omega)^T $$ (60)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Now, we can use the extraction operator

<span id="(61)"> :{| style="width:100%" border="0" $$  \displaystyle \zeta (\mathbf S \boldsymbol \Omega) = e_{ijk}S_p^j \underbrace{\Omega_i^p}_{\frac{1}{2}e_{ipq}\omega_q} \mathbf e_k = \frac{1}{2}\underbrace{e_{ijk}e_{ipq}}_{\delta_{jp}\delta_{kq}-\delta_{jq}\delta_{kp}} S_p^j \omega_q \mathbf e_k = \frac{1}{2}\left[\underbrace{S_j^j}_{=0}\omega_k - \underbrace{S_k^j}_{=S_j^k}\omega_j \right] \mathbf e_k = -\frac{1}{2} \mathbf S \cdot \boldsymbol \omega $$ (61)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

<span id="(62)"> :{| style="width:100%" border="0" $$  \displaystyle \zeta (\mathbf S \boldsymbol \Omega)^T = \underbrace{ e_{ijk}}_{-e_{jik}} S_p^i \Omega_j^p \mathbf e_k = -\frac{1}{2}\underbrace{e_{jik}e_{jpq}}_{\delta_{ip}\delta_{kq}-\delta_{iq}\delta_{kp}} S_p^i \omega_q \mathbf e_k = -\frac{1}{2}\left[\underbrace{S_i^i}_{=0}\omega_k - \underbrace{S_k^i}_{=S_i^k}\omega_i \right] \mathbf e_k = \frac{1}{2} \mathbf S \cdot \boldsymbol \omega $$ (62)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Collecting the terms,

<span id="(63)"> :{| style="width:100%" border="0" $$  \displaystyle \beta = \zeta \mathbf J^2 = \zeta \left[ (\mathbf S \boldsymbol \Omega) - (\mathbf S \boldsymbol \Omega)^T \right] = -\mathbf S \cdot \Omega $$ (63)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

But, $$\displaystyle \mathbf S \cdot \boldsymbol \Omega = \underbrace{(\mathbf S + \boldsymbol \Omega)}_{\mathbf J}\cdot \boldsymbol \omega $$ as $$\displaystyle  \boldsymbol \Omega \boldsymbol \omega = \boldsymbol \omega \times \boldsymbol \omega = \mathbf 0$$. Therefore,

<span id="(64)">
 * {| style="width:100%" border="0"

$$  \displaystyle \beta = -\mathbf J \cdot \omega = -\nabla \mathbf u \cdot \boldsymbol \omega $$ (64)
 * style="width:92%; padding:10px; border:2px solid #8888aa" |
 * style="width:92%; padding:10px; border:2px solid #8888aa" |
 * <p style="text-align:right">
 * }

This completes the proof.

<span id="(Using vector identities)">

Using vector identities
Here, we use the vector calculus identities to derive the vorticity equation. It should be noted that the gradient used in the vector calculus identities is the form used in Eq.(16b), which is the transpose of the gradient, given in Eq.(21), used in this article. The reader is referred to the article Gradient of vector: Two tensor conventions for a clear understanding of the two conventions.

The curl of the Navier-Stokes equation, i.e. Eq.(15), gives

<span id="(65)">
 * {| style="width:100%" border="0"

$$  \displaystyle \nabla \times  \frac {D \mathbf u}  {Dt} = \nabla \times \left( -\frac{1}{\rho} \nabla P \right) + \nabla \times \left(\nu \nabla^2 \mathbf u \right) $$     (65)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Let’s first focus on the left hand side

<span id="(66)">
 * {| style="width:100%" border="0"

$$  \displaystyle \nabla \times  \frac {D \mathbf u}  {Dt} = \nabla \times  \left(\frac{\partial \mathbf u}{\partial  t}\right) + \nabla \times \left(\mathbf u \cdot \overset{\rightarrow}{\nabla} \mathbf u \right) $$     (66)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

But, the partial derivative in Eq.(66) can be written as

<span id="(67)">
 * {| style="width:100%" border="0"

$$  \displaystyle \nabla \times  \left(\frac{\partial \mathbf u}{\partial  t}\right) = \frac{\partial}{\partial t} \underbrace{\left (\nabla \times \mathbf u \right)}_{ \boldsymbol \omega} $$     (67)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Also, using one of the vector dot product identities

<span id="(68)">
 * {| style="width:100%" border="0"

$$  \displaystyle \mathbf{u} \cdot (\overset{\rightarrow}{\nabla} \mathbf{u}) = \frac{1}{2}\nabla \left(\mathbf u \cdot \mathbf u \right) - \mathbf u \times \underbrace{\left (\nabla \times \mathbf u \right)}_{ \boldsymbol \omega} $$     (68)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Now, using one of the vector cross product identities

<span id="(69)">
 * {| style="width:100%" border="0"

$$  \displaystyle \nabla \times \left(\mathbf u \times \boldsymbol \omega \right) = \mathbf u \underbrace{\left( \nabla \cdot \boldsymbol \omega \right)}_{0} - \boldsymbol \omega \underbrace{\left( \nabla \cdot \mathbf u \right)}_{0} + \boldsymbol \omega \cdot \overset{\rightarrow}{\nabla} \mathbf u - \mathbf u \cdot \overset{\rightarrow}{\nabla}\boldsymbol \omega $$     (69)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

The first term on the right of Eq.(69) is zero because it is the divergence of a curl, and the second term is zero because of the incompressibility condition. As the curl of the gradient of a scalar field vanishes

<span id="(70)">
 * {| style="width:100%" border="0"

$$  \displaystyle \nabla \times \left( \frac{1}{2}\nabla \left(\mathbf u \cdot \mathbf u \right) \right) = \mathbf 0 $$     (70)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Grouping all terms on the left hand side

<span id="(71)">
 * {| style="width:100%" border="0"

$$  \displaystyle \nabla \times  \frac{D \mathbf u}{Dt} = \frac{\partial \boldsymbol \omega}{\partial t} - \boldsymbol \omega \cdot \overset{\rightarrow}{\nabla} \mathbf u + \mathbf u \cdot \overset{\rightarrow}{\nabla}\boldsymbol \omega = \frac{D \boldsymbol \omega}{Dt} - \boldsymbol \omega \cdot \overset{\rightarrow}{\nabla} \mathbf u $$ (71)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

The first term on the right hand side of Eq.(65) vanishes as in Eq.(70), i.e.,

<span id="(72)">
 * {| style="width:100%" border="0"

$$  \displaystyle \nabla \times \left( -\frac{1}{\rho} \nabla P \right) = \mathbf 0 $$     (72)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Finally, the last term in Eq.(65) is

<span id="(73)">
 * {| style="width:100%" border="0"

$$  \displaystyle \nabla \times \left(\nu \nabla^2 \mathbf u \right) = \nu \nabla^2 \left(\nabla \times \mathbf u \right) = \nu \nabla^2 \boldsymbol \omega $$     (73)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Combining the left and the right hand sides, we obtain

<span id="(74)">
 * {| style="width:100%" border="0"

$$  \displaystyle \frac{D \boldsymbol \omega}{Dt} - \boldsymbol \omega \cdot \overset{\rightarrow}{\nabla} \mathbf u = \nu \nabla^2 \boldsymbol \omega $$  (74)
 * style="width:92%; padding:10px; border:2px solid #8888aa" |
 * style="width:92%; padding:10px; border:2px solid #8888aa" |
 * <p style="text-align:right">
 * }

<span id="(Meaning of vortex-line stretching)">

Meaning of vortex-line stretching
Vortex-line is a line along the vorticity vector $$\displaystyle \mathbf \omega $$ at a point $$\displaystyle x $$. In a uniform flow field, $$\displaystyle \nabla \mathbf u = \mathbf 0$$, therefore two points on an infinitesimally small vorticity line have the same velocity and the line stays undistorted. However, in a non-uniform field the vortex line is stretched as different velocities are assumed at different points of the line. Therefore, $$\displaystyle \nabla \mathbf u \cdot \boldsymbol \omega $$ represents the vorticity production by vortex-line stretching.

<span id="(Evolution equation for an infinitesimally small line element)">

Evolution equation for an infinitesimally small line element
The mapping at time $$\displaystyle t $$ between the spatial position $$\displaystyle \mathbf x $$ and the material position $$\displaystyle \mathbf X$$ is given with

<span id="(75)">
 * {| style="width:100%" border="0"

$$  \displaystyle \mathbf x = \phi (\mathbf X,t) = \phi_t (\mathbf X) $$ (75)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Then a small displacement with series approximation can be written as,

<span id="(76)">
 * {| style="width:100%" border="0"

$$  \displaystyle \mathbf x + \mathbf h = \phi_t(\mathbf X+ \mathbf H) = \phi_t(\mathbf X) + \frac{\partial \phi(\mathbf X)}{\partial X_k} H_k + h.o.t = \underbrace{ \phi_t(\mathbf X)}_{ \mathbf x} + \underbrace{\left [\rm grad_X \phi_t(\mathbf X) \right] \cdot  \mathbf  H(\mathbf X)}_{\mathbf  h} $$ (76)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

The material time derivative is defined at constant material point as

<span id="(77)">
 * {| style="width:100%" border="0"

$$  \displaystyle
 * style="width:95%" |
 * style="width:95%" |

\frac{D \mathbf h(\mathbf x, t)}{Dt} = \left. \frac{d \mathbf h(\mathbf x, t)}{dt} \right|_{X \ fixed} = \left[ \rm grad_X \underbrace{ \frac{\partial \phi (\mathbf X, t) }{\partial t } }_{\mathbf V(\mathbf X, t)= \mathbf v(\mathbf x, t) } \right] $$     (77)
 * <p style="text-align:right">
 * }

Note that the material description of the velocity field, $$ \displaystyle \mathbf V (\mathbf X, t) $$ is equal to the spatial description of the velocity field, $$ \displaystyle \mathbf v (\mathbf x, t) $$, i.e., $$ \displaystyle \mathbf V (\mathbf X, t) =\mathbf v (\mathbf x, t) $$ and using the chain rule

<span id="(78)">
 * {| style="width:100%" border="0"

$$  \displaystyle \frac{\partial \mathbf v}{\partial x_i} = \frac{\partial \mathbf V}{\partial X_k} \frac{\partial X_k}{\partial x_i} $$     (78)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

To write the gradient of the velocity field with material description in terms of that for spatial description, first the definition for the gradient is introduced as

<span id="(79)">
 * {| style="width:100%" border="0"

$$  \displaystyle \rm grad_X \mathbf V = \frac{\partial V_i}{\partial X_j} \mathbf e_i \otimes \mathbf e_j $$     (79)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Again, using the chain rule,

<span id="(80)">
 * {| style="width:100%" border="0"

$$  \displaystyle \frac{\partial V_i}{\partial X_j} = \underbrace{\frac{\partial V_i}{\partial x_p}}_{\frac{\partial v_i}{\partial x_j}} \underbrace{\frac{\partial x_p}{\partial X_j}}_{F_j^p} $$     (80)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Now, the gradient of the material description of the velocity vector can be expressed in terms of the gradient of the spatial description of the velocity vector as

<span id="(81)">
 * {| style="width:100%" border="0"

$$  \displaystyle \rm grad_X \mathbf V = \rm grad_x \mathbf v \cdot \mathbf F $$ (81)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

where $$ \displaystyle \mathbf F $$ is defined in Eq.(80). Recall that the material or the spatial descriptions of the velocity fields are equal and can be expressed as

<span id="(82)">
 * {| style="width:100%" border="0"

$$  \displaystyle \mathbf V = V_i \mathbf E_i = \mathbf v= v_i \mathbf e_i $$     (82)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

The basis vectors are equal

<span id="(83)">
 * {| style="width:100%" border="0"

$$  \displaystyle \mathbf e_i = \mathbf E_i \, \ i=1,2,3 $$     (83)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Hence, the components of the velocities are equal to each other as

<span id="(84)">
 * {| style="width:100%" border="0"

$$  \displaystyle v_i(\mathbf x,t) = V_i (\mathbf X,t) $$     (84)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

The term $$ \displaystyle \mathbf F $$ in Eq.(80) can be written in component form as

<span id="(85)">
 * {| style="width:100%" border="0"

$$  \displaystyle \mathbf F = F_j^i \mathbf e_i \otimes \mathbf E_j $$     (85)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Combining Eq.(79) and Eq.(83), the gradient of the material description of the velocity field is

<span id="(86)">
 * {| style="width:100%" border="0"

$$  \displaystyle \nabla_X \mathbf V = \frac{\partial V_i}{\partial X_j} \mathbf E_i \otimes \mathbf E_j $$     (86)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

and the gradient of the spatial description of the velocity field is

<span id="(87)">
 * {| style="width:100%" border="0"

$$  \displaystyle \nabla_x \mathbf v = \frac{\partial v_i}{\partial x_j} \underbrace{ \mathbf e_i }_{\mathbf E_i} \otimes \; \mathbf e_j $$     (87)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Then,

<span id="(88)">
 * {| style="width:100%" border="0"

$$  \displaystyle \nabla_X \mathbf V \cdot \mathbf H = \left(\nabla_x \mathbf v \cdot \mathbf F \right) \cdot \mathbf H = \nabla_x(\mathbf v) \cdot \underbrace{\left(\mathbf F \cdot \mathbf H \right) }_{\mathbf h} $$ (88)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

Hence, a small change in the spatial position is expressed as

<span id="(89)">
 * {| style="width:100%" border="0"

$$  \displaystyle \frac{D \mathbf h}{Dt} = \left( \nabla_x \mathbf v \right) \cdot \mathbf h $$ (89)
 * style="width:92%; padding:10px; border:2px solid #8888aa" |
 * style="width:92%; padding:10px; border:2px solid #8888aa" |
 * <p style="text-align:right">
 * }

The corresponding expression in Pope (2000, p.14, Eq.(2.16)) is given as

<span id="(90)">
 * {| style="width:100%" border="0"

$$  \displaystyle \frac{d \mathbf s}{dt} = \mathbf s \cdot \nabla \mathbf u $$ (90)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

and the component form is given as Pope (2000, p.13, Eq.(2.12))

<span id="(91)">
 * {| style="width:100%" border="0"

$$  \displaystyle \underbrace{\mathbf s \cdot \nabla \mathbf u}_{\nabla \mathbf u \cdot \mathbf s \ \rm our \; \rm notation} = s_i \frac{ \partial u_j}{ \partial x_i} \mathbf e_j $$     (91)
 * style="width:95%" |
 * style="width:95%" |
 * <p style="text-align:right">
 * }

= References =