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

3D simulation with blended wingbody. Lift a 1/3 of what it shuld be. Any wrong?

Register Blogs Members List Search Today's Posts Mark Forums Read

Reply
 
LinkBack Thread Tools Display Modes
Old   May 27, 2010, 18:39
Default 3D simulation with blended wingbody. Lift a 1/3 of what it shuld be. Any wrong?
  #1
Member
 
Rickard Hidefjäll
Join Date: Jun 2009
Location: Uppsala, Sweden
Posts: 39
Rep Power: 8
SkunkWorks is on a distinguished road
Hi!
is there any wrong in my simulation ?
To little lift in patchIntegrate, about a 1/3 of what expected.

I`ve made a simulation with simpleFoam, and I let 1000 steps be calculated. dt=1s

The wind speed is 20 m/s and my nu is 1.5*10^-5.

I have 5 degree of angel of attack, and wingtip 3 degree.

Lift = Cl*Wingarea*density*speed²/2=33 N
Cl=0.5 @ Re=5*10⁵ and angel=5 for Naca0008.
Wingarea=0.28m²
Dencity 1.2 kg/m³ @ 20 celcius
Speed 20m/s
__________________________________________________ __________________
Time = 1000
Area vector of patch wing[0] = (-1.55176e-18 7.36813e-19 1.1418e-18)
Area magnitude of patch wing[0] = 0.570385
Reading volScalarField p
Integral of p over vector area of patch wing[0] = (-0.000573922 -3.05068 21.1538)
Integral of p over area magnitude of patch wing[0] = -10.6492

End
__________________________________________________ ___________________
The mesh is about 0.9 million tetrahedra.

checkMesh
__________________________________________________ ___________________
Time = 0

Mesh stats
points: 188033
faces: 1907992
internal faces: 1782516
cells: 922627
boundary patches: 4
point zones: 0
face zones: 0
cell zones: 0

Overall number of cells of each type:
hexahedra: 0
prisms: 0
wedges: 0
pyramids: 0
tet wedges: 0
tetrahedra: 922627
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
wing 82492 41248 ok (closed singly connected)
inlet 5122 2646 ok (non-closed singly connected)
outlet 5210 2690 ok (non-closed singly connected)
channel 32652 16494 ok (non-closed singly connected)

Checking geometry...
Overall domain bounding box (-8.4 -12.5 -8) (8.4 12.5 8)
Mesh (non-empty, non-wedge) directions (1 1 1)
Mesh (non-empty) directions (1 1 1)
Boundary openness (1.00339e-17 -9.39416e-17 -1.34001e-17) OK.
Max cell openness = 2.37738e-16 OK.
Max aspect ratio = 9.32888 OK.
Minumum face area = 6.66656e-09. Maximum face area = 0.377916. Face area magnitudes OK.
Min volume = 2.92226e-13. Max volume = 0.0621314. Total volume = 5276.9. Cell volumes OK.
Mesh non-orthogonality Max: 70.3763 average: 16.8812
*Number of severely non-orthogonal faces: 1.
Non-orthogonality check OK.
<<Writing 1 non-orthogonal faces to set nonOrthoFaces
Face pyramids OK.
Max skewness = 1.89516 OK.

Mesh OK.
__________________________________________________ _______________________

Any one who knows what's wrong? Thankful for answer!
It would be nice to know the centre if lift.

Regards / SkunkWorks
SkunkWorks is offline   Reply With Quote

Old   May 31, 2010, 06:45
Default
  #2
Member
 
Andre Z
Join Date: Dec 2009
Posts: 32
Rep Power: 7
LVDH is on a distinguished road
Hi,
there are two things I see on the first look.

1. I do not think your mesh resolution is sufficient. As far as I can see from the picture the domainsize seems quite small as well.

2. 1000 iterations is probably not enough. You should add the liftdrag function to your controldict and plot your lift over the iterations.
LVDH is offline   Reply With Quote

Old   June 1, 2010, 17:22
Default
  #3
Member
 
Rickard Hidefjäll
Join Date: Jun 2009
Location: Uppsala, Sweden
Posts: 39
Rep Power: 8
SkunkWorks is on a distinguished road
Hi and thanks for answer!

The domainsize is a tube which is 25 meter long and 16.8 meter wide and 16 meter high. The air craft is 0.8 meter long a and equal span.

What size of the tube would be sufficient?
And equally what mesh resolutions would be sufficient?

"You should add the liftdrag function to your controldict and plot your lift over the iterations. "
Thanks a lot for your hint

Best regards / Rickard
SkunkWorks is offline   Reply With Quote

Old   June 2, 2010, 01:10
Default
  #4
Member
 
Andre Z
Join Date: Dec 2009
Posts: 32
Rep Power: 7
LVDH is on a distinguished road
Hi,
in that case your domainsize should be OK.

Now you should first generate the plots like I suggested. You have to rerun the simulation for that.

Right now I do not know how many cells your case will need. But you can find out yourself, as long as you have enough time. After you get converged results with your current mesh you refine it and run another simulation. Then you run a few simulations with even finer meshes. When refining does not significantly alter your results you got the exact mesh size you are looking for. But thats only what the books say. You should have cells on your surfaces which should not exceed a few millimeters and you should have enough fine cells in the regions of interest. That would be everything close to your body, the wake and I would suppose in the wake of the wing tip vortices.
LVDH is offline   Reply With Quote

Old   June 2, 2010, 01:21
Default
  #5
Member
 
Andre Z
Join Date: Dec 2009
Posts: 32
Rep Power: 7
LVDH is on a distinguished road
I do not know what your experience is and what your are exactly doing but maybe you want to take a look at the software xflr5.
LVDH is offline   Reply With Quote

Old   June 2, 2010, 07:43
Default
  #6
Member
 
Rickard Hidefjäll
Join Date: Jun 2009
Location: Uppsala, Sweden
Posts: 39
Rep Power: 8
SkunkWorks is on a distinguished road
Hi!


Then I run " patchIntegrate p wing " I can see that I have a converges the last 80 iterations doesn't change the results at all. And the last 440 iterations change from -10.6123 to -10.6492.

I will tweak liftdrag function and it can take some time. The simulation took little more than a day. I haven't parallelised it yet.

The mesh resolution changes (get finer) with the inverse double integrate i.e 1/f´´(x). There f(x) is the function on the wing profile. That means the more curve bending -> finer mesh.
Then I also put extra refinement the as the mesh closes in to the wing tip, in order to get enough cells on the tip. I does the meshing in Gmsh.
Then I does the 3D mesh in enGrid, and export the mesh for OpenFOAM.
I can make the whole mesh for the aircraft finer with a variable, and another variable for the tube/virtual wind channel.
I can also in the .geo file change the angel of attack and the twist of the wing tip.

The form of the air craft is analytical define from 3 equations. one for the front, one for the back and the last for the wing profile.

I will try to make a finer mesh, and see and the liftdrag function

I intend to bulid it later on with a radio control.
It is of a great importance to get the point of the drag, in order to get it flying.
I does the simulation mostly for fun and that I'm very fascinated of CFD:s and aircraft's . I will also keep my math fresh.

xflr5 looks very interesting!
And thanks again.

Best regards / Rickard
SkunkWorks is offline   Reply With Quote

Old   June 10, 2010, 17:16
Default
  #7
Member
 
Rickard Hidefjäll
Join Date: Jun 2009
Location: Uppsala, Sweden
Posts: 39
Rep Power: 8
SkunkWorks is on a distinguished road
Hi !

I have tried the liftDrag function and got 25.86 N
according to patchIntegrate 10.6492 there is a factor 2,428
look at:
forces function VS patchIntegrate

I think it's a strange difference , has any one any explanation?

25.86 N looks more resonoble to 33 N according to the fact that the wingtip has less attacking angel. But 7 N is pretty much, but on the other hand most lift was generated on the outer parts of the wing.
Does any one have an other opinion? It would be interesting.

I have done a mesh with 2.6 million tetrahedrons witch pass checkMesh with very good results, but the simulation explode after 4 steps.
It isn,t the first time I experience that kind of phenomena. It would be very interesting to know way. Any one who knows? it had a max scrawniness 69 degree and average about 14.

For 2 days ago I experience a disk crash so I may have to come back later with a new case.

Regards Rickard
SkunkWorks is offline   Reply With Quote

Old   April 7, 2011, 10:49
Default
  #8
New Member
 
Silvan Brändli
Join Date: Aug 2009
Posts: 27
Rep Power: 7
s_braendli is on a distinguished road
Hi Rickard

Did you figure out the reason for the difference between simulated and expected lift?

I would be very interested in it, as I am facing the same problem with a 2D-naca case.

Best Regards
Silvan

btw: function forces and patchIntegration give me the same result.
s_braendli is offline   Reply With Quote

Old   April 7, 2011, 15:49
Default
  #9
Member
 
Rickard Hidefjäll
Join Date: Jun 2009
Location: Uppsala, Sweden
Posts: 39
Rep Power: 8
SkunkWorks is on a distinguished road
I have done a new similar case.

No sadly not. The patchIntegrate still give the approx. 1/3 of the calculated lift
And the forcefunction give approx. half 3,479 N.

Lift = Cl*Wingarea*density*speed²/2=7,575 N

Cl=0.5 @ Re=5*10⁵ and angel=5 for Naca0008.
Wingarea=0.2525m²
Dencity 1.2 kg/m³ @ 20 celcius
Speed 10m/s

LiftDragFunction=> 3,479 N

patchIntegrate p wing -time 1000

Create time

Create mesh for time = 1.3

Time = 1.3
Area vector of patch wing[0] = (1.35366e-18 6.55577e-19 -1.26087e-17)
Area magnitude of patch wing[0] = 0.515495
Reading volScalarField p
Integral of p over vector area of patch wing[0] = (-0.000187886 -0.231746 2.94984)
Integral of p over area magnitude of patch wing[0] = -2.01594

End

__________________________________________________ _______________
checkMesh =>


Create polyMesh for time = 0

Time = 0

Mesh stats
points: 147363
faces: 1494142
internal faces: 1395830
cells: 722493
boundary patches: 4
point zones: 0
face zones: 0
cell zones: 0

Overall number of cells of each type:
hexahedra: 0
prisms: 0
wedges: 0
pyramids: 0
tet wedges: 0
tetrahedra: 722493
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
wing 87470 43737 ok (closed singly connected)
inlet 1406 748 ok (non-closed singly connected)
outlet 1388 739 ok (non-closed singly connected)
channel 8048 4112 ok (non-closed singly connected)

Checking geometry...
Overall domain bounding box (-7.5 -10 -7) (7.5 10 7)
Mesh (non-empty, non-wedge) directions (1 1 1)
Mesh (non-empty) directions (1 1 1)
Boundary openness (8.65344e-18 1.78066e-16 1.40199e-17) OK.
Max cell openness = 2.7385e-16 OK.
Max aspect ratio = 8.86279 OK.
Minumum face area = 2.33323e-10. Maximum face area = 0.822625. Face area magnitudes OK.
Min volume = 5.25689e-15. Max volume = 0.266388. Total volume = 3296.44. Cell volumes OK.
Mesh non-orthogonality Max: 73.091 average: 17.2453
*Number of severely non-orthogonal faces: 17.
Non-orthogonality check OK.
<<Writing 17 non-orthogonal faces to set nonOrthoFaces
Face pyramids OK.
Max skewness = 1.29086 OK.

Mesh OK.

End
__________________________________________________ __________

/ Regards Rickard
SkunkWorks is offline   Reply With Quote

Old   April 8, 2011, 04:49
Default
  #10
Senior Member
 
Attesz's Avatar
 
Attesz
Join Date: Mar 2009
Posts: 355
Rep Power: 8
Attesz is an unknown quantity at this point
What about your boundary layers? It is very important to have a fine resolution there to get the right lift. Can you post pictures about is (some cutplanes at different distances from the centerline of the body)
Attesz is offline   Reply With Quote

Old   April 8, 2011, 05:36
Default
  #11
New Member
 
Silvan Brändli
Join Date: Aug 2009
Posts: 27
Rep Power: 7
s_braendli is on a distinguished road
Thanks to both of you for the fast replies.

My case is mainly according to the wingMotion-tutorial from 1.7.1

Following are my boundary layers at the tip, the trailing edge (looks a little bit weird there) and the middle of the airfoil.







The y+ I get is between 5.37 and 48.2. I use the kOmegaSST turbulence-model with the wall funtions omegaWallFunction and kqRWallFunction. Are the layers too thin? In literature I see values of y+=10-100 when using wall functions.

Moreover, I'm thinking about using no turbulence model at all, and, as I am comparing with a potential-method, even using a slip-boundary-condition at the wall. Do you have suggestions on that?

Thanks in advance
Silvan
s_braendli is offline   Reply With Quote

Old   April 11, 2011, 08:55
Default
  #12
Senior Member
 
Attesz's Avatar
 
Attesz
Join Date: Mar 2009
Posts: 355
Rep Power: 8
Attesz is an unknown quantity at this point
I can't see your pictures!

I made a quadrotor helicopter simulation and I get half the thrust with the same settings, however it was made using CFX. I did'n know what was the problem there, but I used coarse mesh because of the PC limitations. y+ may be not problem in first step, but you should always resolve the mesh near the wall if you want to study phenomena there (heat transfer, shear, and lift-drag).

Using no turbulence model will not help you, maybe using the spalart-almaras 1eq model (if any) can be useful (it was developped for airfoil flows) to decrease computational needs. I'm not familiar with the potential-method - i can't help you in this.

You should use the k-w sst and refine the mesh and the boundary layer near the walls - it is the first step.

Last edited by Attesz; April 11, 2011 at 09:57.
Attesz is offline   Reply With Quote

Old   April 15, 2011, 04:15
Default
  #13
Member
 
Rickard Hidefjäll
Join Date: Jun 2009
Location: Uppsala, Sweden
Posts: 39
Rep Power: 8
SkunkWorks is on a distinguished road
Sorry for late response! You got the picture and the complete case in:
http://www.tooslow.net/rattus/CFD_Pi...fa5Tord2Test2/

I have try to make a finer mesh but Engrid and netgen make it only worse.

(I'm still uploading so it can take a lite wile.)

Regards / Rickard
SkunkWorks is offline   Reply With Quote

Old   April 18, 2011, 08:17
Default
  #14
Senior Member
 
Attesz's Avatar
 
Attesz
Join Date: Mar 2009
Posts: 355
Rep Power: 8
Attesz is an unknown quantity at this point
Hi,

I don't see inflation layers on the surfaces! You will not get accurate lift and drag forces unless you build an appropriate boundary layer mesh using at least 10 layers in it!

http://www.telusplanet.net/public/dj...es/naca12b.gif
Attesz is offline   Reply With Quote

Old   April 19, 2011, 03:14
Default
  #15
Member
 
Rickard Hidefjäll
Join Date: Jun 2009
Location: Uppsala, Sweden
Posts: 39
Rep Power: 8
SkunkWorks is on a distinguished road
Hi and thanks a lot!

I was afraid of that. I have tried many time to make a prismatic layer in Engrid, but it doesn't work! I think I had to use snappyhexmesh instead.

Regards / Rickard
SkunkWorks is offline   Reply With Quote

Old   April 19, 2011, 06:54
Default
  #16
Senior Member
 
Attesz's Avatar
 
Attesz
Join Date: Mar 2009
Posts: 355
Rep Power: 8
Attesz is an unknown quantity at this point
I've no experience in the accuracy of the snappy boundary layers. Extruded layers (hexa or prism) is always better. Use expansion factor 1.2 as upper limit, and at least 5 layers. 10 is a good start and 15 would be the best.
Attesz is offline   Reply With Quote

Old   April 19, 2011, 10:37
Default
  #17
New Member
 
Join Date: Mar 2011
Posts: 17
Rep Power: 6
nicolarre is on a distinguished road
Quote:
Originally Posted by s_braendli View Post
Hi Rickard

Did you figure out the reason for the difference between simulated and expected lift?

I would be very interested in it, as I am facing the same problem with a 2D-naca case.

Best Regards
Silvan

btw: function forces and patchIntegration give me the same result.
Im doing some model validation runs with a 2D NACA0009 @ 5º attack angle.

Lift error around 8~9%. Drag, however, nets a whoopn' 370~420% Error.

the Mesh is good. 670.000+ cells and Y+ ranging from 9 to 70

As for the Forces Vs Patchintegrate, the results i get are consistent:

PatchIntegrate * Rho = Forces (Pressure)

Last edited by nicolarre; April 19, 2011 at 11:15.
nicolarre is offline   Reply With Quote

Old   May 4, 2011, 09:20
Default
  #18
Member
 
Rickard Hidefjäll
Join Date: Jun 2009
Location: Uppsala, Sweden
Posts: 39
Rep Power: 8
SkunkWorks is on a distinguished road
Hi !

I've did a last try with enGrid, and did manage to make a prismatic layer.

But I got a bad mesh quality. I did 4 meshes 2 with 10 layer and 2 with 15.
The one with 15 layer bailed out with nNonOrthogonalCorrectors 2;
10 layer ones did converge to cd=0.21 as my small mesh 0.7 million tetrahedrals without any layer.
I did a re-mesh 2 times in enGrid but I've got the same number non-orthogonal faces and Mesh non-orthogonality Max. so I think that belongs to the prismatic layer.
I don't know if a did a bad wrong or the enGrid does a bad mesh when you used prismatic layer. I did also extrude in th mesh.

Does any one have a good clue or advice?

The 4 cases I did is identical to the ones above except for the nNonOrthogonalCorrectors is 2 and not 0.

Here is the checkMesh output. in case of interest.

Best regards / Rickard

rickard@SkunkWorks:/Data/FlygFiler/Projekt20110401/Projekt20110401Alfa5Tord2Layer2$ startFoam171
rickard@SkunkWorks:/Data/FlygFiler/Projekt20110401/Projekt20110401Alfa5Tord2Layer2$ checkMesh
/*---------------------------------------------------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: 1.7.x |
| \\ / A nd | Web: www.OpenFOAM.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
Build : 1.7.x-3776603e4c6c
Exec : checkMesh
Date : May 04 2011
Time : 14:56:59
Host : SkunkWorks
PID : 7660
Case : /Data/FlygFiler/Projekt20110401/Projekt20110401Alfa5Tord2Layer2
nProcs : 1
SigFpe : Enabling floating point exception trapping (FOAM_SIGFPE).

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Create time

Create polyMesh for time = 0

Time = 0

Mesh stats
points: 1404594
faces: 11423843
internal faces: 11290687
cells: 5334445
boundary patches: 4
point zones: 0
face zones: 0
cell zones: 0

Overall number of cells of each type:
hexahedra: 34200
prisms: 1308350
wedges: 0
pyramids: 0
tet wedges: 0
tetrahedra: 3991895
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
wing 124092 62048 ok (closed singly connected)
inlet 1810 1324 ok (non-closed singly connected)
outlet 1806 1322 ok (non-closed singly connected)
channel 5448 4510 ok (non-closed singly connected)

Checking geometry...
Overall domain bounding box (-5.79793 -26.5868 -5.29793) (5.79793 20.6688 5.29793)
Mesh (non-empty, non-wedge) directions (1 1 1)
Mesh (non-empty) directions (1 1 1)
Boundary openness (3.40654e-18 -1.35787e-17 -3.13522e-17) OK.
Max cell openness = 1.30801e-14 OK.
Max aspect ratio = 285.049 OK.
Minumum face area = 3.54516e-11. Maximum face area = 1.57032. Face area magnitudes OK.
Min volume = 1.30018e-16. Max volume = 0.384325. Total volume = 4556.42. Cell volumes OK.
Mesh non-orthogonality Max: 88.6478 average: 14.9831
*Number of severely non-orthogonal faces: 7166.
Non-orthogonality check OK.
<<Writing 7166 non-orthogonal faces to set nonOrthoFaces
Face pyramids OK.
***Max skewness = 23.6909, 1115 highly skew faces detected which may impair the quality of the results
<<Writing 1115 skew faces to set skewFaces

Failed 1 mesh checks.

End
SkunkWorks is offline   Reply With Quote

Old   May 4, 2011, 09:52
Default
  #19
Senior Member
 
Attesz's Avatar
 
Attesz
Join Date: Mar 2009
Posts: 355
Rep Power: 8
Attesz is an unknown quantity at this point
Visually check the mesh, and search for the highly skewed cells, which can be for example near to highly curved faces. Otherwise about 1000 cells compared to 5millions can not cause inaccuracy. if the cells are too bad, the solver fails, otherwise the number of the bad cells are small (of course you have to avoid them).
Other advices:
Check for yplus. it should be under 2 at most of the surfaces of the wing
Also the measurements can be inaccurate.
Give us details about the physical setup, solvers, boundaries. Send a picture about the convergence if you can.
However this is very interesting for me, because I've already met with the same problem. I used ANSYS CFX, but in our case it doesn't matter.

Regards
Attesz is offline   Reply With Quote

Old   May 8, 2011, 16:40
Default
  #20
Member
 
Rickard Hidefjäll
Join Date: Jun 2009
Location: Uppsala, Sweden
Posts: 39
Rep Power: 8
SkunkWorks is on a distinguished road
As I wrote the case is identical except for nNonOrthogonalCorrectors is 2 instead for 0.
and the mesh. The picture is also at the the. http://www.tooslow.net/rattus/CFD_Pi...fa5Tord2Test2/
I also putted the enGrid files: 20110407__50_Exp50_150_NL025lc035_Alfa5Tord2__5_33 4445Mt_Iter3_MultiVingprofile_Short4_Narrow4_Layer 10_Iter15_Extrude.egc and 20110407__50_Exp50_150_NL025lc035_Alfa5Tord2__5_33 4445Mt_Iter3_MultiVingprofile_Short4_Narrow4_Layer 10_Iter15_Extrude.egc.vtu
in the folder so you can run the case by yourself

and

Here you have the graph of the convergens of the simulation.


My experience is that a few bad cells can very quickly damage the simulation.

Best regards / Rickard
SkunkWorks is offline   Reply With Quote

Reply

Thread Tools
Display Modes

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 On
Pingbacks are On
Refbacks are On


Similar Threads
Thread Thread Starter Forum Replies Last Post
How obtain the average of lift over time for a transient simulation? aero ANSYS 0 November 11, 2009 03:00
udf error srihari FLUENT 0 February 9, 2009 10:00
Correct lift but wrong pressure drag - possible? zx Main CFD Forum 4 July 27, 2007 23:38


All times are GMT -4. The time now is 02:34.