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

Perfect fluid implementation - interFoam

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

Reply
 
LinkBack Thread Tools Display Modes
Old   January 25, 2013, 10:02
Question Perfect fluid implementation - interFoam
  #1
New Member
 
Gaetano
Join Date: Jul 2012
Posts: 18
Rep Power: 3
Gaetano is on a distinguished road
Hi all.

This is my first post here, so I'm going to introduce myself: my name's Gaetano and I'm doing a PhD in Chemical Engineering in Naples (Italy).

I've been lurking those forums for at least 6 months and I found answers to almost all my questions. This one, however, seems to be a bit more tricky.

I'm trying to have interFoam working with one of the two phases being a perfect fluid, i.e. having zero density and viscosity (I'd like to simply put "0" in the dictionary transportProperties). In doing so, I need to avoid all the division by density, if any. I found some in class twoPhasesMixture and I created a new class twoPhaseMixturePerfect with a little (tricky) workaround in it. When I tried to solve the damBreack tutorial with my interPerfectFoam I found this error at the very first step:

------------------------------------------------------------------------
Starting time loop

Courant Number mean: 0 max: 0
Interface Courant Number mean: 0 max: 0
deltaT = 0.00119048
Time = 0.00119048

MULES: Solving for alpha1
Phase-1 volume fraction = 0.130194 Min(alpha1) = 0 Max(alpha1) = 1
MULES: Solving for alpha1
Phase-1 volume fraction = 0.130194 Min(alpha1) = 0 Max(alpha1) = 1
#0 Foam::error:rintStack(Foam::Ostream&) in "/share/OpenFOAM/OpenFOAM-2.1.1/platforms/linux64GccDPOpt/lib/libOpenFOAM.so"
#1 Foam::sigFpe::sigHandler(int) in "/share/OpenFOAM/OpenFOAM-2.1.1/platforms/linux64GccDPOpt/lib/libOpenFOAM.so"
#2 __restore_rt at sigaction.c:0
#3 Foam::divide(Foam::Field<double>&, double const&, Foam::UList<double> const&) in "/share/OpenFOAM/OpenFOAM-2.1.1/platforms/linux64GccDPOpt/lib/libOpenFOAM.so"
#4 main in "/local/home/iitcrib3/GaetanoDM/damBreak/interPerfectFoam"
#5 __libc_start_main in "/lib64/libc.so.6"
#6 Foam::regIOobject::writeObject(Foam::IOstream::str eamFormat, Foam::IOstream::versionNumber, Foam::IOstream::compressionType) const in "/local/home/iitcrib3/GaetanoDM/damBreak/interPerfectFoam"
Floating point exception
------------------------------------------------------------------------

It seems that there is a problem with a division (see #3 above) in the p_rgh solver (as the next step should be "DICPCG: Solving for p_rgh ..."), but my insight into the problem ends here.

So, here's my question: what does the above error mean? Is there any division with a density in the denominator anywhere in interFoam?

But, besides the answers, I'm also looking for advice from expert foamers like many here: is there any chance that my attempt to "hack" interFoam could fit in its "architecture"? I mean: does the algorithm require a division by a density? I'd like to modify only the "shell" and not the "core" of interFoam (sorry for the poor analogy: I can't figure out any better way to make my point).

Thanks in advance,
Gaetano
Gaetano is offline   Reply With Quote

Old   January 25, 2013, 10:20
Default
  #2
Member
 
Michiel
Join Date: Oct 2010
Location: Delft, Netherlands
Posts: 96
Rep Power: 5
michielm is on a distinguished road
Are you sure that the density should be 0?! As far as I know (and can find) the definition of a perfect fluid is only that the viscosity is 0, not that the density is 0.
michielm is offline   Reply With Quote

Old   January 25, 2013, 10:33
Default
  #3
New Member
 
Gaetano
Join Date: Jul 2012
Posts: 18
Rep Power: 3
Gaetano is on a distinguished road
Quote:
Originally Posted by michielm View Post
Are you sure that the density should be 0?! As far as I know (and can find) the definition of a perfect fluid is only that the viscosity is 0, not that the density is 0.
I think that would be an Euler fluid, but I might have chosen a wrong name for my fluid. I confirm that I'm looking for something with no inertia and no viscosity.
Gaetano is offline   Reply With Quote

Old   January 25, 2013, 10:44
Default
  #4
Member
 
Michiel
Join Date: Oct 2010
Location: Delft, Netherlands
Posts: 96
Rep Power: 5
michielm is on a distinguished road
Ok, in that case I would suggest using interTrackFoam from the OpenFOAM-extend project. If you define a single mesh for fluid 1 and no mesh for fluid 2, you are basically assuming fluid 2 to have 0 density and 0 viscosity.
michielm is offline   Reply With Quote

Old   January 25, 2013, 12:19
Default
  #5
ngj
Senior Member
 
Niels Gjoel Jacobsen
Join Date: Mar 2009
Location: Rotterdam, The Netherlands
Posts: 1,469
Rep Power: 23
ngj will become famous soon enoughngj will become famous soon enough
And a small detail with respect to the density of the fluid. You will end up dividing by zero, since you create a term, which out of memory, is called rAU in the pEqn.H file. This is essentially an interpolation of your diagonal coefficients in your matrix for the momentum equation. These coefficients scale by rho (density).

The _inverse_ of this field is acting as a diffusion coefficient for the pressure on the Poisson equation, hence you are indeed dividing by zero.

Kind regards,

Niels
ngj is offline   Reply With Quote

Old   January 25, 2013, 14:42
Default
  #6
NJG
Member
 
Nick Gutschow
Join Date: Jan 2013
Posts: 36
Rep Power: 3
NJG is on a distinguished road
I am not too expereinced in OpenFOAM, but can you do the ol'e 0~ some tiny tiny number? That has worked for me in other modeling programs before.

-NG
NJG is offline   Reply With Quote

Old   January 26, 2013, 05:27
Default
  #7
Senior Member
 
akidess's Avatar
 
Anton Kidess
Join Date: May 2009
Location: Delft, Netherlands
Posts: 912
Rep Power: 15
akidess will become famous soon enough
Quote:
Originally Posted by NJG View Post
I am not too expereinced in OpenFOAM, but can you do the ol'e 0~ some tiny tiny number? That has worked for me in other modeling programs before.

-NG
While this is an idea worth trying, I suspect spurious currents to kill you here, as they scale with the density and viscosity ratio between both fluids.
__________________
*On twitter @akidTwit
*Spend as much time formulating your questions as you expect people to spend on their answer.
*Check out the scientific computing exchange http://scicomp.stackexchange.com
akidess is offline   Reply With Quote

Old   January 26, 2013, 12:54
Default
  #8
New Member
 
Gaetano
Join Date: Jul 2012
Posts: 18
Rep Power: 3
Gaetano is on a distinguished road
Well, I've already thought about the possibility of rho=eps, but I'm wondering if it would be possible to have it *exactly* zero.

Here's where my evils turn true (bold mine):

Quote:
Originally Posted by ngj View Post
The _inverse_ of this field (rAU) is acting as a diffusion coefficient for the pressure on the Poisson equation, hence you are indeed dividing by zero.
It seems that I do have to play with the matrix coefficients, and I'd really like to avoid such a task. But still: would it be possible to rewrite the term rAU? And most important: is it worth an attempt? I'm not really a numeric guy...
Gaetano is offline   Reply With Quote

Reply

Tags
floating point exception, interfoam, perfect fluid, zero density, zero viscosity

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


Similar Threads
Thread Thread Starter Forum Replies Last Post
Water subcooled boiling Attesz CFX 7 January 5, 2013 03:32
Locating and observing a transient fluid phenomena Chander CFX 2 September 25, 2011 18:49
error using combination of step function xujjun CFX 1 January 15, 2008 16:46
What is the total energy for incompressible fluid? Harry Dong Main CFD Forum 12 February 4, 2006 00:55
Help With Modeling A Projectile Fluid In Flight. Dzeff G. Main CFD Forum 0 December 10, 1998 21:24


All times are GMT -4. The time now is 00:26.