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

Calculate Lift and Drag Coefficients CL and CD

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

Like Tree3Likes
  • 2 Post By Simon Lapointe
  • 1 Post By desert_1250

Reply
 
LinkBack Thread Tools Display Modes
Old   September 28, 2009, 20:04
Default Calculate Lift and Drag Coefficients CL and CD
  #1
Member
 
Sven Winkler
Join Date: May 2009
Posts: 70
Rep Power: 8
sven is on a distinguished road
Hi,

I am simulating a turbulent flow around a square. For This flow, I want to calculate the lift and the drag coefficients CD and CL. By scanning through posts in this forum, I found that this is possible by adding the following lines to the controlDict file:

Code:
functions
(
forces
{
type forces;
functionObjectLibs ("libforces.so"); //Lib to load
patches (cylinder); // change to your patch name
rhoInf 1.204; //Reference density for fluid
CofR (0 0 0); //Origin for moment calculations
outputControl timeStep;
outputInterval 100;
}

forceCoeffs
{
type forceCoeffs;
functionObjectLibs ("libforces.so");
patches (cylinder); //change to your patch name
rhoInf 1.204;
CofR (0 0 0);
liftDir (0 1 0);
dragDir (1 0 0);
pitchAxis (0 0 0);
magUInf 0.1;
lRef 1;
Aref 1;
outputControl timeStep;
outputInterval 100;
}
);
Unfortunately this does not yield the expected values for CL and CD. For fixing this problem and obtaining the correct CL/CD I first need some information on the code posted above:


  1. rhoInf: What is the reference density? The density in my free stream?
  2. What is magUInf? How can i determine/calculate this value?
  3. The program seems to use "libforces.so". I could find this file, but where can I find the corresponding source code?
  4. Where are the equations for calculating CD and CL defined?
I appreciate your help. Thank you very much
sven is offline   Reply With Quote

Old   September 28, 2009, 22:17
Default
  #2
Member
 
Simon Lapointe
Join Date: May 2009
Location: Qubec, Qc, Canada
Posts: 33
Rep Power: 8
Simon Lapointe is on a distinguished road
Hi,

1. rhoInf is the freestream density.
2. magUinf is the magnitude of the freestream velocity.
3. The source code for the forces and forceCoeffs are located, respectively, in OpenFOAM-1.6/src/postProcessing/forces/ and OpenFOAM-1.6/src/postProcessing/forceCoeffs/
4. The equations of CL and CD are defined in OpenFOAM-1.6/src/postProcessing/forceCoeffs/forceCoeffs.C

Hope this helps
ruby_nuaa and FrankFlow like this.
Simon Lapointe is offline   Reply With Quote

Old   September 28, 2009, 22:54
Default
  #3
Member
 
Sven Winkler
Join Date: May 2009
Posts: 70
Rep Power: 8
sven is on a distinguished road
Thanks for your answer Simon. I still have a problem with magUInf, the magnitude of the free stream velocity. What exactly is the magnitude of my free stream velocity? My free stream velocity is 0.322, so the magnitude is 0.1? I tried it with that value (0.1), but it does not yield the correct results! Thanks again for your help!
sven is offline   Reply With Quote

Old   September 29, 2009, 02:53
Default
  #4
New Member
 
bigred's Avatar
 
Matthew Philpott
Join Date: Aug 2009
Location: Belgium
Posts: 24
Rep Power: 8
bigred is on a distinguished road
If your free stream velocity is 0.322, then I'm guessing the magnitude is 0.322. If that is the only velocity component ie, value (0.322 0 0). Magnitude is calculated using typical vector mathematics v=vx+vy+vz. The square of the velocity magnitude is the sum of the squares of your velocity components.
__________________
CAELinux 2009 + OF1.5
Ubuntu 9.04 x64 (jaunty jackalope) + OF1.6
bigred is offline   Reply With Quote

Old   September 29, 2009, 21:01
Default
  #5
Member
 
Sven Winkler
Join Date: May 2009
Posts: 70
Rep Power: 8
sven is on a distinguished road
Ok..now I got it! Thank you all very much for your help. Greetings from California

Sven
sven is offline   Reply With Quote

Old   October 2, 2009, 21:19
Default
  #6
Member
 
Sven Winkler
Join Date: May 2009
Posts: 70
Rep Power: 8
sven is on a distinguished road
Well, I looked up the definition for CD and CL in the source Code. First of all, the right location for the files is:
/home/user/OpenFOAM/OpenFOAM-1.6/src/postProcessing/functionObjects/forces

There are two folders (forces and forceCoeffs) where you can find the appropriate .C and .H files.

I first had a look at the forceCoeffs.C file, here is part of that:

Code:
void Foam::forceCoeffs::write()
{
    if (active_)
    {
        // Create the forces file if not already created
        makeFile();

        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);

        if (Pstream::master())
        {
            forcesFilePtr_()
                << obr_.time().value() << tab
                << Cd << tab << Cl << tab << Cm << endl;

            if (log_)
            {
                Info<< "forceCoeffs output:" << nl
                    << "    Cd = " << Cd << nl
                    << "    Cl = " << Cl << nl
                    << "    Cm = " << Cm << nl
                    << endl;
            }
        }
    }
}
OpenFOAM first calculates a total force, by multiplying it with the directions of each component, this force is divided into its components (lift/drag). These components are then divided by the referenceDensity multiplied with the refArea which leads to CL and CD.
Unfortunately I do not understand, how OpenFOAM calculates the total force mentioned above. I tried to look up the definition of first() and second(), which are defined in Tuple2.H (OpenFoam/InInclude/Tuple2.H). Here is the definition of these functions:

Code:
//- Return first
        inline const Type1& first() const
        {
            return f_;
        }

        //- Return first
        inline Type1& first()
        {
            return f_;
        }

        //- Return second
        inline const Type2& second() const
        {
            return s_;
        }

        //- Return second
        inline Type2& second()
        {
            return s_;
        }
I don't understand the definitions of these functions. Can anyone help me with that?
sven is offline   Reply With Quote

Old   October 2, 2009, 21:21
Default
  #7
Member
 
Sven Winkler
Join Date: May 2009
Posts: 70
Rep Power: 8
sven is on a distinguished road
In the forces.C file, I also found this:


Code:
Foam::forces::forcesMoments Foam::forces::calcForcesMoment() const
{
    forcesMoments fm
    (
        pressureViscous(vector::zero, vector::zero),
        pressureViscous(vector::zero, vector::zero)
    );

    if (directForceDensity_)
    {
        const volVectorField& fD = obr_.lookupObject<volVectorField>(fDName_);

        const fvMesh& mesh = fD.mesh();

        const surfaceVectorField::GeometricBoundaryField& Sfb =
            mesh.Sf().boundaryField();

        forAllConstIter(labelHashSet, patchSet_, iter)
        {
            label patchi = iter.key();

            vectorField Md = mesh.C().boundaryField()[patchi] - CofR_;

            scalarField sA = mag(Sfb[patchi]);

            // Normal force = surfaceUnitNormal * (surfaceNormal & forceDensity)
            vectorField fN =
                Sfb[patchi]/sA
               *(
                    Sfb[patchi] & fD.boundaryField()[patchi]
                );

            fm.first().first() += sum(fN);
            fm.second().first() += sum(Md ^ fN);

            // Tangential force (total force minus normal fN)
            vectorField fT = sA*fD.boundaryField()[patchi] - fN;

            fm.first().second() += sum(fT);
            fm.second().second() += sum(Md ^ fT);
        }
    }
    else
    {
        const volVectorField& U = obr_.lookupObject<volVectorField>(UName_);
        const volScalarField& p = obr_.lookupObject<volScalarField>(pName_);

        const fvMesh& mesh = U.mesh();

        const surfaceVectorField::GeometricBoundaryField& Sfb =
            mesh.Sf().boundaryField();

        tmp<volSymmTensorField> tdevRhoReff = devRhoReff();
        const volSymmTensorField::GeometricBoundaryField& devRhoReffb
            = tdevRhoReff().boundaryField();

        forAllConstIter(labelHashSet, patchSet_, iter)
        {
            label patchi = iter.key();

            vectorField Md = mesh.C().boundaryField()[patchi] - CofR_;

            vectorField pf = Sfb[patchi]*p.boundaryField()[patchi];

            fm.first().first() += rho(p)*sum(pf);
            fm.second().first() += rho(p)*sum(Md ^ pf);

            vectorField vf = Sfb[patchi] & devRhoReffb[patchi];

            fm.first().second() += sum(vf);
            fm.second().second() += sum(Md ^ vf);
        }
    }

    reduce(fm, sumOp());

    return fm;
}
Can anyone by this understand, how OpenFOAM calculates these forces?

Thanks a lot for all your help!
sven is offline   Reply With Quote

Old   August 10, 2011, 05:21
Default
  #8
Member
 
s.rasoul_varedi
Join Date: Feb 2010
Posts: 82
Rep Power: 6
desert_1250 is an unknown quantity at this point
Send a message via Yahoo to desert_1250
Quote:
Originally Posted by sven View Post
In the forces.C file, I also found this:


Code:
Foam::forces::forcesMoments Foam::forces::calcForcesMoment() const
{
    forcesMoments fm
    (
        pressureViscous(vector::zero, vector::zero),
        pressureViscous(vector::zero, vector::zero)
    );

    if (directForceDensity_)
    {
        const volVectorField& fD = obr_.lookupObject<volVectorField>(fDName_);

        const fvMesh& mesh = fD.mesh();

        const surfaceVectorField::GeometricBoundaryField& Sfb =
            mesh.Sf().boundaryField();

        forAllConstIter(labelHashSet, patchSet_, iter)
        {
            label patchi = iter.key();

            vectorField Md = mesh.C().boundaryField()[patchi] - CofR_;

            scalarField sA = mag(Sfb[patchi]);

            // Normal force = surfaceUnitNormal * (surfaceNormal & forceDensity)
            vectorField fN =
                Sfb[patchi]/sA
               *(
                    Sfb[patchi] & fD.boundaryField()[patchi]
                );

            fm.first().first() += sum(fN);
            fm.second().first() += sum(Md ^ fN);

            // Tangential force (total force minus normal fN)
            vectorField fT = sA*fD.boundaryField()[patchi] - fN;

            fm.first().second() += sum(fT);
            fm.second().second() += sum(Md ^ fT);
        }
    }
    else
    {
        const volVectorField& U = obr_.lookupObject<volVectorField>(UName_);
        const volScalarField& p = obr_.lookupObject<volScalarField>(pName_);

        const fvMesh& mesh = U.mesh();

        const surfaceVectorField::GeometricBoundaryField& Sfb =
            mesh.Sf().boundaryField();

        tmp<volSymmTensorField> tdevRhoReff = devRhoReff();
        const volSymmTensorField::GeometricBoundaryField& devRhoReffb
            = tdevRhoReff().boundaryField();

        forAllConstIter(labelHashSet, patchSet_, iter)
        {
            label patchi = iter.key();

            vectorField Md = mesh.C().boundaryField()[patchi] - CofR_;

            vectorField pf = Sfb[patchi]*p.boundaryField()[patchi];

            fm.first().first() += rho(p)*sum(pf);
            fm.second().first() += rho(p)*sum(Md ^ pf);

            vectorField vf = Sfb[patchi] & devRhoReffb[patchi];

            fm.first().second() += sum(vf);
            fm.second().second() += sum(Md ^ vf);
        }
    }

    reduce(fm, sumOp());

    return fm;
}
Can anyone by this understand, how OpenFOAM calculates these forces?

Thanks a lot for all your help!
Hi, for calculating total Forces in OpenFOAM, integrated the pressure Field over the surfaces and projection of this forces in drag or lift direction Division Pdyn*Aref gives the cd and cl Coeff

_____
Rasoul
ehsan likes this.
desert_1250 is offline   Reply With Quote

Old   June 13, 2012, 06:20
Default
  #9
New Member
 
Malhar Malushte
Join Date: May 2012
Posts: 16
Rep Power: 5
malhar is on a distinguished road
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.
malhar is offline   Reply With Quote

Old   December 17, 2012, 05:53
Default
  #10
Member
 
Malik
Join Date: Dec 2012
Location: Austin, USA
Posts: 52
Rep Power: 4
malaboss is on a distinguished road
Quote:
Originally Posted by desert_1250 View Post
Hi, for calculating total Forces in OpenFOAM, integrated the pressure Field over the surfaces and projection of this forces in drag or lift direction Division Pdyn*Aref gives the cd and cl Coeff
Hi,
by saying this, do you mean that the pressure value used is the one on the boundary of the cells which touch the surface on which drag and lift force are applied ?

I am a new Foamer but i think that we forget the viscous force if OpenFOAM only does this. In the code, there is a tangential force calculated, but i cannot say if it concerns the viscous force or not. Does anyone have an idea ?
malaboss is offline   Reply With Quote

Old   August 1, 2014, 01:29
Default
  #11
Senior Member
 
Tushar Chourushi
Join Date: Jul 2009
Location: IIT-Indore, India
Posts: 319
Blog Entries: 1
Rep Power: 9
Tushar@cfd is on a distinguished road
Understanding of the terms...

I am actually looking out for the term mentioned below in the previous posts on this thread:

const volVectorField& fD = obr_.lookupObject<volVectorField>(fDName_);

Where is the term located in the case directory?

Correct me if I am wrong

Thanks in Advance!
Tushar@cfd is offline   Reply With Quote

Old   August 1, 2014, 01:38
Default
  #12
Senior Member
 
Tushar Chourushi
Join Date: Jul 2009
Location: IIT-Indore, India
Posts: 319
Blog Entries: 1
Rep Power: 9
Tushar@cfd is on a distinguished road
Oh, it's clear the term is for compressible flows.

Please, excuse for previous post. I did mistake it for incompressible flows.

Last edited by Tushar@cfd; August 2, 2014 at 01:24.
Tushar@cfd 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



All times are GMT -4. The time now is 16:07.