Hello,

for compressible viscous computations I intend to use/modify the sonicFoam solver, but a few things I don't understand:

- the viscous term is given as a laplacian (the case forwardStep is discussed in sec. 3.3 of the P-guide) and I can not see why. Only for incompressible flows the stress reduces to laplacian(mu, U). Is there a reason for it?

I tried to incorporate the full stress tensor, but this is only possible in explicit form because the function grad(U) is only available in the fvc class and now I model it by:

volTensorField gradU = fvc::grad(U);

mu*(( gradU + gradU.T() ) - 2.0/3.0*fvc::div(U)*tensor::I);

but the solver gets less stable, especially for larger mu.
Is there a reason why I can't use grad(U) in implicit form (why it is not available in the fvm class)? Is there another way to do this implicitly?

-Furthermore I see from the results that the choice of outflow boundary conditions (zerogradient on p, T, U) is not ideal since it generates an artificial boundary layer on the outflow, see also figure 3.8 in the P-guide). For my work a correct boundary treatment is essential.
Is there a specific reason for this choice of boundary conditions? It seems to mee that it is not the usual case for Navier-Stokes outflow conditions.

regards, Jelmer

[hjasak]

Hi,

1) Laplacian
The viscous term is a laplacian for numerical reasons. If you have a variable density flow (like you do) or wish to incorporate some external stress tensor, it is advisable to split it into the laplacian + correction terms. This is because the laplacian is critical for stability, as your fully explicit implementation shows. You will need to read up on discretisation to get more details on precisely why your algorithm fails.

2) gradient
It's quite easy to implicitly discretise the gradient, but that's only the beginning of your problems. As you know, the diagonal contribution from the discretised gradient is zero (!) and this is only useful in cross-equation coupling, eg. thr grad p term in the momentum equation. BTW, how would you then discretise a divergence from the discretised gradient and what happens to your computational molecule? As for another way, see above: this comes from controlled assembly of the matrix for consistency, stability and iterative solvers.

3) outflow boundary
There is no "artificial boundary layer" at the outlet: you've messed up your boundary conditions (I bet you've got a subsonic outlet and the b.c.-s you list are supersonic.

4) usual case for Navier-Stokes
What do you mean? The b.c.-s depend on the problem and physics and there's no such thing as "usual"

Hrv

Thank you for the reply. I am not familiar with the laplacian + correction term for the Navier-Stokes equations, but I will have a look in literature.

Cosidering the outflow: it is supersonic (velocity above 1) and I copied the sonicFoam forwardStep case 1 to 1 into my directory, including the boundary conditions. However, I was wrong about the zero gradient (sorry, mistakes happen), but from the files I have:
pressure:
{
type pressureTransmissive;
pInf 1;
lInf 3;
}

temperature:
{
type inletOutlet;
inletValue uniform 1;
value uniform 1;
}

velocity:
{
type inletOutlet;
inletValue uniform (3 0 0);
value uniform (3 0 0);
}

Considering the pressureTransmissive condition:
the Usersguide on page 143 state that only pInf is required but in the file it also gives a lInf variable. The code does not run without it, but what is it's meaning and purpose?

I assume the inletOutlet condition on the velocity only 'turns on' with a value of 3 when there is inflow.

The results (with the original sonicFoam solver) I get are equal to the one from fig 3.8 in which also the 'artificial boundary layer' on the outlet is visible in the pressure field. Playing the intermediate results sequentially also shows the layer exists at all timesteps. In the picture below it shows the layer for the pressure. A similar phenomenon happens for the temperature.

Since the outflow is supersonic the conditions on temperature and velocity are 'not in use' so only the pressure is of importance. Could it be the discretisation of the pressureTransmissive condition that causes the effect?



regards, Jelmer

Dear prof Jasak,

I really can't find the cause of this layer. I checked other (inviscid) compressible solvers with their corresponding tutorials and there's no such problem with them (although the rhoSonicFoam shows a terrible decoupling in the solution with the tutorial settings). In spite the supersonic outflow where NO boundary condition should be given, there's zeroGradient for the p, T and U in the tutorial bc files. I don't see the reason why. In previous Euler solvers I wrote myself I never put any bc on supersonic outflow, all characteristics are outgoing...

I also tried upwind schemes for the sonicFoam forwardStep case, because in other cases the discretisation may use the boundary value for p given on the outlet (as for the pressureTransmissive condition), which might have caused the artificial layer. But with upwind schemes the same layer occurs.

Furthermore I had a look at the PISO correction loop since it contains a laplacian on the pressure, from which I expected it to require a boundary condition for the elliptic term. But neither does it have influence on the layer.

So far, I am out of ideas to solve this issue.
What is your opinion about this, you have any suggestions?

regards, Jelmer