
[Sponsors] 
January 6, 2014, 11:22 
Too large negative power coefficient for VAWT

#1 
Member
Join Date: Nov 2012
Posts: 62
Rep Power: 5 
Hello,
I am currently using OpenFOAM for 3D simulation of Vertical Axis Wind Turbine. Actually I moved on to 3D simulation after 2D simulation. My current problem is I am getting too large negative Cp. Infact sometimes it may go beyond 1. Here is the result of my CheckMesh utility: /**\  =========    \\ / F ield  OpenFOAM: The Open Source CFD Toolbox   \\ / O peration  Version: 2.2.2   \\ / A nd  Web: www.OpenFOAM.org   \\/ M anipulation   \**/ Build : 2.2.29739c53ec43f Exec : checkMesh Date : Jan 06 2014 Time : 23:47:53 Host : "hpcmaster01" PID : 17376 Case : /home/Turbine_BUET_0_exp. nProcs : 1 sigFpe : Floating point exception trapping  not supported on this platform fileModificationChecking : Monitoring runtime modified files using timeStampMaster allowSystemOperations : Disallowing usersupplied system call operations // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // Create time Create polyMesh for time = 0 Time = 0 Mesh stats points: 1001176 faces: 6019344 internal faces: 5930454 cells: 2665263 faces per cell: 4.48353 boundary patches: 10 point zones: 0 face zones: 0 cell zones: 2 Overall number of cells of each type: hexahedra: 603260 prisms: 0 wedges: 0 pyramids: 82226 tet wedges: 0 tetrahedra: 1979777 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: 2 The mesh has multiple regions which are not connected by any face. <<Writing region information to "0/cellToRegion" Checking patch topology for multiply connected surfaces... Patch Faces Points Surface topology blade2 31703 30935 ok (closed singly connected) Rotor 5144 2692 ok (nonclosed singly connected) blade1 31703 30935 ok (closed singly connected) RotorAMI 5782 5900 ok (nonclosed singly connected) topAndBottom 2678 1447 ok (nonclosed singly connected) Stator 5144 2692 ok (nonclosed singly connected) inlet 171 200 ok (nonclosed singly connected) sides 612 700 ok (nonclosed singly connected) outlet 171 200 ok (nonclosed singly connected) StatorAMI 5782 5900 ok (nonclosed singly connected) Checking geometry... Overall domain bounding box (15 10 2.5) (25 10 3) Mesh (nonempty, nonwedge) directions (1 1 1) Mesh (nonempty) directions (1 1 1) Boundary openness (1.48176e18 7.67635e17 1.5866e18) OK. Max cell openness = 3.39725e16 OK. Max aspect ratio = 61.3245 OK. Minimum face area = 2.18162e08. Maximum face area = 1.57259. Face area magnitudes OK. Min volume = 1.75721e12. Max volume = 0.535229. Total volume = 4400. Cell volumes OK. Mesh nonorthogonality Max: 70.075 average: 19.2599 *Number of severely nonorthogonal (> 70 degrees) faces: 2. Nonorthogonality check OK. <<Writing 2 nonorthogonal faces to set nonOrthoFaces Face pyramids OK. Max skewness = 1.66154 OK. Coupled point location match (average 0) OK. Mesh OK. End As you may see from the checkMesh utility the mesh is fine. But there is only one thing that may be a reason for abnormality and that is: "The mesh has multiple regions which are not connected by any face." What does it mean? I got the same message while 2D simulation too. But in case of 2D I got good results. I have modelled my whole problem after the propeller tutorial. In my opinion the problem may arise from the fvSchemes or fvSolutions too. I am posting it here. FoamFile { version 2.0; format ascii; class dictionary; location "system"; object fvSchemes; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // ddtSchemes { default Euler; } gradSchemes { default Gauss linear; grad(p) Gauss linear; grad(U) cellLimited Gauss linear 1; } divSchemes { default none; // div(phi,U) Gauss upwind; div(phi,U) Gauss linearUpwind grad(U); div(phi,k) Gauss upwind; div(phi,omega) Gauss upwind; div((nuEff*dev(T(grad(U))))) Gauss linear; } laplacianSchemes { default Gauss linear corrected; } interpolationSchemes { default linear; } snGradSchemes { default corrected; } fluxRequired { default no; pcorr ; p ; } // solvers { pcorr { solver GAMG; tolerance 1e7; relTol 0.005; smoother GaussSeidel; nPreSweeps 0; nPostSweeps 2; cacheAgglomeration off; agglomerator faceAreaPair; nCellsInCoarsestLevel 20; mergeLevels 1; maxIter 100; minIter 1; } p { $pcorr; tolerance 1e5; relTol 0.01; } pFinal { $p; tolerance 1e6; relTol 0; } "(UkepsilonomeganuTilda)" { solver smoothSolver; smoother GaussSeidel; nSweeps 1; tolerance 1e07; relTol 0.05; maxIter 100; minIter 1; }; "(UkepsilonomeganuTilda)Final" { solver smoothSolver; smoother GaussSeidel; nSweeps 2; tolerance 1e07; relTol 0.05; maxIter 100; minIter 1; } } PISO { nCorrectors 2; // PISO loop // not used in turbSimpleDyMFoam nOuterCorrectors 2; // PIMPLE loop nNonOrthogonalCorrectors 0; correctPhi true; } PIMPLE { nOuterCorrectors 2; nCorrectors 2; nNonOrthogonalCorrectors 0; correctPhi true; } SIMPLE { nNonOrthogonalCorrectors 0; pRefCell 0; pRefValue 0; } relaxationFactors { p 0.3; U 0.7; k 0.6; omega 0.6; epsilon 0.6; } I am at lost at the moment. I would really appreciate any help from you.
__________________
Happy Foaming 

January 6, 2014, 12:09 

#2 
Member
Join Date: Nov 2012
Posts: 62
Rep Power: 5 
For your convenience I am attaching some pictures. It is the picture of the overall computational domain:
Capture.jpg These are the AMI PATCHes: Capture2.JPG Thanks
__________________
Happy Foaming 

January 8, 2014, 01:57 

#3 
Member
Join Date: Nov 2009
Posts: 62
Rep Power: 8 
Hi Naruto
I also simulate micro VAWT. Are your simulations successful for stopped or very slow rotation cases ? In addition, plz. check in your constant/controlDict for axis definition of moment and air density etc. Best regards, waku2005 

January 9, 2014, 08:33 

#4 
Member
Join Date: Nov 2012
Posts: 62
Rep Power: 5 
Thanks for your reply. Actually I am simulating a high speed turbine. The Reynolds number is quiet high. But I am getting negative power coefficient. Although I was able to solve the problem a little. Previously I used to get power coeff around 0.5.
Now it is around 0.2. I achieved this change by changing fvSolution. solvers { pcorr { solver GAMG; tolerance 1e7; relTol 0.005; smoother GaussSeidel; nPreSweeps 0; nPostSweeps 2; cacheAgglomeration off; agglomerator faceAreaPair; nCellsInCoarsestLevel 10; mergeLevels 1; maxIter 100; minIter 1; } p { $pcorr; tolerance 1e5; relTol 0.01; } pFinal { $p; tolerance 1e6; relTol 0; } "(UkepsilonomeganuTilda)" { solver smoothSolver; smoother GaussSeidel; nSweeps 1; tolerance 1e07; relTol 0.05; maxIter 100; minIter 1; }; "(UkepsilonomeganuTilda)Final" { solver PBiCG; preconditioner DILU; nSweeps 2; tolerance 1e07; relTol 0.05; maxIter 100; minIter 1; } } PISO { nCorrectors 2; // PISO loop // not used in turbSimpleDyMFoam nOuterCorrectors 2; // PIMPLE loop nNonOrthogonalCorrectors 0; correctPhi true; } PIMPLE { nOuterCorrectors 2; nCorrectors 2; nNonOrthogonalCorrectors 0; correctPhi no; } relaxationFactors { p 0.3; U 0.7; k 0.6; omega 0.6; epsilon 0.6; } cache { grad(U); } It is my changed fvSolution. My rotational axis definition is ok. Could you suggest me anything else?
__________________
Happy Foaming 

January 10, 2014, 05:51 

#5 
Member
Join Date: Nov 2009
Posts: 62
Rep Power: 8 
Hi Naruto
Well, I think fvSolution is not the reason of negative Cp. I'd like to know your TSR and free stream velocity. If you set too high TSR, it may bring the negative Cp. How about set a relatively slow TSR such as around 1. Best regards. waku2005 

January 10, 2014, 07:30 

#6 
Member
Join Date: Nov 2012
Posts: 62
Rep Power: 5 
Hi,
Thanks for your interaction. My freestream velocity is 7m/s. The TSR is only 3.0. I really could not figure out the problem. For your conveinence I am attaching my case files. If you could manage time please take a look. I would be very glad https://dl.dropboxusercontent.com/u/...T.S.R.3.0.rar
__________________
Happy Foaming 

January 11, 2014, 01:31 

#7 
Member
Join Date: Nov 2009
Posts: 62
Rep Power: 8 
Hi Naruto
Thanks for your case files. I made a run of the case without any changes on my desktop (OF 2.2.x/pimpleDyMFoam), and got Cp and Torque variations until t=0.6 as attached figures. Cp is relatively small and negative in some period, but I think that is due to its wing section and solidity( NACA0024 ?) , and these results are seemed to be reasonable. # Cp will vary and some averaging process would be required. For the above Cp calculation, I used below values.
Best regards, waku2005 

January 11, 2014, 12:50 

#8 
Member
Join Date: Nov 2012
Posts: 62
Rep Power: 5 
Dear Wakku,
I think, I have partially solved the problem. Actually I am trying to validate with an experimental work. According to the experiment, the average power coefficient is 0.15. Anyways, actually my boundary layer meshing was incorrect. After increasing the boundary layer thickness a little bit, I got a little better result. Although, It is still not sufficient. So far I had been able to get an average of 0.08. That's horrible. I have found one very interesting thing though. If you use upwind scheme for div(phi,omega) and div(phi,k), you would get better result. From my very limited knowledge of CFD, I am guessing that at the current configuration the flow is fully turbulent. That is why it is more stable with upwind scheme instead of linear upwind one which is more suitable in viscous sublayer. I have two questions. Like for setting the value of k I use the following formula k=1.5(U*I)*(UI) where I is the turbulent intesity. As I am using SST komega model for omega my equation is omega=sqrt(k)/((turbulent length scale)*Cmu^0.25). I got this value after studying the code of OpenFOAM. Am I correct? I am normally choosing a turbulence intensity of 1% and turbulent length scale of like 0.25 at the moment. Could you suggest anything else to improve my result? Thanks.
__________________
Happy Foaming 

January 11, 2014, 22:20 

#9 
Member
Join Date: Nov 2009
Posts: 62
Rep Power: 8 
Hi naruto
I think turb. length of 0.25 is somewhat long. How about 0.1 or 0.05? And if you use more refined boundary layers, please consider to use komegaSST_LowRe for low Re effect. Sincerery yours waku2005 

January 12, 2014, 04:19 

#10 
Member
Join Date: Nov 2012
Posts: 62
Rep Power: 5 
Dear wakku,
I heard about low Re komega SST model. But my current RE_no is aroudn 1.1e05. Would this turbulence model be efficient for this high Re number? And I will use length of 0.05 and 0.1. Thank you
__________________
Happy Foaming 

January 12, 2014, 20:34 

#11 
Member
Join Date: Nov 2009
Posts: 62
Rep Power: 8 
Hi Naruto
RASmodels such as komegaSST_lowRe includes low Re number effects near wall. Those models are not for low Re cases, and also appricable for hi Re cases. It might improve numerical results but you will require more finer boundary mesh and change wall treatments of wings. I hope this will help you. Best regards waku2005 

January 12, 2014, 23:55 

#12 
Member
Join Date: Nov 2012
Posts: 62
Rep Power: 5 
Dear wakku,
Thank you for your suggestion. Could you please give me a little hint on how to give boundary condition using kklOmega model? Thank you.
__________________
Happy Foaming 

January 13, 2014, 19:14 

#13 
Member
Join Date: Nov 2009
Posts: 62
Rep Power: 8 
Hi Naruto
Sorry, I've never used kklOmega model for low Re RANS simulation. Please refer the below thread: [OpenFOAM2.1.0] kklOmega RAS Turbulence Model (low Re) Regards waku2005 

Thread Tools  
Display Modes  


Similar Threads  
Thread  Thread Starter  Forum  Replies  Last Post 
nonorthogonal faces and incorrect orientation?  nennbs  OpenFOAM Native Meshers: blockMesh  7  April 17, 2013 05:42 
error message with modeling a cube with a hold at the center  hsingtzu  OpenFOAM Native Meshers: blockMesh  2  March 14, 2012 10:56 
mesh airfoil NACA0012  anand_30  OpenFOAM Meshing & Mesh Conversion  12  December 12, 2011 05:16 
BlockMesh FOAM warning  gaottino  OpenFOAM Native Meshers: blockMesh  7  July 19, 2010 14:11 
Axisymmetrical mesh  Rasmus Gjesing (Gjesing)  OpenFOAM Native Meshers: blockMesh  10  April 2, 2007 14:00 