CFD Online Discussion Forums (http://www.cfd-online.com/Forums/)
-   OpenFOAM (http://www.cfd-online.com/Forums/openfoam/)
-   -   NACA0012: Residual convergence issue (http://www.cfd-online.com/Forums/openfoam/91033-naca0012-residual-convergence-issue.html)

 Ofoam200_nwb July 28, 2011 07:42

NACA0012: Residual convergence issue

4 Attachment(s)
Hi everyone,

I am just about getting started with OpenFOAM, and I was trying to simulate a NACA0012 2D section at 5 degrees angle of attack in a Reynolds number of 5 million (medium: water. So rho = 1e3 and nu = 1e-6). However, when I run the simulation, the residuals, especially pressure one, refuse to converge. The residuals for Ux and Uy hover around 1E5 and for p around 1E3!
1. Could this be corrected by the use of a higher order numerical scheme?
Or
2. Could this be because of the grid? It has some high aspect ratio cells, but all of them I suppose are away from the airfoil. Around the airfoil I tried to maintain a y+ of 30 and obtained an average y+ of 26, so I think that is fine. But then again, the solution itself is faulty so I don't know how reliable this value is (I get a Cd of .05 and Cl of .29, when, from JavaFoil, I was expecting Cl = .6 and Cd = .01!)

Turbulence Model - Spallart Allmaras

For your reference: blockMeshDict, controlDict, transportProperties and fvSolution are attached.

Thanks a lot.
--
M. Butc

 Ofoam200_nwb July 28, 2011 09:40

Any sort of suggestions are most welcome.
Thank you!

 FelixL July 28, 2011 10:06

Hello, Mitch,

first of all: residual convergence and solution convergence basically are two different things. If the residuals are below a certain level, the solution might have converged, but this is not neccessarily the case, so its usually better to check for convergence via integral parameters like drag coefficient, pressure drop and so on.
Anyway, it's the same vice versa: the solution might reach a convergent state even when the residuals (especially pressure) seem to stagnate or won't fall below a certain level. Here, too, it is important to check the integral parameters and not the residuals.

Going back to your problem: Have you checked if your force coefficients are converged? It shouldn't be too much of an issue if the pressure residual doesn't drop anymore, as long as other indicators like the continuity errors are reasonably low, everything should be alright.

I had the experience that my pressure residuals stagnate or don't stagnate depending on the setup - still the results were nearly the same. E.g. when starting the simulation with a good initial field (for example mapped from another case), the pressure residual often stagnates quickly, but the results aren't affected.

Regarding your questions: yes you could try another scheme (which did you use, anyway?) - but residual convergence is not a matter of scheme accuracy. It's more about how the residuals are defined and how the schemes work. I made the experience that using a TVD scheme (like limitedLinear) lead to stagnating pressure residuals and using a different scheme (like linearUpwind) provided pressure residual convergence. The results weren't affected then, too, though.

Poor quality grids also affect the residual convergence, yes. And I suspect there's something with your mesh. My suspicion isn't because of the convergence issues, though, but for the poor results (especially the lift is way too off!). Try increasing the domain size, refine the grid (esp. near the wall), check the boundary conditions etc.pp.
There are a lot of threads about naca 0012 airfoils here in the forums, I think you'll get some ideas there.

Greetings,
Felix.

 Ofoam200_nwb July 28, 2011 12:30

Felix,

Thanks a lot. I get what yo're saying about convergence, and I do suspect that my grid is maybe not large enough (I did make sure that it was sufficiently fine close to the airfoil.) Anyway, I'll monitor force coefficients, make my grid slightly larger, and let you know if it works.
About the scheme, I think it is set to GaussSeidel as of now.
Thanks again.

--
M

 Ofoam200_nwb July 29, 2011 09:47

4 Attachment(s)
Felix,

I enlarged my domain, waited for the integral quantities to converge and the results I finally obtained were:
Cl = 0.52; Cd = 0.045! (large by a factor of 5!)

I have attached my y+ plot over the airfoil surface. (Nose - Top - Bottom)
The only place where the y+ is less than 30 is the nose, but that is such a small region that I feel like it can't be the cause of this large an error (or can it?).

I've also attached three pictures of my mesh close to the foil. Noticeable is the slightly abrupt jumps in cell sizes from one block to the other, however these are far from the airfoil surface. Close to the surface the variation is mostly smooth. Could you please take a quick look and see if the mesh seems to be problematic?

Thanks a lot, again. Look forward to hearing from you soon.
--
M

 Ofoam200_nwb July 29, 2011 09:52

I should probably clarify the yplus plot. The order in which I defined the airfoil faces are:
1. nose 2. top (nose-end to trailing edge) 3. bottom (nose-end to trailing edge)
As a result of this, the first dip is at the nose, then comes the top part (from where the nose ends till the trailing edge), and then the bottom part (from nose-end to trailing edge), and hence the apparent discontinuity.

Hello,
Have you seen this discussion: http://www.cfd-online.com/Forums/ope...igh-drag.html? I found it really useful...

 FelixL July 29, 2011 13:49

Hello, Mitch,

the lift coefficient looks much better now, good work.

Regarding GaussSeidel - this is an algorithm for solving linear systems of equations, not an interpolation scheme. Could you attach your fvSchemes file? That might help getting an idea about the rest of your numerical setup.

What boundary conditions are you using? Especially for nut? In case you're using blended wall functions, your y+ should be fine. If you're using HighRe-Wallfunctions, your y+ values could be bigger (between 30 and 100). And in the case of no Wall functions, your near wall grid is much too coarse (should be y+=1...5).

Oh, and I nearly forgot: these cell size jumps should be avoided. Your case seems to converge and the results don't look too bad, but still - I always get a bad feeling in my numerical guts when looking at such cell size ratios.

Greetings,
Felix

 Martin Hegedus July 29, 2011 15:42

In regards to residual not converging:

A short while ago we had a topic which included the airfoil2D which is under the simpleFoam tutorial directory. http://www.cfd-online.com/Forums/ope...oam-v-v-2.html At zero degrees angle of attack the problem would not converge in regards to getting the residuals to machine zero because of the use of freestream boundary condition. Using other boundary conditions such as combinations of fixedvalue and zerogradient worked. Here is felixL's post for BCs regarding BCs. http://www.cfd-online.com/Forums/ope...tml#post313369. The problem with freestream bc may be because the flow at zero degrees angle of attack is parallel to the top and bottom faces of the box which defines the outer boundary.

The solution for a flat plate using the freeStream BC did not converge to machine zero either (residuals). The top boundary condition of the outer edge for this case is also parallel to the far field velocity vector.

 Ofoam200_nwb July 30, 2011 05:48

Thanks, everyone! Lots of good stuff out there that I've seemed to completely miss. Sorry about that.

Mad: Thanks. This is certainly germane. I will try switching models keeping in mind their y+ restrictions and see if it helps.

Felix: As I said, I'm just a beginner, so the bcs and schemes I'm using for nut and everything else are the same as for the simpleFoam airFoil2D case. I do suspect however that I may be using a too high inlet turbulent viscosity. Anyway, the relevant files are attached along with this message. I am using the nutUSpaldingWallFunction, which frankly I can't even find in the UserGuide, but for a first simulation I just copied this from the tutorial. What you say about the cell jumps is worth listening to. Well, I don't have meshing experience so that is something I am working on too.

nut BC
Code:

```/*--------------------------------*- C++ -*----------------------------------*\ | =========                |                                                | | \\      /  F ield        | OpenFOAM: The Open Source CFD Toolbox          | |  \\    /  O peration    | Version:  2.0.0                                | |  \\  /    A nd          | Web:      www.OpenFOAM.com                      | |    \\/    M anipulation  |                                                | \*---------------------------------------------------------------------------*/ FoamFile {     version    2.0;     format      ascii;     class      volScalarField;     object      nut; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // dimensions      [0 2 -1 0 0 0 0]; internalField  uniform 0.14; boundaryField {     inlet     {         type            freestream;         freestreamValue uniform 0.14;     }     upperLower     {         type            freestream;         freestreamValue uniform 0.14;     }     outlet     {         type            freestream;         freestreamValue uniform 0.14;     }     airFoil     {         type            nutUSpaldingWallFunction;         value          uniform 0;     }     frontAndBack     {         type            empty;     } } // ************************************************************************* //```
p BC
Code:

```/*--------------------------------*- C++ -*----------------------------------*\ | =========                |                                                | | \\      /  F ield        | OpenFOAM: The Open Source CFD Toolbox          | |  \\    /  O peration    | Version:  2.0.0                                | |  \\  /    A nd          | Web:      www.OpenFOAM.com                      | |    \\/    M anipulation  |                                                | \*---------------------------------------------------------------------------*/ FoamFile {     version    2.0;     format      ascii;     class      volScalarField;     object      p; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // dimensions      [0 2 -2 0 0 0 0]; internalField  uniform 0; boundaryField {     inlet     {         type            freestreamPressure;     }         upperLower     {         type            freestreamPressure;     }     outlet     {         type            freestreamPressure;     }     airFoil     {         type            zeroGradient;     }     frontAndBack     {         type            empty;     } } // ************************************************************************* //```
fvSchemes
Code:

```/*--------------------------------*- C++ -*----------------------------------*\ | =========                |                                                | | \\      /  F ield        | OpenFOAM: The Open Source CFD Toolbox          | |  \\    /  O peration    | Version:  2.0.0                                | |  \\  /    A nd          | Web:      www.OpenFOAM.com                      | |    \\/    M anipulation  |                                                | \*---------------------------------------------------------------------------*/ FoamFile {     version    2.0;     format      ascii;     class      dictionary;     location    "system";     object      fvSchemes; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // ddtSchemes {     default        steadyState; } gradSchemes {     default        Gauss linear;     grad(p)        Gauss linear;     grad(U)        Gauss linear; } divSchemes {     default        none;     div(phi,U)      Gauss linearUpwind grad(U);     div(phi,nuTilda) Gauss linearUpwind grad(nuTilda);     div((nuEff*dev(T(grad(U))))) Gauss linear; } laplacianSchemes {     default        none;     laplacian(nuEff,U) Gauss linear corrected;     laplacian((1|A(U)),p) Gauss linear corrected;     laplacian(DnuTildaEff,nuTilda) Gauss linear corrected;     laplacian(1,p)  Gauss linear corrected; } interpolationSchemes {     default        linear;     interpolate(U)  linear; } snGradSchemes {     default        corrected; } fluxRequired {     default        no;     p              ; } // ************************************************************************* //```
Martin: Thanks for the links. I am in fact using the freestream boundary conditions, and for the angle of attack, I meshed my geometry with a rotated Foil and hence the velocity is indeed parallel to the top and bottom faces. I'll try playing around with the bcs then.

Thanks, again!
--
M

 Ofoam200_nwb July 30, 2011 17:38

3 Attachment(s)
Hey everyone,

To start with the improvements you all suggested, I refined my grid a bit from 80k cells to 120k cells, reduced the abruptness of cell size jumps a bit, and reduced freestream turbulent viscosity from 0.14 (a value I had copied from the turotrial file) to 5e-7 (so no incoming turbulence, I guess). I have what may seem like good news, as well as bad news.

1. I got a lift coefficient within 1% of the 2pi slope value (0.548; JavaFoil gives a much higher value of 0.596! Anyone know why this is?), and a drag coefficient value within 1% of the JavaFoil value (.0107). My velocity and pressure residuals all reduced to orders of magnitude of 1e-07.

2. The snag is that my domain doesn't seem to be large enough, and my yPlus values are out of whack. My yPlus values lie in the buffer region (>5, <30), and hence any of the wall function models shouldn't be applicable (right? I use the Spalding wall function.). And for my domain, there seems to be not enough space for the wake region to decay! (Plots attached.) The turbulent viscosity also doesn't decay, but I remember reading in "A First Course in Turbulence" (by Tennekes and Lumley) that a plane wake, once turbulent, remains turbulent. So that seems fine. But the yPlus values are definitely problematic, right?

I am attaching plots of U,p and nut along with this message. Look forward to hearing from you all soon!

--
M

 Martin Hegedus July 30, 2011 22:06

Lift curve slopes for airfoils will be higher than 2pi due to boundary layer thickness. I view it as the boundary layer creating an airfoil with a longer chord than the original. OK, not entirely correct. But, that's how I view it. Can't comment on JavaFoil. But at an Re of 5e6 your results should probably be a couple percent higher than 2pi.

I wouldn't be too concerned about the wake leaving your outer boundary. However, you do want the outer boundary far enough so that the wake and the outer boundary condition don't fight each other too much.

A bigger concern is to have the outer boundary far enough so that the error of not including the airfoil vortex at the outer boundary is negligible. Usually, for an airfoil, I have my outer boundary at 100 to 150 times the chord. Best approach is to expand your outer boundary and see how much your results change.

I'm not a fan of wall functions, so I do not use them. Therefore, I usually have y+=1 at about the 10% chord location.

 Ofoam200_nwb July 31, 2011 06:42

Hi, Martin,
JavaFoil is kind of like XFoil but less rigourous in its approach. Yes, I remember reading about effective-lengthening of geometry somewhere.
Anyway, regarding the wake, in the pictures that I posted above, the outlet was 15 chord-lengths away. I increased that to 20 chord-lengths and got virtually no difference (a difference of the order of magnitude 1e-4 for the lift coefficient), so I think the domain may be fine now. The lift coefficient is <1% below the 2pi gradient value, and the drag seems to be spot on. Now I've got to try and see what I did right. And check the pressure distribution as well.

Thanks for all your help, people. These forums are awesome.
--
M

 FelixL July 31, 2011 09:24

Hello, everyone,

actually, Mitch, your lift results seem pretty good to me. My experimental data available reads an cl of about 0.54 so I think JavaFoil is a bit off due to modeling errors. Xfoil yields ~0.55.

The nutUSpaldingWallFunction you're using is a continuous wall function. This means the wall function is applicable to the viscous sublayer, the buffer layer and the log layer. Your y+-values should be fine then, so no worries there!

Glad we could be of help.

Greetings,
Felix

Hi FelixL and Foamers,
I am simulating axial flow fan, I am using MRFSimpleFoam.
The settings are almost like the tutorial case (mixer2DVessle).

The fan section is with shroud, my chord length is 0.186m, I am using 3m inlet section and 7m outlet section. my simulation is not at all getting converged.
mesh size is --> around 4.5 millions
the y+ value near to the fan blade --> 3
other than fan blade, shroud and other walls y+ value is around 50
turbulence model --> kOmegaSST

fvScheme looks like:

ddtSchemes
{
}

{
default Gauss linear;
}

divSchemes
{
default none;
div(phi,U) Gauss linear;
div(phi,k) Gauss upwind;
div(phi,omega) Gauss upwind;
div(phi,epsilon) Gauss upwind;
}

laplacianSchemes
{
default Gauss linear corrected;
}

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

{
default corrected;
}

fluxRequired
{
default no;
p ;
}

for wall function:
k --> kqRWallFunction;
epsilon --> epsilonWallFunction;
nut --> nutkWallFunction;

can you guys help me to improve my simulation?
any sort of hit and tips are appreciated.