SPINS Science Files

SPINS comes with a variety of Science files which help with simulation data analysis. Below is a list of currently available functions with a brief description of their purpose and the assumptions behind their implementation. These Science files are in spins/src/Science.

Unless otherwise mentioned, the functions below work for both unmapped and mapped domains where the mapping is along the bottom of the domain (where the bottom is at the 0-index in the z-dimension, as opposed to at the $N_{z}-1$ index). Another way to say this is that z increases with the index.

Units

The units of each of the functions have been built under the assumption that physical units are being used in the simulation. These computations will work for nondimensional simulations, but more care must be taken to ensure the scalings are appropriate. Physical units are typically:

Length	m
Time	s
Density	kg/m^3

The reference density should be 1000 kg/m^3 for the stress and dissipation calculations to be physically correct. The density field can be (and often is) non-dimensionalized to $\rho '=\rho /\rho _{0}-1$. If nondimensionalized, compute_Background_PE.cpp and compute_BPE_from_internal.cpp assume this form of scaling. These two functions also need to know if the density field is dimensional or not.

bottom_slope.cpp

Compute the vector of the bottom slope at each x-grid point. This is used for calculating the stresses along the bottom boundary (cf. bottom_stress_x.cpp and bottom_stress_y.cpp).

bottom_stress_x.cpp

Compute the 2D vector of the x-component of stress on the bottom boundary. The z boundary condition must be no slip.

If the case is unmapped:

$t_{x}=\mu u_{z}$

where $\mu =\nu \rho _{0}$ is the dynamic viscosity.

If the case is mapped:

$t_{x}={\frac {\mu }{\alpha ^{2}(x)}}\left(2h'(x)(w_{z}-u_{x})+(u_{z}+w_{x})(1-h'^{2}(x))\right)$

where $\alpha (x)={\sqrt {1+h'^{2}(x)}}$ and where $h'(x)$ is the bottom slope.

bottom_stress_y.cpp

Compute the 2D vector of the y-component of stress on the bottom boundary. The z boundary condition must be no slip.

If the case is unmapped:

$t_{y}=\mu v_{z}$

where $\mu =\nu \rho _{0}$ is the dynamic viscosity.

If the case is mapped:

$t_{y}={\frac {\mu }{\alpha (x)}}\left(v_{z}-h'(x)v_{x}\right)$

where $\alpha (x)={\sqrt {1+h'^{2}(x)}}$ and where $h'(x)$ is the bottom slope.

compute_Background_PE.cpp

Compute the value of the background potential energy,

$E_{b}=g\int _{V}\rho _{*}z_{*}dV$

where $z_{*}$ is the depth of an element of the sorted density field which minimized the potential energy (ie. $\rho _{*}$). For more information see Winters et al. (1995).

compute_BPE_from_internal.cpp

Compute the rate of energy transfer from internal energy into background potential energy.

$\phi _{i}=-g\kappa \int _{A}\rho (z=z_{u})-\rho (z=z_{l})dA$

where $z_{u}$ and $z_{l}$ are the depths of the upper and lower boundaries, respectively. The integral is taken over the horizontal extent. For more information see Winters et al. (1995).

compute_quadweights.cpp

Compute the quadrature weights for each dimension.

If the boundary condition is no-slip, use the Clenshaw-Curtis quadrature weights:

$w=...$

If the boundary condition is free-slip, use the trapezoid rule:

$w=L_{x}/N_{x}=\Delta x$

compute_vorticity.cpp

Compute all vorticity components.

compute_vort_x.cpp

Compute the x-component of vorticity: $\omega _{x}=w_{y}-v_{z}$

compute_vort_y.cpp

Compute the y-component of vorticity: $\omega _{y}=u_{z}-w_{x}$

compute_vort_z.cpp

Compute the z-component of vorticity: $\omega _{z}=v_{x}-u_{y}$

dissipation.cpp

Compute the viscous dissipation $\epsilon =2\mu e_{ij}e_{ij}$ where

$e_{ij}={\frac {1}{2}}({\frac {\partial u_{i}}{\partial x_{j}}}+{\frac {\partial u_{j}}{\partial x_{i}}})$ is the strain rate tensor.

enstrophy_density.cpp

Compute the enstrophy density $\Omega ={\frac {1}{2}}(\omega _{x}^{2}+\omega _{y}^{2}+\omega _{z}^{2})$

enstrophy_stretch_production.cpp

Compute the enstrophy production due to vortex tilting/stretching $\sum _{ij}\omega _{i}\omega _{j}e_{ij}=\sum _{ij}\omega _{i}\omega _{j}{\frac {\partial u_{i}}{\partial x_{j}}}$ where,

$e_{ij}={\frac {1}{2}}({\frac {\partial u_{i}}{\partial x_{j}}}+{\frac {\partial u_{j}}{\partial x_{i}}})$ is the strain rate tensor.

cf. The Enstrophy equation.