# src/navier-stokes/swirl.h

# Azimuthal velocity for axisymmetric flows

The centered Navier–Stokes solver can be combined with the axisymmetric metric but assumes zero azimuthal velocity (“swirl”). This file adds this azimuthal velocity component: the $w$ field.

Assuming that $x$ is the axial direction and $y$ the radial direction, as in axi.h, the incompressible, variable-density and viscosity, axisymmetric Navier–Stokes equations (with swirl) can be written $${\partial}_{x}{u}_{x}+{\partial}_{y}{u}_{y}+\frac{{u}_{y}}{y}=0$$ $${\partial}_{t}{u}_{x}+{u}_{x}{\partial}_{x}{u}_{x}+{u}_{y}{\partial}_{y}{u}_{x}=-\frac{1}{\rho}{\partial}_{x}p+\frac{1}{\rho y}\nabla \cdot (\mu y\nabla {u}_{x})$$ $${\partial}_{t}{u}_{y}+{u}_{x}{\partial}_{x}{u}_{y}+{u}_{y}{\partial}_{y}{u}_{y}-\frac{{w}^{2}}{y}=-\frac{1}{\rho}{\partial}_{y}p+\frac{1}{\rho y}(\nabla \cdot (\mu y\nabla {u}_{y})-2\mu \frac{{u}_{y}}{y})$$ $${\partial}_{t}w+{u}_{x}{\partial}_{x}w+{u}_{y}{\partial}_{y}w+\frac{{u}_{y}w}{y}=\frac{1}{\rho y}[\nabla \cdot (\mu y\nabla w)-w(\frac{\mu}{y}+{\partial}_{y}\mu )]$$ We will thus need to solve and advection-diffusion equation for $w$.

```
#include "tracer.h"
#include "diffusion.h"
scalar w[], * tracers = {w};
```

The azimuthal velocity is zero on the axis of symmetry ($y=0$).

`w[bottom] = dirichlet(0);`

We will need to add the acceleration term ${w}^{2}/y$ in the evolution equation for ${u}_{y}$. If the acceleration field is not allocated yet, we do so.

```
event defaults (i = 0)
{
if (is_constant(a.x)) {
a = new face vector;
foreach_face()
a.x[] = 0.;
boundary ((scalar *){a});
}
}
```

The equation for ${u}_{y}$ is solved by the centered Navier–Stokes solver combined with the axisymmetric metric, but the acceleration term ${w}^{2}/y$ is missing. We add it here, taking care of the division by zero on the axis, and averaging $w$ from cell center to cell face.

```
event acceleration (i++)
{
face vector av = a;
foreach_face (y)
av.y[] += y > 0. ? sq(w[] + w[0,-1])/(4.*y) : 0.;
}
```

The advection of $w$ is done by the tracer solver, but we need to add diffusion. Using the diffusion solver, we solve $$\theta {\partial}_{t}w=\nabla \cdot (D\nabla w)+\beta w$$ Identifying with the diffusion part of the equation for $w$ above, we have $$\begin{array}{ccc}\hfill \theta & \hfill =\hfill & \rho y\hfill \\ \hfill D& \hfill =\hfill & \mu y\hfill \\ \hfill \beta & \hfill =\hfill & -(\rho {u}_{y}+\frac{\mu}{y}+{\partial}_{y}\mu )\hfill \end{array}$$ Note that the *rho* and *mu* fields (defined by the Navier–Stokes solver) already include the metric (i.e. are $\rho y$ and $\mu y$), which explains the divisions by $y$ in the code below.

```
event tracer_diffusion (i++)
{
scalar β[], θ[];
foreach() {
θ[] = ρ[];
double muc = (μ.x[] + μ.x[1] + μ.y[] + μ.y[0,1])/4.;
double dymu = (μ.y[0,1]/fm.y[0,1] - μ.y[]/fm.y[])/Δ;
β[] = - (ρ[]*u.y[] + muc/y)/y - dymu;
}
diffusion (w, dt, μ, θ = θ, β = β);
}
```