CFD Online Discussion Forums

CFD Online Discussion Forums (http://www.cfd-online.com/Forums/)
-   OpenFOAM Post-Processing (http://www.cfd-online.com/Forums/openfoam-post-processing/)
-   -   Forces on cylinder in Openfoam (http://www.cfd-online.com/Forums/openfoam-post-processing/98644-forces-cylinder-openfoam.html)

krsp March 15, 2012 08:50

Forces on cylinder in Openfoam
 
Hi everyone,

I am trying to determine the forces and the force coefficients on a cylinder in OpenFoam 1.6-ext. I have typed the following in my controlDict file:

libs ("libmyIncompressibleRASModels.so" "libOpenFOAM.so" "libgroovyBC.so");
application osciTurbFoam;
startFrom latestTime;
deltaT 2.7e-4;
writeControl adjustableRunTime;
startTime 0;
stopAt endTime;
endTime 200.0;
writeInterval 1;
purgeWrite 0;
writeFormat ascii;
writePrecision 6;
writeCompression uncompressed;
timeFormat general;
timePrecision 6;
runTimeModifiable yes;
adjustTimeStep yes;
maxCo 0.5;
maxDeltaT 0.01;

functions
(
forces
{
type forces;
functionObjectLibs ("libforces.so");
patches (fixedWall);
pName p;
UName U;
rhoName rhoInf;
rhoInf 1000; // Reference density, fluid
CofR (0 0 0); // Origin for moment calculations
outputControl timeStep;
outputInterval 1;
}
forceCoeffs
{
type forceCoeffs;
functionObjectLibs ("libforces.so");
patches (fixedWall);
// pName p;
// UName U;
rhoName rhoInf;
rhoInf 1000;
CofR (0 0 0);
liftDir (0 0 1);
dragDir (1 0 0);
pitchAxis (0 0 0);
magUInf 0.05; // Free stream velocity
lRef 0.3; // Diameter of cylinder?
Aref 0.0707; // Ref. Area = cross sectional area?
outputControl timeStep;
outputInterval 1;
}
);


My problem is, that in the two directories "forces" and "forceCoeffs" I only get a 0 directory (only one .dat file for each). Does someone know, why the other timesteps do not generate files for forces and force coefficients?

Thank you in advance :)

niaz March 17, 2012 09:54

it is correct. the folder name shows the first time step. but it writes timesteps in the files.

November March 20, 2012 21:32

Hi

If you look at the forceCoeffs.C file, the meaning of magUInf, lRef, Aref are quite clear:

Directions for lift and drag forces, and pitch moment
00073 dict.lookup("liftDir") >> liftDir_;
00074 dict.lookup("dragDir") >> dragDir_;
00075 dict.lookup("pitchAxis") >> pitchAxis_;
00076 00077 // Free stream velocity magnitude
00078 dict.lookup("magUInf") >> magUInf_;
00079 00080 // Reference length and area scales
00081 dict.lookup("lRef") >> lRef_;
00082 dict.lookup("Aref") >> Aref_;
...
00118 scalar pDyn = 0.5*rhoRef_*magUInf_*magUInf_;
...
00127 scalar Cl = liftForce/(Aref_*pDyn);
00128 scalar Cd = dragForce/(Aref_*pDyn);
00129 scalar Cm = pitchMoment/(Aref_*lRef_*pDyn);

in your post, I am wondering if the rhoInf should be 1 rather than 1000 for water??
Please point out if I am wrong...

malhar June 13, 2012 06:25

hello every1,
i am also new to openfoam. i have done similar sim ulation for flow acrodss a cylinder. i get drag n lift forces correctly, i.e accor ding to vortex shedding i am getting variation in lift forces. but the coeff of drag n lift Cd n Cl, i am getting them constant throughout. i am unable ti digest this contrasting behavour.
please guide me thr this.

thanks n regards malhar.

krsp June 15, 2012 05:00

Hi Malhar,

I'll put my answer to your question here, just in case it could be usefull to others :) I am also quite new to CFD, but I will try to help as much as I can :)

Just to understand (please correct me if I'm wrong):
You get oscillating drag and lift forces, corresponding to the expected vortex shedding. But your drag and lift coefficients does not oscillate.

Since CD = FD/(1/2*rho*g*D*U^2) and CL = FL/(1/2*rho*g*D*U^2), a quick guess could be, that you have puttet the reference values for either rho, D or U equal to 0 in the controlDict file. (D corresponds to Aref).

And to your other question, you're right, that the total force on the cylinder comes from both pressure differences around the cylinder and from friction. The drag force (also called the in-line force) acts in the direction of the flow, while the lift force acts in the transverse direction.

Regards..

malhar June 15, 2012 05:04

dear krsp,
thanks for the reply. well u got my problem correctly. i wil post the control dict file which i have employed for this case. then u can easily find out the fault in it.

thanks n regards,
malhar

malhar June 15, 2012 05:04

controldict
 
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: 2.1.0 |
| \\ / A nd | Web: www.OpenFOAM.org |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
location "system";
object controlDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

application icoFoam;

startFrom startTime;

startTime 0;

stopAt endTime;

endTime 200;

deltaT 0.01;

writeControl timeStep;

writeInterval 50;

purgeWrite 0;

writeFormat ascii;

writePrecision 6;

writeCompression off;

timeFormat general;

timePrecision 6;

runTimeModifiable true;

functions
(
forces
{
type forces;
functionObjectLibs ("libforces.so");
patches (cylinder);
pName p;
UName U;
rhoName rhoInf;
rhoInf 1000; // Reference density, fluid
CofR (-2.5 0 0); // Origin for moment calculations
outputControl timeStep;
outputInterval 50;
}
forceCoeffs
{
type forceCoeffs;
functionObjectLibs ("libforces.so");
patches (cylinder);
// pName p;
// UName U;
rhoName rhoInf;
rhoInf 1000;
CofR (-2.5 0 0);
liftDir (-2.5 0.5 0);
dragDir (-2 0 0);
pitchAxis (-2.5 0 0);
magUInf 1; // Free stream velocity
lRef 1; // Diameter of cylinder?
Aref 3.142; // Ref. Area = cross sectional area?
outputControl timeStep;
outputInterval 50;
}
);

// ************************************************** *********************** //

malhar June 15, 2012 05:10

forcesMoments fm = forces::calcForcesMoment();

scalar pDyn = 0.5*rhoRef_*magUInf_*magUInf_;

vector totForce = fm.first().first() + fm.first().second();
vector totMoment = fm.second().first() + fm.second().second();

scalar liftForce = totForce & liftDir_;
scalar dragForce = totForce & dragDir_;
scalar pitchMoment = totMoment & pitchAxis_;

scalar Cl = liftForce/(Aref_*pDyn);
scalar Cd = dragForce/(Aref_*pDyn);
scalar Cm = pitchMoment/(Aref_*lRef_*pDyn);



above is the code written in the forceCoeffs.C file for calculating the Cl, Cd,Cm. i am not understanding what are time(), first(), second()

malhar

krsp June 15, 2012 05:49

I don't understand the directions you typed for drag and lift. If your steady flow runs parallel to the x-axis, I think they should be something like:

liftDir (0 0 1); (OR (0 1 0); if the y-axis is your vertical axis)
dragDir (1 0 0);

Also I think your outputInterval is quite high. If your time step is 0.01s, you only get the results for every 0.5s. Just make sure, that it is often enough to capture the vortex shedding frequency.
To your reference area I would put it equal to the cylinder diameter. (Of course that depends on which formula you refer to, but that is what is done in the book "Hydrodynamics around cylindrical structures" by B. Mutlu Sumer and Jørgen Fredsøe)

Finally as posted above, I think that if you solve an incompressible problem rhoInf should be 1 instead of 1000, because the value is already included in the pressure..

Hopefully this will help a bit :)

malhar June 15, 2012 05:56

dear krsp.....
it may be co-incidental... there u found out the fault n here me too:-)...i got dat prob. i have set the code to run. wil msg u the results once i get them.

thank you very very much.

wishes,
malhar

amin66 July 2, 2012 18:03

i've just done the simulation around cylinder and use above controlDict. in Re=150 experiment data says Cd=1.5 but i get Cd=1.6 that is not accurate. i changed my mesh(gambit) better(more resolution) and used smaller time step(0.01) but the result became worse! Cd=2.7 ??!!!:eek:
i dont know what did happen? with a better mesh and smaller time step and similar other conditions my result became worse!
please help what maybe happened that i got this bad result!??

malhar July 3, 2012 00:35

dear amin,
wat is ur domain size?? is ur cylinder in a confined channel??? wat inlet boiundary condition have you taken???

regards Malhar.

amin66 July 3, 2012 10:11

4 Attachment(s)
dear malhar

1. for curved mesh and time step=0.1 taken this result for Cd:
fluent:1.484 , Openfoam:1.6 , experiment:1.5

2. for rectangular domain mesh and timestep=0.001 taken this result for Cd:
fluent:1.32 , Openfoam:2.55 , experiment:1.5

velocity inlet condition for inlet

malhar July 3, 2012 13:15

dear Amin
i wanted to know your type of boundary conditin. have you taken fully developed flow condition i.e parabolic for your prob at inlet???? if nt den its necessary if ur upstream length is nt enough long for the flow to get developed. most probably the incorrectness in Cd Cl is because of dat.
regards,
Malhar

niaz July 3, 2012 13:15

Dear amin
your solution does not converge change the number of correction in cvsolution text.
but i am sure that the openfoam predict drag more than real value.
it may occure because of pressure calculation

amin66 July 3, 2012 14:24

1 Attachment(s)
dear malhar, i attached my case study.

dear niaz, how can i know that my solution has converged or not? and how can plot the residuals like fluent?

malhar July 3, 2012 14:42

dear Amin,
i am too new to OF...even i dont know how to know that solution has converged or nt and to plot residuals.....in fact i am nt understanding wat u meant by convergence in ur case n how it may be related to ur prob.

niaz July 3, 2012 15:18

Dear amin
I check it by saving it in a text file. and check it by eye.
secondly, It is a point which I catch it by experiment. if you import meshes from other sources, you need to do more loops for convergance of solution. try it :)

Quote:

Originally Posted by amin66 (Post 369581)
dear malhar, i attached my case study.

dear niaz, how can i know that my solution has converged or not? and how can plot the residuals like fluent?


niaz July 3, 2012 15:22

other points, your mesh shows that your blockage ratio is not enough. so don`t compare your result with one which has wider blockage ratios.
it is important point. you can find this point in a paper from sanjey mittal.
only search it. it is in his site :rolleyes:

amin66 July 3, 2012 16:53

Dear A_R, tnx for your attention. i just say in rectangular domain mesh i have wider blockage ratio and more node numbers than the curved domain mesh.(as you can see up the pictures) but the result became worse.


All times are GMT -4. The time now is 18:21.