hi @thc
A few comments:
1. NRR
I can’t read Chinese, but I used Google translate:
(注: 你截的图里省略了这段，但我帮你加上了，并已在 USE rxns 里声明了 NRR 和 ReactionRates，确保 ParaView 能出图)
The USE rxns, only : NRR, ReactionRates does not by itself “ensure ParaView can plot it” — you also need to set NRR
In this case, you should set it to 2 because you have 2 reactions:
(this control is in the “Model setup” pane, but perhaps it should be in "Chemistry’)
Since NRR is 0 (not set in your file), these lines:
usr_rates.f: IF (CH4_Combustion <= NRR) &
usr_rates_des.f: IF (Fe2O3_Reduction <= NRR) &
are preventing anything from happening.
2. Array indices in ReactionRates:
@(RXNS)
CH4_combustion {
chem_eq = "CH4 + 2*O2 --> CO2 + 2*H2O"
}
@(END)
@(DES_RXNS)
Fe2o3_reduction {
chem_eq = "12*Fe2O3 + CH4 --> CO2 + 8*Fe3O4 + 2*H2O"
}
@(DES_END)
x280-mfix-git:) cat species.inc
!This file was generated by /home/cgw/Work/NETL/mfix/model/rxn_preproc.py
!from /tmp/thc/4.8.mfx on Fri Apr 17 10:17:04 2026
!Edit at your own risk!!!
INTEGER, PARAMETER :: CH4 = 1
INTEGER, PARAMETER :: O2 = 2
INTEGER, PARAMETER :: CO2 = 3
INTEGER, PARAMETER :: H2O = 4
INTEGER, PARAMETER :: N2 = 5
INTEGER, PARAMETER :: Fe2O3 = 1
INTEGER, PARAMETER :: Fe3O4 = 2
INTEGER, PARAMETER :: CH4_combustion = 1
INTEGER, PARAMETER :: Fe2o3_reduction = 1
x280-mfix-git:) grep -i ch4_comb *f
usr_rates.f: RATES(CH4_Combustion) = EP_g(IJK) * 1.58489d13 * exp(-24356.5d0/T_g(IJK)) * &
usr_rates.f: RATES(CH4_Combustion) = 0.0d0
usr_rates.f: IF (CH4_Combustion <= NRR) &
usr_rates.f: ReactionRates(IJK, CH4_Combustion) = RATES(CH4_Combustion)
x280-mfix-git:) grep -i Fe2o3_reduction *f
usr_rates_des.f: DES_RATES(Fe2O3_Reduction) = Rate_kmol_s
usr_rates_des.f: DES_RATES(Fe2O3_Reduction) = 0.0d0
usr_rates_des.f: IF (Fe2O3_Reduction <= NRR) &
usr_rates_des.f: ReactionRates(IJK, Fe2O3_Reduction) = ReactionRates(IJK, Fe2O3_Reduction) + DES_RATES(Fe2O3_Reduction)
`RXNS` and `DES_RXNS` are numbered separately, so both `CH4_Combustion` and `Fe2O3_Reduction` are set to `1` in `species.inc`. This is OK when they are used as indices into `RATES` and `DES_RATES`, but when you use them as indices into `ReactionRates`, they collide. Column 2 in that array is never being set. You need to manage the indices in `ReactionRates` yourself.
In `usr_rates.f`:
```fortran
INTEGER, PARAMETER :: NRR_GAS_CH4 = 1
...
IF (NRR_GAS_CH4 <= NRR) &
ReactionRates(IJK, NRR_GAS_CH4) = RATES(CH4_Combustion)
In usr_rates_des.f:
INTEGER, PARAMETER :: NRR_DES_FE2O3 = 2
...
IF (NRR_DES_FE2O3 <= NRR) &
ReactionRates(IJK, NRR_DES_FE2O3) = ReactionRates(IJK, NRR_DES_FE2O3) + DES_RATES(Fe2O3_Reduction)
3. SI vs CGS
It looks like your usr_rates.f may be using CGS units while MFIX only supports SI:
units = 'SI'
If this is true, you need to convert your code to use SI:
In usr_rates.f:
DOUBLE PRECISION, PARAMETER :: SI_TO_CGS_CONC = 1.0d-3 ! kmol/m^3 -> mol/cm^3
DOUBLE PRECISION, PARAMETER :: CGS_TO_SI_RATE = 1.0d3 ! mol/cm^3/s -> kmol/m^3/s
c_O2_cgs = (RO_g(IJK) * X_g(IJK,O2) / MW_g(O2)) * SI_TO_CGS_CONC
c_CH4_cgs = (RO_g(IJK) * X_g(IJK,CH4) / MW_g(CH4)) * SI_TO_CGS_CONC
RATES(CH4_Combustion) = EP_g(IJK) * 1.58489d13 * exp(-24356.5d0/T_g(IJK)) * &
(c_CH4_cgs**0.7d0) * (c_O2_cgs**0.8d0) * CGS_TO_SI_RATE
While you’re in there:
I think your 1.0d-10 cutoff was written assuming CGS (mol/cm³). In SI (kmol/m³), 1e-10 kmol/m³ = 1e-13 mol/cm³, which is effectively never trips. In practice that’s fine — the rate goes to zero smoothly anyway — but if you fix units and keep the threshold, it’ll now activate at a very different physical concentration than you intended. Consider 1.0d-13 (SI) or just drop the guards once the expression is correct.
I think usr_rates_des.f is OK because mass fractions are dimensionless, and PMASS/MW gives kmol in SI or mol in CGS, so the expressions are the same.
Let us know if this gives you better results!
– Charles
