Skip to content
GitLab
Closed
Issue created by Robin Kamenicky@Robin

Turbulence models (kEpsilon, kOmegaSSTSato) TwoPhaseSystem - preconditioner

Summary

When kEpsolin, kOmegaSST, kOmegaSSTSato is used and the alpha fields are initialized via utility such as setField to value 0 or 1 the preconditioner gives floating point exception during usage of preconditioner class. This happens during first iteration.

Steps to reproduce

Use the tutorial $FOAM_TUTORIALS/heatTransfer/chtMultiRegionTwoPhaseEulerFoam/solidQuenching2D. Change the constant/water/turbulenceProperties.liquid from

simulationType RAS; RAS { RASModel mixtureKEpsilon;//kEpsilon; turbulence on; printCoeffs on; }

to

simulationType RAS; RAS { RASModel kEpsilon; turbulence on; printCoeffs on; }

and change file constant/water/turbulenceProperties.gas from

simulationType RAS; RAS { RASModel mixtureKEpsilon; turbulence on; printCoeffs on; }

to

simulationType RAS; RAS { RASModel kEpsilon; turbulence on; printCoeffs on; }

The same happens if we use for instance a combination for gas-liquid kEpsilon-laminar or kOmegaSST-kOmegaSSTSato.

The error is

#1  Foam::sigFpe::sigHandler(int) at ~/OpenFOAM/OpenFOAM-v1906-debug/src/OSspecific/POSIX/signals/sigFpe.C:130                                       
#2  ? in /lib/x86_64-linux-gnu/libc.so.6                                                                                                     
#3  Foam::DILUPreconditioner::calcReciprocalD(Foam::Field<double>&, Foam::lduMatrix const&) at ~/OpenFOAM/OpenFOAM-v1906-debug/src/OpenFOAM/matrices/lduMatrix/preconditioners/DILUPreconditioner/DILUPreconditioner.C:80 (discriminator 2)

It happens already at first construction of the class. It seems that division by 0 occurs on line 80 in DILUPreconditioner.C. Change of preconditioner or solver does not help. The only option is not to use preconditioner at all. When that is done, the simulation runs smoothly and the results seem to be fine.

Example case

Described in the steps to reproduce.

What is the current bug behaviour?

SIGFPE when preconditioner used with kEpsilon, kOmegaSST, kOmegaSSSato and internalField of phase volume fractions initialized to 0 or 1.

What is the expected correct behavior?

No SIGFPE.

Relevant logs and/or images

Environment information

  • OpenFOAM version : v1906
  • Operating system : Ubuntu 16.04.6 LTS
  • Hardware info :
  • Compiler : gcc

Possible fixes

Bound the phase fraction in $FOAM_SRC/phaseSystemModels/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseSystem/twoPhaseSystem.C:342

from

alpha1.clip(0, 1);

to

alpha1.clip(SMALL, 1 - SMALL);

the same will probably happen also for twoPhaseEulerFoam. I have not tested that. the fix would be in $FOAM_APP/solvers/multiphase/twoPhaseEulerFoam/twoPhaseSystem/twoPhaseSystem.C

where lines 532, 533 would need to be changed from

alpha1.max(0); alpha1.min(1);

to

alpha1.max(SMALL); alpha1.min(1 - SMALL);

Edited by Robin Kamenicky