reactingTwoPhaseEulerFoam添加空化模型，发生相变时焓方程或者迭代求温度发散

[zwl]

相方程和压力方程的源相均是类似interPhaseChangeFoam的源相添加方式添加，能量方程是按照质量传递乘以相应的焓，但是在求解过程中，一旦发生相变，焓方程或者迭代求温度就会发散。（我并没有按照本求解器原有的蒸发模型添加，因为相变原理不同）要以什么样的数值方式处理才能解决这个问题呢？有经验的专家能给些建议吗？@东岳
代码如下：

E1Eqn =
            (
                E1Eqn
             ==
               *heatTransfer[phase1.name()]
              +(fluid.mDotAlphal()[0]*limitedAlpha2+fluid.mDotAlphal()[1]*limitedAlpha1)*phase2.thermo().he()
              + alpha1*rho1*(U1&g)
              + fvOptions(alpha1, rho1, phase1.thermo().he())
            );

            E1Eqn->relax();
            fvOptions.constrain(E1Eqn.ref());
            E1Eqn->solve();

        }
    }

    {
        tmp<fvScalarMatrix> E2Eqn(phase2.heEqn());

        if (E2Eqn.valid())
        {
            E2Eqn =
            (
                E2Eqn
             ==
               *heatTransfer[phase2.name()]
              -fvm::Sp((fluid.mDotAlphal()[0]*limitedAlpha2+fluid.mDotAlphal()[1]*limitedAlpha1)*2, phase2.thermo().he())
              +(fluid.mDotAlphal()[0]*limitedAlpha2+fluid.mDotAlphal()[1]*limitedAlpha1) * phase1.thermo().he()
              + alpha2*rho2*(U2&g)
              + fvOptions(alpha2, rho2, phase2.thermo().he())
            );

            E2Eqn->relax();elax();