
[Sponsors] 
November 24, 2007, 11:42 
I'm working on a NACA0012 in 2

#1 
Senior Member
Emanuele
Join Date: Mar 2009
Posts: 110
Rep Power: 9 
I'm working on a NACA0012 in 2D and i must calculate the lift and drag coefficent. Can someone post here a complete liftDrag tool with all files needed to compile and install it? Not only the files under directory src/postProcessing/incompressible/liftDrag but also the files under application/utilities/postProcessing/miscellaneous/liftDrag
.Can someone explain me how the tool calculate the lift and drag coefficents? In the forum there are many thread on this argument but no one is complete. Thanks a lot Emanuele 

November 25, 2007, 13:54 
Everything that you've asked f

#2 
Senior Member
Srinath Madhavan (a.k.a pUl)
Join Date: Mar 2009
Location: Edmonton, AB, Canada
Posts: 703
Rep Power: 13 
Everything that you've asked for (and more) is in this board already. Please use the search utility to the left to find it.


December 3, 2007, 11:45 
I think the question can be re

#3 
Senior Member
Alberto Passalacqua
Join Date: Mar 2009
Location: Ames, Iowa, United States
Posts: 1,907
Rep Power: 27 
I think the question can be reformulated as "How can I compile liftDrag on OpenFOAM 1.4 and following?" ;)
Regards, Alberto
__________________
Alberto Passalacqua GeekoCFD  A free distribution based on openSUSE 64 bit with CFD tools, including OpenFOAM. Available as in both physical and virtual formats. OpenQBMM  An opensource implementation of quadraturebased moment methods. To obtain more accurate answers, please specify the version of OpenFOAM you are using. 

December 3, 2007, 17:34 
Dear Emanuele,
I have been

#4 
Member
Alessandro Spadoni
Join Date: Mar 2009
Location: Atlanta, GA
Posts: 65
Rep Power: 9 
Dear Emanuele,
I have been working on a naca0012 as well. These are the things I found out so far: 1. The solution is very very sensitive to your mesh quality. At the moment, I am using a structured mesh with maximum nonorthogonality of 9 deg, average 4 deg. Maximum aspect ratio 808. Max skewness 0.28. You can find out this information by running checkMesh. 2. You should initialize your case with potential foam to get good initial conditions for velocity. Copy potentialFoam/yourcase/0/U into simpleFoam/yourCase/0 3. At this point you should start your simulation with turbulence off (simpleFoam/yourCase/constant/tubulenceProperties) Let the pressure and velocity initial residuals drop to say 1e3, and then turn the turbulence on. Remember that to be able to have your code see changes as it is running you need to have runTimeModifiable yes in simpleFoam/yourcase/system/controlDict. 4. You should converge your results in simpleFoam with very conservative schemes, i.e.: file simpleFoam/yourCase/system/fvSchemes divSchemes { default none; div(phi,U) Gauss upwind; div(phi,k) Gauss upwind; div(phi,epsilon) Gauss upwind; div(phi,omega) Gauss limitedLinear 1; div(phi,R) Gauss upwind; div(R) Gauss linear; div(phi,nuTilda) Gauss upwind; div((nuEff*dev(grad(U).T()))) Gauss linear; } 5. After your solution has converged to 1e6 (residual for all variables), you can go more aggressive with the discretization schemes. for example: file simpleFoam/yourCase/system/fvSchemes divSchemes { default none; div(phi,U) Gauss limitedLinearV 1; div(phi,k) Gauss limitedLinear 1; div(phi,epsilon) Gauss limitedLinear 1; div(phi,omega) Gauss limitedLinear 1; div(phi,R) Gauss linear; div(R) Gauss linear; div(phi,nuTilda) Gauss linear; div((nuEff*dev(grad(U).T()))) Gauss linear; } If you compute your drag at steps 4 and 5, you may notice a difference of 20%. This is due to the better results you get with the more aggressive discretization schemes. 6. Your drag estimate is going to be very sensitive to your boundary conditions at the airfoil also. At the moment I am having very satisfactory results using nutWallFunctions for the airfoil with a y+ of about 7. You can check your y+ during your solution or when the solution is converged issuing checkYplus. I hope this is a good start. Keep me posted about your findings, Alessandro 

December 5, 2007, 07:34 
Good morning Alessandro! I hav

#5 
Senior Member
Emanuele
Join Date: Mar 2009
Posts: 110
Rep Power: 9 
Good morning Alessandro! I have followed your advices but i have a problem with the boundary conditions: where should i set the nutWallFunctions ?
Thanks in advance Emanuele 

December 5, 2007, 07:58 
what turbulence model do you u

#6 
Senior Member
Emanuele
Join Date: Mar 2009
Posts: 110
Rep Power: 9 
what turbulence model do you use?


December 5, 2007, 09:41 
Emanuele,
At the moment, I

#7 
Member
Alessandro Spadoni
Join Date: Mar 2009
Location: Atlanta, GA
Posts: 65
Rep Power: 9 
Emanuele,
At the moment, I am experimenting with: 1. kOmegaSST (HighRe) meaning you use wall functions. 2. LaunderSharmaKe (LowRe) meaning you do not use wall functions. To set wall functions, edit /simpleFoam/yourCase/constant/polymesh/boundary. For case 1: frontAndBack { type empty; ... ... } inlet { type patch; physicalType inlet; ... ... } outlet { type patch physicalType outlet; } airfoil { type wall; physicalType nutWallFunctions; (or physicalType wallFunctions); ... ... } /simpleFoam/yourCase/0/p frontAnfback { type empty; } inlet { type zeroGradient; } outlet { type fixedValue; value uniform 0; } airfoil { type zeroGradient; } /simpleFoam/yourCase/0/U frontAnfback { type empty; } inlet { type fixedValue; value uniform (50 0 0) } outlet { type zeroGradient; } airfoil { type fixedValue; value uniform (0 0 0) } /simpleFoam/yourCase/0/k frontAnfback { type empty; } inlet { type fixedValue; value uniform 0.375; } outlet { type zeroGradient; } airfoil { type zeroGradient; } /simpleFoam/yourCase/0/omega frontAnfback { type empty; } inlet { type fixedValue; value uniform 12800; } outlet { type zeroGradient; } airfoil { type zeroGradient; } For case 2, I am not so sure about how to enforce no wall functions but I think it is like this: frontAndBack { type empty; ... ... } inlet { type patch; physicalType inlet; ... ... } outlet { type patch physicalType outlet; } airfoil { type wall; physicalType wall; ... ... } /simpleFoam/yourCase/0/p frontAnfback { type empty; } inlet { type zeroGradient; } outlet { type fixedValue; value uniform 0; } airfoil { type zeroGradient; } /simpleFoam/yourCase/0/U frontAnfback { type empty; } inlet { type fixedValue; value uniform (50 0 0) } outlet { type zeroGradient; } airfoil { type fixedValue; value uniform (0 0 0) } /simpleFoam/yourCase/0/k frontAnfback { type empty; } inlet { type fixedValue; value uniform 0.375; } outlet { type zeroGradient; } airfoil { type fixedValue; value uniform 1.0e10; } /simpleFoam/yourCase/0/omega or epsilon (since I use LaunderSharmaKE) frontAnfback { type empty; } inlet { type fixedValue; value uniform 12800; Make sure this value is consistent with your analysis } outlet { type zeroGradient; } airfoil { type fixedValue; value uniform 1.0e10; } I hope this helps, Alessandro 

December 5, 2007, 12:35 
Hello,
I thought physicalTy

#8 
Member
Rosario Russo
Join Date: Mar 2009
Location: Trieste, Italy
Posts: 56
Rep Power: 9 
Hello,
I thought physicalType was only used in FoamX and that it didn't have any effect on solvers. Wall functions should be used only if an high Reynolds turbulent model is employed otherwise they are not used, regardless of what physicalType" one is using. Please correct me if I'm wrong or something has been changed in openfoam about this. Ciao. 

December 6, 2007, 10:11 
this is the output of checkYPl

#9 
Senior Member
Emanuele
Join Date: Mar 2009
Posts: 110
Rep Power: 9 
this is the output of checkYPlus:
Create mesh for time = 1000 Time = 1000 Reading field p Reading field U Reading/calculating face flux field phi Selecting incompressible transport model Newtonian Selecting turbulence model kOmegaSST Patch 3 named ala y+ : min: 0.0658652 max: 102.546 average: 65.7859 End 

December 6, 2007, 11:24 
Emanuele,
The inlet velocit

#10 
Member
Alessandro Spadoni
Join Date: Mar 2009
Location: Atlanta, GA
Posts: 65
Rep Power: 9 
Emanuele,
The inlet velocity seems very low. What Re are you trying to obtain? What is your kinematic viscosity? Also, remember that in liftDrag a value of Aref is computed, which is equivalent to the wetted area of your wing, but your are interested in your sectional drag coefficient I assume. This is what I do: 1. Set Aref=1 in src/postProcessing/incompressible/liftDrag/liftDrag.C recomplile your liftdrag tool by issuing ./Allwmake. Then divide the drag coefficient that you obtain this way by your depth in the zdirection. Make sure that in applications/postProcessing/miscellaneous/liftDrag/liftdrag.C turbDragCoeffcient is called, instead of the default dragCoefficient. If you change the function call I mentioned, recompile by issuing wclean && wmake. As an after thought, does your velocity field look reasonable? What about your pressure field? From "Theory of wing sections", you get the maximum velocity over a naca0012 at 0deg aoa as v_local/Vinf = 1.188. Does your velocity field agree with this? Let me know, Alessandro 

December 6, 2007, 12:58 
Hi Alessandro,
in transportPr

#11 
Senior Member
Emanuele
Join Date: Mar 2009
Posts: 110
Rep Power: 9 
Hi Alessandro,
in transportProperties file there is this line: nu nu [0 2 1 0 0 0 0] 1.4e06; I'm trying to obtain Re = 1e+06 If my depth in the zdirection is 1, can i unchange the liftDrag.C file? The velocity field looks reasonable in paraView also the pressure field. 

December 6, 2007, 13:15 
Emanuele,
What is the chord

#12 
Member
Alessandro Spadoni
Join Date: Mar 2009
Location: Atlanta, GA
Posts: 65
Rep Power: 9 
Emanuele,
What is the chord length of your airfoil? To get Re=1e6, given your nu, it looks like the chord is ~ 0.72. So, the aspect ratio of your faces in the zdirection is very high. If you search the postings, you will find a posting where Dr. Jasak mentioned the procedure that is used to compute the drag coefficient; he also mentioned that while making your zdepth = 1 is appealing, since you don't have do recompute your Cd (by normalizing about zdirection, if Aref = 1), this leads to significant errors. I think you should reduce your zdimension to something like 1%10% of your chord, to be safe. Furthermore, I do not know the details of the Navierstokes formulation in OpenFOAM, but I remember reading that solving the NavierStokes equations at really low speed may be difficult (I think there are special formulations for this, but look this up I am not a CFD expert). So I would suggest changing the inlet velocity to something reasonable like 50, and playing with nu to get your desired Re. Alessandro 

December 6, 2007, 13:18 
Emanuele,
Sorry about the c

#13 
Member
Alessandro Spadoni
Join Date: Mar 2009
Location: Atlanta, GA
Posts: 65
Rep Power: 9 
Emanuele,
Sorry about the chord length; This is what I mean Re=chord*U/nu; I messed that up in my previous posting. I see your chord must be 1, right? Alessandro 

December 6, 2007, 14:01 
i have change my zdepth to 1%

#14 
Senior Member
Emanuele
Join Date: Mar 2009
Posts: 110
Rep Power: 9 
i have change my zdepth to 1% of chord (my chord's lenght is 1 > zdepth 0.01) and i've adjusted nu to obtain Re = 1e6 with a velocity inlet of (50*cos(18°), 50*sin(18°), 0) supposing an aoa = 18 degree but the lift coefficent is too low.


December 6, 2007, 14:13 
Emanuele,
With aoa = 18 deg

#15 
Member
Alessandro Spadoni
Join Date: Mar 2009
Location: Atlanta, GA
Posts: 65
Rep Power: 9 
Emanuele,
With aoa = 18 deg, you have massive separation. SimpleFoam, being a steady solver, is not applicable to this case. You should use turbFoam. For a naca0012, I think you are safe using simpleFoam up to 89 deg maximum. Alessandro 

December 10, 2007, 14:48 
Alessandro,
with simpleFoam

#16 
Senior Member
Emanuele
Join Date: Mar 2009
Posts: 110
Rep Power: 9 
Alessandro,
with simpleFoam now i obtain a reasonable value of lift coefficent from 0 to 6 deg. I'm traying to use turbFoam for aoa > 6 deg but i have always the floating point error. In what way can i solve this error?? Thanks a lot Emanuele 

December 10, 2007, 18:06 
Emanuele,
Could you be more

#17 
Member
Alessandro Spadoni
Join Date: Mar 2009
Location: Atlanta, GA
Posts: 65
Rep Power: 9 
Emanuele,
Could you be more specific about your turbFoam case: 1. turbulence model 2. checkMesh output and y+ 3. Boundary conditions 4. Initial conditions 5. deltaT and Courant number you observe during the simulation. 6. fvSchemes and fvSolution specifications. Alessandro 

December 11, 2007, 22:19 
Emanuele,
I don't know whet

#18 
Member
Alessandro Spadoni
Join Date: Mar 2009
Location: Atlanta, GA
Posts: 65
Rep Power: 9 
Emanuele,
I don't know whether or not the spalartAllmaras uses wall functions. If it does, you should set zeroGradient for all turbulence variables at the boundary patch ala like you have done. If spalartAllmaras does not use wall functions, you should set all turbulence variables to 1e10 at the patch ala. I think there is something wrong since your y+ is 0 everywhere. How do you specify nuTilda? I think you need nuTilda for the spalartAllmaras turbulence model. In order to find out which variables you need to specify for the SpalartAllmaras turbulence model edit /src/turbulenceModels/incompressible/SpalartAllmaras/SpalartAllmaras.C At line 125 you can see that this turbulence model definitely wants to read nuTilda, so you need to specify it. for external flows, (please double check what I am saying) nuTilda = beta * nu, where 1 <= beta <= 10. For more details see: http://www.cfdonline.com/OpenFOAM_D...es/1/3476.html fvSchemes and fvSolutions look good. Later on when you get your model working you might what to tighten up convergence tolerances to at least 1e6. I am no expert in CFD and much less in turbulence modeling, but are you sure you can use spalartAllmaras for separated flows? (I assume you are still trying to analyze your naca0012 at 18deg aoa). If you are still trying to simulate a naca0012 at 18deg aoa, I think you would be better off with kOmegaSST. I assume you farfield1 patch is the inlet correct? I hope this helps, Alessandro 

December 12, 2007, 12:31 
Hello,
yPlus are equal to

#19 
Member
Rosario Russo
Join Date: Mar 2009
Location: Trieste, Italy
Posts: 56
Rep Power: 9 
Hello,
yPlus are equal to yPlus_* =Cmu^0.25 sqrt(k) yp / nu so for Spalart Allmaras the field is always zero since in that model k is not used and its value is set to zero. On the other hand I think it could be useful to have the "other" yPlus (without *): yPlus= u_* yp / nu yPlus_* and yPlus should be comparable in equilibrium turbulent boundary layers. I'm attaching here an application which computes both yPlus* and yPlus and writes them on files so that the fields can be visualized by paraFoam. I hope it can be useful to anybody, and/or that anybody will correct it if I made some mistakes ;) . computeYplus.tgz Just one more thing: in the liftDrag utility the computation of the turbulent force is mesleading it is better to compute laminar and turbulent forces at the same time: laminarTurbForce += gSum ( turbulence>nuEff()().boundaryField()[patchI]* U.boundaryField()[patchI].snGrad()* mesh.magSf().boundaryField()[patchI] ); This explains why in the automotive test case thread I mentioned strange values obtained for turbulent forces. Ciao. 

December 13, 2007, 00:50 
Rosario,
Thank you for your

#20 
Member
Alessandro Spadoni
Join Date: Mar 2009
Location: Atlanta, GA
Posts: 65
Rep Power: 9 
Rosario,
Thank you for your elucidations on the SpalartAllmaras model, and the associated method to compute y+. Now I have a question for you. I have been going crazy with trying to understand wallFunctions. I now understand that there are HiRe models that do use wallFunctions such as kOmega, kEpsilon on so forth. I experimented with LaunderSharma (no wallFunctions) but I get a really bad looking velocity field at the leading edge of my naca0012. I went through all threads discussing wallFunctions, and the only thing I saw was people enforcing the use of them with physicalType .... As you have pointed out though, physicalType specifications are only used by FoamX. Do hiRe turbulence models automatically use wallFunctions, no matter what you specify? what about lowRe turbulence models? Do you set you wall patch as type wall? Thank you, Alessandro 

Thread Tools  
Display Modes  


Similar Threads  
Thread  Thread Starter  Forum  Replies  Last Post 
Question to liftDrag  hoochie  OpenFOAM PostProcessing  29  September 19, 2014 03:38 
LiftDrag for 141  ryan_m  OpenFOAM Running, Solving & CFD  2  August 24, 2009 21:26 
Liftdrag calculation  marco  OpenFOAM PostProcessing  10  March 6, 2009 10:51 
LiftDrag coefficient in LES  fabian_korn  OpenFOAM PostProcessing  1  September 22, 2008 02:34 
LiftDrag utility not available  guggi  OpenFOAM Running, Solving & CFD  1  August 2, 2006 12:36 