CFD Online Logo CFD Online URL
www.cfd-online.com
[Sponsors]
Home > Forums > Software User Forums > OpenFOAM

second order schemes

Register Blogs Community New Posts Updated Threads Search

Like Tree58Likes

Reply
 
LinkBack Thread Tools Search this Thread Display Modes
Old   October 19, 2010, 13:45
Default
  #21
Senior Member
 
Alberto Passalacqua
Join Date: Mar 2009
Location: Ames, Iowa, United States
Posts: 1,912
Rep Power: 36
alberto will become famous soon enoughalberto will become famous soon enough
Quote:
Originally Posted by vkrastev View Post
div(phi, U) Gauss linearUpwindV cellMDLimited Gauss linear 1;
div(phi, k) Gauss linearUpwindV cellMDLimited Gauss linear 1;
div(phi, epsilon) Gauss linearUpwindV cellMDLimited Gauss linear 1;
linearUpwind not linearUpwindV for k and eps.

Quote:
PS-I'm sorry if I'm repetitive, but I'll be very glad if you can direct me to some more informations about the theoretical basis of the linearUpwind scheme in its bounded formulation
Check the implementation in the code. I am not aware of any specific references (it is pretty standard stuff).

Best,
dickcruz likes this.
__________________
Alberto Passalacqua

GeekoCFD - A free distribution based on openSUSE 64 bit with CFD tools, including OpenFOAM. Available as in both physical and virtual formats (current status: http://albertopassalacqua.com/?p=1541)
OpenQBMM - An open-source implementation of quadrature-based moment methods.

To obtain more accurate answers, please specify the version of OpenFOAM you are using.
alberto is offline   Reply With Quote

Old   October 19, 2010, 13:52
Default
  #22
Senior Member
 
Vesselin Krastev
Join Date: Jan 2010
Location: University of Tor Vergata, Rome
Posts: 368
Rep Power: 20
vkrastev is on a distinguished road
Quote:
Originally Posted by vkrastev View Post
First of all, thanks a lot for your suggestions. About the nonOrth correctors, the 8-option was only a trial to see if there is some significant dependence of such parameter in my case, but actually the residuals' graphs that you have seen are all referred to simulations with nonOrth correctors set to 3. About your other indications, If I understand properly from your previous posts the right setting in the fvSchemes dictionary should be:

div(phi, U) Gauss linearUpwindV cellMDLimited Gauss linear 1;
div(phi, k) Gauss linearUpwindV cellMDLimited Gauss linear 1;
div(phi, epsilon) Gauss linearUpwindV cellMDLimited Gauss linear 1;

am I right?
Apart from this, I'll try to follow your advices and then I'll let you know what happens.
Thank you once again

Regards

V.

PS-I'm sorry if I'm repetitive, but I'll be very glad if you can direct me to some more informations about the theoretical basis of the linearUpwind scheme in its bounded formulation
Of course linearUpwindV should become linearUpwind for the scalars...sorry for the cut-paste error....
Ramzy1990 likes this.
vkrastev is offline   Reply With Quote

Old   October 19, 2010, 13:54
Default
  #23
Senior Member
 
Vesselin Krastev
Join Date: Jan 2010
Location: University of Tor Vergata, Rome
Posts: 368
Rep Power: 20
vkrastev is on a distinguished road
Quote:
Originally Posted by alberto View Post
linearUpwind not linearUpwindV for k and eps.

Check the implementation in the code. I am not aware of any specific references (it is pretty standard stuff).

Best,
Ok (for both the observations), thanks
vkrastev is offline   Reply With Quote

Old   October 21, 2010, 13:20
Default
  #24
Senior Member
 
Vesselin Krastev
Join Date: Jan 2010
Location: University of Tor Vergata, Rome
Posts: 368
Rep Power: 20
vkrastev is on a distinguished road
Quote:
Originally Posted by alberto View Post
I have a few comments on what you posted.

- First your case barely converges also with upwind (residuals of the pressure are the highest (yellow line)). Reduce the tolerances on the linear solvers to something closer to your machine precision (10^-12 for p, 10^-10 for the rest: yes it will take more iterations, it does not matter).

- Relax the turbulent quantities more than the velocity. Typically an URF = 0.4 works well.

- As schemes, linearUpwindV for div(phi, U) and linearUpwind for the rest, with cellLimited modifier for gradients should work just fine.

- Your mesh does not suffer of strong non-orthogonality, so I would not push the correctors too much (surely not to 8!).

- Stay away from SFCD, QUICK, UMIST. They won't give you any significant advantage. All the love for QUICK comes from the fact that it is formally third-order accurate, but its dissipation error is still high, and its stability is not good.

I hope this helps.

Best,
Hi Alberto,
as I said I've tried to follow your suggestions, but unfortunately the results were not satisfactory...In particular, I've tried the following new settings:

- URF factors for k and epsilon reduced to 0.4
- Tolerances of the solver reduced to 10^-10 for p and 10^-08 for the other quantities
- Relative tolerance for p reduced to 10^-03 (from 10^-02)
- divSchemes set to "Gauss linearUpwindV cellMDLimited Gauss linear" for the velocities and the same without "V" for k and epsilon

Additionally, before switching to linearUpwind I let the simulation go on with all the divSchemes set to the basic upwind scheme for 250 iterations, hoping this would improve the convergence of the run...But, as I said before, the results were not good. In particular, after a few iterations the residuals of k fall suddenly to a very low value, so in practice for the rest of the run the k-equation is no more solved. The other equations instead, keep to be solved and the relative residuals reduce slowly to lower values (with the exception of epsilon, whose residuals remain fairly constant). Apart from this, the solution after 1500 iterations is definetly not in accordance with what one should expect from the Ahmed body case: in particular, the tke distribution in the symmetry plane is completely unphysical (but this should be expected from a practically not resolved field), and the drag coefficient associated with the nose part of the body is heavily underestimated comparing with both experimental and other numerical data (including my previous all-upwind run). For more clearness, I will attache the images concerning the convergence patterns, the tke distribution in the symmetry plane and the nose drag coefficients calculated in two different runs (the linearUpwind one and a previous all-upwind one). Of course any help would be really appreciated...

Regards

V.
Attached Images
File Type: png residuals.png (34.0 KB, 535 views)
File Type: jpg tke_symmetry_plane.jpg (15.9 KB, 399 views)
File Type: png noseDrag.png (15.0 KB, 365 views)
Ramzy1990 likes this.
vkrastev is offline   Reply With Quote

Old   October 21, 2010, 13:53
Default
  #25
Senior Member
 
Alberto Passalacqua
Join Date: Mar 2009
Location: Ames, Iowa, United States
Posts: 1,912
Rep Power: 36
alberto will become famous soon enoughalberto will become famous soon enough
Hello,

the solution is not converged. You are doing a steady state simulation: residuals should be significantly lower. Set the convergence criterion at 10^-10 and wait until that value is reached or the residuals become completely flat for a number of iterations.

I would not start with upwind and then switch scheme. It is of no use. If the problem is set up correctly, it will converge also with the second-order scheme.

Since k is not being solved anymore, the under-relaxation of 0.4 is probably too strong (assuming the physical setup is correct). You might want to increase the URF again.

Best,
__________________
Alberto Passalacqua

GeekoCFD - A free distribution based on openSUSE 64 bit with CFD tools, including OpenFOAM. Available as in both physical and virtual formats (current status: http://albertopassalacqua.com/?p=1541)
OpenQBMM - An open-source implementation of quadrature-based moment methods.

To obtain more accurate answers, please specify the version of OpenFOAM you are using.
alberto is offline   Reply With Quote

Old   October 21, 2010, 14:24
Default
  #26
Senior Member
 
Vesselin Krastev
Join Date: Jan 2010
Location: University of Tor Vergata, Rome
Posts: 368
Rep Power: 20
vkrastev is on a distinguished road
Quote:
Originally Posted by alberto View Post
Hello,

the solution is not converged. You are doing a steady state simulation: residuals should be significantly lower. Set the convergence criterion at 10^-10 and wait until that value is reached or the residuals become completely flat for a number of iterations.

I would not start with upwind and then switch scheme. It is of no use. If the problem is set up correctly, it will converge also with the second-order scheme.

Since k is not being solved anymore, the under-relaxation of 0.4 is probably too strong (assuming the physical setup is correct). You might want to increase the URF again.

Best,
Hi,
ok, I'll try with these new settings but let me add some more observations:

- The same behaviour of tke's residuals can be viewed with the URFs set to 0.7: does it mean that even 0.7 is a too strong URF for k and epsilon?

- Let's say that lowering further the tolerance I will come to some satisfactory solution: it seems to me (but maybe I'm wrong about it) that this is a quite restrictive condition... should I amplify many times the total runtime in any situation that would require a lower numerical diffusion? Or maybe, as you said, this is simply a typical behaviour of the steady-state solvers (in this last case, It seems like it would be much more convenient to run a transient simulation and wait until it will converge to a steady-state configuration)?

However, thanks once again for your replies

regards

V.
vkrastev is offline   Reply With Quote

Old   October 21, 2010, 15:39
Default
  #27
Senior Member
 
Alberto Passalacqua
Join Date: Mar 2009
Location: Ames, Iowa, United States
Posts: 1,912
Rep Power: 36
alberto will become famous soon enoughalberto will become famous soon enough
Quote:
Originally Posted by vkrastev View Post
- The same behaviour of tke's residuals can be viewed with the URFs set to 0.7: does it mean that even 0.7 is a too strong URF for k and epsilon?
No, I don't think 0.7 is too strong. Try to figure out why k is not solved anymore, while all the rest keeps changing (if there is a reason different from being it converged).

Quote:
- Let's say that lowering further the tolerance I will come to some satisfactory solution: it seems to me (but maybe I'm wrong about it) that this is a quite restrictive condition... should I amplify many times the total runtime in any situation that would require a lower numerical diffusion? Or maybe, as you said, this is simply a typical behaviour of the steady-state solvers (in this last case, It seems like it would be much more convenient to run a transient simulation and wait until it will converge to a steady-state configuration)?
It is really common practice, independently from the code you use, to require convergence of residuals at machine precision for steady-state solutions (or very low). It has nothing to do with numerical diffusion but with accuracy.

Running an unsteady case would require significantly longer in general, since obtaining those residual values is fast with steady solvers if the mesh is good and the setup is correct.

Best,
__________________
Alberto Passalacqua

GeekoCFD - A free distribution based on openSUSE 64 bit with CFD tools, including OpenFOAM. Available as in both physical and virtual formats (current status: http://albertopassalacqua.com/?p=1541)
OpenQBMM - An open-source implementation of quadrature-based moment methods.

To obtain more accurate answers, please specify the version of OpenFOAM you are using.
alberto is offline   Reply With Quote

Old   October 22, 2010, 10:21
Default
  #28
Senior Member
 
Vesselin Krastev
Join Date: Jan 2010
Location: University of Tor Vergata, Rome
Posts: 368
Rep Power: 20
vkrastev is on a distinguished road
Quote:
It is really common practice, independently from the code you use, to require convergence of residuals at machine precision for steady-state solutions (or very low). It has nothing to do with numerical diffusion but with accuracy.
I got your point, but I didn't want to say that lowering the residuals means to have a less diffusive solution. I was talking about the fact that if one have to reach a really good overall accuracy there has to be two subsequential steps to follow: first, a higher order scheme is needed, because if you lowered a lot the residuals with a first order scheme you will have at most the best accuracy permitted by the scheme itself which, in some cases, maybe would be not satisfactory; secondly, as you explained above, the residuals should fall below a really small value. So, at the end, having in mind that convergence seems to be quite faster with more diffusive schemes, one have to enlarge significantly the runtime to reach real accuracy.

Apart from this, I lowered the tolerance of the solvers to 10^-11 for p and to 10^-10 for the other quantities and (as it could be seen from the residuals pattern atached) there are no substantial differences, because tke's residual still falls down below the imposed tolerance after a few iterations, and after that all the other residuals does not converge at all...I made a review of all the settings of the case, but honestly I cannot find a reason for which the k-equation behaves like this (apart from the numerical scheme chosen for the convection term). I can post again my fvSchemes and fvSolution dictionary, maybe they could add some more information which probably I'm missing...

Regards

V.

fvSchemes:
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

ddtSchemes
{
default steadyState;
}

gradSchemes
{
default Gauss linear;
grad(p) Gauss linear;
grad(U) Gauss linear;
}

divSchemes
{
default none;
div(phi,U) Gauss linearUpwindV cellMDLimited Gauss linear 1;
div(phi,k) Gauss linearUpwind cellMDLimited Gauss linear 1;
div(phi,epsilon) Gauss linearUpwind cellMDLimited Gauss linear 1;
div(R) Gauss linear;
div((nuEff*dev(grad(U).T()))) Gauss linear;
}

laplacianSchemes
{
default none;
laplacian(nuEff,U) Gauss linear corrected;
laplacian((1|A(U)),p) Gauss linear corrected;
laplacian(DkEff,k) Gauss linear corrected;
laplacian(DepsilonEff,epsilon) Gauss linear corrected;
}

interpolationSchemes
{
default linear;
interpolate(U) linear;
}

snGradSchemes
{
default corrected;
}

fluxRequired
{
default no;
p ;
}


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

fvSolution:
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

solvers
{
p
{
solver PCG;
preconditioner DIC;
tolerance 1e-11;
relTol 0.001;
}

U
{
solver PBiCG;
preconditioner DILU;
tolerance 1e-10;
relTol 0.1;
}

k
{
solver PBiCG;
preconditioner DILU;
tolerance 1e-10;
relTol 0.1;
}

epsilon
{
solver PBiCG;
preconditioner DILU;
tolerance 1e-10;
relTol 0.1;
}

}

SIMPLE
{
nNonOrthogonalCorrectors 3;
}

relaxationFactors
{
p 0.3;
U 0.7;
k 0.7;
epsilon 0.7;
}


// ************************************************** *********************** //
Attached Images
File Type: png residuals2.png (26.3 KB, 274 views)
vkrastev is offline   Reply With Quote

Old   October 22, 2010, 15:51
Default
  #29
Senior Member
 
Alberto Passalacqua
Join Date: Mar 2009
Location: Ames, Iowa, United States
Posts: 1,912
Rep Power: 36
alberto will become famous soon enoughalberto will become famous soon enough
Quote:
Originally Posted by vkrastev View Post
I got your point, but I didn't want to say that lowering the residuals means to have a less diffusive solution. I was talking about the fact that if one have to reach a really good overall accuracy there has to be two subsequential steps to follow: first, a higher order scheme is needed, because if you lowered a lot the residuals with a first order scheme you will have at most the best accuracy permitted by the scheme itself which, in some cases, maybe would be not satisfactory; secondly, as you explained above, the residuals should fall below a really small value. So, at the end, having in mind that convergence seems to be quite faster with more diffusive schemes, one have to enlarge significantly the runtime to reach real accuracy.
Well, finite volume schemes are globally at best second order accurate, and the increase in computational time is really negligible. ;-)

Quote:
Apart from this, I lowered the tolerance of the solvers to 10^-11 for p and to 10^-10 for the other quantities and (as it could be seen from the residuals pattern atached) there are no substantial differences, because tke's residual still falls down below the imposed tolerance after a few iterations, and after that all the other residuals does not converge at all...I made a review of all the settings of the case, but honestly I cannot find a reason for which the k-equation behaves like this (apart from the numerical scheme chosen for the convection term). I can post again my fvSchemes and fvSolution dictionary, maybe they could add some more information which probably I'm missing...
I would start looking more carefully at other aspects of the case setup then. The numerical setup seems fine.

Best,
__________________
Alberto Passalacqua

GeekoCFD - A free distribution based on openSUSE 64 bit with CFD tools, including OpenFOAM. Available as in both physical and virtual formats (current status: http://albertopassalacqua.com/?p=1541)
OpenQBMM - An open-source implementation of quadrature-based moment methods.

To obtain more accurate answers, please specify the version of OpenFOAM you are using.
alberto is offline   Reply With Quote

Old   October 23, 2010, 05:56
Default
  #30
Senior Member
 
Markus Rehm
Join Date: Mar 2009
Location: Erlangen (Germany)
Posts: 184
Rep Power: 17
markusrehm is on a distinguished road
Hello Vesselin,

have you tried to set k and epsilon tolerance to real machine precision (1e-16)? Does make a difference?

Regards, Markus.
markusrehm is offline   Reply With Quote

Old   October 25, 2010, 05:41
Default
  #31
Senior Member
 
Vesselin Krastev
Join Date: Jan 2010
Location: University of Tor Vergata, Rome
Posts: 368
Rep Power: 20
vkrastev is on a distinguished road
Quote:
Originally Posted by markusrehm View Post
Hello Vesselin,

have you tried to set k and epsilon tolerance to real machine precision (1e-16)? Does make a difference?

Regards, Markus.

Hi, Markus
if you're talking about the absolute tolerances I didn't try such a setting (at most I have lowered the tolerances to 10^-10, 10^-11), but i don't think that it could make some difference, because the run starts to have some problems with the k-equation after a few iterations and after this point the residuals stop to converge at all (you can check this behaviour in the image attached to my previous post), remaining at values of about 10^-03/10^-04.

Regards

V.
vkrastev is offline   Reply With Quote

Old   October 25, 2010, 12:17
Default
  #32
Senior Member
 
Vesselin Krastev
Join Date: Jan 2010
Location: University of Tor Vergata, Rome
Posts: 368
Rep Power: 20
vkrastev is on a distinguished road
However, my case is a 3D external aerodynamics case, with a symmetric domain. The domain itself is composed by inlet and outlet sections, the symmetry plane, a floor and top and side patches. I've tried the following settings:

setting 1:
-top and side treated as generic patches, with a zeroGradient condition for all the quantities involved in the calculation

Of course the inlet and outlet are generic patches and the symmetry plane is symmetryPlane. Here there is the content of my 0/ folder for this setting.

U:

internalField uniform (40 0 0);

boundaryField
{
inlet
{
type fixedValue;
value uniform (40 0 0);
}

outlet
{
type zeroGradient;
}

body
{
type fixedValue;
value uniform (0 0 0);
}

nose
{
type fixedValue;
value uniform (0 0 0);
}

slant
{
type fixedValue;
value uniform (0 0 0);
}

back
{
type fixedValue;
value uniform (0 0 0);
}

floor
{
type fixedValue;
value uniform (0 0 0);
}

side
{
type zeroGradient;
}

top
{
type zeroGradient;
}

symmetry
{
type symmetryPlane;
}
}


p:

internalField uniform 0;

boundaryField
{
inlet
{
type zeroGradient;
}

outlet
{
type fixedValue;
value uniform 0;
}

body
{
type zeroGradient;
}

nose
{
type zeroGradient;
}

slant
{
type zeroGradient;
}

back
{
type zeroGradient;
}

floor
{
type zeroGradient;
}

symmetry
{
type symmetryPlane;
}

side
{
type zeroGradient;
}

top
{
type zeroGradient;
}
}


k:

internalField uniform 0.0096;

boundaryField
{
inlet
{
type fixedValue;
value uniform 0.0096;
}
outlet
{
type zeroGradient;
}
body
{
type kqRWallFunction;
value uniform 0.0096;
}

nose
{
type kqRWallFunction;
value uniform 0.0096;
}

slant
{
type kqRWallFunction;
value uniform 0.0096;
}

back
{
type kqRWallFunction;
value uniform 0.0096;
}

floor
{
type kqRWallFunction;
value uniform 0.0096;
}

side
{
type zeroGradient;
}

top
{
type zeroGradient;
}

symmetry
{
type symmetryPlane;
}
}


epsilon:

internalField uniform 0.155;

boundaryField
{
inlet
{
type fixedValue;
value uniform 0.00155;
}
outlet
{
type zeroGradient;
}
body
{
type epsilonWallFunction;
value uniform 0.00155;
}
nose
{
type epsilonWallFunction;
value uniform 0.00155;
}
slant
{
type epsilonWallFunction;
value uniform 0.00155;
}
back
{
type epsilonWallFunction;
value uniform 0.00155;
}
floor
{
type epsilonWallFunction;
value uniform 0.00155;
}

side
{
type zeroGradient;
}

top
{
type zeroGradient;
}

symmetry
{
type symmetryPlane;
}
}


nut:

internalField uniform 0;

boundaryField
{
inlet
{
type calculated;
value uniform 0;
}
outlet
{
type calculated;
value uniform 0;
}
body
{
type nutWallFunction;
value uniform 0;
}
nose
{
type nutWallFunction;
value uniform 0;
}
slant
{
type nutWallFunction;
value uniform 0;
}
back
{
type nutWallFunction;
value uniform 0;
}
floor
{
type nutWallFunction;
value uniform 0;
}

side
{
type zeroGradient;
value uniform 0;
}

top
{
type zeroGradient;
value uniform 0;
}

symmetry
{
type symmetryPlane;
}
}

k and epsilon values are calculated imposing a turbulence inlet intensity of 0.2% (reffered to the inlet bulk velocity).


setting 2: is very close to setting 1, but the top and side patches are imposed to be symmetry planes. All the entries in the 0/ folder are changed consequently

As you already know, both of the settigns does not work with higher order convection schemes, even with a mesh with quite good quality parameters (you can check all of them in my previous posts).
Any additional ideas?

V.
vkrastev is offline   Reply With Quote

Old   October 27, 2010, 07:04
Default
  #33
Senior Member
 
Vesselin Krastev
Join Date: Jan 2010
Location: University of Tor Vergata, Rome
Posts: 368
Rep Power: 20
vkrastev is on a distinguished road
Does anybody think that changing the linear solver of some of the quantities (for instance, switching to GAMG for p) would help?

V.
vkrastev is offline   Reply With Quote

Old   November 2, 2010, 09:34
Default
  #34
Senior Member
 
Vesselin Krastev
Join Date: Jan 2010
Location: University of Tor Vergata, Rome
Posts: 368
Rep Power: 20
vkrastev is on a distinguished road
Quote:
Originally Posted by alberto View Post
Well, finite volume schemes are globally at best second order accurate, and the increase in computational time is really negligible. ;-)



I would start looking more carefully at other aspects of the case setup then. The numerical setup seems fine.

Best,
Ok, after about 50 runs I have to admit that I'm starting to get a bit frustrated...Could you please post the fvSchemes and fvSolution dictionaries of a well working case (I mean a case that you have already ran successfully) concerning a steady-state simulation on a fully tetrahedral 3D mesh, with a higher order convection scheme configuration? I will be also very glad if I can see the checkMesh report of such a case. Just for information I will post down below the fvSchemes, fvSolution, checkMesh report and residuals pattern of one of my last attempts (I tried to run the case with a linearUpwindV scheme imposed only on the velocity, hoping that the upwind choice for at least the turbulent quantities would improve convergence and stability...when I tried the linearUpwind schemes for k and epsilon the simulation simply blows up, also with URFs set to 0.5). In addition (I don't know if I have mentioned it before), I use the 1.6 version.

Thank you in advance

Best Regards

V.

fvSchemes:
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

ddtSchemes
{
default steadyState;
}

gradSchemes
{
default leastSquares;
grad(p) leastSquares;
grad(U) leastSquares;
}

divSchemes
{
default none;
div(phi,U) Gauss linearUpwindV faceMDLimited leastSquares 1;
div(phi,k) Gauss upwind;
div(phi,epsilon) Gauss upwind;
div(R) Gauss linear;
div((nuEff*dev(grad(U).T()))) Gauss linear;
}

laplacianSchemes
{
default Gauss linear limited 0.5;
// laplacian(nuEff,U) Gauss linear corrected;
// laplacian((1|A(U)),p) Gauss linear corrected;
// laplacian(DkEff,k) Gauss linear corrected;
// laplacian(DepsilonEff,epsilon) Gauss linear corrected;
}

interpolationSchemes
{
default linear;
interpolate(U) linear;
}

snGradSchemes
{
default limited 0.5;
}

fluxRequired
{
default no;
p ;
}


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

fvSolution:
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

solvers
{
// p
// {
// solver PCG;
// preconditioner DIC;
// tolerance 1e-11;
// relTol 0.001;
// }

p
{
solver GAMG;
tolerance 1e-11;
relTol 0.0001;
smoother GaussSeidel;
cacheAgglomeration true;
nCellsInCoarsestLevel 1000;
agglomerator faceAreaPair;
mergeLevels 1;
}

U
{
solver PBiCG;
preconditioner DILU;
tolerance 1e-10;
relTol 0.1;
}

k
{
solver PBiCG;
preconditioner DILU;
tolerance 1e-10;
relTol 0.1;
}

epsilon
{
solver PBiCG;
preconditioner DILU;
tolerance 1e-10;
relTol 0.1;
}

}

SIMPLE
{
nNonOrthogonalCorrectors 3;
}

relaxationFactors
{
p 0.3;
U 0.7;
k 0.5;
epsilon 0.5;
}


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

checkMesh report:
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Create time

Create polyMesh for time = 0

Time = 0

Mesh stats
points: 357845
faces: 3218365
internal faces: 3105296
cells: 1527609
boundary patches: 10
point zones: 0
face zones: 0
cell zones: 0

Overall number of cells of each type:
hexahedra: 0
prisms: 213225
wedges: 0
pyramids: 0
tet wedges: 0
tetrahedra: 1314384
polyhedra: 0

Checking topology...
Boundary definition OK.
Point usage OK.
Upper triangular ordering OK.
Face vertices OK.
Number of regions: 1 (OK).

Checking patch topology for multiply connected surfaces ...
Patch Faces Points Surface topology
symmetry 22600 12584 ok (non-closed singly connected)
floor 14660 7609 ok (non-closed singly connected)
inlet 426 242 ok (non-closed singly connected)
outlet 432 245 ok (non-closed singly connected)
side 2228 1204 ok (non-closed singly connected)
top 1648 910 ok (non-closed singly connected)
body 54899 27768 ok (non-closed singly connected)
nose 8704 4461 ok (non-closed singly connected)
slant 3948 2058 ok (non-closed singly connected)
back 3524 1841 ok (non-closed singly connected)

Checking geometry...
Overall domain bounding box (-1.26 0 0.1945) (6.27 1.44 1.1945)
Mesh (non-empty, non-wedge) directions (1 1 1)
Mesh (non-empty) directions (1 1 1)
Boundary openness (8.32191e-20 -1.52667e-19 7.68048e-20) OK.
Max cell openness = 2.7637e-16 OK.
Max aspect ratio = 5.97752 OK.
Minumum face area = 1.71072e-06. Maximum face area = 0.0159908. Face area magnitudes OK.
Min volume = 1.10751e-09. Max volume = 0.000575106. Total volume = 10.7878. Cell volumes OK.
Mesh non-orthogonality Max: 54.8744 average: 14.5203
Non-orthogonality check OK.
Face pyramids OK.
Max skewness = 1.98187 OK.

Mesh OK.

End
Attached Images
File Type: png residualsRun1.png (31.9 KB, 223 views)
vkrastev is offline   Reply With Quote

Old   November 2, 2010, 09:52
Default
  #35
Senior Member
 
Vesselin Krastev
Join Date: Jan 2010
Location: University of Tor Vergata, Rome
Posts: 368
Rep Power: 20
vkrastev is on a distinguished road
Quote:
Originally Posted by vkrastev View Post
Ok, after about 50 runs I have to admit that I'm starting to get a bit frustrated...Could you please post the fvSchemes and fvSolution dictionaries of a well working case (I mean a case that you have already ran successfully) concerning a steady-state simulation on a fully tetrahedral 3D mesh, with a higher order convection scheme configuration? I will be also very glad if I can see the checkMesh report of such a case. Just for information I will post down below the fvSchemes, fvSolution, checkMesh report and residuals pattern of one of my last attempts (I tried to run the case with a linearUpwindV scheme imposed only on the velocity, hoping that the upwind choice for at least the turbulent quantities would improve convergence and stability...when I tried the linearUpwind schemes for k and epsilon the simulation simply blows up, also with URFs set to 0.5). In addition (I don't know if I have mentioned it before), I use the 1.6 version.

Thank you in advance

Best Regards

V.

fvSchemes:
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

ddtSchemes
{
default steadyState;
}

gradSchemes
{
default leastSquares;
grad(p) leastSquares;
grad(U) leastSquares;
}

divSchemes
{
default none;
div(phi,U) Gauss linearUpwindV faceMDLimited leastSquares 1;
div(phi,k) Gauss upwind;
div(phi,epsilon) Gauss upwind;
div(R) Gauss linear;
div((nuEff*dev(grad(U).T()))) Gauss linear;
}

laplacianSchemes
{
default Gauss linear limited 0.5;
// laplacian(nuEff,U) Gauss linear corrected;
// laplacian((1|A(U)),p) Gauss linear corrected;
// laplacian(DkEff,k) Gauss linear corrected;
// laplacian(DepsilonEff,epsilon) Gauss linear corrected;
}

interpolationSchemes
{
default linear;
interpolate(U) linear;
}

snGradSchemes
{
default limited 0.5;
}

fluxRequired
{
default no;
p ;
}


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

fvSolution:
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

solvers
{
// p
// {
// solver PCG;
// preconditioner DIC;
// tolerance 1e-11;
// relTol 0.001;
// }

p
{
solver GAMG;
tolerance 1e-11;
relTol 0.0001;
smoother GaussSeidel;
cacheAgglomeration true;
nCellsInCoarsestLevel 1000;
agglomerator faceAreaPair;
mergeLevels 1;
}

U
{
solver PBiCG;
preconditioner DILU;
tolerance 1e-10;
relTol 0.1;
}

k
{
solver PBiCG;
preconditioner DILU;
tolerance 1e-10;
relTol 0.1;
}

epsilon
{
solver PBiCG;
preconditioner DILU;
tolerance 1e-10;
relTol 0.1;
}

}

SIMPLE
{
nNonOrthogonalCorrectors 3;
}

relaxationFactors
{
p 0.3;
U 0.7;
k 0.5;
epsilon 0.5;
}


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

checkMesh report:
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Create time

Create polyMesh for time = 0

Time = 0

Mesh stats
points: 357845
faces: 3218365
internal faces: 3105296
cells: 1527609
boundary patches: 10
point zones: 0
face zones: 0
cell zones: 0

Overall number of cells of each type:
hexahedra: 0
prisms: 213225
wedges: 0
pyramids: 0
tet wedges: 0
tetrahedra: 1314384
polyhedra: 0

Checking topology...
Boundary definition OK.
Point usage OK.
Upper triangular ordering OK.
Face vertices OK.
Number of regions: 1 (OK).

Checking patch topology for multiply connected surfaces ...
Patch Faces Points Surface topology
symmetry 22600 12584 ok (non-closed singly connected)
floor 14660 7609 ok (non-closed singly connected)
inlet 426 242 ok (non-closed singly connected)
outlet 432 245 ok (non-closed singly connected)
side 2228 1204 ok (non-closed singly connected)
top 1648 910 ok (non-closed singly connected)
body 54899 27768 ok (non-closed singly connected)
nose 8704 4461 ok (non-closed singly connected)
slant 3948 2058 ok (non-closed singly connected)
back 3524 1841 ok (non-closed singly connected)

Checking geometry...
Overall domain bounding box (-1.26 0 0.1945) (6.27 1.44 1.1945)
Mesh (non-empty, non-wedge) directions (1 1 1)
Mesh (non-empty) directions (1 1 1)
Boundary openness (8.32191e-20 -1.52667e-19 7.68048e-20) OK.
Max cell openness = 2.7637e-16 OK.
Max aspect ratio = 5.97752 OK.
Minumum face area = 1.71072e-06. Maximum face area = 0.0159908. Face area magnitudes OK.
Min volume = 1.10751e-09. Max volume = 0.000575106. Total volume = 10.7878. Cell volumes OK.
Mesh non-orthogonality Max: 54.8744 average: 14.5203
Non-orthogonality check OK.
Face pyramids OK.
Max skewness = 1.98187 OK.

Mesh OK.

End
PS - I use ANSA as mesh generator
vkrastev is offline   Reply With Quote

Old   November 2, 2010, 11:41
Default
  #36
Senior Member
 
Alberto Passalacqua
Join Date: Mar 2009
Location: Ames, Iowa, United States
Posts: 1,912
Rep Power: 36
alberto will become famous soon enoughalberto will become famous soon enough
Hi,

I have used the setup I told you successfully in almost all my runs, so I would double-check the boundary conditions. Your checkMesh does not have anything wrong.

Out of curiosity, what happens to the residuals if you let the case run more than 2000 iterations? Do they keep increasing or do they oscillate?

Best,
__________________
Alberto Passalacqua

GeekoCFD - A free distribution based on openSUSE 64 bit with CFD tools, including OpenFOAM. Available as in both physical and virtual formats (current status: http://albertopassalacqua.com/?p=1541)
OpenQBMM - An open-source implementation of quadrature-based moment methods.

To obtain more accurate answers, please specify the version of OpenFOAM you are using.
alberto is offline   Reply With Quote

Old   November 2, 2010, 12:07
Default
  #37
Senior Member
 
Vesselin Krastev
Join Date: Jan 2010
Location: University of Tor Vergata, Rome
Posts: 368
Rep Power: 20
vkrastev is on a distinguished road
Quote:
Originally Posted by alberto View Post
Hi,

I have used the setup I told you successfully in almost all my runs, so I would double-check the boundary conditions. Your checkMesh does not have anything wrong.

Out of curiosity, what happens to the residuals if you let the case run more than 2000 iterations? Do they keep increasing or do they oscillate?

Best,
The boundary conditions are described in details in one of my last posts in this tread but, briefly talking, they are:
-Symmetrycal domain composed by a half Ahmedd model and the surrounding wind tunnel
-The symmetry plane is (obviously) set as symmetryPlane bc
-The inlet is set to fixed values for U, k and epsilon and to zeroGradient for p
-The outlet is set to zeroGradient for all the quantities except p, which is fixed to 0.
-The floor and the body are set as walls, so fixedValue (0 0 0) for U, zeroGradient for p and wallFunctions for k, epsilon and nut (I'm using the high-Re realizable k-epsilon model).
-For the side and top patches I decided to use the zeroGradient condition for all the quantities. I also tried to switch to symmetryPlane, without (seemingly) any benefits.

To my knowledge, and having a look around for the literature on the Ahmed body (incompressible) simulations, such a BC setting seems to my quite standard.

About the second issue, the same run is now going on further: I could tell you more about it maybe this evening or tomorrow morning

Best Regards

V.
vkrastev is offline   Reply With Quote

Old   November 3, 2010, 06:29
Default
  #38
Senior Member
 
Vesselin Krastev
Join Date: Jan 2010
Location: University of Tor Vergata, Rome
Posts: 368
Rep Power: 20
vkrastev is on a distinguished road
Quote:
Originally Posted by alberto View Post
Out of curiosity, what happens to the residuals if you let the case run more than 2000 iterations? Do they keep increasing or do they oscillate?
This is the residuals pattern after 2000 iterations: indeed, they do oscillate over practically constant values.

V.
Attached Images
File Type: png run1_2000-5500.png (30.6 KB, 326 views)

Last edited by vkrastev; November 3, 2010 at 12:05.
vkrastev is offline   Reply With Quote

Old   November 4, 2010, 04:46
Default
  #39
Senior Member
 
Matthias Voß
Join Date: Mar 2009
Location: Berlin, Germany
Posts: 449
Rep Power: 20
mvoss is on a distinguished road
hi,

since i´m not very familiar with all of the impacts of the schemes and everything this is maybe smth. totally wrong: if i look at the flow and keep in mind that my residuals aren't converging tight enough... isn't that supposing a transient effect (back flow etc.) in my flow field leading to the "poor" convergence?
--> check for the high residuals within the domain.. if they are at the back of the body this could be an hint on transient effects.

neewbie
mvoss is offline   Reply With Quote

Old   November 4, 2010, 06:21
Default
  #40
Senior Member
 
Vesselin Krastev
Join Date: Jan 2010
Location: University of Tor Vergata, Rome
Posts: 368
Rep Power: 20
vkrastev is on a distinguished road
Quote:
Originally Posted by neewbie View Post
hi,

since i´m not very familiar with all of the impacts of the schemes and everything this is maybe smth. totally wrong: if i look at the flow and keep in mind that my residuals aren't converging tight enough... isn't that supposing a transient effect (back flow etc.) in my flow field leading to the "poor" convergence?
--> check for the high residuals within the domain.. if they are at the back of the body this could be an hint on transient effects.

neewbie
Hi neewbie,
about the backflow I don't see any trace of it during the runs (the outlet patch is at 5 body lenghts downstream of the body itself: at this distance the flow is again quite unidirectional, so having a backflow there would be a really uncommon situation...).
About to check where the residuals are higher in the domain, can you give me a hint on how it could be done? Thank you

V.
vkrastev is offline   Reply With Quote

Reply


Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off
Trackbacks are Off
Pingbacks are On
Refbacks are On


Similar Threads
Thread Thread Starter Forum Replies Last Post
4th order schemes in channelOodles maka OpenFOAM Bugs 9 January 19, 2009 12:58
Help for fourth order accurate convective schemes Z.C.Wang Main CFD Forum 0 January 15, 2009 07:53
2nd order conservative schemes taw Main CFD Forum 1 September 16, 2008 08:05
CFL condition for higher order schemes Shyam Main CFD Forum 2 February 14, 2008 15:24
High order compact finite difference schemes Mikhail Main CFD Forum 6 August 5, 2003 11:36


All times are GMT -4. The time now is 16:23.