[Oliver Meng]

Hi lovely CFD People,


i'm running a compressible simulation with rhoSimpleFoam but it does not converge. Even worse, the residual of p is really really weird.


Anyone happens to know the reason?


I basically started with Gauss Upwind to get the stability (but honestly, the p is not stable at all) and switched to Gauss Upwind linear at 100. Iteration.


Thx for any answers in advance!

[JulioPieri]

Hey! This can be caused by a lot of reasons. Poor mesh would be my first guess, either from low quality elements or by overly coarse mesh. You may also have inconsistent BCs. Lastly, your problem might be so complex that the solver can't find a solution. Breaking it down into chunks of the domain, or gradually accelerate it (actual velocity or thermal stuff) may give you some insights.

Try solving a smaller problem with similar set of properties and BCs to easily debug and test for solutions.

[Oliver Meng]

Hi Julio,

thanks for your reply! I found out that it was a plotting mistake from me in my gnuplot script. It shoudl look like this one in the attachment.

But it is still not satisfying because it didnt converge to my criterien.
I changed the schemes step by step during the simulation. Basically from loose to tight. Do you have an idea why the residuals didnt go lower?

/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  v2106                                 |
|   \\  /    A nd           | Website:  www.openfoam.com                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    object      fvSchemes;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

ddtSchemes
{
    default         steadyState;
}

gradSchemes
{
    // default         cellMDLimited Gauss linear 0; // super loose
    // grad(U)         cellMDLimited Gauss linear 1; // super loose
     
    // default         leastSquare; // loose
    // default         cellMDLimited Gauss linear 1; // loose
    default         Gauss linear; // tight
    grad(U)         cellLimited<cubic> 1.5 Gauss linear 1; // in Einsatz ab 1100 Iter.
}

divSchemes
{
    default         none;

    // div(phi,U)      bounded Gauss upwind; // loose
    div(phi,U)      Gauss linearUpwind grad(U); // tight
    
    // energy          bounded Gauss upwind; // loose
    energy          Gauss linearUpwind default; // tight

    div(phi,e)      $energy;
    div(phi,h)      $energy;
    div(phi,Ekp)    $energy;

    turbulence      Gauss upwind;
    // turbulence      bounded Gauss upwind;
    // turbulence      Gauss linearUpwind default;

    div(phi,k)      $turbulence;
    div(phi,omega)  $turbulence;

    div(((rho*nuEff)*dev2(T(grad(U))))) Gauss linear;
}

laplacianSchemes
{
    default         Gauss linear corrected;
    // default         Gauss linear corrected; 
    // default         Gauss linear limited 0.33; 
    // default         Gauss linear limited 0.777; 
    // default         Gauss linear orthogonal; 
    // default         Gauss linear uncorrected;
}

interpolationSchemes
{
    default         linear;
}

snGradSchemes
{
   default         corrected;
   // default         limited 0.777; // tight
   // default        limited 0.5; // loose
   // default limited 0.33; // super loose
}

wallDist
{
method meshWave;
}



// ************************************************************************* //

/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  v2106                                 |
|   \\  /    A nd           | Website:  www.openfoam.com                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    object      fvSolution;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

solvers
{
    p
    {
        solver          GAMG;
        tolerance       1e-08; // loose=1e-5, tight=1e-8
        relTol          0.01;  // super loose=0.5, tight=0.01
        smoother        GaussSeidel;
        nCellsInCoarsestLevel 100; //100
        minIter         3;
 
    }

     "(U|k|omega|e)"
    {
        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-8; // loose=1e-6, tight=1e-8
        relTol          0.0;
    }
}

SIMPLE // semi-implicit method for pressure-linked euqations
{
    consistent          yes; // SIMPLEC formulation for more robust and faster convergence
    nNonOrthogonalCorrectors 2; // loose=3, tight=2
    pMinFactor          0.4; //0.4
    pMaxFactor          1.5; //1.5

    residualControl
    {
        p               1e-8; //eventuell noch kleiner
        U               1e-7;
        "(k|omega|e)"   1e-7;
    }
}

relaxationFactors
{
    fields
    {
        p               0.4; // loose 0.1
        rho             0.4; //loose 0.01
    }
    equations
    {
        U               0.7;
        "(k|omega)"     0.7;
        e               0.7;
    }
}


// ************************************************************************* //

Originally Posted by JulioPieri

Hey! This can be caused by a lot of reasons. Poor mesh would be my first guess, either from low quality elements or overly coarse mesh. You may also have inconsistent BCs. Lastly, your problem might be so complex that the solver can't find a solution. Breaking it down into chunks of the domain, or gradually accelerate it (actual velocity or thermal stuff) may give you some insights.

Try solving a smaller problem with similar set of properties and BCs to easily debug and test for solutions.

[JulioPieri]

Hmm, have you tried changing the relaxation factors back to 1 (gradually)? They affect the solution as they can make residuals appear to have stabilized (even if they are oscillating like yours).

The residuals are not the best way to check convergence. It's important to add some physical quantity as well, like force, heat flux, etc. If those have not converged, then your problem really didn't converge. However, if they are stable (and you don't have super low relaxation factors), it might be a stronger indicative you've reached a solution than assessing the residuals.

This wiki might also be helpful, though I think you should first tweek your relaxation.

https://www.cfd-online.com/Wiki/Ansy...gence_criteria

[Oliver Meng]

Hi Julio,

thanks again for you idea. I am still doing some optimising with my mesh. And i will try a new simulation tonight with your advice of changing the relaxation step by step to 1. (My new work station with 16 cores just arrived this afternoon)

And this is my force profile, it looks quite stable.

Originally Posted by JulioPieri

Hmm, have you tried changing the relaxation factors back to 1 (gradually)? They affect the solution as they can make residuals appear to have stabilized (even if they are oscillating like yours).

The residuals are not the best way to check convergence. It's important to add some physical quantity as well, like force, heat flux, etc. If those have not converged, then your problem really didn't converge. However, if they are stable (and you don't have super low relaxation factors), it might be a stronger indicative you've reached a solution than assessing the residuals.

This wiki might also be helpful, though I think you should first tweek your relaxation.

https://www.cfd-online.com/Wiki/Ansy...gence_criteria

[JulioPieri]

Congrats on the workstation!

As for the plot, it seems converged to me. I've seem some rare cases that changing the relax factor would change that profile, though that shouldn't happen. Run the simulation with increased relax factor and post back!

[mAlletto]

Hello, I found that running rhoPimpleFoam with a local Euler sometimes helps to get a steady solution. If there are recirculation regions then it is very hard to get a converged steady solution