Cantera vs. CHEMKIN: Different Laminar Flame Speed Sensitivity Results with Extended Mechanism

[Loraine Sanson]

Hello colleagues,

I hope I can describe my question as clearly as possible to minimize any unnecessary confusion.

I’m currently working on combustion inhibition research involving hydrocarbon fuels and Halon replacements. Recently, I’ve transitioned from the CHEMKIN platform to Cantera. The reaction mechanisms I use include the USC II mechanism along with decomposition and combustion reactions of halogenated inhibitors.

However, I’ve encountered a puzzling issue: the sensitivity data for laminar flame speed computed by Cantera is inconsistent with the results from CHEMKIN — not just in magnitude, but even in the ranking order of the top contributing elementary reactions (based on absolute values).

To explore the issue, I performed the following comparison:

Platforms:

• CHEMKIN

• Cantera

Mechanisms:

• Mechanism A: USC II

• Mechanism B: USC II + halogenated inhibitor reactions
  (based on J.L. Pagliaro, N. Bouvet, G.T. Linteris, Premixed flame inhibition by CF₃Br and C₃H₂F₃Br (2-BTP), Combustion and Flame, 169 (2016) 272–286. https://doi.org/10.1016/j.combustflame.2016.04.017)

Computed Parameters:

• Laminar flame speed

• Peak mole fraction of active radicals (H, O, OH)

• Sensitivity of laminar flame speed

Observations:

• CHEMKIN+A, CHEMKIN+B, and Cantera+A produced highly consistent results for all three parameters and required relatively short simulation times.

• Cantera+B produced notably different sensitivity results. While the laminar flame speed and radical peaks (H, O, OH) were still consistent with the other three cases, the sensitivity rankings and values differed substantially — both in magnitude and absolute order.

Below are two plots showing the top 10 sensitivities calculated by Cantera for Mechanism A and Mechanism B. The discrepancies in Mechanism B are what concern me the most.

Interpretation:

• Since CHEMKIN+A and CHEMKIN+B agree well, I believe the mechanisms themselves (A and B) are valid and correctly implemented on the CHEMKIN platform.

• Since Cantera+A agrees with both CHEMKIN results, I believe my Cantera code is functioning properly.

• Since Cantera+B gives consistent flame speed and radical profiles, I believe the converted Mechanism B YAML file is likely correct, despite the discrepancy in sensitivity data.

This comparison was my attempt to address the sensitivity anomaly, as it's quite challenging to directly verify the correctness of a converted mechanism with 181 species and 1513 reactions.

I’m also aware that both CHEMKIN and Cantera output normalized sensitivity coefficients by default. However, it seems Cantera might not exclude irrelevant reactions — even under pure hydrogen conditions, it still evaluates the influence of all reactions, including those involving fluorinated species. For instance, in Cantera+B, some fluorine-containing reactions have larger sensitivity values than hydrogen-related reactions, which seems unexpected. In contrast, CHEMKIN might handle this exclusion more effectively.

This might also explain why Cantera+B simulations take much longer — even 6× longer. Perhaps unintended fluorinated reactions are being evaluated? Interestingly, the flame speed remains unaffected.

I’ve consulted some colleagues, but most of them focus on simpler fuels and do not encounter merged mechanisms like Mechanism B. I’m eager for this issue to be discussed and hopefully resolved — not just for myself, but for others who may face similar problems when transitioning to Cantera with complex, merged combustion-inhibition mechanisms.

Thanks in advance for your insights and suggestions!

—————————————————————————————————————————————————
—————————————————————————————————————————————————

Here is the core part of my Cantera calculation code:

import cantera as ct

import numpy as np
import pandas as pd

# Suppress Cantera warning messages
warnings.filterwarnings('ignore')

# ===== Basic initial parameters =====
T0 = 293.0  # Initial temperature [K]
p0 = 101325.0  # Initial pressure [Pa]
air_ratio = 3.76  # Molar ratio of N2/O2 in air
volume_L = 0.5  # Volume of the container [L]
volume_m3 = volume_L / 1000  # Convert to m³, required by Cantera
mech_file = "chem.yaml"  # Replace with the mechanism file path you are using (CTI/YAML)

# ===== Simulation settings =====
width = 0.05  # Domain width [m], must be large enough to contain the flame structure

# ===== Operating condition variables =====
# This section of the code is omitted #

# === Create gas object ===
gas = ct.Solution(mech_file)
gas.TPX = T0, p0, X

# Set solver refinement criteria
flame.set_refine_criteria(ratio=3.0, slope=0.06, curve=0.12)
flame.max_grid_points = 3000

# === Solve the flame ===
flame.solve(loglevel=1, auto=True)

# === Extract results ===
# 1. Laminar flame speed
flame_speed = flame.velocity[0]  # m/s

# 2. Flame thickness (based on maximum temperature gradient)
temp_grad = np.gradient(flame.T, flame.grid)
flame_thickness = (max(flame.T) - min(flame.T)) / max(abs(temp_grad)) * 1000  # Convert to mm

# 3. Temperature and pressure profiles
distance = flame.grid * 1000  # Convert to mm
temperature = flame.T  # K
pressure = flame.P / 1000  # Convert to kPa

# 4. Elemental reaction heat release rates ×××
heat_release_rates = flame.heat_release_rate / 1e6  # Convert to MW/m³
reaction_names = gas.reaction_equations()

# Select top 30 reactions by absolute heat release rate
n_reactions = min(30, len(heat_release_rates))
sorted_indices = np.argsort(np.abs(heat_release_rates))[::-1][:n_reactions]

# 5. Sensitivity analysis (core section) ×××
# Create a sensitivity DataFrame containing reaction index, equation, and sensitivity value
sens = pd.DataFrame(index=gas.reaction_equations(), columns=["sensitivity"])

# Get sensitivity data
sens.sensitivity = flame.get_flame_speed_reaction_sensitivities()
sens.head(10)  # Show top 10 sensitivities

# Generate reaction IDs (e.g., R1, R2, ..., Rn) and add to the DataFrame
reaction_ids = [f"R{i+1}" for i in range(gas.n_reactions)]
sens["reaction_id"] = reaction_ids
sens["reaction_eqn"] = gas.reaction_equations()

# Sort by absolute value of sensitivity in descending order
sens_sorted = sens.iloc[(-sens['sensitivity'].abs()).argsort()]
print(sens_sorted.head(10))

# Extract top 50 reactions by absolute sensitivity
top_sens = sens_sorted.head(50)[["reaction_id", "reaction_eqn", "sensitivity"]]
top_sens.to_csv("Top50_Sensitivities.csv", index=False)

# 6. ROP (Rate of Production) analysis ×××

# 7. Key species mole fractions √√√
species_mole_fractions = {}
for species in key_species:
      if species in gas.species_names:
          spec_idx = gas.species_index(species)
          mole_fractions = flame.X[spec_idx, :]  # Get mole fraction
          species_mole_fractions[species] = mole_fractions

[Loraine Sanson]

Sorry, I should have added one more detail for better clarity:
The calculations above were performed using H₂ as the fuel and an oxidizer mixture of O₂:N₂ = 1:3.76, with an equivalence ratio of 0.6.

[Ray Speth]

Hi,

Thanks for reporting this unexpected behavior. It appears that it can be traced back to the very common issue of excessively high reverse rate constants. These can be identified with a simple check:

kr = gas.reverse_rate_constants for i,R in enumerate(gas.reactions()): if kr[i] > 1e18: print(i, R.equation, kr[i])

which, at an initial temperature of 300 K, gives the following list:

327 OH + aC3H5 <=> C2H3CHO + 2 H 5.935380791512869e+20 363 C3H6 + O <=> C2H3CHO + 2 H 2.8724419912602785e+22 784 HF + M <=> F + H + M 3.0336962187667686e+33 797 CF4 + M <=> CF3 + F + M 1.2161333198330624e+19 1370 CF3CHO (+M) <=> CF2CO + HF (+M) 5.854179315298665e+21 1455 BR + CH2BR <=> BR2 + CH2 3.888297162403887e+22

I found that just disabling the one reaction HF + M <=> F + H + M was sufficient to change the behavior and give more normal-looking sensitivity coefficients. My best guess is that this rate constant is having an adverse effect on the condition number of the Jacobian matrix, and that this is then affecting the resulting sensitivities. While there are other reactions with very high reverse rate constants, they don’t seem to throw off the sensitivity calculation.

You can also verify that this reaction is the problem by adding just this one reaction to an existing, well-behaved mechanism. I’ll create an issue on our tracker documenting this, since we don’t want to return incorrect values to the user even if there is a problem with the input data.

Regards,
Ray