doubt about axisymmetry source for SST

[bigfooted]

Starting point is the SST model, specifically the production term, here are some equations:

(#1551)
The final production term is:

I will focus on the axisymmetry contribution of this term only. For axisymmetry, the strain does not contribute because it is invariant under coordinate transformations.
[EDIT] The tensor is invariant under rotations, the magnitude in cylindrical coordinates with axisymmetry assumption has an additional component $S_{\theta \theta} = \frac{v}{r}$

So the terms that are left are the divergence terms.

Divergence in cylindrical coordinates is $$\nabla\cdot u = (\frac{\partial u}{\partial x} + \frac{\partial v}{\partial r}) + \frac{v}{r}$$

The second divergence term written in cylindrical coordinates leads to an additional term $$-\frac{2}{3}\rho k\frac{v}{r}$$.

The squared divergence term written in cylindrical coordinates leads to an additional axisymmetric term $$-\frac{2}{3}\mu_t((\frac{v}{r})^2 +2\frac{v}{r}(\frac{\partial u}{\partial x} + \frac{\partial v}{\partial r}))$$

with $$\mu_t = \frac{\rho k}{\zeta}$$, we can write it as:

$$ -\frac{2}{3}\frac{\rho v k}{r \zeta}((\frac{v}{r}) +2(\frac{\partial u}{\partial x} + \frac{\partial v}{\partial r}))$$

The production term when axisymmetry is active for the SST-model is implemented as:

const su2double pk_axi = max(
        0.0, 2.0 / 3.0 * rhov * k * ((2.0 * yinv * V_i[idx.Velocity() + 1] - PrimVar_Grad_i[idx.Velocity()+1][1] - PrimVar_Grad_i[idx.Velocity()][0]) / zeta - 1.0));

This does not match with the derivation. The final "-1" at the end is the contribution of the second divergence term, this term matches with the derivation. The other part does not match. It looks like the factor of 2.0 should be -1.0, and a factor 2.0 should be in front of the gradients.

Could somebody recheck my derivation?

[bigfooted]

By the way, I'd like to add that the SST_xxxx_m versions neglect the divergence terms, so these terms should be 0 for the SST_xxxx_m models!

[pcarruscag]

What are the steps necessary for convection and diffusion?

[bigfooted]

If the goal is to get rid of the inline function entirely, we could simply directly add the remaining source terms to ComputeResidual, I mean after all we are already calling ResidualAxisymmetric there anyway (and nowhere else for that matter). So would you say it's prefered to remove the inline function completely and just have an if (axisymmetric) {} block directly in ComputeResidual? It's only 10 lines an the inline is not used elsewhere so maybe this is the most clear?

[pcarruscag]

No, the goal is to implement something correctly, that doesn't leave someone checking the terms every now and again.
Therefore, if for the source term all we need is to add v / r to the "cartesian divergence" that we already compute in ComputeResidual we simply add one line there. However, if it is not that simple, because of convection and diffusion for example, we keep the function.
So, what do we need for the convection and diffusion terms?

(We need more, but smaller, functions, not bigger ones)

[bigfooted]

What we do now is we split cylindrical terms like this: $$\frac{1}{r}\frac{\partial r \phi}{\partial x} = \frac{\phi}{r} + \frac{\partial \phi}{\partial x}$$ and then we can treat Cartesian and cylindrical setups in the same way, except for the occurrence of the additional source term. We could deal with the convection and diffusion terms directly, in their respective schemes, but this looks to me more difficult than working with an additional source term. The convection and diffusion source terms are not difficult so it is easy to check their implementation.

So then my proposal is to only move the divergence source terms to the production routine, and leave the rest of the axisymmetry implementation the same.

I'll create a PR now and then we can have a look at it.

[bigfooted]

I've created a WIP PR: #1784