CFD Online Discussion Forums - Add extra body force (e.g. gravity) in simpleFoam: cavity example

CFD Online Discussion Forums (https://www.cfd-online.com/Forums/)

- OpenFOAM Running, Solving & CFD (https://www.cfd-online.com/Forums/openfoam-solving/)

- - Add extra body force (e.g. gravity) in simpleFoam: cavity example (https://www.cfd-online.com/Forums/openfoam-solving/145285-add-extra-body-force-e-g-gravity-simplefoam-cavity-example.html)

Add extra body force (e.g. gravity) in simpleFoam: cavity example

Hi everyone,

I heard lots of praises about OpenFOAM so I'm trying to get into it, but I have to admit it isn't as straightforward as I thought to set up a new model ;)

I need a steady-state solver for incompressible turbulent flows that are subject to an additional non-uniform body force. As a first step, I just want to use the cavity example (from icoFoam) with laminar flow, no initial pressure gradient, no initial velocity field and no lid movement, but taking the constant uniform field of gravity (0, -9.81, 0) into account. With this very simple example I expect to observe a nice pressure gradient from top to bottom of the cavity (constant in the x-direction).

Here is what I did so far. Following this thread, I created a new solver "myFoam" based on simpleFoam. In order to read a constant force field I added this code

Code:

    Info<< "Reading field F\n" << endl;

    volVectorField F

    (

        IOobject

        (

            "F",

            runTime.constant(),

            mesh,

            IOobject::MUST_READ,

            IOobject::NO_WRITE

        ),

        mesh

    );

in the file myFoam/createFields.H. Then I modified the file myFoam/UEqn.H like this

Code:

solve(UEqn() == -fvc::grad(p) + F );

and added the momentum corrector

Code:

U -= rAU*(fvc::grad(p) - F );

in the file myFoam/pEqn.H. I assume that the field F is already divided by rho so that the dimensions match. Compilation of the new "myFoam" application worked perfectly.

Next, I took the cavity example and added the file "F" describing the gravity (internalField uniform (0 -9.81 0);) in the constant directory. In the file 0/U I also set the initial velocity of the lid to zero. Since the parameters of the cavity example are intended for icoFoam and not simpleFoam, I then moved the files

- constant/RASProperties
- constant/transportProperties
- system/controlDict
- system/fvSchemes
- system/fvSolution

from the simpleFoam/pitzDaily example to my new cavity example. I set the RASmodel to laminar and call myFoam instead of simpleFoam in the controlDict file. Now when I run myFoam I get a "Floating point exception (core dumped)" after 50 steps, which isn't very nice... but what bothers me even more is the pressure field I obtain after 10 steps (see attached figure), which is far from resembling the one I was expecting :D

Am I using the wrong approach or doing something wrong? Or has anyone another idea as to how I could achieve my objective?

Thanks a lot in advance for your help! Best regards,

Julien

Have a look at the buoyantBoussinesqSimpleFoam and buoyantSimpleFoam solvers.

They include buoyancy forces. If you just want to see the effect of gravity on the pressure field (set the thermal expansion coefficient to 0 in the transportProperties file of the buoyantBoussinesqSimpleFoam solver).

Tutorials are in $FOAM_TUTORIALS/heatTransfer/ directory

Hi J Herb,

I have ran two similar cases of buoyantBoussinesqSimpleFoam solver, the only difference is that in one case gravity is enabled and in another without gravity. In both the cases I found the velocity profiles to be same. By the way as you suggested, I have set beta to be zero and there are no temperature gradients in the cases [to simulate simpleFoam kind of case with gravity included].
I am wondering that, the effect of gravity works in buoyantBoussinesqSimpleFoam when there is a difference in density of the fluid due to temperature?

Kindly clarify my doubt.

thanks

Quote:

Originally Posted by rakesh927 (Post 628618)

I am wondering that, the effect of gravity works in buoyantBoussinesqSimpleFoam when there is a difference in density of the fluid due to temperature?
thanks

If you look at the source code of it (https://github.com/OpenFOAM/OpenFOAM...nesqSimpleFoam) you see, that the density difference is used here:
https://github.com/OpenFOAM/OpenFOAM...oam/UEqn.H#L28
https://github.com/OpenFOAM/OpenFOAM...Foam/pEqn.H#L8

The density is calculated from the the temperature field:
https://github.com/OpenFOAM/OpenFOAM...oam/TEqn.H#L26

So if you set beta to zero, there is no density gradient and therefore no difference between simpleFoam and buoyantBoussinesqSimpleFoam. Does this answer your question?

Quote:

Originally Posted by jherb (Post 629418)

So if you set beta to zero, there is no density gradient and therefore no difference between simpleFoam and buoyantBoussinesqSimpleFoam. Does this answer your question?

Thanks Jherb for your clarification.

Why you add the force in the momentum corrector?

Sent from my GT-I8190L using CFD Online Forum mobile app

Do you know this: https://openfoamwiki.net/index.php/B...sinesqPisoFoam

Of course it's the PISO and not the SIMPLE algorithm but it might be a start: Check for the differences due to the gravitational force with the normal algorithm and then apply similar changes to the SIMPLE algorithm.