twoPhaseMixture类修改，编译通过，运行出错

Rachel0096

各位老师，下午好
我目前尝试在compressibleInterFoam添加组分方程，将twoPhaseMixture类做了以下修改：

thermo1_ = rhoThermo::New(U.mesh(), phase1Name());
thermo2_ = rhoReactionThermo::New(U.mesh(), phase2Name());

只修改thermo2的原因是假设液相为单组分系统，气相为多组分系统。
目前编译正常通过，在运行算例时会遇到下面的问题：

#0  Foam::error::printStack(Foam::Ostream&) at ??:?
#1  Foam::sigFpe::sigHandler(int) at ??:?
#2  ? in /lib/x86_64-linux-gnu/libpthread.so.0
#3  Foam::multiComponentMixture<Foam::sutherlandTransport<Foam::species::thermo<Foam::janafThermo<Foam::perfectGas<Foam::specie> >, Foam::sensibleEnthalpy> > >::patchFaceVolMixture(double, double, int, int) const at ??:?
#4  Foam::heRhoThermo<Foam::rhoReactionThermo, Foam::SpecieMixture<Foam::singleStepReactingMixture<Foam::sutherlandTransport<Foam::species::thermo<Foam::janafThermo<Foam::perfectGas<Foam::specie> >, Foam::sensibleEnthalpy> > > > >::calculate(Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh> const&, Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh>&, Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh>&, Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh>&, Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh>&, Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh>&, Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh>&, bool) at ??:?
#5  Foam::rhoReactionThermo::addfvMeshConstructorToTable<Foam::heRhoThermo<Foam::rhoReactionThermo, Foam::SpecieMixture<Foam::singleStepReactingMixture<Foam::sutherlandTransport<Foam::species::thermo<Foam::janafThermo<Foam::perfectGas<Foam::specie> >, Foam::sensibleEnthalpy> > > > > >::New(Foam::fvMesh const&, Foam::word const&) at ??:?
#6  Foam::autoPtr<Foam::rhoReactionThermo> Foam::basicThermo::New<Foam::rhoReactionThermo>(Foam::fvMesh const&, Foam::word const&) at ??:?
#7  Foam::rhoReactionThermo::New(Foam::fvMesh const&, Foam::word const&) at ??:?
#8  Foam::twoPhaseMixtureThermo::twoPhaseMixtureThermo(Foam::GeometricField<Foam::Vector<double>, Foam::fvPatchField, Foam::volMesh> const&, Foam::GeometricField<double, Foam::fvsPatchField, Foam::surfaceMesh> const&) at ??:?
#9  ? in ~/OpenFOAM/dyfluid-v2012/platforms/linux64GccDPInt32Opt/bin/atomizationFoam
#10  __libc_start_main in /lib/x86_64-linux-gnu/libc.so.6
#11  ? in ~/OpenFOAM/dyfluid-v2012/platforms/linux64GccDPInt32Opt/bin/atomizationFoam
Floating point exception (core dumped)

我通过逐行插入输出语句发现运行卡顿在下面这一句：

    Info<<"     twoPhaseMixtureThermo::debug2222222"<<endl;
    thermo1_ = rhoThermo::New(U.mesh(), phase1Name());
    thermo2_ = rhoReactionThermo::New(U.mesh(), phase2Name());   
    Info<<"     twoPhaseMixtureThermo::debug3333333"<<endl;

输出如下：

debug0000000
Constructing twoPhaseMixtureThermo
Selecting surfaceTensionModel constant
     twoPhaseMixtureThermo::debug00000000
         phase1Name:water
     twoPhaseMixtureThermo::debug111111111
         phase2Name:air
     twoPhaseMixtureThermo::debug2222222
use rhoThermo::New->basicThermo::New
      into basicThermo::New
      into lookupThermo
        actually into lookupThermo
Selecting thermodynamics package 
{
    type            heRhoThermo;
    mixture         pureMixture;
    transport       const;
    thermo          hConst;
    equationOfState rhoConst;
    specie          specie;
    energy          sensibleEnthalpy;
}
    debug::into rhoReactionThermo1111111111
      into basicThermo::New
      into lookupThermo
        actually into lookupThermo
Selecting thermodynamics package 
{
    type            heRhoThermo;
    mixture         singleStepReactingMixture;
    transport       sutherland;
    thermo          janaf;
    energy          sensibleEnthalpy;
    equationOfState perfectGas;
    specie          specie;
}
Selecting chemistryReader foamChemistryReader
chemistryReader::debug111111111
    foamChemistryReader::debug22222222222
    foamChemistryReader::debug22222222222
    foamChemistryReader::debug22222222222
    foamChemistryReader::debug22222222222
    foamChemistryReader::debug22222222222
Fuel heat of combustion :46447948
stoichiometric air-fuel ratio :15.602043
stoichiometric oxygen-fuel ratio :3.635401
Maximum products mass concentrations:
    H2O: 0.098613587
    CO2: 0.18067909
    N2: 0.72070733

而后再对multiComponentMixture进行逐行输出，发现卡顿在下面这段循环语句：
想问一下，有没有老师遇到过类似的情况，或者出现这个错误可能出现的原因？
此外，液体入口处是通过snapyyHexMesh生成的，如图所示，这会不会也对case的运行有影响呢？

wangfei9088

@Rachel0096 在 twoPhaseMixture类修改，编译通过，运行出错中说：

而后再对multiComponentMixture进行逐行输出，发现卡顿在下面这段循环语句：

好像碰到过这个问题。临时的解决的办法是初始化密度的倒数为极小数，不要初始化为0。
scalar rhoInv = VSMALL;
先试试这个，具体的问题可能还要想其他办法。也有可能你的问题不是这个原因。

Rachel0096

@wangfei9088 您是指在createField.H中初始化为极小数嘛？

Rachel0096

@wangfei9088
您好，我尝试将其修改为scalar rhoInv = VSMALL，但还其还是卡在了下面的for循环当中，不知道您还有什么建议嘛？

输出：

celli:516314.............
into cellVolMixture
   into forAll(speciesData_, i)
    rhoInv:0.71703454
    rhoInv:0.71703454
    rhoInv:0.71703454
    rhoInv:0.71703454
    rhoInv:3.4130544
      into mixtureVol_
      into for(label n=1; n<Y_.size(); n++)
celli:516315.............
into cellVolMixture
   into forAll(speciesData_, i)
    rhoInv:0.71703454
    rhoInv:0.71703454
    rhoInv:0.71703454
    rhoInv:0.71703454
    rhoInv:3.4130544
      into mixtureVol_
      into for(label n=1; n<Y_.size(); n++)

wangfei9088

@Rachel0096 在 twoPhaseMixture类修改，编译通过，运行出错中说：

还是卡在了下面的for循环当中

没瞧出有啥问题。看输出的信息就是遍历网格，计算每个网格内的密度倒数。
如果猜测是卡在for循环里了，这个cellVolMixture函数在heRhoThermo.C文件里也有调用，那里面也有一个for循环，瞅瞅是不是那里边有问题。

CFD中文网

twoPhaseMixture类修改，编译通过，运行出错