4.4.12.4. Reaction Blocks

Table of Contents

Stoichiometry
Reaction Rates
Heat Release
Reaction Phases
Reversible Reactions

Both CHEMEQ and General Chemistry use a common syntax for specifying reactions. Each reaction block defines the reaction stoichiometry, the different rate-controlling functions, and optional parameters like the heat of reaction and temperature or pressure phases (not valid for CHEMEQ).

Begin Reaction [REACTION NAME]
  Reaction is A + 1.2B -> 3C
  Multiplier Function = [Optional]
  Rate Function = [Required]
  Pressure Function = [Optional]
  Concentration Function = [Optional]
  Heat of Reaction = [Value]
  Temperature Phase = [GAS_PHASE, SOLID_PHASE, NO_MATERIAL_PHASE (default)]
  Pressure Phase = [GAS_PHASE, SOLID_PHASE, NO_MATERIAL_PHASE (default)]
  Origin Phase Heat Fraction = [Value]
End

The syntax rules for specifying reactions using these blocks is shown in the following sections.

4.4.12.4.1. Stoichiometry

The reaction stoichiometry is defined by the Reaction is line. The left hand side and right hand side must be separated by either “->” or “=>” and any numeric constants before a species name is interpreted as a stoichiometric coefficient. When using mass-based units (density or mass fraction) it is required that the stoichiometric coefficient of reactants and products be equal to conserve mass. For molar units (mole fractions or molar concentration) there is no such requirement.

Reaction is A + 2B -> C    # not valid for mass-based concentrations
Reaction is A + 2B -> 3C   # valid for mass-based concentrations

4.4.12.4.2. Reaction Rates

The overall rate of each reaction is given by

$\omega_r = f_m(t,T,P,area,Y) f_r(t,T,P,area,Y) f_p(P,T) f_c(Y)$

where $f_m$ is the multiplier function, $f_r$ is the rate function, $f_p$ is the pressure function, and $f_c$ is the concentration function. The available options for each of these functions are outlined in the following section. They can be selected in the input file by adding as many reaction blocks in the parent chemistry description block.

Multiplier Functions

Unity

The default multiplier function is unity,

f_m(t,T,P,area,Y) = 1

and can be used by either omitting the Multiplier Function line or using

Multiplier Function = Unity

Tanh Time

The tanh_time multiplier function uses a hyperbolic tangent function of time to increase the reaction rate at a specified time.

$f_m(t,T,P,area,Y) = \frac{1}{2} \left[ 1 + \tanh \left( x_{span} (t-x_{ref}) \right) \right]\left(f_{end} - f_{start}\right) + f_{start}$

The two inputs for this are $x_{ref}$ and $x_{span}$ . If you do not specify an input for the span, it uses $x_{span} = 1/x_{ref}$ , unless $x_{ref}$ was set as 0, in which case an error is thrown. The default start and end points ( $f_{start}$ and $f_{end}$ ) of the tanh blending function are 0 and 1, respectively. These can be overridden by user-specified start and end values if desired. Example inputs are shown below.

Multiplier Function = Tanh_Time x_ref = 54
Multiplier Function = Tanh_Time x_ref = 2.5 span = 0.1
Multiplier Function = Tanh_Time x_ref = 2.5 span = 0.1 start = 1 end = 100

Tanh Temperature

The tanh_temperature multiplier function uses a hyperbolic tangent function of temperature to increase the reaction rate at a specified temperature.

$f_m(t,T,P,area,Y) = \frac{1}{2} \left[ 1 + \tanh \left( x_{span} (T-x_{ref}) \right) \right]\left(f_{end} - f_{start}\right) + f_{start}$

The arguments for this function as the same format as for the tanh_time function.

Multiplier Function = Tanh_Temperature x_ref = 450
Multiplier Function = Tanh_Temperature x_ref = 600 span = 10
Multiplier Function = Tanh_Temperature x_ref = 600 span = 10 start = 100 end = 5

Thermodynamic Equilibrium

Note

Thermodynamic_Equilibrium is not supported with CHEMEQ

The thermodynamic_equilibrium multiplier function can be used to calculate the inverse equilibrium constant $\frac{1}{K_c}$ using thermodynamic properties. When multiplied by the forward rate function, this will provide a reverse rate function. The equation for $K_c$ satisfies the Van’t Hoff equation and uses Gibbs free energy.

$f_m(t,T,P,area,Y) = \exp{\left(-\sum_{i=1}^{n_s} \frac{\nu_i}{RT} (T S^{o}_{i} - H^{o}_{i})\right)} \left(\frac{RT}{P_{atm}}\right)^{\sum_{i=1}^{n_s} \nu_i}$

In this equation, $P_{atm}$ is the atmospheric pressure, $R$ is the universal gas constant, $n_s$ is the number of species, $\nu_i$ is the stoichiometric coefficient of the $i^{th}$ species, and $S^{o}_{i}$ and $H^{o}_{i}$ are the standard state molar entropy and molar enthalpy of the $i^{th}$ species respectively.

The two optional inputs for this are atmospheric pressure, $P_{atm}$ , and the universal gas constant, $R$ . The default values for these are 101325 Pa and 8.314 J/mol-K respectively. Heat capacity evaluators must be provided with the species properties in order to use this function (the enthalpy and entropy evaluators will be derived from the provided heat capacity evaluators). Example inputs are shown below.

Multiplier Function = Thermodynamic_Equilibrium
Multiplier Function = Thermodynamic_Equilibrium  Patm = 101325 R = 8.314

Simplified Equilibrium

Note

Simplified_Equilibrium is not supported with CHEMEQ

The simplified_equilibrium multiplier function can be used to calculate the inverse equilibrium constant $\frac{1}{K_c}$ using a known value of the equilibrium constant (in terms of partial pressures) at a specific temperature and the heat of reaction. When multiplied by the forward rate function, this will provide a reverse rate function. The equation for $K_c$ is the integrated Van’t Hoff equation.

$f_m(t,T,P,area,Y) = \left(\frac{1}{K_{eq}}\right) \exp{\left(\frac{\Delta H^o}{R}\left(\frac{1}{T} - \frac{1}{T_{eq}}\right)\right)} \left(\frac{RT}{P_{atm}}\right)^{\sum\limits_{i=1}^{n_s} \nu_i}$

In this equation, $P_{atm}$ is the atmospheric pressure, $R$ is the universal gas constant, $n_s$ is the number of species, $\nu_i$ is the stoichiometric coefficient of the $i^{th}$ species, $K_{eq}$ is the known equilibrium constant (in terms of partial pressures) at temperature $T_{eq}$ , and $\Delta H^o$ is the heat of reaction.

The required input for this is $T_{eq}$ . The optional inputs include $K_{eq}$ , $P_{atm}$ , and $R$ defaulted to 1.0, 101325 Pa, and 8.314 J/mol-K respectively. Either the heat of reaction or the enthalpies of formation must be provided with the reactions in the input file to use this equation accurately. Example inputs are shown below.

Multiplier Function = Simplified_Equilibrium  Teq=800
Multiplier Function = Simplified_Equilibrium  Teq=800 Keq=0.6
Multiplier Function = Simplified_Equilibrium  Teq=800 Keq=0.6 Patm=101325 R=8.314

General

The general multiplier function uses a function string supplied in the input file to evaluate the rate generally. Allowable variables in this string include time (time), T (temperature), P (pressure), area, and any valid species name from the current model. The value used for a species is in whatever unit system the general chemistry model is using. These functions can use all the functions available in string function expressions, but cannot access arbitrary expression values.

An example input for this option is

Multiplier Function = General  Function = 0.5*(1+tanh((T-350)^2))*max(T,300)

Rate Functions

Unity

The rate function is the only required function, so it has no default like the others. However, you can still select Unity for it if desired, using

Rate Function = Unity

f_r(t,T,P,area,Y) = 1

Arrhenius

The Arrhenius rate function is the most common, and uses

$f_r(t,T,P,area,Y) = A T^{\beta} \exp \left( \frac{-E_a}{R T} \right)$

The pre-exponential constant, $A$ , is the only required input. Optional inputs include $E_a$ , $\beta$ , and $R$ . If omitted, the defaults for these are 0 for $E_a$ and $\beta$ and 8.314 for $R$ . Some example lines to select this rate function is

Rate Function = Arrhenius  A = 100 Ea = 20000
Rate Function = Arrhenius  A = 100 Ea = 350 R = 1
Rate Function = Arrhenius  A = 100 beta = 2 Ea = 350 R = 1

Distributed Arrhenius

The Distributed Arrhenius rate function uses a cumulative normal distribution to vary the activation energy with respect to a selected progress species.

$f_r(t,T,P,area,Y) = A T^{\beta} \exp \left( \frac{-(E_a + \sigma \zeta)}{R T} \right)$

$\zeta = CDF^{-1}(1 - Y_j/Y_0)$

where $CDF^{-1}$ is the clipped inverse cumulative normal distribution function, $Y_j$ is the concentration of the selected progress species, $Y_0$ is the specified reference concentration, and $\sigma$ is the activation energy variation. In addition to these inputs, the same input rules as for the Arrhenius reaction apply (only $A$ is required of those again). The name specified for the progress species must be a valid species name.

Some example inputs for this rate type are (note that these are shown on multiple lines to fit in the manual, but must be on the same line in your input file):

Rate Function = Distributed_Arrhenius           \$
    A = 100 Ea = 20000                          \$
    sigma = 1000 Y0 = 1.2 ProgressSpecies = N2

Rate Function = Distributed_Arrhenius          \$
    A = 100 beta = 2 Ea = 350 R = 1            \$
    sigma = 100 Y0 = 1.2 ProgressSpecies = N2

General

The general multiplier function uses a function string supplied in the input file to evaluate the rate generally. Allowable variables in this string include time (time), T (temperature), P (pressure), area, and any valid species name from the current model. The value used for a species is in whatever unit system the general chemistry model is using. These functions can use all the functions available in string function expressions, but cannot access arbitrary expression values.

An example input for this option is:

Rate Function = General Function = 2.1*exp(-350/T)*N2^3.2*(1-CO2)*max(CO2,0.1)
Rate Function = General Function = (T>300 ? 5.0 : 0.01)*CO2*exp(-350/T)

Pressure Functions

Unity

The default pressure function is unity,

f_p(P,T) = 1

and can be used by either omitting the Pressure Function line or using

Pressure Function = Unity

Exponential

The exponential pressure function is evaluated as

$f_p(P,T) = \left(\frac{P}{P_{ref}}\right)^n,$

where $P_{ref}$ must be positive. An example input for this option is

Pressure Function = Exponential  P_ref = 1 n = 1.5

Concentration Functions

Unity

The default concentration function is unity,

f_c(Y) = 1

and can be used by either omitting the Concentration Function line or using

Concentration Function = Unity

Standard

The standard concentration function is evaluated as

$f_c(Y) = \prod_{i=1}^{N_s} Y_i^{\mu_{i}}$

The values of $\mu$ can be specified manually using a list of entries for $\mu$ . This is done by entering a list of name-value pairs for each entry with a non-zero $\mu$ value. Alternately, you can provide a list of $\mu$ values for all species.

If your reaction rate only depends on the concentration of the reactants, and the value of $\mu$ is the same as the stoichiometric coefficient, then you can have the values set automatically from the reaction stoichiometry by using the keyword Automatic.

For example, if your reaction is 2A + B -> C the automatic option would do the same thing as specifying mu = A 2 B 1. If you have repeated species as both reactants and products, only the reaction portion of the coefficient is considered. For example, A + B -> 2A + C would result in automatic exponents of mu = A 1 B 1.

Concentration Function = Standard  mu = A 2 B 1
Concentration Function = Standard  mu = 2 1 0
Concentration Function = Standard  mu = Automatic

Prout Tompkins

The Prout Tompkins concentration function can only be used with a reaction that has a single reactant and a single product, and is evaluated as

$f_c(Y) = \alpha^n \left(1 - \left(1 - 10^{-p} \right) \alpha \right)^m$

where $\alpha$ is the concentration of the reactant. An example input for this model is

Concentration Function = Prout_Tompkins  n = 1.2 m = 2.5 p = 9

Evaporation

Note

Evaporation is not supported with CHEMEQ

The evaporation concentration function is evaluated as

$f_c(Y) = A \sqrt{\frac{2}{\pi M R T}}\frac{1}{2-C_c}\left( C_e p_v(T) x_{liq} - C_c p x_{vap} \right)$

$p_v(T) = p_0 \exp \left( \frac{L}{R} \left( \frac{1}{T_b} - \frac{1}{T} \right)\right)$

where $C_c$ and $C_e$ are the condensation and evaporation coefficients, $L$ is the molar enthalpy of vaporization, $T_b$ is the boiling temperature at $p_0$ , $p$ is the gas phase pressure, $A$ is the surface area, $x_{liq}$ is the mole fraction of the liquid that is evaporating, and $x_{vap}$ is the gas phase mole fraction of the evaporating species. If $L$ is provided in mass units, $R$ must be converted appropriately so their units match. If you provide only $C_c$ or only $C_e$ , then the coefficients are assumed to be the same.

For reactions that do not include evaporation, species that are defined in the material but not involved in any reactions may be omitted from the chemistry definition block. When using evaporation reactions, this is not allowed since it will alter the value of the mole fractions in the equations above. You must also use a species variable name of species when using evaporation since this requires molar concentrations to get mole fractions.

The only rate, pressure, and multiplier functions allowed with an Evaporation concentration function are the Unity functions.

Concentration Function = Evaporation  R=8.314  M=18e-3  Cc=0.01  Tb=373  L=4e4

General Evaporation

Note

General Evaporation is not supported with CHEMEQ

This is a more general case of the evaporation function specified above. With general evaporation, different forms of the saturation pressure may be specified through a fugacity functor. The general evaporation function is evaluated as

$f_c(Y) = A \sqrt{\frac{2}{\pi M R T}}\frac{1}{2-C_c}\left( C_e f_l - C_c f_v \right)$

where $f_l$ and $f_v$ are the fugacity of the liquid and vapor respectively. Choices for the liquid fugacity include:

HENRY_CONSTANT: returns $x_{liq} H:math:$ for a specified Henry’s law constant $H$ .

ANTOINE: returns $x_{liq} 10^{\left(A - \frac{B}{C + T}\right)}$ for specified Antoine coefficients $A$ , $B$ , and $C$ .

VANT_HOFF_BOILING: returns $x_{liq} P_{atm} \exp \left( \frac{L}{R} \left( \frac{1}{T_b} - \frac{1}{T} \right)\right)$ for required parameters of molar enthalpy of vaporization $L$ and boiling temperature $T_b$ at $P_{atm}$ ; and optional parameters of atmospheric pressure $P_{atm}$ and universal gas constant $R$ with default values of 101325 Pa and 8.314 J/mol-K respectively.

Choices for the vapor fugacity include:

IDEAL_GAS: returns $x_{vap} P$ (no parameters to specify).

In the equation for $f_c(Y)$ , inputs $R$ , $M$ , $C_c$ , and $C_e$ work the same as they would for the Evaporation case above. The rules defined above for using Evaporation also apply to General_Evaporation.

For an example reaction Liq -> Gas, general evaporation may be specified as

Concentration Function = General_Evaporation  R=8.314  M=18e-3  Cc=0.01
Begin Fugacity For Liq
    Type = VANT_HOFF_BOILING
    Parameter Patm = 101325
    Parameter R = 8.314
    Parameter Tb = 373
    Parameter L = 4e4
End
Begin Fugacity For Gas
    Type = IDEAL_GAS
End

which would give the same result as the Evaporation example given above. Other options for the liquid fugacity might look like

Begin Fugacity For Liq
    Type = ANTOINE
    Parameter A = 8
    Parameter B = 1000
    Parameter C = 200
End

Begin Fugacity For Liq
    Type = HENRY_CONSTANT
    Parameter H = 1
End

Condensation

Note

Condensation is not supported with CHEMEQ

The condensation concentration function is evaluated as

$f_c(Y) = \frac{3 D M (f_v - f_l)}{r^2 R T}$

where $f_l$ and $f_v$ are the fugacity of the liquid and vapor respectively. Choices for the fugacity are given in the General Evaporation description above.

In the equation for $f_c(Y)$ , $D$ is the diffusion coefficient, $M$ is the molar mass of the condensing species, $r$ is the droplet radius, $R$ is the universal gas constant, and $T$ is temperature. The required parameter for condensation is $D$ , while the optional parameter of $R$ has the default value 8.314 J/mol-K.

Similar rules apply for using Condensation as using General_Evaporation defined above. However, you must use a mass basis for $Y$ when using condensation and provide the molecular weights of all species present in order to calculate the mole fractions in the fugacity functions.

For an example reaction Gas -> Liq, condensation may be specified as

Concentration Function = Condensation  R=8.314  D=1e-4
Begin Fugacity For Gas
    Type = IDEAL_GAS
End
Begin Fugacity For Liq
    Type = VANT_HOFF_BOILING
    Parameter Patm = 101325
    Parameter R = 8.314
    Parameter Tb = 373
    Parameter L = 4e4
End

4.4.12.4.3. Heat Release

Reaction heat release can be specified by setting per-species formation enthalpies or by supplying a heat of reaction in each reaction block, using the Heat of Reaction = X command. The units of the provided value should be energy per mass (e.g. J/kg) when using mass-based units and energy per mol when using molar units.

Note

The sign convention for heat of reaction in the chemistry blocks is that endothermic reactions have a positive heat of reaction and exothermic reactions have a negative heat of reaction.

When supplying heat of reaction values, formation enthalpies will be calculated automatically. If you supply a set of heats of reaction for which formation enthalpies cannot be calculated you will get an error. For example, in the following scenario there are no formation enthalpies that could be constructed that satisfy the given heats of reaction.

Begin Reaction R1
  Reaction is A -> B
  Heat of Reaction = 100
End

Begin Reaction R2
  Reaction is B -> A
  Heat of Reaction = 100
End

Similarly, using the same reaction stoichiometry twice with different heat release values will generate a similar error.

Begin Reaction R1
  Reaction is A -> B
  Heat of Reaction = 100
End

Begin Reaction R2
  Reaction is A -> B
  Heat of Reaction = 200
End

Sometimes one may want to use two identical reaction stoichiometries to represent two reaction pathways that occur at different temperatures. The best way to accomplish that is to create different product species.

Begin Reaction R1
  Reaction is A -> Blt
  Heat of Reaction = 100
End

Begin Reaction R2
  Reaction is A -> Bht
  Heat of Reaction = 200
End

Alternately, instead of specifying heats of reaction per reaction you can specify enthalpies of formation per species. You should not provide both heats of reaction and heats of formation however. Enthalpies of formation are supplied per species outside the reaction blocks, using a single line for CHEMEQ

species names are A B C D
Enthalpies of Formation = 0 3.5564e8 -7.1126e8 -99.5786e8

and multiple lines for General Chemistry

species names = A B C D
Enthalpy of Formation for A = 0
Enthalpy of Formation for B = 3.5564e8
Enthalpy of Formation for C = -7.1126e8
Enthalpy of Formation for D = -99.5786e8

4.4.12.4.4. Reaction Phases

Note

Reaction phase arguments are not supported with CHEMEQ

If the pressure phase argument is omitted but temperature phase is provided, the pressure phase defaults to the same phase as the temperature. By default, all of the reaction heating is applied in the temperature phase (e.g. solid phase). However, if you specify a value for the “Origin Phase Heat Fraction” between 0 and 1 you can divide that heat between phases. The fraction specified controls how much energy goes into the origin phase (a.k.a. temperature phase). The remaining energy is split evenly between the other phases that participate in the chemistry model.

4.4.12.4.5. Reversible Reactions

Note

Reversible reactions are not supported with CHEMEQ

If a reversible reaction is desired, the forward and reverse reaction must be specified in two separate reaction blocks. However, there is an option to say Reaction is Reverse of [REACTION NAME] giving the name of the forward reaction block. This will copy over the reaction string with the opposite signs for the stoichiometric coefficients, the rate function, the temperature and pressure phases, and the equal and opposite heat of reaction. If the forward reaction block contains a concentration function that is not Unity, this statement will also create a Standard concentration function for the reverse reaction with automatic coefficients. These default rate and concentration functions can be overwritten by specifying a new function inside the reverse reaction block. An example of how this option might be used for an equilibrium model is shown below.

Begin Reaction [FORWARD REACTION NAME]
    Reaction is A -> B
    Rate Function = [Required]
    Concentration Function = [Optional]
    Heat of Reaction = [Value]
End
Begin Reaction [REVERSE REACTION NAME]
    Reaction is Reverse of [FORWARD REACTION NAME]
    Multiplier Function = [Thermodynamic_Equilibrium, Simplified_Equilibrium] [args]
End