Segfaults when creating Solution from scratch

[ischoegl]

Problem description

When two Solution objects share Reaction objects, there is a risk of segfaults.

Steps to reproduce

Run the following script:

import cantera as ct
from math import exp

gas0 = ct.Solution('gri30.yaml')
species = gas0.species()
reactions = gas0.reactions()

custom_reactions = gas0.reactions()  # problematic
# custom_reactions = [_ for _ in reactions]  # safe alternative that does not segfault
custom_reactions[2] = ct.Reaction(
    equation='H2 + O <=> H + OH',
    rate=lambda T: 38.7 * T**2.7 * exp(-3150.1542797022735/T))

gas1 = ct.Solution(thermo='ideal-gas', kinetics='gas',
                   species=species, reactions=custom_reactions)

# custom_reactions = gas0.reactions()  # problematic
# custom_reactions = [_ for _ in reactions]  # safe alternative that does not segfault
custom_reactions[2] = ct.Reaction(
    equation='H2 + O <=> H + OH',
    rate=lambda T: 38.7 * T**2.7 * exp(-3150.1542797022735/T))

gas2 = ct.Solution(thermo='ideal-gas', kinetics='gas',
                   species=species, reactions=custom_reactions)

gas2.net_production_rates  # works
# gas0.net_production_rates  # works
gas1.net_production_rates  # segfaults

Behavior

The variants

custom_reactions = [_ for _ in reactions]

are safe, whereas

custom_reactions = gas0.reactions()

cause segfaults (or bus errors resulting in crashes).

System information

• Cantera version: latest main
• OS: macOS
• Python/MATLAB/other software versions: Python 3.11.5

Other remarks

• This issue surfaced in the context of the custom_reactions.py sample.
• May be related to Eliminate traps associated with accessing ThermoPhase objects in use by Reactors enhancements#201 ... it appears that instantiation of the gas2 object corrupts gas1?

[bryanwweber]

I wonder if we're not keeping a reference to the reaction objects so they get garbage collected, where the explicit list comprehension creates that reference so things work?

[ischoegl]

I wonder if we're not keeping a reference to the reaction objects so they get garbage collected, where the explicit list comprehension creates that reference so things work?

Could be. An easy fix would be to move the list comprehension into the Solution constructor, but that may just mask a deeper underlying issue. What's surprising to me is that the C++ layer should create shared_ptrs internally, which presumably persist even if an external Python object gets destroyed.

[ischoegl]

@bryanwweber ... on second thought, I believe that it may be the Python lambda function

rate=lambda T: 38.7 * T**2.7 * exp(-3150.1542797022735/T)

that gets garbage collected. C++ will have all pointers in place, but the callback is gone.

[speth]

I think it's the CustomRate object that's getting garbage collected, rather than the lambda. Here's a simplified example that demonstrates this:

import cantera as ct
from math import exp

gas0 = ct.Solution('gri30.yaml')
species = gas0.species()
reactions = gas0.reactions()

custom_reactions = gas0.reactions()

L = lambda T: 38.7 * T**2.7 * exp(-3150.15/T)
rate1 = ct.CustomRate(L)
custom_reactions[2] = ct.Reaction(
    equation='H2 + O <=> H + OH',
    rate=rate1)

gas1 = ct.Solution(thermo='ideal-gas', kinetics='gas',
                   species=species, reactions=custom_reactions)

# Removing either of these 'del' statements avoids segfault
del custom_reactions
del rate1
gas1.net_production_rates  # segfaults

Interestingly, I also find that there's no segfault if the rate is a different type, for example ct.ArrheniusRate(1000, 1, 101), or if you use gas.add_reaction rather than inserting it into the array of reactions used to initially construct the mechanism.

[ischoegl]

@speth ... thanks for shedding more light on this issue. The Python CustomRate wrapper, together with Python Func1 hold objects that should prevent garbage collection, but that also means that once a CustomRate object is gone, there's nothing to call from C++. It appears that we need to hold references to anything that involves Python code within Python Solution objects?

I haven't investigated whether this also affects ExtensibleRate ... although this may involve manually deleting class definitions from scope/workspace?

[speth]

Oh, well, here's the answer. We do explicitly hold onto these objects in add_reaction:

cantera/interfaces/cython/cantera/kinetics.pyx

Lines 197 to 202 in 9d0e54a

	def add_reaction(self, Reaction rxn):
	""" Add a new reaction to this phase. """
	self.kinetics.addReaction(rxn._reaction)
	if isinstance(rxn.rate, (CustomRate, ExtensibleRate)):
	# prevent premature garbage collection
	self._custom_rates.append(rxn.rate)

It's just the pathway of constructing a Solution object from a list of reactions calls the C++ method directly instead of using add_reaction (in _SolutionBase._init_parts):

cantera/interfaces/cython/cantera/solutionbase.pyx

Lines 219 to 230 in 9d0e54a

	if isinstance(self, Kinetics):
	self.base.setKinetics(newKinetics(stringify(kinetics)))
	self.kinetics = self.base.kinetics().get()
	self.kinetics.addThermo(self.base.thermo())
	for phase in adjacent:
	# adjacent bulk phases for a surface phase
	self.kinetics.addThermo(phase.base.thermo())
	self.kinetics.init()
	self.kinetics.skipUndeclaredThirdBodies(True)
	for reaction in reactions:
	self.kinetics.addReaction(reaction._reaction, False)
	self.kinetics.resizeReactions()

And yes, this affects ExtensibleRate as well as CustomRate.