CFD Online Discussion Forums

CFD Online Discussion Forums (https://www.cfd-online.com/Forums/)
-   OpenFOAM Running, Solving & CFD (https://www.cfd-online.com/Forums/openfoam-solving/)
-   -   VOF interface (https://www.cfd-online.com/Forums/openfoam-solving/59335-vof-interface.html)

michele May 26, 2005 09:38

I'm just starting to use the O
 
I'm just starting to use the OpenFOAM. I'm interested in VOF methods for ship hydrodynamics simulations. I had previously used other CFD tools when I worked at University.
The first test I'm making is a very simple undisturbed free surface 2D case. An uniform velocity field enters from the left and, after having crossed the domain, leaves the domain through an outlet boundary condition. I'm using the rasInterFoam, as I'm interested in getting a model more and more sophisticated (the complete simulation will require a 3D full rigid body dynamics for a boat... I would like to make roll damping calculations and performance analyses).

But I must confess I hadn't obtained good results even for this simple problem. The problem is well initialised and starts well. But after about 2.6 seconds the free surface is distorted (I was expecting a trivial solution of steady undisturbed free surface). The link to a picture of gamma is the following:
http://www.lem3.it/tests/gamma_2.57102.png

The process becomes unstable some time after, as can be seen in the following picture:
http://www.lem3.it/tests/gamma_3.15095.png

The time step is limited by a Courant number of 0.2, the fields are initialised accordingly with the boundary conditions...

I was considering what I was missing... Thanks in advance.
Michele.

hjasak May 26, 2005 13:36

The second picture looks like
 
The second picture looks like a problem of a wave trying to reflect back from an inlet boundary...

Hrv

eugene May 26, 2005 16:05

What are your left and right b
 
What are your left and right boundary conditions on gamma? They should be fixed value on the left and zeroGradient on the outlet.

I assume in terms of numerics you used all the default settings from the rasInterFoam test case? Try running it with both
div(phi,U) Gauss upwind;
div(rho*phi,U) Gauss upwind;
in your fvSchemes dictionary. If this doesn't solve your problem (and even if it does), it is more than likely an issue with the step change between the air and liquid at the inlet.

Could you please pack up the case and post it to this forum? I would like to take a look at it some time.

michele May 27, 2005 03:33

Thanks for your kind replies.
 
Thanks for your kind replies. I think I imposed correctly the boundary conditions. The numerical scheme used for the momentum equation is upwind (at the first trials I realised that the gamma schemes where more unstable).
Herebelow, as Eugene requested to me, I attach the whole case.

By he way, its name, "waves", is misleading... it is a simple steady case (when the current issues will be solved, I will implement a time/space varying inlet boundary condition accordingly with potential wave theory).
I'm however still having problems with the case attached... any suggestions?

Thanks,
Michele.

michele May 27, 2005 03:39

Thanks for your kind replies.
 
Thanks for your kind replies. I think I imposed correctly the boundary conditions. The numerical scheme used for the momentum equation is upwind (at the first trials I realised that the gamma schemes where more unstable).
Herebelow, as Eugene requested to me, is a link to the case (it's too big to be posted to this forum):
www.lem3.it/tests/waves.tar.gz

By he way, its name, "waves", is misleading... it is a simple steady case (when the current issues will be solved, I will implement a time/space varying inlet boundary condition accordingly with potential wave theory).
I'm however still having problems with the case attached... any suggestions?

Thanks,
Michele.

chris May 27, 2005 06:33

It might be worth running with
 
It might be worth running with interFoam (rather than rasInterFoam) with upwind on the momentum equation. See if you hit the same problem.

michele May 27, 2005 09:44

I tried running the case even
 
I tried running the case even with interFoam (with upwind discretisation for the momentum equations), but I experienced the same issues, i.e.
1. during the first part of the simulation the free surface continues going down (instead of keeping the original height)
2. after this starting phase, instabilities near the inlet boundary are observable

Michele.

henry May 27, 2005 09:53

The conservation problem you h
 
The conservation problem you have sounds like your mesh is incorrect, have you tried running checkMesh on it? Have you tried simpler meshes?

hjasak May 27, 2005 10:01

Hello Michele, I've got an
 
Hello Michele,

I've got an idea for you, but before I can think of a fix, I need you to try something out for me. Could you please repeat the same run but set the gravity vector to zero:

constant/environmentalProperties

g g [0 1 -2 0 0 0 0] (0 0 0);

Now you should get a perfect undisturbed surface, right? If not, something else is badly wrong, but if you do, the problem will be easier to solve.

Please let me know,

Hrv

henry May 27, 2005 10:32

Check the number of faces of t
 
Check the number of faces of the front and back "empty" patches, these should be equal to the number of cells otherwise the case is not 2D. Are you sure the mesh is 1-cell thick everywhere and is properly connected? How did you generate the mesh?

henry May 30, 2005 07:03

> The checkMesh tool suggests
 
> The checkMesh tool suggests using zipUpMesh, but nothing happens

Are you sure nothing happens? Doesn't it write out the corrected mesh? If you run it a second time what does it tell you? If you run checkMesh after running zipUpMesh do you still get errors/warnings?

> Another strange issue: the mesh non-orthogonality should be zero

Does your mesh contain refined regions? If so the mesh will not be perfectly orthogonal.

michele May 30, 2005 08:46

Thank you, Henry. Yes, indeed
 
Thank you, Henry.
Yes, indeed the mesh contains refined regions.
So I assume that the non-orthogonality and skewness parameters are non zero due to the refined faces. That's OK.

Indeed, if I run the zipUpMesh tool twice, it always outputs the following message:

-------------------------------------------------------------------------------
/*---------------------------------------------------------------------------*\
| ========= | |
| \ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \ / O peration | Version: 1.1 |
| \ / A nd | Web: http://www.openfoam.org |
| \/ M anipulation | |
\*---------------------------------------------------------------------------*/

Exec : zipUpMesh . waves
Date : May 30 2005
Time : 14:34:03
Host : binah
PID : 6479
Root : /home/michele/OpenFOAM/michele-1.1/run/simulations/rasInterFoam
Case : waves
Nprocs : 1
Create time

Reading polyMesh
Cycle 1 changed 3600 faces.
Cycle 2 changed 0 faces.
-------------------------------------------------------------------------------

No one of the files into the polyMesh directory is modified (I verified to have write access to these files). There is the
polyMesh/sets/zipUpCells
file, generated during the checkMesh execution.
I don't understand what is going wrong with this mesh.

Michele

henry May 30, 2005 08:53

I am not sure why zipUpMesh do
 
I am not sure why zipUpMesh doesn't write a new mesh, we will investigate. However the test checkMesh fails on and requests running zipUpMesh is only important is you want to run finite-elements on the mesh, i.e. mesh-motion, because it relates to edge connectivity which is only used for that purpose.

eugene May 31, 2005 08:48

I haven't had a chance to look
 
I haven't had a chance to look at this case yet. Has anything changed since last week?

michele June 1, 2005 04:38

Dear Henry, thanks for the sup
 
Dear Henry, thanks for the support.
I've just uploaded on my webspace the whole test case (www.lem3.it/tests/waves_050601/waves.tar.gz). Note that in the constant/polyMesh/star_mesh directory reside the exported star mesh (the simple mesh is "star" and the mesh containing embedded regions is "star2").

I however assume that for this simulation the mesh "star2" (the one that fails the checkMesh tests) is correct.
Regarding the rasInterfoam simulation made with this mesh, the following observations came out:

1. The visualiser shows skewed elements in the refinement interface region (as it can be seen two refinements were applied)
(see fig: http://www.lem3.it/tests/waves_05060...30s/1_mesh.png)
The mesh instead is made of hexaedra with split faces on refinement. (Is it due to this the fact that checkMesh reveals skewed elements?)
Is it a problem of starToFoam conversion or only of visualisation? I think the latter, as I manually verified the vertex positions (into the imported vertices files) and they seemed to be correct.

2. One small issue: the runFoamX applications, when saving the simulation files, forgets adding the
div((nuEff*dev(grad(U).T()))) divScheme in the fvSchemes dictionary... only a little bug.


Regarting the simulations results, at time 30s, I have to make the following observations:

a. Gamma shows a wavy pattern (see fig: http://www.lem3.it/tests/waves_05060...0s/2_gamma.png).
I think that this type of instability has no hydrodynamic sense. What's happening?

b. Turbulence grows when the wavy instability starts as a coupled effect (see fig: http://www.lem3.it/tests/waves_050601/res_30s/3_k.png and http://www.lem3.it/tests/waves_05060.../4_epsilon.png).

c. pressure shows some spikes (of large magnitude, up to pd= 426Pa/density, i.e. the pressure goes up to 4.26 bar see fig: http://www.lem3.it/tests/waves_050601/res_30s/5_pd.png)

d. The graphs of velocity also show these kind of instabilities. (see Ux in fig: http://www.lem3.it/tests/waves_050601/res_30s/6_Ux.png and Uy in fig: http://www.lem3.it/tests/waves_050601/res_30s/7_Uy.png).

The fvSchemes employed in this simulation are of high order. No numerical instability problems arose.

I will appreciate any comment.
Yours sincerely,
Michele.

mattijs June 1, 2005 05:43

Hi Michele, we found a bug
 
Hi Michele,

we found a bug in the writing of zipUpMesh. It does do its work but doesn't write the new mesh!

Attached a version of zipUpMesh which does. I have however no files to test it with so let me know if it actually helps.

Compile as usual:
- unpack
- wclean
- wmake

http://www.cfd-online.com/OpenFOAM_D...hment_icon.gif zipUpMesh.tgz

henry June 1, 2005 05:43

1) The split-hex cells are de
 
1) The split-hex cells are decomposed into pyramids for ParaView because it cannot directly handle cells "higher" than hex.

2) Yes this is an omission in the FoamX configuration files which is easy to fix.

Have you tried:

making the walls slip?
lower-order schemes?
or any other schemes?
compression coefficient of 1?

henry June 1, 2005 06:28

I also notice you have a relat
 
I also notice you have a relative tolerance set for pd, I think this is unwise although it what is set for the tutorial case. For the next release I will reset and test all the settings for tutorial cases of all the interFoam-based codes. I suggest you try:

Setting relative tolerance on pd to 0

cGamma to 1

the coefficient on the Gamma201 to 1 on both terms

nNonOrthogonalCorrectors to 0

eugene June 1, 2005 06:29

I had a quick look at this cas
 
I had a quick look at this case.

Neither upwind, setting cGamma=0, changing surface tension to zero nor any combination of these factors changes the occurrance of the initial perturbation near the inlet (lower order schemes do of course damp out the larger perturbations further downstream).

However, switching off gravity makes the problem disappear very rapidly.

I'm not quite sure why this happens, but it is noteworthy that a fixedValue outlet as used in this case is not valid. Pressure should increase with depth, not remain constant.

It is not clear whether this erronous BC influences the deccelelation of the interfacial gas phase near the inlet. I think it somewhat unlikely, but I dont have the time to pinpoint the origin at the moment.

henry June 1, 2005 06:34

> but it is noteworthy that a
 
> but it is noteworthy that a fixedValue outlet as used in this case is not valid. Pressure should increase with depth, not remain constant.

Why do you think there is a problem with the oulet pressure BC? Remember we are solving for pd which is p - rho*g*h, i.e. is constant with depth at the outlet.

eugene June 1, 2005 06:46

Aah, my mistake. I was wonderi
 
Aah, my mistake. I was wondering where g was hiding.

Still, the pressure does behave strangely toward the bottom of the outlet and the removal of gravity makes the problem as described above go away. I will have another look at the formulation when I have some time.

eugene June 1, 2005 06:55

Belay that, after a second loo
 
Belay that, after a second look I dont believe there is anything wrong at or near the outlet.

michele June 1, 2005 07:07

Dear Mattijs, thanks for the
 
Dear Mattijs,
thanks for the help. I tried to compile the new zipUpMesh, but I experienced problems (the compilation of the old zipUpMesh executable was instead straightforward).
I obtained the following errors after issuing wmake:
Making dependency list for source file zipUpCells.C
Making dependency list for source file zipUpMesh.C

SOURCE_DIR=.
SOURCE=zipUpCells.C ; g++ -m64 -DlinuxAMD64 -Wall -W -Wno-unused-parameter -march=opteron -O3 -ffast-math -fno-gcse -DNoRepository -ftemplate-depth-30 -I/home/michele/OpenFOAM/OpenFOAM-1.1/src/OpenFOAM/lnInclude -IlnInclude -I. -fPIC -c $SOURCE -o Make/linuxAMD64Opt/zipUpCells.o
zipUpCells.C: In member function `void Foam::polyMesh::zipUpCells()':
zipUpCells.C:134: error: `WarningIn' undeclared (first use this function)
zipUpCells.C:134: error: (Each undeclared identifier is reported only once for each function it appears in.)
make: *** [Make/linuxAMD64Opt/zipUpCells.o] Error 1

However, if you are interested on a test mesh, the following link is that of the case under discussion:
www.lem3.it/tests/waves_050601/waves.tar.gz

Michele.
PS. Henry, I've just started a simulation with low order approximations, walls whose condition was set to 2 m/s (in rasInterFoam is not coded the slipWall condition, am I right?), and compression coefficient for gamma of 1.0. I will obtain results in about 3 hours...

mattijs June 1, 2005 07:29

Hi Michele, Sorry, didn't t
 
Hi Michele,

Sorry, didn't try it under 1.1.

Replace in zupUpCells.C

WarningIn("XXXX")

with

Warning<< "XXXX"

michele June 1, 2005 08:38

Hi Mattijs, thank you. Now
 
Hi Mattijs,
thank you.
Now the utility zipUpMesh works correctly. The corrected mesh passes all the tests of the checkMesh.

Thanks,
Michele.

michele June 1, 2005 09:10

Dear Henry: >> Why do you t
 
Dear Henry:

>> Why do you think there is a problem with the oulet pressure BC? Remember we are solving for pd which is p - rho*g*h, i.e. is constant with depth at the outlet.

Excuse me, Henry, isn't pd = p - g*h ?
I expect for pd the dimensions [m^2/s^2], am I correct?
Another question about the pressure. I'm interested in computing the pressure around floating bodies (and integrate pressure and viscous/turbulent stress over surfaces in order to obtain overall forces).
It is rather complicated to obtain pressure: if I use the expression p = pd + gh, h is not a constant as it depends on the wave elevation above each point.
That is to say that
p = pd + int(-gamma*g, ds)
an integral to be computed along the g direction (gamma is zero in air and one in water, introduced in order to capture the surface).
But below a floating body it is impossible to integrate up to a free surface (this space is occupied by the floating body itself!).
Have you any suggestion?

Ragarding the simulation, I'm making the run with the parameters you suggested.
However my impression is the following: it seems to me that boundary conditions don't pose any problem on the solution: all the instability grows well after the flow enters the domain. It resembles more a turbulence generated instability...
By the way, I would like to get the eddy viscosity among the output results... is it automatically possible or does this require coding?

henry June 1, 2005 09:49

> Excuse me, Henry, isn't pd =
 
> Excuse me, Henry, isn't pd = p - g*h ?

No pd = p - rho*(g.h) and you can calculate p from pd using this expression. Notice the dimension of pd:

dimensions [1 -1 -2 0 0 0 0];

h is the the position vectors of the cell centres and g is the gavitational force vector.

You can get the eddy viscosity from the turbulence model and write it out.

nico765 April 2, 2006 09:24

Hello, I also have a simila
 
Hello,

I also have a similar problem at the inlet of my 3d domain. As shown in the following pictures (with isosurface at gamma=0.5) (t=0.03), soon after the start of the computation, a 'wave' is developing near the inlet, and gets bigger and bigger.

meshCheck is ok.

my boundaries
pd:
boundaryField
{
inlet
{
type zeroGradient;
}

inlet:002
{
type zeroGradient;
}

outlet
{
type zeroGradient;
}

bottom
{
type zeroGradient;
}

top
{
type totalPressure;
p0 uniform 0;
value uniform 0;
}


right
{
type zeroGradient;
}
hull
{
type zeroGradient;
}
sym
{
type zeroGradient;
}
}

U:
boundaryField
{
inlet
{
type fixedValue;
value uniform (0.5 0 0);
}
inlet:002
{
type fixedValue;
value uniform (2 0 0);
}
outlet
{
type zeroGradient;
}
top
{
type pressureInletOutletVelocity;
phi phi;
rho rho;
value uniform (0 0 0);
}
bottom
{
type fixedValue;
value uniform (0 0 0);
}
sym
{
type fixedValue;
value uniform (0 0 0);
}
right
{
type fixedValue;
value uniform (0 0 0);
}
hull
{
type fixedValue;
value uniform (0 0 0);
}
}

http://www.cfd-online.com/OpenFOAM_D...ges/1/2101.jpg
http://www.cfd-online.com/OpenFOAM_D...ges/1/2103.jpg
http://www.cfd-online.com/OpenFOAM_D...ges/1/2102.jpg

liu April 2, 2006 14:38

As I said before, OF is not re
 
As I said before, OF is not ready for open channel yet. At least the inlet and outlet boundary conditions will cause problem.

kumar2 April 2, 2006 20:56

Hi liu,nicolas i did numeri
 
Hi liu,nicolas

i did numerical simulations of a 2d hydrofoil in a channel. my preliminary results are 100% ok. i also had issues with the correct BC. but the VOF works for open channel flows. nicolas , please see my posts rasInterFoam - STRANGE RESULTS AT BOUNDARY . after i incorporated Hrv's suggestions my simulations went smooth

regards

kumar

liu April 3, 2006 11:47

Hi, kumar, Which post are you
 
Hi, kumar,
Which post are you refering to? I can't find it. Can you upload your case? I am so interested to see where I did wrong.

Thank you.

nico765 April 3, 2006 14:22

It's here http://www.cfd-onli
 
It's here
http://www.cfd-online.com/cgi-bin/Op...show.cgi?1/656

nico765 April 3, 2006 14:25

oops, actually it s this one
 
oops, actually it s this one

http://www.cfd-online.com/OpenFOAM_D...es/1/2014.html

nico765 April 4, 2006 04:07

I have checked the thread, but
 
I have checked the thread, but i do something similar in my case. My problem is at the inlet and not at the outlet.

eugene July 4, 2006 10:06

Can anyone supply or point me
 
Can anyone supply or point me toward an implementation of outlet boundary conditions for free surface wave transmission?

khleitz November 5, 2007 10:27

Hallo, I have the following
 
Hallo,
I have the following problem: I have a solid surrounded by air that moves from the left to the right with a constant velocity. However I get a very high pressure wave on the interface solid-air at the inlet and at the outlet. I have the feeling, that something with the boundary conditions is wrong. Can anybody give me a hit what the problem could be?
Regards,
Karl-Heinz
P.S. I would like to add a picture of my problem. Can anybody tell me how to do that?

khleitz November 5, 2007 10:35

Hallo, concerning my problem
 
Hallo,
concerning my problem with the strange pressure patterns on the interface solid air. Here is a picture of my current model.
http://www.cfd-online.com/OpenFOAM_D...ges/1/5878.jpg
I use a modified interFoam solver. However these pressure waves make my silulation unstable.
Best regards,
Karl-Heinz

santiagomarquezd September 29, 2009 18:31

Hi all, I'm dealing with surface problems in interFoam, the post is:

http://www.cfd-online.com/Forums/ope...-sloshing.html

Any comments are welcome. Regards.

ehsan January 7, 2012 22:43

VOF Method
 
Dear All

Could I ask you whether you could tell me which VOF technique is used in OpenFOAM? I mean there are different VOF models like Hirt and Nicols, Youngs, and so on. Which one is used in OpenFOAM?

Regards
Ehsan

nimasam January 8, 2012 09:45

Dear ehsan
1) please ask your question just in one post, not in several posts!
2) look this paper:
Quote:

BRACKBILL, J. U., KOTHE, D. B. & ZEMACH, C. 1992. A Continuum Method for Modeling Surface Tension. Journal of Computational Physics, 100, 335-354.
this method reduces need for interface reconstruction, this method interface is smeared among 2 or 3 cells and it uses a high order schemes
and interface compression method to keep interface sharp
3) if you uses interFoam solver look at here:
Quote:

BERBEROVIĆ, E., HINSBERG, N. P. V., JAKIRLIĆ, S., ROISMAN, I. V. & TROPEA, C. 2009. Drop impact onto a liquid layer of finite thickness: Dynamics of the cavity evolution. The American Physical Society, 79, 036306 (15).
4)if you look for more, look at here:
Quote:

Rusche, H., Computational Fluid Dynamics of Dispersed Two-Phase Flows at High Phase Fractions,
in Department of Mechanical Engineering. 2002, Imperial college of Science: University of London


All times are GMT -4. The time now is 12:09.