Updates to make HP queries return values instead of failing for higher pressures

[wandadars]

I was running a test case for a 1D counterflow diffusion flame at very high pressures using the Peng-Robinson equation of state, and the initial condition that used the HP equilibrate call return EoS failures. Codex 5.3 identified this as a possible reason for why the real/ideal gas mixed PREOS might be failing, and it did resolve the EoS that was happening for high pressure H-P equilibrate calls.

AI Statement (required)

• Extensive use of generative AI.
  Significant portions of code or documentation were generated with AI, including
  logic and implementation decisions. All generated code and documentation were
  reviewed and understood by the contributor. Examples: Output from agentic coding
  tools and/or substantial refactoring by LLMs (web-based or local).

For additional information on Cantera's AI policy, see
https://github.com/Cantera/cantera/blob/main/CONTRIBUTING.md#restrictions-on-generative-ai-usage -->

[codecov]

Codecov Report

❌ Patch coverage is 59.67742% with 50 lines in your changes missing coverage. Please review.
✅ Project coverage is 77.49%. Comparing base (05de7fe) to head (2901fd9).

Files with missing lines	Patch %	Lines
src/thermo/MixtureFugacityTP.cpp	51.35%	17 Missing and 1 partial ⚠️
src/thermo/PengRobinson.cpp	63.63%	10 Missing and 6 partials ⚠️
src/thermo/RedlichKwongMFTP.cpp	62.79%	10 Missing and 6 partials ⚠️

Additional details and impacted files

@@            Coverage Diff             @@
##             main    #2087      +/-   ##
==========================================
- Coverage   77.53%   77.49%   -0.05%     
==========================================
  Files         451      451              
  Lines       52810    52928     +118     
  Branches     8828     8859      +31     
==========================================
+ Hits        40947    41016      +69     
- Misses       8888     8924      +36     
- Partials     2975     2988      +13     

☔ View full report in Codecov by Sentry.
📢 Have feedback on the report? Share it here.

🚀 New features to boost your workflow:

• ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

[speth]

I think this may be partially on the right track, but I'd suggest focusing on a narrower definition of the problem. The equilibrium solver here is acting as a way of "finding" conditions that cause trouble for the equation of state solver, but isn't itself part of the problem.

The equations of the P-R model do not like the case where $a=b=0$. There are a lot of places where this results in 0/0, such that all the resulting thermodynamic properties are nan. For example:

import cantera as ct

def partial_lookup_peng_robinson_phase():
    """Build a mixed-parameter Peng-Robinson phase used by the original test."""
    species = ct.Species.list_from_file("test/data/co2_PR_example.yaml")
    phase_species = []
    for sp in species:
        data = dict(sp.input_data)
        if data["name"] in {"CH4", "O2"}:
            data.pop("equation-of-state", None)
        else:
            data["equation-of-state"] = {
                "model": "Peng-Robinson",
                "a": 0,
                "b": 0,
                "acentric-factor": 0.0,
            }
        phase_species.append(ct.Species.from_dict(data))

    return ct.Solution(
        thermo="Peng-Robinson",
        kinetics="bulk",
        species=phase_species,
        reactions=[],
        name="partial-lookup-pr",
    )

gas = partial_lookup_peng_robinson_phase()
gas.TPX = 500, 1e6, 'CO2:0.333, H2O:0.666'
gas()

which outputs:

  partial-lookup-pr:

       temperature   500 K
          pressure   1e+06 Pa
           density   6.4177 kg/m^3
  mean mol. weight   26.68 kg/kmol
   phase of matter   unspecified

                          1 kg             1 kmol     
                     ---------------   ---------------
          enthalpy   nan               nan              J
   internal energy   nan               nan              J
           entropy   nan               nan              J/K
    Gibbs function   nan               nan              J
 heat capacity c_p   nan               nan              J/K
 heat capacity c_v   nan               nan              J/K

                      mass frac. Y      mole frac. X     chem. pot. / RT
                     ---------------   ---------------   ---------------
               CO2           0.54984           0.33333   nan            
               H2O           0.45016           0.66667   nan            
     [   +5 minor]                 0                 0  

So it's pretty clear why the equilibrium solver can't do anything with such a mixture.

If you want to pick out specific cases where the cubic solver is causing trouble, I'd put a try/catch around the place in the equilibrium solver where it sets the mixture state and print out the full state that it's trying to set if that fails. You probably also want to do any exploration using a specific equilibrium solver (vcs or gibbs), as the automatic switching behavior may have the effect of hiding some errors that arise in the cubic equation of state solver.

[wandadars]

Good catch @speth . We have checks that watch for a=0, b=0 and use ideal-gas fallbacks for those quantities. Also those fallback were moved up to the base class since this is a common issue in both PREOS and RK cubic equations of state.