[shipman]

Hi to everyone,

I am trying to simulate cavitation inside nozzle using interPhaseChangeFoam. I am using kOmega SST turbulence model. When i run the case, after 40 or 50 time steps always giving following floating point error:

smoothSolver:  Solving for alpha.water, Initial residual = 4.1184e-06, Final residual = 3.80347e-09, No Iterations 1
Phase-1 volume fraction = 0.98612  Min(alpha1) = 0.000310125  Max(alpha1) = 1
MULES: Correcting alpha.water
MULES: Correcting alpha.water
Liquid phase volume fraction = 0.98612  Min(alpha1) = 0.000310125  Max(alpha1) = 1
smoothSolver:  Solving for Ux, Initial residual = 1.98514e-05, Final residual = 1.02782e-09, No Iterations 2
smoothSolver:  Solving for Uy, Initial residual = 5.44074e-06, Final residual = 1.94705e-10, No Iterations 2
smoothSolver:  Solving for Uz, Initial residual = 8.8507e-05, Final residual = 3.80435e-09, No Iterations 2
GAMG:  Solving for p_rgh, Initial residual = 2.79212e-08, Final residual = 2.04895e-10, No Iterations 1
GAMG:  Solving for p_rgh, Initial residual = 6.59936e-09, Final residual = 6.59936e-09, No Iterations 0
GAMG:  Solving for p_rgh, Initial residual = 6.61102e-09, Final residual = 6.61102e-09, No Iterations 0
GAMGPCG:  Solving for p_rgh, Initial residual = 6.61138e-09, Final residual = 6.61138e-09, No Iterations 0
smoothSolver:  Solving for omega, Initial residual = 5.4048e-06, Final residual = 6.75835e-10, No Iterations 2
smoothSolver:  Solving for k, Initial residual = 9.03591e-06, Final residual = 2.43949e-09, No Iterations 2
ExecutionTime = 97153.3 s  ClockTime = 98393 s

Max pressure: 323086
Min pressure: -35072.5
Max velocity: 28.9586

Min velocity: 0

Courant Number mean: 0.00779459 max: 0.0996774
deltaT = 3.50877e-08
Time = 0.01051926

smoothSolver:  Solving for alpha.water, Initial residual = 4.11863e-06, Final residual = 3.80331e-09, No Iterations 1
Phase-1 volume fraction = 0.986119  Min(alpha1) = 0.000310159  Max(alpha1) = 1
MULES: Correcting alpha.water
MULES: Correcting alpha.water
Liquid phase volume fraction = 0.986119  Min(alpha1) = 0.000310159  Max(alpha1) = 1
smoothSolver:  Solving for Ux, Initial residual = 1.98021e-05, Final residual = 1.02634e-09, No Iterations 2
smoothSolver:  Solving for Uy, Initial residual = 5.44207e-06, Final residual = 1.94729e-10, No Iterations 2
smoothSolver:  Solving for Uz, Initial residual = 8.84228e-05, Final residual = 3.80449e-09, No Iterations 2
GAMG:  Solving for p_rgh, Initial residual = 2.58283e-08, Final residual = 1.86337e-10, No Iterations 1
[5] #0  Foam::error::printStack(Foam::Ostream&) at ??:?
[5] #1  Foam::sigFpe::sigHandler(int) at ??:?
[5] #2   in "/lib/x86_64-linux-gnu/libc.so.6"
[5] #3  Foam::divide(Foam::Field<double>&, Foam::UList<double> const&, Foam::UList<double> const&) at ??:?
[5] #4  Foam::tmp<Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh> > Foam::operator/<Foam::fvPatchField, Foam::volMesh>(Foam::tmp<Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh> > const&, Foam::tmp<Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh> > const&) at ??:?
[5] #5  Foam::phaseChangeTwoPhaseMixtures::SchnerrSauer::pCoeff(Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh> const&) const at ??:?
[5] #6  Foam::phaseChangeTwoPhaseMixtures::SchnerrSauer::mDotP() const at ??:?
[5] #7  Foam::phaseChangeTwoPhaseMixture::vDotP() const at ??:?
[5] #8
[5]  at ??:?
[5] #9  __libc_start_main in "/lib/x86_64-linux-gnu/libc.so.6"
[5] #10
[5]  at ??:?
[baris-desktop:10415] *** Process received signal ***
[baris-desktop:10415] Signal: Floating point exception (8)
[baris-desktop:10415] Signal code:  (-6)
[baris-desktop:10415] Failing at address: 0x3e8000028af
[baris-desktop:10415] [ 0] /lib/x86_64-linux-gnu/libc.so.6(+0x36d40) [0x7ff21608cd40]
[baris-desktop:10415] [ 1] /lib/x86_64-linux-gnu/libc.so.6(gsignal+0x39) [0x7ff21608ccc9]
[baris-desktop:10415] [ 2] /lib/x86_64-linux-gnu/libc.so.6(+0x36d40) [0x7ff21608cd40]
[baris-desktop:10415] [ 3] /home/baris/OpenFOAM/OpenFOAM-2.3.1/platforms/linux64GccDPOpt/lib/libOpenFOAM.so(_ZN4Foam6divideERNS_5FieldIdEERKNS_5UListIdEES6_+0xd5) [0x7ff21734cc05]
[baris-desktop:10415] [ 4] /home/baris/OpenFOAM/baris-2.3.1/platforms/linux64GccDPOpt/lib/libphaseChangeTwoPhaseMixtures.so(_ZN4FoamdvINS_12fvPatchFieldENS_7volMeshEEENS_3tmpINS_14GeometricFieldIdT_T0_EEEERKS8_SA_+0x19b) [0x7ff21ae2aaab]
[baris-desktop:10415] [ 5] /home/baris/OpenFOAM/baris-2.3.1/platforms/linux64GccDPOpt/lib/libphaseChangeTwoPhaseMixtures.so(_ZNK4Foam27phaseChangeTwoPhaseMixtures12SchnerrSauer6pCoeffERKNS_14GeometricFieldIdNS_12fvPatchFieldENS_7volMeshEEE+0x33c) [0x7ff21ae2f3cc]
[baris-desktop:10415] [ 6] /home/baris/OpenFOAM/baris-2.3.1/platforms/linux64GccDPOpt/lib/libphaseChangeTwoPhaseMixtures.so(_ZNK4Foam27phaseChangeTwoPhaseMixtures12SchnerrSauer5mDotPEv+0xff) [0x7ff21ae2ff8f]
[baris-desktop:10415] [ 7] /home/baris/OpenFOAM/baris-2.3.1/platforms/linux64GccDPOpt/lib/libphaseChangeTwoPhaseMixtures.so(_ZNK4Foam26phaseChangeTwoPhaseMixture5vDotPEv+0xa2) [0x7ff21ae18242]
[baris-desktop:10415] [ 8] MRinterPhaseChangeFoam() [0x431350]
[baris-desktop:10415] [ 9] /lib/x86_64-linux-gnu/libc.so.6(__libc_start_main+0xf5) [0x7ff216077ec5]
[baris-desktop:10415] [10] MRinterPhaseChangeFoam() [0x433792]
[baris-desktop:10415] *** End of error message ***

When i check the results Pressure, Velocity and alpha distributions seem OK. At the beginning I thought that it is due to mesh problem. Although gradually i have increased the mesh number, still i am getting same ERROR. Finally, i recognized that I am getting this floating error since Omega value is excessively high.

My initial condition for k and omega are as follows:

==>k:

internalField   uniform 0.0234;

boundaryField
{
    inlet
    {
         type            fixedValue;
        value           $internalField;
    }
    outlet
    {
       type            inletOutlet;
        inletValue      $internalField;
        value           $internalField;
    }
    wall
    {
        type            kqRWallFunction;
        value           uniform 0.0234;
    }
 
}

==>Omega:

internalField   uniform 174.69;

boundaryField
{
    inlet
    {
         type            fixedValue;
        value           $internalField;
    }
    outlet
    {
        type            inletOutlet;
        inletValue      $internalField;
        value           $internalField;
    }
    wall
    {
        type            epsilonWallFunction;
        value           uniform 0.315;
    }

Any comment and help will be really appreciated.

Thanks in advance.

[alexeym]

Hi,

Are you sure that the problem is caused by kOmegaSST? Maybe divergence of omega equation is a consequence of certain other reasons?

As to your error, it happens in SchnerrSauer:Coeff method and the error is division by zero, the only division is in the return statement of the method:

Foam::tmp<Foam::volScalarField>
Foam::phaseChangeTwoPhaseMixtures::SchnerrSauer::pCoeff
(
    const volScalarField& p
) const
{
    ...
    return
        (3*rho1()*rho2())*sqrt(2/(3*rho1()))
       *rRb(limitedAlpha1)/(rho*sqrt(mag(p - pSat()) + 0.01*pSat()));
}

Guess rho1 can not be zero, so this expression

rho*sqrt(mag(p - pSat()) + 0.01*pSat())

is zero. And it can be zero only if rho is zero or your pSat is negative.

[shipman]

Hi Alexey,

Thank you for prompt answer. Yes you are right. I also recognized it and therefore i deleted 0.01*pSat() term in order avoid negative value inside of square root.

So i used like

Foam::tmp<Foam::volScalarField>
Foam::phaseChangeTwoPhaseMixtures::SchnerrSauer::pCoeff
(
    const volScalarField& p
) const
{
    ...
    return
        (3*rho1()*rho2())*sqrt(2/(3*rho1()))
       *rRb(limitedAlpha1)/(rho*sqrt(mag(p - pSat()));
}

In this way inside of square root should be always positive even pSat is negative. However, still i am getting same error.

Do you have any idea?

Thank you

[alexeym]

Hi,

And if p == pSat()?

As I don't know even order of magnitudes for the pressure, my proposition will be quite naive, rewrite denominator

(rho*sqrt(mag(p - pSat()) + 0.01*pSat()))

as

(rho*sqrt(mag(p - pSat()) + 0.01*pSat()) + <small value>)

where <small value> is a dimensionedScalar with "rho*sqrt(mag(p - pSat()) + 0.01*pSat())" units. As rho >= 0, sqrt(...) >= 0, <small value> is positive, so denominator should be always positive.

[shipman]

Hi Alexey,

Thank you so much.

Answer is quite clear.

BR.

[esujby]

Hello,

I have a similar issue, the case was initially running, but now i get this error and i can't figure out what exactly happened. please help..

Region: inlet Courant Number mean: 3.31401e-08 max: 6.96144e-07
Region: outlet Courant Number mean: 0.912229 max: 4418.7
Region: insulator Diffusion Number mean: 1.90139e-11 max: 1.53073e-10
Region: s1 Diffusion Number mean: 2.86584e-08 max: 6.24037e-08
Region: s2 Diffusion Number mean: 2.87506e-08 max: 6.11088e-08
Region: s3 Diffusion Number mean: 2.90832e-08 max: 5.7235e-08
Region: s4 Diffusion Number mean: 2.92442e-08 max: 5.69394e-08
Region: s5 Diffusion Number mean: 2.90227e-08 max: 5.37904e-08
Region: s6 Diffusion Number mean: 2.90552e-08 max: 5.36443e-08
Region: s7 Diffusion Number mean: 2.89756e-08 max: 5.78914e-08
Region: s8 Diffusion Number mean: 2.91473e-08 max: 6.40512e-08
Region: s9 Diffusion Number mean: 2.86969e-08 max: 6.01415e-08
Region: s10 Diffusion Number mean: 2.92766e-08 max: 5.74572e-08
Region: s11 Diffusion Number mean: 1.95274e-08 max: 3.83532e-08
Region: s12 Diffusion Number mean: 1.96676e-08 max: 3.79916e-08
Region: s13 Diffusion Number mean: 1.96613e-08 max: 4.33353e-08
Region: s14 Diffusion Number mean: 1.95072e-08 max: 4.04921e-08
Region: s15 Diffusion Number mean: 1.9569e-08 max: 4.03973e-08
Region: lens Diffusion Number mean: 5.26653e-12 max: 2.47215e-11
deltaT = 6.63239e-14
Time = 2.42058e-05


Solving for fluid region inlet
diagonal:  Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
DILUPBiCG:  Solving for Ux, Initial residual = 0.01466, Final residual = 9.35381e-11, No Iterations 2
DILUPBiCG:  Solving for Uy, Initial residual = 0.013276, Final residual = 8.35837e-11, No Iterations 2
DILUPBiCG:  Solving for Uz, Initial residual = 0.0133429, Final residual = 1.20532e-10, No Iterations 2
DILUPBiCG:  Solving for h, Initial residual = 0.894301, Final residual = 9.70702e-09, No Iterations 2
Min/max T:300 300
GAMG:  Solving for p_rgh, Initial residual = 1.16314e-05, Final residual = 1.34599e-10, No Iterations 1
diagonal:  Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
time step continuity errors (inlet): sum local = 1.50502e-16, global = -3.36401e-19, cumulative = -1.3015e-10
GAMG:  Solving for p_rgh, Initial residual = 1.43257e-09, Final residual = 1.43257e-09, No Iterations 0
diagonal:  Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
time step continuity errors (inlet): sum local = 1.50614e-16, global = -2.90421e-19, cumulative = -1.3015e-10

Solving for fluid region outlet
diagonal:  Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
DILUPBiCG:  Solving for Ux, Initial residual = 0.00636657, Final residual = 7.72967e-09, No Iterations 2
DILUPBiCG:  Solving for Uy, Initial residual = 0.00697353, Final residual = 4.10965e-08, No Iterations 2
DILUPBiCG:  Solving for Uz, Initial residual = 0.0064832, Final residual = 1.03779e-09, No Iterations 2
DILUPBiCG:  Solving for h, Initial residual = 0.858225, Final residual = 2.4626e-09, No Iterations 3
[13] 
[13] 
[13] --> FOAM FATAL ERROR: 
[13] Maximum number of iterations exceeded
[13] 
[13]     From function Foam::scalar Foam::species::thermo<Thermo, Type>::T(Foam::scalar, Foam::scalar, Foam::scalar, Foam::scalar (Foam::species::thermo<Thermo, Type>::*)(Foam::scalar, Foam::scalar)const, Foam::scalar (Foam::species::thermo<Thermo, Type>::*)(Foam::scalar, Foam::scalar)const, Foam::scalar (Foam::species::thermo<Thermo, Type>::*)(Foam::scalar)const) const [with Thermo = Foam::hConstThermo<Foam::perfectGas<Foam::specie> >; Type = Foam::sensibleEnthalpy; Foam::scalar = double; Foam::species::thermo<Thermo, Type> = Foam::species::thermo<Foam::hConstThermo<Foam::perfectGas<Foam::specie> >, Foam::sensibleEnthalpy>]
[13]     in file /home/parallels/OpenFOAM/OpenFOAM-3.0.x/src/thermophysicalModels/specie/lnInclude/thermoI.H at line 66.
[13] 
FOAM parallel run aborting
[13] 
[13] #0  Foam::error::printStack(Foam::Ostream&) at ??:?
[13] #1  Foam::error::abort() at ??:?
[13] #2  Foam::species::thermo<Foam::hConstThermo<Foam::perfectGas<Foam::specie> >, Foam::sensibleEnthalpy>::THs(double, double, double) const at ??:?
[13] #3  Foam::heRhoThermo<Foam::rhoThermo, Foam::pureMixture<Foam::constTransport<Foam::species::thermo<Foam::hConstThermo<Foam::perfectGas<Foam::specie> >, Foam::sensibleEnthalpy> > > >::calculate() at ??:?
[13] #4  Foam::heRhoThermo<Foam::rhoThermo, Foam::pureMixture<Foam::constTransport<Foam::species::thermo<Foam::hConstThermo<Foam::perfectGas<Foam::specie> >, Foam::sensibleEnthalpy> > > >::correct() at ??:?
[13] #5  ? at ??:?
[13] #6  __libc_start_main in "/lib/x86_64-linux-gnu/libc.so.6"
[13] #7  ?--------------------------------------------------------------------------
MPI_ABORT was invoked on rank 13 in communicator MPI_COMM_WORLD 
with errorcode 1.

NOTE: invoking MPI_ABORT causes Open MPI to kill all MPI processes.
You may or may not see output from other processes, depending on
exactly when Open MPI kills them.
--------------------------------------------------------------------------
 at ??:?
--------------------------------------------------------------------------
mpirun has exited due to process rank 13 with PID 3457 on
node ubuntu exiting improperly. There are two reasons this could occur:

1. this process did not call "init" before exiting, but others in
the job did. This can cause a job to hang indefinitely while it waits
for all processes to call "init". By rule, if one process calls "init",
then ALL processes must call "init" prior to termination.

2. this process called "init", but exited without calling "finalize".
By rule, all processes that call "init" MUST call "finalize" prior to
exiting or it will be considered an "abnormal termination"

This may have caused other processes in the application to be
terminated by signals sent by mpirun (as reported here).
--------------------------------------------------------------------------

kind regards