[jp3g12]

Hello,

I am new to OpenFOAM (I have been using it for about two weeks now) and I decided to start learning the software by trying to simulate a simple case such as a NACA 0012 airfoil at an angle of attack of zero degrees. The specification of the problem can be found here: http://turbmodels.larc.nasa.gov/naca0012_val.html

I used the mesh they provide (I converted it to OpenFOAM format using pointwise).

I have been working for a number of days on this problem without much success. I have attached the residuals and the force coefficients. As it can be seen, the lift is pretty much converged after 6000 iterations. However, I have been experiencing difficulties with drag convergence. In essence, it is going down but it is taking forever to do so! Moreover, the residuals are not that great either (ignore the xlabel, the number of iterations should be 6000).

I am not sure what is wrong with my set up. I stumbled upon this thread http://www.cfd-online.com/Forums/ope...t-airfoil.html and I modified some of the files on the system folder accordingly but it did not help with convergence of results.

I suspect there may be a mistake or something wrong with the settings on the 0 folder. As I said, I am new to OpenFOAM and I am not used to this environment at all. I am using the Spalart-Allmaras. I have attached the system folder and the 0 folder, and perhaps someone can take a look and tell me possible mistakes or things you would do differently. On the nuTilda folder I am using a value that is suggested by Spalart and Rumsey for external aerodynamics calculations based on the velocity that I am using (52 m/s -- Reynolds number = 6,000,000). Spalart and Rumsey suggest a value of nuTilda/nut = 3.

Any help is very much appreciated, especially if someone has attempted the same NASA problem before!

Thank you very much for the help,

JP

[msuaeronautics]

As you may know, it's generally quite difficult regardless to accurately capture drag in a stable manner. The NASA grids look fine to me- though my eyes aren't quite trained- so I cannot comment on that.

I have run a 0012 sim and noticed that this parameter definitely depends on how you set nut and nuTilda for the SA model. Thus, I also suggest you try other models, out of curiosity, to observe the turbulence init cond. effects. Otherwise, it's worth waiting for.

[immortality]

hi
Send me your case to test it.

[jp3g12]

I have tried also the SST with similar results (convergence wise in drag and residuals). I am assuming that the grids NASA provides are good; as far as I could see, the mesh looked fine to me, although I don't know if the conversion from pointwise to OpenFOAM can cause any problems.

[immortality]

Hi
as I replied you the file is too big to download for me.delete all files in polyMesh folder except to blockMeshDict(they can be created later by blockMesh command),and also all time folders except to 0 folder so that i can run a test.

[FelixL]

Hello James,


your boundary conditions seem incorrect. This may cause the errorenous results.

Two critical things:

a.) your simulation is incompressible - therefore using dirichlet bcs (fixedValue) for both velocity and pressure at the same boundary is inconsistent. You should set

u -> fixedValue
p -> zeroGradient

at the inlet and

u -> zeroGradient
p -> fixedValue

at the outlet.

b.) nut is derived from the turbulent quantity nuTida, so setting it to zero at all boundaries may lead to inconsistencies. Especially at the outlet this setting is wrong, because in the wake (which usually is turbulent), nut is clearly non-zero. You can solve this by setting nut to calculated everywhere.


I hope this helps,


Greetings
Felix

[jp3g12]

Felix,

Thank you very much for your input. I will try your suggestions and I will let you know the outcome. I have another question: I have observed that the meshes NASA provides have very large AR in the wake of the airfoil. Should I be concerned about this? can this lead to trouble?

Thank you very much for the help,

James

[FelixL]

Hello James,


grid aspect ratio is not really critical in regions, where the flow is highly anisotropic (e.g. boundary layers or the wake). Ideally, the high aspect ratio grid region is aligned to the local streamlines - but this is usually the case in boundary layer meshes.

I recommend you to perform studies on a simple case (flat plate, for example) or have a look in the literature if you want to find out, when high aspect ratios become important.


Greetings
Felix

[immortality]

Hi Felix
could you explain a bit that when a flow is anisotropic?
I knew that but I have forgotten the subject.thanks.

[jp3g12]

Hello Felix,

So I have managed to do a simulation with your suggestions and here are the results (find attached a plot of the Residuals and the Force Coefficients). Looking at the residuals I'd say the solution is pretty much converged as confirmed by the force coefficient plot. However, as it can be seen, the drag value is very very large (about 1.39764) and the lift coefficient is (0.00033) which is higher than what I would expect from a NACA 0012. I am forwarding you a link where I have posted my case (0, constant and system folder) in case you want to take a look at my setup.

https://www.dropbox.com/sh/ylxws1mx6w3ki9f/2mDt6r0CHw

To give you some background as to what I am doing, I am testing the NACA 0012 for a Reynolds number of 6,000,000 and zero AOA using the SA turbulence model. I am using the grid with 2 x 449 x 129 (2 x 257 points on airfoil surface) mesh NASA provides (the middle one). I am using standard values for rho and nut and I calculated my speed based on the Re number. The nutilda value I obtained from a paper from Spalart and Rumsey where they advise using nutilda/nut = 3 to 5. I think I am using in this particular case 5.

On the system folder I am using two different methods for the fvSchemes. First I ran it with the divSchemes as:

 div(phi,U)      Gauss upwind;
    div(phi,nuTilda) Gauss upwind;

Once the residuals have dropped significantly I switch to

     div(phi,U)      Gauss linearUpwind grad(U);
    div(phi,nuTilda) Gauss linearUpwind grad(nuTilda);

Can you spot anything fundamentally wrong that would cause the solution to be so far from the values that can be expected?? Can it be that the domain is GIGANTIC?

Any help is appreciated!!! Thank you very much!

James

[FelixL]

Hello James,


for further comments, can you please provide us with your checkMesh results?


Greetings
Felix

[jp3g12]

Hi Felix,

Here is the checkMesh results... As I said, it complains about high aspect ratio cells. They are mostly located on the wake of the airfoil towards the outlet. Also high aspect ratio cells can be found near the airfoil surface but I suppose that is unavoidable.

/*---------------------------------------------------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  2.1.1                                 |
|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
Build  : 2.1.1-221db2718bbb
Exec   : checkMesh
Date   : Jun 04 2013
Time   : 12:51:58
Host   : "lyceum2"
PID    : 1883
Case   : /home/jp3g12/OpenFOAM/jp3g12-2.1.1/run/sa/airfoil-med
nProcs : 1
sigFpe : Enabling floating point exception trapping (FOAM_SIGFPE).
fileModificationChecking : Monitoring run-time modified files using timeStampMaster
allowSystemOperations : Disallowing user-supplied system call operations
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Create time
Create polyMesh for time = 0
Time = 0
Mesh stats
    points:           461056
    internal points:  0
    faces:            918656
    internal faces:   457600
    cells:            229376
    boundary patches: 4
    point zones:      0
    face zones:       0
    cell zones:       0
Overall number of cells of each type:
    hexahedra:     229376
    prisms:        0
    wedges:        0
    pyramids:      0
    tet wedges:    0
    tetrahedra:    0
    polyhedra:     0
Checking topology...
    Boundary definition OK.
    Cell to face addressing 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                  
    frontAndBack        458752   461056   ok (non-closed singly connected)  
    inlet               896      1794     ok (non-closed singly connected)  
    outlet              512      1026     ok (non-closed singly connected)  
    wall                896      1792     ok (non-closed singly connected)  
Checking geometry...
    Overall domain bounding box (-484.457 -1 -507.807) (501 0 507.807)
    Mesh (non-empty, non-wedge) directions (1 0 1)
    Mesh (non-empty) directions (1 0 1)
    All edges aligned with or perpendicular to non-empty directions.
    Boundary openness (-1.61899e-20 -4.24013e-16 -2.5902e-19) OK.
 ***High aspect ratio cells found, Max aspect ratio: 3.17847e+07, number of cells 29164
  <<Writing 29164 cells with high aspect ratio to set highAspectRatioCells
    Minumum face area = 9.5796e-11. Maximum face area = 843.777.  Face area magnitudes OK.
    Min volume = 9.5796e-11. Max volume = 843.777.  Total volume = 875558.  Cell volumes OK.
    Mesh non-orthogonality Max: 19.7377 average: 1.63527
    Non-orthogonality check OK.
    Face pyramids OK.
    Max skewness = 0.201756 OK.
    Coupled point location match (average 0) OK.
Failed 1 mesh checks.
End

Thanks for the Help!!

[FelixL]

Hello James,


your boundary Definition is wrong. Your mesh is not stitched together in the wake behind the airfoil. You have some kind of long flat plate with a tiny airfoil at the leading edge. I found out about this, because the boundary file says 896 faces for patch wall instead of 256, as the nada description of the mesh suggests.

Please correct this and everything should be fine now, I double checked everything else.

A hint for the future: always take a look at the mesh and boundaries using paraFoam before running a case. In 50% of all problematic cases its the mesh or boundary definition, in 30% its the boundary conditions, 15% numerical settings and 5% its OpenFOAM


Greetings
Felix

[jp3g12]

Hi Felix,

OH DAMN!!! Yep you are very much correct . I just inspected the mesh in pointwise (everything seemed to be Ok in paraview) and it looks like the airfoil has a flat plate coming out of the trailing edge and it extends all the way to the outlet.

Even though this question is a bit off topic, how the heck am I supposed to define the boundary condition at the airfoil if it is attached to this flat plate?? It appears as though the airfoil and that "flat plate extension" are just one part so that whenever you select the airfoil to prescribe the wall BC you are also selecting the flat plate. Clearly setting wall BC does not work. how did you do it when you attempted this problem?

Nice catch!!! Cheers,

James

[FelixL]

Hello James,


you can see you individual boundary patches in paraView by selecting and deselecting them in the patch-tree (don't know how it's called right now) in paraView. Then you should see the mishap with the boundary patch definition.

As for resolving the problem: that's up to you or someone, who is familiar with pointwise. I used the OF Tools plot3dToFoam and autoPatch to convert the NASA meshes, but the procedure is difficult to explain. Just play with the OF provided tools and you'll find your own procedure. Since you are a beginner with OF, it doesn't hurt to get to know the mesh utilities.


Greetings
Felix

[jp3g12]

Thanks a lot for the help Felix!!!

I think I have figured out how to isolate the airfoil geometry from the flat plate extension in pointwise (crossing fingers)... I am running a simulation right now, and I will post the outcomes of that once I am done!! Hopefully the issues are resolved now !

Cheers,

James

[jp3g12]

Felix,

All the changes you suggested WORKED perfectly (at least for the zero angle of attack case)!!! Now I get values for the lift and drag that are very close to what NASA gives!!!

Thank you very much for taking the time to help me out! Now, on to attempt the same problem with different turbulence models !!

Thanks a lot

James

[Combas]

Hello James,

I am also trying to validate OpenFoam on an airfoil, and my first results are not very good on the drag coefficient.

Since you say that you got good results with the NACA0012, would it be possible that:
- you give the Cd you get at AoA = 0°
- you share the input files of your case (folders "0", "system" and "constant")

Thank you in advance for your help!
Laurent

[K31th]

Hi James,

I faced the same issue on your last part, where the airfoil has a flat plate coming out of the trailing edge and it extends all the way to the outlet.
wonder how did you resolve the problem ?

[K31th]

Hi Felix,

I did the same way as you ( plot3dToFoam and autoPatch on NASA mesh), however i face the problem of the airfoil has a flat plate coming out of the trailing edge and it extends all the way to the outlet.

I wonder if i miss out any step in between ?