CFD Online Discussion Forums

CFD Online Discussion Forums (
-   OpenFOAM Running, Solving & CFD (
-   -   Understanding interFoam (

sega July 10, 2008 09:35

Hello. I'm trying to get a

I'm trying to get a rather complete overview about what interFoam is doing and how it works.

What do you think is the best way of understanding a solver?

I had a quick look into the Doxygen documentation.
Starting from interFoam.C I looked up the first header-file: fvCFD.H. There I find 14 other .H-files. Some of them linking even more dependencies.

So, at first glance it seems to be impossible to read this whole crowd of codes.

Where can I find the essentials of the underlying finite volume method, like how surface tension is implemented, how the interface is tracked/reconstructed, or the divergence and laplacian schemes are treated just to give some examples?

Maybe you can give me a hint where to start?

jaswi July 11, 2008 02:17

Hi Sebastian I have some p
Hi Sebastian

I have some personal notes on interFOAM which were made for learning purposes. I will send them over to you after Tuesday. Please just send me an e-mail next week . Regarding the interface tracking and reconstruction , interFOAM uses compressive schemes so there is no need of reconstruction. The gamma scalar transport equation is solvd with compressive schemes which not only keeps the interface sharp but also compresses it. I will give you some more details on that next week .

I started with that same approach as you are trying now but it did not work. All you need is to pick some major classes and study them.

My advice would be leave the OpenFOAM/primitives folder for now and try taking a look at :

primitiveMesh, polyMesh, fvMesh and classes related to patches (there is a whole hierarchy for that) . Once you understand the structure of mesh and patches then you will understand the fvPacthes and fvPatchFields as it is from where the boundary conditions come into play. If you will be planning to go for a deeper study then we can divide the stuff among us and try to create some easy to understand guide or documentation.I have done some code reading and can share my experience with you.

I will be actively pursuing this casue after next week. Let me know if you want to join me :-).

With Best Regards

sega July 21, 2008 09:44

Having a deeper look into the
Having a deeper look into the PhD-thesis of Henrik Rusche I get the feeling, that this is a rather complete overview what interFoam is doing.

Is interFoam a result from this work?
And what is (still) applying to the current solver?

sega August 13, 2008 05:43

Hello again! After looking
Hello again!

After looking into the UEqn.H file of the interFoam solver I gathered some more knowledge about the solver. But of course there are some further questions.

The UEqn.H defines a (left-hand-side) equation and solves the variables in it until they satisfy another (right-hand-side) equation?

Regarding the first equation:

fvm::ddt(rho, U)
+ fvm::div(rhoPhi, U)
- fvm::laplacian(muf, U)
- (fvc::grad(U) & fvc::grad(muf))
//- fvc::div(muf*(fvc::interpolate(dev(fvc::grad(U))) & mesh.Sf()))

muf: This must be the "average" dynamic viscosity from the two fluids?
rhoPhi: Don't know what this is. Phi is a cell face motion flux (per definition). What does rhoPhi mean?

What about these logical operators (&) in an algebraic equation? What kind of operation is this?

What is mesh.Sf() ?

Second equation:

( fvc::interpolate(interface.sigmaK())*fvc::snGrad(g amma)
- ghf*fvc::snGrad(rho)
- fvc::snGrad(pd)
) * mesh.magSf()

The hydrostatic phenomena is hidden behind snGrad(pd)? But why is is a surface normal gradient and not a "simple" gradient?

What are ghf and mesh.magSf()

I would appreciate any answers.
Thanks so far!

andersking August 13, 2008 09:48

Hi, In regards to '&' (dot pro
Hi, In regards to '&' (dot product) and mesh.Sf() (face surface area vectors) the Programmer's Guide gives lots of info. (in ~/OpenFOAM/OpenFOAM-1.?/doc/Guides-a4/ )

they also give a lot of other vector/tensor operators.

As for the C++ code dependencies the online Doxygen docs are useful. ( )


jaswi August 14, 2008 04:35

Hi Sebastian I tried to con
Hi Sebastian

I tried to contact you on related to interfoam , but the mail bounced back . Could you please confirm the address or please mail me your another email address.

With Best Regards

sega August 15, 2008 08:56

Try this one: sebastian.gatzka
Try this one:

sega August 20, 2008 11:24

I had a look into the Programm
I had a look into the Programmers Guide.
It looks like the code is getting clearer to me, now.

Just to make sure I understand it right:

In the case of an incompressible fluid with constant properties the expressions grad(U) is zero and so

fvc::grad(U) & fvc::grad(muf)


fvc::div(muf*(fvc::interpolate(dev(fvc::grad(U))) & mesh.Sf()))

will be equal to zero?

With constant fluid properties snGrad(rho) will be zero as well?

Is this correct so far?

eugene August 20, 2008 11:50

grad(U) isn't zero div(U) is
grad(U) isn't zero
div(U) is zero

There is a big difference.

sega August 21, 2008 05:16

Yes, of course. Sorry. So m
Yes, of course. Sorry.

So my problem persists:
There are three components in the equation I don't understand. I haven't seen these terms in any incompressible Navier-Stokes so far.

1) fvc::grad(U) & fvc::grad(muf)
2) fvc::div(muf*(fvc::interpolate(dev(fvc::grad(U))) & mesh.Sf()))
3) ghf*fvc::snGrad(rho)


1) Of course grad(U) is NOT zero. What about grad(muf)? Is there a gradient in the neighbourhood of the interface?

2) Whats this whole expression? I haven't seen anything similar in the works of Rusche and Ubbink or any other textbook about CFD so far. Any hints, where to look?

3) Maybe the same like grad(muf). Is snGrad(rho) nonzero in the sourrounding of the interface?

Thanks so far.
Greetings from sunny Germany.

david_h August 21, 2008 08:40

Gatzka, 1) this term is wha

1) this term is what is left of term (2) if the differential operators commute and div(U)=0

2) this term is the part of the stress tensor which cannot be represented by the laplacian

(I assume (1) was selected over (2) to improve the stability)

3) this term is to cancel the density gradient which is introduced due to using the gradient of dyanmic pressure, pd = p - rho * (gh), in the momentum equation.

hope this helps

sega August 24, 2008 04:22

Thank you Dave. Finally, I
Thank you Dave.

Finally, I can see where this will lead, and now I recognize the analogy to the equations in Ubbinks PhD-thesis. Obviously I have to be very carefully with the differential operators in the momentum equation.

I gave it a try to recalculate the interFoam UEqn.
On the way I tried to treat all operators right.
Maybe you want to have a look at the point where I am stuck:

Obviously there a still many differences between my equation (which is not jet finished), equation (2.15) in Ubbinks work and the final UEqn used with interFoam.

1) First of all Ubbink is losing the transposed gradients of the velocity-field somewhere between the lines. I don't know how.

2) My equation is stating a "real" laplacian operator. Ubbinks has non, just a "pre"-laplacian (?) div(mu*grad(U)). In interFoam there is a laplacian which is including the muf (face interpolated dynamic viscosity).

3) interFoam has these deviatoric components. Dave said they are components left over from the stress tensor. Maybe they are some representation from the transposed gradients?

You see, I want to recap the development of the interFoam solver my way. Please tell me if I'm completely wrong. Regarding the deviatoric components I have read the section in the programmers guide but don't know jet how to use these information in my equations.

Regards. S.

david_h August 28, 2008 09:47

Sebastian, 1) The 4th term

1) The 4th term on right hand side of the last equation is zero, if:
- the gradient and divergence operators commute
- and the velocity is divergence free

In the indicial notation:
mu ( d_j ( d_i u_j ) ) = mu ( d_i (d_j u_j) )

2) If you combine terms 2 and 3 of RHS you will get the modified laplacian in interFoam.

Also, the density should be outside the D/Dt operator on the left hand side of the first equation. If combined with the continuity equation you will be able to get your equation to look like eq. 2.15 of the Thesis and interFoam

sega August 30, 2008 07:37

Thanks! I made the correcti

I made the corrections to the PDF-file. Have a look.

Further unclarities are marked with the question mark. These are

- in the material derivative: Is there an identity U*grad(rho*U) = grad(rho*U*U) ?
- is the commuted divergence operator correct (righ hand side 4)
- is the "combination" of RHS 2 + 3 correct this way? If not how does it work? Whats behind this operation?

sega August 30, 2008 07:38

Sorry, I have missed one impor
Sorry, I have missed one important question.

What about the transposed gradient?
How do I treat this one?

david_h September 8, 2008 11:43

I think what happens is you co
I think what happens is you converted

div( mu grad(U)^T ) = grad(U)^T . grad(mu) (using div(U) == 0)

as opposed to

div( mu grad(U)^T ) = grad(mu) . grad(U)^T = grad(U) . grad(mu)

sega September 13, 2008 10:53

Ok, so div( mu grad(U)^T )
Ok, so

div( mu grad(U)^T ) = grad(mu) . grad(U)^T = grad(U) . grad(mu)

is correct?

I converted like this

div( mu grad(U)^T ) = mu div( grad(U)^T ) + grad(U)^T grad(mu)

and that is wrong?

Would you mind having a look at the other questions, too?
Obviously my knowledge about differential operators is limited ...

sega September 24, 2008 08:50

Hello. I'm back again and m

I'm back again and made some changes to the file above with the new information. Would anyone mind having a look at the equations?

There are still some questions (some of them already mentioned above without answer so far).

1) Expression div(mu(grad(u)) must be the source of to the interFoam expression fvc::div(muf*(fvc::interpolate(dev(fvc::grad(U)))

How can this be done? Maybe someone can tell me in operator syntax?

2) I didn't get the thing about dynamic pressure.
If I try to use dynamic pressure instead of static pressure pd = p - rho*g*h and transform to the interFoam expression ghf*fvc::snGrad(rho) and fvc::snGrad(pd) there is still g*rho left from the original equation. What happens to it?

Thanks so far.

tomislav_maric September 30, 2008 09:30

Jaswinder: I'm trying to lea
I'm trying to learn interFoam also, so I would really benefit from Your notes. Can You please send them to:

Thank You!

isabel November 12, 2009 05:28

There is something in interFoam I don't understand. At createFields.H:

volScalarField rho
gamma*rho1 + (scalar(1) - gamma)*rho2,

Does anybody know what "gamma.boundaryfield().types()" and "rho.oldTime();" means?

All times are GMT -4. The time now is 10:49.