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

Floating point exceptions with own modified k-omega SST with Reboud's correction

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

Reply
 
LinkBack Thread Tools Display Modes
Old   September 23, 2013, 09:37
Default Floating point exceptions with own modified k-omega SST with Reboud's correction
  #1
Member
 
Join Date: Jul 2011
Posts: 33
Rep Power: 5
A_Pete is on a distinguished road
Hi everyone!

I have implemented the correction of turbulent viscosity by Reboud in interface regions for cavitating flows. Physically it does everything I want pretty well, but from time to time my solutions will abort due to floating point exception errors.
For me this is an indication that I didn't prevent the solver from dividing by zero somewhere. So I found a spot where this might happen during the solution, altough it shouldn't and implemented a small prevention.
But that still doesn't get the job done, because I will still get the same FPEs and have no idea where to look for. This is often accompanied by omega getting very large first. But since this doesn't happen with other turbulence models in my simulations something has to be wrong with my code.

What I have done is take the standard kOmegaSST turbulence model of OpenFOAM and modify it the way I need it.

The error says this:
Code:
[0] #0  Foam::error::printStack(Foam::Ostream&) in
"/home/peters/OpenFOAM/OpenFOAM-2.2.x/platforms/linux64GccDPOpt/lib/libOpenFOAM.so"
[0] #1  Foam::sigFpe::sigHandler(int) in
"/home/peters/OpenFOAM/OpenFOAM-2.2.x/platforms/linux64GccDPOpt/lib/libOpenFOAM.so"
[0] #2   in "/lib/x86_64-linux-gnu/libc.so.6"
[0] #3  double Foam::sumProd<double>(Foam::UList<double> const&,
Foam::UList<double> const&) in
"/home/peters/OpenFOAM/OpenFOAM-2.2.x/platforms/linux64GccDPOpt/lib/libOpenFOAM.so"
[0] #4  Foam::PBiCG::solve(Foam::Field<double>&, Foam::Field<double> const&,
unsigned char) const in
"/home/peters/OpenFOAM/OpenFOAM-2.2.x/platforms/linux64GccDPOpt/lib/libOpenFOAM.so"
[0] #5  Foam::fvMatrix<double>::solveSegregated(Foam::dictionary const&) in
"/home/peters/OpenFOAM/OpenFOAM-2.2.x/platforms/linux64GccDPOpt/lib/libfiniteVolume.so"
[0] #6  Foam::fvMatrix<double>::solve(Foam::dictionary const&) in
"/home/peters/OpenFOAM/OpenFOAM-2.2.x/platforms/linux64GccDPOpt/bin/interPhaseChangeDyMFoam"
[0] #7  Foam::fvMatrix<double>::solve() in
"/home/peters/OpenFOAM/OpenFOAM-2.2.x/platforms/linux64GccDPOpt/bin/interPhaseChangeDyMFoam"
[0] #8  Foam::SolverPerformance<double>
Foam::solve<double>(Foam::tmp<Foam::fvMatrix<double> > const&) in
"/home/peters/OpenFOAM/OpenFOAM-2.2.x/platforms/linux64GccDPOpt/bin/interPhaseChangeDyMFoam"
[0] #9  Foam::incompressible::RASModels::modKOmegaSST::correct() in
"/home/peters/OpenFOAM/OpenFOAM-2.2.x/platforms/linux64GccDPOpt/lib/libincompressibleRASModels.so"
[0] #10  
[0]  in
"/home/peters/OpenFOAM/OpenFOAM-2.2.x/platforms/linux64GccDPOpt/bin/interPhaseChangeDyMFoam"
[0] #11  __libc_start_main in "/lib/x86_64-linux-gnu/libc.so.6"
[0] #12  
[0]  in
"/home/peters/OpenFOAM/OpenFOAM-2.2.x/platforms/linux64GccDPOpt/bin/interPhaseChangeDyMFoam"
[prop-04:28417] *** Process received signal ***
[prop-04:28417] Signal: Floating point exception (8)
...
Something seems to go wrong at #9 in the correct() method of modKOmegaSST, but I can't find it. Below are what I suppose to be the interesting parts of my code regarding this issue.


In the header file I added this:
Code:
// under member functions
//- Calculate Reboud's correction and return the normalized density function

virtual volScalarField normDensFunc();
And in the source file the interesting parts are the following:
Code:
// definition of the calculated function under member functions

volScalarField modKOmegaSST::normDensFunc()
{
    //- Reading fields for Reboud's correction
    const volScalarField& alpha1 = mesh_.lookupObject<volScalarField>("alpha1");
    const volScalarField& alpha2 = mesh_.lookupObject<volScalarField>("alpha2");

    //- Lookup transportProperties dictionary to recieve fluid properties
    const dictionary& transportProperties =
        mesh_.lookupObject<IOdictionary>("transportProperties");
    const dictionary& phase1 = transportProperties.subDict("phase1");
    const dictionary& phase2 = transportProperties.subDict("phase2");
    const dimensionedScalar& rho1 = phase1.lookup("rho");
    const dimensionedScalar& rho2 = phase2.lookup("rho");

    //- Calculation of mixture density
    const volScalarField rho_m = alpha1*rho1 + alpha2 * rho2;

    //- Calculation of density function
    const volScalarField f_rho =
        pow(( rho_m - rho2 ), nReb_) / pow((rho1-rho2), (nReb_-1)) + rho2;
    
    //- Calculation of ratio of density function and mixture density
    const volScalarField fRhoNorm = 
        f_rho / max(rho_m, rho1);

    //const volScalarField& test = mesh_.lookupObject<volScalarField>("alpha1");
    return fRhoNorm; 
}

...

void modKOmegaSST::correct()
{
...

// function is called and used for calculation of nut_ below

    //- Get normalized density function for Reboud's correction
    const volScalarField& fRhoNorm_(normDensFunc());

    // Re-calculate viscosity; standard calculation is the same, just w/o fRhoNorm
    nut_ = a1_*k_*fRhoNorm_/max(a1_*omega_, b1_*F23()*sqrt(S2));
    nut_.correctBoundaryConditions();

}
So if anybody has any idea what I could try out I would be very happy and thankful.


Thanks in advance.

Best regards,
Andreas
A_Pete is offline   Reply With Quote

Old   September 24, 2013, 10:12
Default
  #2
Member
 
Join Date: Jul 2011
Posts: 33
Rep Power: 5
A_Pete is on a distinguished road
Any hints?
A_Pete is offline   Reply With Quote

Old   September 28, 2013, 11:49
Default
  #3
New Member
 
Join Date: Apr 2013
Posts: 24
Rep Power: 4
GPesch is on a distinguished road
I dont actually have a suggestion of how to avoid the FPE but I just want to give you a hint on how to find at which part of the code there is something going wrong:

a. Either get a debug version of OpenFOAM in which you can directly find out which line is causing the error.

b. Use Info<<"..."<<endl; to write some lines into the terminal at relevant parts of the code. From that you can see up to which point your code is still running. From that you can narrow done the line of code which produces the error..
That was helping me always a lot (since I didnt had a debug version of OF)
GPesch is offline   Reply With Quote

Reply

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
icoFoam - Floating point exception (core dumped) File OpenFOAM Running, Solving & CFD 22 August 4, 2013 12:04
a reconstructPar issue immortality OpenFOAM Post-Processing 8 June 16, 2013 11:25
rhoPimpleFoam Floating point exception after mapFields yossi OpenFOAM Running, Solving & CFD 0 June 10, 2013 09:03
an odd(at least for me!) reconstructPar error on a field immortality OpenFOAM Running, Solving & CFD 3 June 3, 2013 22:36
CFX4.3 -build analysis form Chie Min CFX 5 July 12, 2001 23:19


All times are GMT -4. The time now is 18:23.