[seb_210]

Hello everyone,

I am coming back after fixing the problem I used to have here : http://www.cfd-online.com/Forums/ope...ation-2-a.html

Now my simulation converges well with first order schemes. The thing is that now, I would like to use Second Order schemes (linear or QUICK) in the divergence schemes to increase the model's precision. I use the result of the first order simulation as initial condition for the second order but it crashes sooner or later with the use of such schemes.

Pressure and velocity explode and reach unrealistic values.

Is there anything doable to avoid this numerical behaviour ?

Thanks in advance for your help,

Seb


PS: this problem is referenced in the following thread as well and nobody answered

http://www.cfd-online.com/Forums/ope...r-schemes.html

[seb_210]

Still nobody ?

Seb

[GDTech]

Hi,

If your Mach number is above 0.75 in some places, I suggest you to activate transonic flag in fvSolution :

SIMPLE
{
transonic true;
.
.
.
}

From my experience, you need a really good mesh (essentially hex) to be able to use second order schemes with rhoSimpleFoam and high Mach number.

Good luck.

[seb_210]

Hello

Thanks for the tip. The thing is when I use "transonic true;", rhoSimpleFoam rises :

Starting time loop

Time = 3001

smoothSolver:  Solving for Ux, Initial residual = 2.32034e-06, Final residual = 1.84847e-08, No Iterations 2
smoothSolver:  Solving for Uy, Initial residual = 0.00324492, Final residual = 2.34545e-05, No Iterations 2
smoothSolver:  Solving for Uz, Initial residual = 0.000158151, Final residual = 1.15094e-06, No Iterations 2
DILUPBiCG:  Solving for e, Initial residual = 2.32848e-05, Final residual = 3.35775e-07, No Iterations 1


--> FOAM FATAL IO ERROR: 
Unknown asymmetric matrix smoother DICGaussSeidel

Valid asymmetric matrix smoothers are :

5
(
DILU
DILUGaussSeidel
GaussSeidel
nonBlockingGaussSeidel
symGaussSeidel
)


file: /scratch/*****.rep/system/fvSolution.solvers.p from line 22 to line 29.

    From function lduMatrix::smoother::New
    in file matrices/lduMatrix/lduMatrix/lduMatrixSmoother.C at line 126.

FOAM exiting

I've never experienced such an error. Do you know what might be going on here ?

[GDTech]

yes, pressure equation is different now and you have to use an asymmetric smoother.
Replace DICGaussSeidel by DILUGaussSeidel for pressure equation in fvSolution.

[seb_210]

Hello Laurent,

I followed your advice and until now, the solver seems to run. But it takes a incredibly long time to go through "p" equation, iterating it 1000 times for every new iteration. Plus, the residuals seem to oscillate a lot.

I tried lowering the tolerances for "p" solver in the fvSolution dictionnary but sooner or later these 1000 iterations end up coming back.

Do you have an idea of what is happening here ?

Thanks for your help

Seb

[GDTech]

Hi Seb,

p equation going to the maximum default number of iterations (1000) without convergence means that something is wrong in your setup (i.e. boundary conditions or mesh).

Residual oscillations is another story and probably comes from the unsteadiness in your flow or pressure wave reflections.

Maybe I would be able to help you more if you post some pictures of your mesh and flow field.

Regards.

[seb_210]

Hi Laurent,

Here are the pictures you asked for.

Here is the result of checkMesh:

Checking geometry...
    Overall domain bounding box (-5.16 -5.16 -5.16) (5.16 5.16 5.16)
    Mesh (non-empty, non-wedge) directions (1 1 1)
    Mesh (non-empty) directions (1 1 1)
    Boundary openness (-2.76836e-18 4.16889e-17 -1.94542e-17) OK.
    Max cell openness = 5.28473e-15 OK.
    Max aspect ratio = 106.421 OK.
    Minimum face area = 1.69335e-09. Maximum face area = 4.26113.  Face area magnitudes OK.
    Min volume = 2.9706e-12. Max volume = 8.79548.  Total volume = 1099.1.  Cell volumes OK.
    Mesh non-orthogonality Max: 73.5736 average: 9.737
   *Number of severely non-orthogonal faces: 160.
    Non-orthogonality check OK.
  <<Writing 160 non-orthogonal faces to set nonOrthoFaces
    Face pyramids OK.
    Max skewness = 3.27693 OK.
    Coupled point location match (average 0) OK.

Mesh OK.

End

Concerning the boundary/initial conditions, here is what I have for Spalart-Allmaras model:

- p field

internalField uniform 1.01325e+05

boundaryField
domain limits: freeStreamPressure
object surface: zeroGradient

- U field

internalField uniform ( 172 0 15) // Mach0.5, incidence 5°

boundaryField
domain limits: freestream
freestreamValue uniform ( 172 0 15)
object surface: fixedValue
value uniform (0 0 0) // no slip

- nuTilda field:

internalField uniform 0.001 

boundaryField
domain limits: freestream
freestreamValue $internalField
object surface: fixedValue
value uniform 0

- mut field:

internalField uniform 0

boundaryField
domain limits: calculated
value $internalField
object surface: mutUWallFunction
value $internalField

Any ideas ?

Thanks a lot

Seb

[GDTech]

Hi,

your settings seem to be ok.
Be sure your freestream boundaries are far away from your object (at least 20x object size).
You should refine your mesh in the wake since size transition occurs too close to the object.
If it does not solve your problem, try with more relaxation.

Good luck !

[seb_210]

Okay, so, I refined the wake a little more, and the boundaries are far enough (more than 20 times the object length).

Unfortunately, the second order schemes still doesn't manage to converge

I use to do :
All upwind and first order schemes with low relaxation factors (0.1, sometimes lower) for the first 500 iterations.
Then I rise the relaxation factors (up to 0.5 - 0.6, step by step if necessary) for the 500 or more next iterations.
Then I finally switch to 2nd order schemes (linear/QUICK/linearUpwind for fivSchemes) with lower relaxation factor (0.01) and run for 1000 iterations at least. But even at this point, the solver crashes ...

Is there something I'm doing wrong in this method ? Should I use a coarser mesh first ? Should I lower the relaxation factors even more ? I'm running out of solutions ^^

Plus, I feel like the output Aerodynamic Coefficients don't seem right. With positive incidence, I get a negative lift coefficient and wrong values in general. Doesn't sound right to me... Might it be related with the convergence problem ?

Here is how I call them:

functions
 {
  coefficientsAero
   {
 // type               forces;
    type               forceCoeffs;
    functionObjectLibs ("libforces.so");
 // outputControl      outputTime;
    outputControl      timeStep;
    log                yes;
    patches            ( "geometrie_.*" ); // I have many patches
    rhoInf               1.17659 ;
    CofR               ( 0 0 0 ); 
    lRef                 0.037 ;       
    Aref                 0.0010752 ;     
    dragDir            ( 1 0 0 ) ;
    liftDir            ( 0 0 1 ) ;
    pitchAxis          ( 0 1 0 ) ;
    magUInf              173.612 ; 
   }
 }

PS: even with first order schemes where I reach a good convergence, the coefficient appear to be this way and they do not evolve to physically correct values (i.e. positive lift)

If you have any clue about this, I'm interested !

Thanks again for your help and patience

Seb

[GDTech]

when I can't make rhoSimpleFoam running fine, I switch to sonicFoam !
sonicFoam is often slower (unstready solver -> limited by courant number) but it is really more stable than rhoSimpleFoam for transonic flows.

[seb_210]

Hello Laurent,

Thanks for your help again, I'll see what I can do with this

Seb