[Tarak]

Hii,

I am running the case of flow over a square cylinder using LES by ICEM generated grid. But only after a few time steps, the code diverges, with Courant number crossing 4 or 5. I tried to use Upwind Schemes as well as first order time stepping to tackle the problem, but still the issue wasn't solved. At last, I increased the number of PISO non-orthogonal correctors to 5 (since increasing to 3 didn't help), and now its converging, but its painstakingly slow.

Can you please let me know of a way by which I can run OF with tethrahedral grids smoothly.

Thanks,
Tarak

[maddalena]

I am also facing difficulties on getting a 3D tetra mesh converging with OF 1.6.x.
The case is quite standard: some pipes with two fans that put the fluid in motion. The geometry is complex and I have no possibility to use a hexa mesh. The case is steady state and turbulent (simpleFoam); BC are as tested on similar (working) cases.
fvSchemes are:

gradSchemes: faceMDLimited Gauss linear 0.5;
divSchemes: Gauss linearUpwind cellLimited Gauss linear 1;
laplacianSchemes: Gauss linear limited 0.5;

fvSolution:

p
    {
        solver          GAMG;
        tolerance       1e-12;
        relTol          0;
        smoother        GaussSeidel;
        nPreSweeps      0;
        nPostSweeps     2;
        cacheAgglomeration true;
        nCellsInCoarsestLevel 10;
        agglomerator    faceAreaPair;
        mergeLevels     1;
    }

    U - epsilon - k
    {
        solver          smoothSolver;
        smoother        GaussSeidel;
        tolerance       1e-10;
        relTol          0;
    }

SIMPLE
{
    nNonOrthogonalCorrectors 2;
    pRefCell        0;
    pRefValue       0;
}

relaxationFactors  ///really low in the beginning to let the simulation start.
{
    default         0;
    p               0.05;
    U               0.2;
    k               0.2;
    epsilon            0.2;
    nuTilda            0.3;
}

The solution converges for the first 150 steps nicely, but after on it start to get unstable and finally crashes. I really do not know what to do for keep it running...
Suggestions?

mad

[vkrastev]

Originally Posted by maddalena

I am also facing difficulties on getting a 3D tetra mesh converging with OF 1.6.x.
The case is quite standard: some pipes with two fans that put the fluid in motion. The geometry is complex and I have no possibility to use a hexa mesh. The case is steady state and turbulent (simpleFoam); BC are as tested on similar (working) cases.
fvSchemes are:

gradSchemes: faceMDLimited Gauss linear 0.5;
divSchemes: Gauss linearUpwind cellLimited Gauss linear 1;
laplacianSchemes: Gauss linear limited 0.5;

fvSolution:

p
    {
        solver          GAMG;
        tolerance       1e-12;
        relTol          0;
        smoother        GaussSeidel;
        nPreSweeps      0;
        nPostSweeps     2;
        cacheAgglomeration true;
        nCellsInCoarsestLevel 10;
        agglomerator    faceAreaPair;
        mergeLevels     1;
    }

    U - epsilon - k
    {
        solver          smoothSolver;
        smoother        GaussSeidel;
        tolerance       1e-10;
        relTol          0;
    }

SIMPLE
{
    nNonOrthogonalCorrectors 2;
    pRefCell        0;
    pRefValue       0;
}

relaxationFactors  ///really low in the beginning to let the simulation start.
{
    default         0;
    p               0.05;
    U               0.2;
    k               0.2;
    epsilon            0.2;
    nuTilda            0.3;
}

The solution converges for the first 150 steps nicely, but after on it start to get unstable and finally crashes. I really do not know what to do for keep it running...
Suggestions?

mad

Hi Maddalena,
if you are using linearUpwind for all the divSchemes, try instead to set div(phi,U) to Gauss linearUpwindV faceMDLimited Gauss linear 1, and div(phi,k/epsilon) to Gauss upwind. By the way, why so strong underrelaxation factors? Have you tried to leave the URF for p at the 0.3 standard value, for U at 0.7 and to set the URF's for k and epsilon to 0.4-0.5?

Best Regards

Vesselin

[vkrastev]

Originally Posted by Tarak

Hii,

I am running the case of flow over a square cylinder using LES by ICEM generated grid. But only after a few time steps, the code diverges, with Courant number crossing 4 or 5. I tried to use Upwind Schemes as well as first order time stepping to tackle the problem, but still the issue wasn't solved. At last, I increased the number of PISO non-orthogonal correctors to 5 (since increasing to 3 didn't help), and now its converging, but its painstakingly slow.

Can you please let me know of a way by which I can run OF with tethrahedral grids smoothly.

Thanks,
Tarak

Hi Tarak,
two hints: 1) use hexa-meshes for LES calculations; 2) you cannot use first order schemes for LES, because they are too diffusive: you'll have to use at least a limited second order scheme, such as limitedLinear.

Best regards

Vesselin

[vkrastev]

Originally Posted by vkrastev

By the way, why so strong underrelaxation factors? Have you tried to leave the URF for p at the 0.3 standard value, for U at 0.7 and to set the URF's for k and epsilon to 0.4-0.5?

Oops, I haven't noticed the comment in the fvSolution dictionary So
, let me change the question: what are your URF's after the first n-steps in which you leave so low values?

[maddalena]

Originally Posted by vkrastev

what are your URF's after the first n-steps in which you leave so low values?

Unfortunately, I cannot get too far to have the time to increase them! After 120 steps the solution starts to be unstable, thus increase them sounds as not a good idea...

Originally Posted by vkrastev

div(phi,k/epsilon) to Gauss upwind

you want to increase diffusion of turbulence quantities, do not you?

BTW... sounds a bit weird that, whenever I change my fvSchemes & fvSolution, I cannot get too far with my simulation. May there be something more subtle than this? May be due to the 2 fan BC I imposed?


mad

PS: seems like in the last day we are the only active people in the forum

[vkrastev]

Originally Posted by maddalena

you want to increase diffusion of turbulence quantities, do not you?

Yes, It's a matter of trade-off between accuracy and stability: I've noticed in some of my runs that changing only the div(phi,U) scheme to something more accurate than the first order upwind, gives very significant improvements to the results comparing with an all-upwind setting for the divSchemes. On the other hand, till now, with tetra-dominant meshes I was not able to reach satisfactory convergemce patterns for k and epsilon with higher order schemes...So, better a quite accurate and converging result, than a very accurate but no-converging-at-all result...

Originally Posted by maddalena

BTW... sounds a bit weird that, whenever I change my fvSchemes & fvSolution, I cannot get too far with my simulation. May there be something more subtle than this? May be due to the 2 fan BC I imposed?

Yes, the problem could be easily far from the fvSchemes/Solution tuning (what is your checkMesh "opinion" about it?)...unfortunately I have no experience with fan BC's, so I cannot give you any useful advice...

Originally Posted by maddalena

PS: seems like in the last day we are the only active people in the forum

So it seems...

Vesselin

[maddalena]

Originally Posted by vkrastev

better a quite accurate and converging result, than a very accurate but no-converging-at-all result...

Especially if you have some close-by deadlines...

Originally Posted by vkrastev

what is your checkMesh "opinion" about it?

    Overall domain bounding box (-12 -1.804 -7.19731e-17) (14.905 1.804 10.0839)
    Mesh (non-empty, non-wedge) directions (1 1 1)
    Mesh (non-empty) directions (1 1 1)
    Boundary openness (-3.66027e-18 -7.23166e-17 -1.6178e-15) OK.
    Max cell openness = 2.94713e-16 OK.
    Max aspect ratio = 7.21422 OK.
    Minumum face area = 2.38932e-06. Maximum face area = 0.734457.  Face area magnitudes OK.
    Min volume = 2.87124e-09. Max volume = 0.183334.  Total volume = 327.312.  Cell volumes OK.
    Mesh non-orthogonality Max: 65.408 average: 20.4121
    Non-orthogonality check OK.
    Face pyramids OK.
    Max skewness = 0.763815 OK.

checkMesh says that the problem lies somewhere else...
I will give a try to upwind on k and epsilon. Thank you in the meanwhile..

mad

[vkrastev]

Originally Posted by maddalena

Especially if you have some close-by deadlines...

That's for sure!

Originally Posted by maddalena

    Overall domain bounding box (-12 -1.804 -7.19731e-17) (14.905 1.804 10.0839)
    Mesh (non-empty, non-wedge) directions (1 1 1)
    Mesh (non-empty) directions (1 1 1)
    Boundary openness (-3.66027e-18 -7.23166e-17 -1.6178e-15) OK.
    Max cell openness = 2.94713e-16 OK.
    Max aspect ratio = 7.21422 OK.
    Minumum face area = 2.38932e-06. Maximum face area = 0.734457.  Face area magnitudes OK.
    Min volume = 2.87124e-09. Max volume = 0.183334.  Total volume = 327.312.  Cell volumes OK.
    Mesh non-orthogonality Max: 65.408 average: 20.4121
    Non-orthogonality check OK.
    Face pyramids OK.
    Max skewness = 0.763815 OK.

checkMesh says that the problem lies somewhere else...

I agree with "him"

Originally Posted by maddalena

I will give a try to upwind on k and epsilon. Thank you in the meanwhile..

Good luck!

V

[santiagomarquezd]

Hi Tarak, Did you try renumberMesh after importing the mesh to FOAM? I've found that meshes imported from Fluent are not very efficient in the numeration and this impacts strongly in iterative solvers.

Regards.

[maddalena]

Originally Posted by santiagomarquezd

Hi Tarak, Did you try renumberMesh after importing the mesh to FOAM? I've found that meshes imported from Fluent are not very efficient in the numeration and this impacts strongly in iterative solvers.

NB: this will not work if you have cyclic patch! When renumbering mesh, the cyclics will be lost!

mad

[maddalena]

Hi everybody,
seems like all my problems are still there, on a case that is similar to what posted above. I wrote on a different thread: http://www.cfd-online.com/Forums/ope...tml#post293910
Maybe you can take a look to it and suggest something more on the subject...

mad

[tcarrigan]

This is the thread I meant to post on. I started clicking your links and posted to the wrong one . I'll repeat my post:

----------------------------------------------------------------
Just curious, have you tried using leastSquares for the gradScheme?

I did some 2D calculations for a NACA airfoil using both structured and unstructured grids. I too suffered convergence issues when running the calculation for the unstructured case. However, switching the gradScheme to a cellLimited leastSquares happened to solve the problem.

Let me know if this works.
----------------------------------------------------------------

And to answer your question, there may not be an advantage to using leastSquares, it's just an observation made from several simulations I've run for unstructured grids. Most of the time I have issues obtaining convergence when using tet grids and the Gauss linear gradScheme. Just and observation, was curious to see if it had any effect on your calculation.

[maddalena]

Originally Posted by tcarrigan

And to answer your question, there may not be an advantage to using leastSquares, it's just an observation made from several simulations I've run for unstructured grids. Most of the time I have issues obtaining convergence when using tet grids and the Gauss linear gradScheme. Just and observation, was curious to see if it had any effect on your calculation.

Ok, I'll give it a try...

mad

[Tarak]

Hii,

I would recommend using Gamma Differencing Scheme with appropriate blending factor for obtaining good convergence/stability for tetrahedral grids. This is how I managed to get the code run.

Thanks,
Tarak

[maddalena]

Originally Posted by tcarrigan

Just curious, have you tried using leastSquares for the gradScheme?

I did some 2D calculations for a NACA airfoil using both structured and unstructured grids. I too suffered convergence issues when running the calculation for the unstructured case. However, switching the gradScheme to a cellLimited leastSquares happened to solve the problem.

So, here an answer to your observations: http://www.cfd-online.com/Forums/ope...tml#post294338

mad

[hjasak]

In 1.6-ext, we've got special discretisation scheme for tet meshes: use it on convection in the momentum equation and on interpolate(U). The scheme is called reconCentral, and the settings for a tet mesh look like this:

gradSchemes
{
default cellLimited leastSquares 1.0;
}

divSchemes
{
div(phi,U) Gauss reconCentral cellLimited leastSquares 1.0;
//...
}

laplacianSchemes
{
default Gauss linear limited 0.5;
}

interpolationSchemes
{
default linear;
interpolate(HbyA) linear;
interpolate(U) reconCentral phi cellLimited leastSquares 1.0;
}

Enjoy,

Hrv

[Mehrez]

Hi Foamers,

I am facing the same problem using multiphaseEulerFoam (OpenFoam 2.3.0) on Ubuntu 16.
In fact, my geometry is discretized by a tetrahedral mesh and when I run checkMesh, I get the Mesh OK message.
My problem is that the solver craches because the simulations become quickly unstable and the Courant number diverges. Using adaptable time steps didn't help.
I have tried to modify my fvSchemes as described above but this didn't help. Here is my fvSchemes :


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

ddtSchemes
{
default Euler;
}

gradSchemes
{
default cellLimited leastSquares 1.0;
}

divSchemes
{
default none;

"div\(phi,alpha.*\)" Gauss reconCentral cellLimited leastSquares 1.0;
"div\(phir,alpha.*,alpha.*\)" Gauss vanLeer;

"div\(phiAlpha.*,U.*\)" Gauss limitedLinearV 1;
div(Rc) Gauss linear;
"div\(phi.*,U.*\)" Gauss limitedLinearV 1;
}

laplacianSchemes
{
default Gauss linear corrected;
}

interpolationSchemes
{
default linear;
interpolate(HbyA) linear;
interpolate(U) reconCentral phi cellLimited leastSquares 1.0;
}

snGradSchemes
{
default corrected;
}

fluxRequired
{
default no;
p ;
pcorr ;
}

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

I saw different threads exposing this problem (unstable OpenFoam solvers with tetrahedral meshes) but didn't find a solution. Your help will be appreciated.
Thank you.

Mhrz

Originally Posted by hjasak

In 1.6-ext, we've got special discretisation scheme for tet meshes: use it on convection in the momentum equation and on interpolate(U). The scheme is called reconCentral, and the settings for a tet mesh look like this:

gradSchemes
{
default cellLimited leastSquares 1.0;
}

divSchemes
{
div(phi,U) Gauss reconCentral cellLimited leastSquares 1.0;
//...
}

laplacianSchemes
{
default Gauss linear limited 0.5;
}

interpolationSchemes
{
default linear;
interpolate(HbyA) linear;
interpolate(U) reconCentral phi cellLimited leastSquares 1.0;
}

Enjoy,

Hrv

[tcarrigan]

Mehrez,

I run hybrid prism, pyramid, tet meshes using OpenFOAM pretty frequently and have come up with a standard fvSchemes dictionary that I start with every time. I also have a fvSolution dictionary that I've created that works equally well coupled with the changes I made to fvSchemes. Both are listed below for your reference.





fvSchemes

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

ddtSchemes
{
default steadyState;
}
gradSchemes
{
default faceMDLimited Gauss linear 0.333;
}

divSchemes
{
default none;
div(phi,U) bounded Gauss linearUpwindV grad(U);
div(phi,U) bounded Gauss upwind;
div(phi,nuTilda) bounded Gauss upwind;
div((nuEff*dev(T(grad(U))))) Gauss linear;
}

laplacianSchemes
{
default Gauss linear limited 0.333;
}

interpolationSchemes
{
default linear;
}

snGradSchemes
{
default limited 0.333;
}

fluxRequired
{
default no;
pcorr ;
p ;
}


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





fvSolution

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

solvers
{
p
{
solver GAMG;
tolerance 1e-6;
relTol 0.01;
smoother FDIC;
nPreSweeps 1;
nPostSweeps 2;
nFinestSweeps 2;
cacheAgglomeration on;
agglomerator faceAreaPair;
nCellsInCoarsestLevel 1000;
mergeLevels 1;
}

U
{
type coupled;
solver PBiCCCG;
preconditioner DILU;
tolerance (1e-6 1e-6 1e-6);
relTol (0.01 0.01 0.01);
}

"(k|omega|nuTilda)"
{
solver PBiCG;
preconditioner DILU;
tolerance 1e-6;
relTol 0.01;
}

}

SIMPLE
{
nNonOrthogonalCorrectors 3;
}

potentialFlow
{
nNonOrthogonalCorrectors 10;
}

relaxationFactors
{
fields
{
p 0.2;
}
equations
{
U 0.5;
k 0.5;
omega 0.5;
nuTilda 0.5;
}
}

cache
{
grad(U);
}

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

[Mehrez]

Dear tcarrigan,

Thank you for your answer and for sharing your files.
I have edited my fvSchemes and fvSolution files as follows:

fvSchemes
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
ddtSchemes
{
default steadyState;//Euler;
}

gradSchemes
{
default faceMDLimited Gauss linear 0.333;//Gauss linear;
}

divSchemes
{
default none;

"div\(phi,alpha.*\)" bounded Gauss upwind;//Gauss vanLeer;
"div\(phir,alpha.*,alpha.*\)" bounded Gauss upwind;//Gauss vanLeer;

"div\(phiAlpha.*,U.*\)" bounded Gauss linearUpwindV grad(U);//Gauss limitedLinearV 1;
"div\(phiAlpha.*,U.*\)" bounded Gauss upwind;//Gauss limitedLinearV 1;
div(Rc) Gauss linear;
"div\(phi.*,U.*\)" bounded Gauss linearUpwindV grad(U);//Gauss limitedLinearV 1;
"div\(phi.*,U.*\)" bounded Gauss upwind;//Gauss limitedLinearV 1;
}

laplacianSchemes
{
default Gauss linear limited 0.333;//Gauss linear corrected;
}

interpolationSchemes
{
default linear;
}

snGradSchemes
{
default limited 0.333;//corrected;
}

fluxRequired
{
default no;
p ;
pcorr ;
}
// ************************************************** *********************** //


fvSolution
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
solvers
{
alpha
{
nAlphaSubCycles 2;
}

p
{
solver GAMG;
tolerance 1e-6;
relTol 0.01;
smoother FDIC;
nPreSweeps 1;
nPostSweeps 2;
nFinestSweeps 2;
cacheAgglomeration on;
agglomerator faceAreaPair;
nCellsInCoarsestLevel 1000;
mergeLevels 1;
}

pFinal
{
$p;
tolerance 1e-14;
relTol 0;
}

pcorr
{
$p;
tolerance 1e-14;
relTol 0;
}

"U.*"
{
type coupled;
solver PBiCCCG;
preconditioner DILU;
tolerance (1e-6 1e-6 1e-6);
relTol (0.01 0.01 0.01);
}

}

PIMPLE
{
nOuterCorrectors 1;
nCorrectors 2;
nNonOrthogonalCorrectors 3;
}

relaxationFactors
{
fields
{
p 0.2;
}
equations
{
U 0.5;
}
}

cache
{
grad(U);
}
// ************************************************** *********************** //


I have tried to run a simulation using these files but it goes quickly unstable. Here is a copy of my log file where the deltaT is going down:

Courant Number mean: 0.0021156010062783 max: 0.049408733271849
deltaT = 7.6566763768132e-19
Time = 6.18252270356846e-06

MULES: Solving for alpha.gas
gas volume fraction, min, max = 0.30019748795334 0.23897356139625 0.41705604707805
MULES: Solving for alpha.liquid
liquid volume fraction, min, max = 0.70087632594298 0.3 0.86565308392893
Phase-sum volume fraction, min, max = 1.0010738138963 0.59999999999999 1.1699030447477
MULES: Solving for alpha.gas
gas volume fraction, min, max = 0.30019798328777 0.2386702783917 0.41704447896555
MULES: Solving for alpha.liquid
liquid volume fraction, min, max = 0.70087752844775 0.3 0.86618197518773
Phase-sum volume fraction, min, max = 1.0010755117355 0.59999999999999 1.1703801282208
GAMG: Solving for p, Initial residual = 0.23419353824923, Final residual = 0.0023325036575349, No Iterations 2
GAMG: Solving for p, Initial residual = 0.052855980175197, Final residual = 0.0003817033201444, No Iterations 2
GAMG: Solving for p, Initial residual = 0.014025336427512, Final residual = 8.8046943721164e-05, No Iterations 2
GAMG: Solving for p, Initial residual = 0.0048650690815121, Final residual = 4.0782603314908e-05, No Iterations 2
time step continuity errors : sum local = 9.3649763390591e-09, global = 7.9493516767848e-10, cumulative = 8.2571867820764e-08
GAMG: Solving for p, Initial residual = 0.20867553926075, Final residual = 0.00073649732328669, No Iterations 3
GAMG: Solving for p, Initial residual = 0.050332050313325, Final residual = 0.00029712401985861, No Iterations 2
GAMG: Solving for p, Initial residual = 0.013334743089887, Final residual = 8.741559351822e-05, No Iterations 2
GAMG: Solving for p, Initial residual = 0.0048312307962256, Final residual = 7.7390130904333e-15, No Iterations 26
time step continuity errors : sum local = 1.8252376105909e-18, global = -3.8452675495394e-20, cumulative = 8.2571867820725e-08
ExecutionTime = 424.27 s


Expression Re_l : average=2.7438854159252e+20 min=6.7983965132898e+18 max=4.2263684536994e+21
Expression Re_g : average=1.2645680808052e+20 min=3.1815443817709e+18 max=2.7037109361158e+21
Expression rho_avg : average=489.19817329068 min=447.12039199452 max=600.27512277584
Expression head : average=9022.1840703653
Expression alpha_g on outlet_gas: average=0.30283498588815 min=0.2386702783917 max=0.36663713347938
Expression alpha_g on outlet_liquid: average=0.3006588178025 min=0.2974219141048 max=0.30429082494955
Expression alpha_g on inlet: average=0.3 min=0.3 max=0.3
Expression Q_all on outlet_gas: sum=-24210659648582 min=-4793324780552.7 max=247183322085.2
Expression Q_all on outlet_liquid: sum=17816623940489 min=193138733172.99 max=521444731230.64
Expression Q_all on inlet: sum=6394035708093.2 min=-1215092696559.3 max=597712717234.6
Expression Q_l on outlet_gas: sum=-17262638871291 min=-3156101143802.7 max=235004594945.67
Expression Q_l on outlet_liquid: sum=12558463445253 min=137931428364.81 max=379354434355.94
Expression Q_l on inlet: sum=1580810441283.7 min=-596686394611.29 max=299229574724.19
Expression Q_g on outlet_gas: sum=-6948020777290.6 min=-1637223636750 max=794642130398.71
Expression Q_g on outlet_liquid: sum=5258160495235.8 min=51791175935.031 max=161116494634.18
Expression Q_g on inlet: sum=4813225266809.5 min=-618406301947.99 max=305405324047.52
Expression p_avg on outlet_gas: average=0 min=0 max=0
Expression p_avg on outlet_liquid: average=7000 min=7000 max=7000
Expression p_avg on inlet: average=9003.4061139999 min=9003.4061139999 max=9003.4061139999
Expression recycle on outlet_liquid: sum=2.7864442355142
Expression p_inlet on inlet: average=9003.4061139999
Expression p_inlet_proporionalTerm on inlet: average=0



Thank you for your help,

Mhrz