NumbericalProblem is thrown in NcpFlash, when no water initially with a dead oil dead gas three phase system

[GitPaean]

The message is as follows,

Program threw an exception: [/home/kaib/OPM-master-test/debug/opm-material/opm/material/constraintsolvers/NcpFlash.hpp:239] NcpFlash solver failed: {c_alpha^kappa} = {4865.86 0 0}, T = 15.555

Any clues about that? I can share the example.

[GitPaean]

(gdb) backtra
#0  0x00007ffff5acbf1d in __cxa_throw () from /usr/lib/x86_64-linux-gnu/libstdc++.so.6
#1  0x0000000000700e39 in Opm::NcpFlash<double, Opm::FluidSystems::BlackOil<double> >::solve<Opm::EclMultiplexerMaterial<Opm::ThreePhaseMaterialTraits<double, 0, 1, 2>, Opm::EclHysteresisTwoPhaseLaw<Opm::EclEpsTwoPhaseLaw<Opm::PiecewiseLinearTwoPhaseMaterial<Opm::TwoPhaseMaterialTraits<double, 1, 2>, Opm::PiecewiseLinearTwoPhaseMaterialParams<Opm::TwoPhaseMaterialTraits<double, 1, 2> > >, Opm::EclEpsTwoPhaseLawParams<Opm::PiecewiseLinearTwoPhaseMaterial<Opm::TwoPhaseMaterialTraits<double, 1, 2>, Opm::PiecewiseLinearTwoPhaseMaterialParams<Opm::TwoPhaseMaterialTraits<double, 1, 2> > > > >, Opm::EclHysteresisTwoPhaseLawParams<Opm::EclEpsTwoPhaseLaw<Opm::PiecewiseLinearTwoPhaseMaterial<Opm::TwoPhaseMaterialTraits<double, 1, 2>, Opm::PiecewiseLinearTwoPhaseMaterialParams<Opm::TwoPhaseMaterialTraits<double, 1, 2> > >, Opm::EclEpsTwoPhaseLawParams<Opm::PiecewiseLinearTwoPhaseMaterial<Opm::TwoPhaseMaterialTraits<double, 1, 2>, Opm::PiecewiseLinearTwoPhaseMaterialParams<Opm::TwoPhaseMaterialTraits<double, 1, 2> > > > > > >, Opm::EclHysteresisTwoPhaseLaw<Opm::EclEpsTwoPhaseLaw<Opm::PiecewiseLinearTwoPhaseMaterial<Opm::TwoPhaseMaterialTraits<double, 0, 1>, Opm::PiecewiseLinearTwoPhaseMaterialParams<Opm::TwoPhaseMaterialTraits<double, 0, 1> > >, Opm::EclEpsTwoPhaseLawParams<Opm::PiecewiseLinearTwoPhaseMaterial<Opm::TwoPhaseMaterialTraits<double, 0, 1>, Opm::PiecewiseLinearTwoPhaseMaterialParams<Opm::TwoPhaseMaterialTraits<double, 0, 1> > > > >, Opm::EclHysteresisTwoPhaseLawParams<Opm::EclEpsTwoPhaseLaw<Opm::PiecewiseLinearTwoPhaseMaterial<Opm::TwoPhaseMaterialTraits<double, 0, 1>, Opm::PiecewiseLinearTwoPhaseMaterialParams<Opm::TwoPhaseMaterialTraits<double, 0, 1> > >, Opm::EclEpsTwoPhaseLawParams<Opm::PiecewiseLinearTwoPhaseMaterial<Opm::TwoPhaseMaterialTraits<double, 0, 1>, Opm::PiecewiseLinearTwoPhaseMaterialParams<Opm::TwoPhaseMaterialTraits<double, 0, 1> > > > > > >, Opm::EclMultiplexerMaterialParams<Opm::ThreePhaseMaterialTraits<double, 0, 1, 2>, Opm::EclHysteresisTwoPhaseLaw<Opm::EclEpsTwoPhaseLaw<Opm::PiecewiseLinearTwoPhaseMaterial<Opm::TwoPhaseMaterialTraits<double, 1, 2>, Opm::PiecewiseLinearTwoPhaseMaterialParams<Opm::TwoPhaseMaterialTraits<double, 1, 2> > >, Opm::EclEpsTwoPhaseLawParams<Opm::PiecewiseLinearTwoPhaseMaterial<Opm::TwoPhaseMaterialTraits<double, 1, 2>, Opm::PiecewiseLinearTwoPhaseMaterialParams<Opm::TwoPhaseMaterialTraits<double, 1, 2> > > > >, Opm::EclHysteresisTwoPhaseLawParams<Opm::EclEpsTwoPhaseLaw<Opm::PiecewiseLinearTwoPhaseMaterial<Opm::TwoPhaseMaterialTraits<double, 1, 2>, Opm::PiecewiseLinearTwoPhaseMaterialParams<Opm::TwoPhaseMaterialTraits<double, 1, 2> > >, Opm::EclEpsTwoPhaseLawParams<Opm::PiecewiseLinearTwoPhaseMaterial<Opm::TwoPhaseMaterialTraits<double, 1, 2>, Opm::PiecewiseLinearTwoPhaseMaterialParams<Opm::TwoPhaseMaterialTraits<double, 1, 2> > > > > > >, Opm::EclHysteresisTwoPhaseLaw<Opm::EclEpsTwoPhaseLaw<Opm::PiecewiseLinearTwoPhaseMaterial<Opm::TwoPhaseMaterialTraits<double, 0, 1>, Opm::PiecewiseLinearTwoPhaseMaterialParams<Opm::TwoPhaseMaterialTraits<double, 0, 1> > >, Opm::EclEpsTwoPhaseLawParams<Opm::PiecewiseLinearTwoPhaseMaterial<Opm::TwoPhaseMaterialTraits<double, 0, 1>, Opm::PiecewiseLinearTwoPhaseMaterialParams<Opm::TwoPhaseMaterialTraits<double, 0, 1> > > > >, Opm::EclHysteresisTwoPhaseLawParams<Opm::EclEpsTwoPhaseLaw<Opm::PiecewiseLinearTwoPhaseMaterial<Opm::TwoPhaseMaterialTraits<double, 0, 1>, Opm::PiecewiseLinearTwoPhaseMaterialParams<Opm::TwoPhaseMaterialTraits<double, 0, 1> > >, Opm::EclEpsTwoPhaseLawParams<Opm::PiecewiseLinearTwoPhaseMaterial<Opm::TwoPhaseMaterialTraits<double, 0, 1>, Opm::PiecewiseLinearTwoPhaseMaterialParams<Opm::TwoPhaseMaterialTraits<double, 0, 1> > > > > > > > >, Opm::CompositionalFluidState<double, Opm::FluidSystems::BlackOil<double>, true> > (fluidState=..., matParams=..., paramCache=..., globalMolarities=..., tolerance=2.2204460492503131e-08, tolerance@entry=-1)
    at /home/kaib/OPM-master-test/debug/opm-material/opm/material/constraintsolvers/NcpFlash.hpp:236
#2  0x0000000000701f4e in Ewoms::BlackOilPrimaryVariables<Ewoms::Properties::TTag::EclFlowProblem>::assignMassConservative<Opm::CompositionalFluidState<double, Opm::FluidSystems::BlackOil<double>, true> > (isInEquilibrium=false, matParams=..., fluidState=..., this=0x100f980)
    at /home/kaib/OPM-master-test/debug/ewoms/ewoms/models/blackoil/blackoilprimaryvariables.hh:237
#3  Ewoms::EclProblem<Ewoms::Properties::TTag::EclFlowProblem>::initial<Ewoms::FvBaseElementContext<Ewoms::Properties::TTag::EclFlowProblem> > (timeIdx=0, spaceIdx=0, context=..., values=..., this=<optimized out>) at /home/kaib/OPM-master-test/debug/ewoms/ebos/eclproblem.hh:944
#4  Ewoms::FvBaseDiscretization<Ewoms::Properties::TTag::EclFlowProblem>::applyInitialSolution (this=0xf59990) at /home/kaib/OPM-master-test/debug/ewoms/ewoms/disc/common/fvbasediscretization.hh:548
#5  0x00000000007d58c8 in Opm::FlowMainEbos<Ewoms::Properties::TTag::EclFlowProblem>::setupEbosSimulator (this=this@entry=0x7fffffffdbc0) at /home/kaib/OPM-master-test/debug/opm-simulators/opm/autodiff/FlowMainEbos.hpp:422
#6  0x00000000007d71a8 in Opm::FlowMainEbos<Ewoms::Properties::TTag::EclFlowProblem>::execute (this=this@entry=0x7fffffffdbc0, argc=<optimized out>, argv=<optimized out>) at /home/kaib/OPM-master-test/debug/opm-simulators/opm/autodiff/FlowMainEbos.hpp:106
#7  0x00000000005d21a9 in Opm::flowEbosBlackoilMain (argc=2, argv=0x7fffffffde28) at /home/kaib/OPM-master-test/debug/opm-simulators/opm/simulators/flow_ebos_blackoil.cpp:59
#8  0x00007ffff547b830 in __libc_start_main (main=0x5a9590 <main(int, char**)>, argc=2, argv=0x7fffffffde28, init=<optimized out>, fini=<optimized out>, rtld_fini=<optimized out>, stack_end=0x7fffffffde18) at ../csu/libc-start.c:291
#9  0x00000000005cd0e9 in _start ()

[GitPaean]

Initially, there is no water and gas.

[andlaus]

it seems like the flash solver cannot find a thermodynamically consistent solution for the condition that you specified. For example, this happens if all fluids that are present are incompressible.

[GitPaean]

The oil data looks originally like the following. You are probably right, when I increase the comprehensibility, it works.

200 PVDO
201     20.0 1.0100001   10.000
202     900. 1.01        20.000
203 /

It was an example I wrote when we used flow_legacy. It works with flow_legacy, should we try to make it work with the current flow, or you think it probably does not matter.

[andlaus]

It is not a big problem to make the code work with this: just always use values.assignNaive() in EclProblem::initial(). That said, because it is extremely easy to specify thermodynamically inconsistent initial conditions, I'm not sure if this is a good idea (i.e., this would basically mean that we fall back to flow_legacy's approach of "undefined behaviour by default"?)

[GitPaean]

Okay. I do not have an opinion about this and just out of curiosity. Thanks for the feedback. I am closing it now.