Fail to converge for zeroD simulation

[BangShiuh]

Problem description
This problem usually occurs for a reaction mechanism with reverse reactions. The reverse reaction rates can become very large for some conditions which might be caused by using the mechanism outside of the temperature range.

Steps to reproduce
error from running:

import cantera as ct
import numpy as np

ct.suppress_thermo_warnings()
gas = ct.Solution("https://raw.githubusercontent.com/BangShiuh/reaction-mechanisms/main/input/A2NOx.yaml", transport_model=None)

solutions = ct.SolutionArray(gas, extra=["t"])
air_composition = 'N2:3.76, O2:1.0'
fuel_composition = "POSF10325:1.0"

phi = 0.1
gas.set_equivalence_ratio(phi, fuel_composition, air_composition)
gas.TP = 300, ct.one_atm
inlet = ct.Reservoir(gas)

# initiate ideal gas reactor
gas.equilibrate('HP', estimate_equil=-1)
reactor = ct.IdealGasReactor(gas)
reactor.volume = 1.0

def mdot(t):
    return reactor.volume * reactor.density / 1e-3

inlet_mfc = ct.MassFlowController(inlet, reactor, mdot=mdot)
exhaust = ct.Reservoir(gas)
ct.PressureController(reactor, exhaust, master=inlet_mfc, K=1e-5)
sim = ct.ReactorNet([reactor])

rt = 0.001
solutions.append(reactor.thermo.state, t=sim.time)

while sim.time < rt * 10:
    # try:
        sim.step()
        solutions.append(reactor.thermo.state, t=sim.time)
    # except:
    #     break

# solutions.write_hdf("result-issue.h5", group='phi')

To see the large reverse reaction constant:

import cantera as ct
import matplotlib.pyplot as plt

ct.suppress_thermo_warnings()
gas = ct.Solution("/home/bschen/reaction-mechanisms/input/A2NOx.yaml")
gas.TP = 300, ct.one_atm

plt.plot(gas.reverse_rate_constants, ".")

print(gas.reaction(1153))
print(gas.reverse_rate_constants[1153])

After making this reaction irreversible, the simulation converged.

Behavior

Traceback (most recent call last):
File "phi-issue.py", line 34, in
sim.step()
File "build/python/cantera/reactor.pyx", line 1158, in cantera.reactor.ReactorNet.step
cantera._utils.CanteraError:

---

CanteraError thrown by CVodesIntegrator::step:
CVodes error encountered. Error code: -3

At t = 0.0012841 and h = 2.1008e-11, the error test failed repeatedly or with |h| = hmin.
Components with largest weighted error estimates:
196: -215.00762894680796
197: 210.9919998140404
21: 5.7857926262174235
2: -0.008801723314520293
23: 0.0030090242342010846
4: 1.0254549364550695e-06
6: 7.229548501420058e-07
14: 5.002437405334156e-07
12: 4.414005225560132e-07
8: 3.0764574908068557e-07

---

System information

• Cantera version: 3.0.0

Attachments

Additional context

[BangShiuh]

In GasKinetics.cpp,

void GasKinetics::updateKc()
{
    thermo().getStandardChemPotentials(m_grt.data());
    fill(m_delta_gibbs0.begin(), m_delta_gibbs0.end(), 0.0);

    // compute Delta G^0 for all reversible reactions
    getRevReactionDelta(m_grt.data(), m_delta_gibbs0.data());

    double rrt = 1.0 / thermo().RT();
    for (size_t i = 0; i < m_revindex.size(); i++) {
        size_t irxn = m_revindex[i];
        m_rkcn[irxn] = std::min(
            exp(m_delta_gibbs0[irxn] * rrt - m_dn[irxn] * m_logStandConc),
            BigNumber);
    }

    for (size_t i = 0; i != m_irrev.size(); ++i) {
        m_rkcn[ m_irrev[i] ] = 0.0;
    }
}

Can we use a "smaller" BigNumber to solve this problem? Or maybe there is a more physical way to set the limitation to the equilibrium constant.

[BangShiuh]

The reaction mechanism can be downloaded at https://raw.githubusercontent.com/BangShiuh/reaction-mechanisms/main/input/A2NOx.yaml