Computing Reaction Rates and Transport Properties in C++

Learn how to use Kinetics and Transport objects to calculate reaction rates and transport properties for a phase.

The following program demonstrates the general method for accessing the following object types from a Solution object:

  • ThermoPhase: Represents the thermodynamic properties of mixtures containing one or more species. Accessed using the thermo() method on the Solution object.

  • Kinetics: Represents a kinetic mechanism involving one or more phases. Accessed using the kinetics() method on the Solution object.

  • Transport: Computes transport properties for a ThermoPhase. Accessed using the transport() method on the Solution object.

#include "cantera/core.h"
#include "cantera/kinetics/Reaction.h"
#include <iostream>

using namespace Cantera;

// The actual code is put into a function that can be called from the main program.
void simple_demo2()
{
    // Create a new 'Solution' object that provides access to ThermoPhase, Kinetics and
    // Transport objects.
    auto sol = newSolution("gri30.yaml", "gri30");

    // Access the ThermoPhase object. Based on the phase definition in the input file,
    // this will reference an IdealGasPhase object.
    auto gas = sol->thermo();

    // Access the Kinetics object. Based on the phase definition in the input file,
    // this will reference a GasKinetics object.
    auto kin = sol->kinetics();

    // Set an "interesting" state where we will observe non-zero reaction rates.
    gas->setState_TPX(500.0, 2.0*OneAtm, "CH4:1.0, O2:1.0, N2:3.76");
    gas->equilibrate("HP");
    gas->setState_TP(gas->temperature() - 100, gas->pressure());

    // Get the net reaction rates.
    vector_fp wdot(kin->nReactions());
    kin->getNetRatesOfProgress(wdot.data());

    writelog("Net reaction rates for reactions involving CO2\n");
    size_t kCO2 = gas->speciesIndex("CO2");
    for (size_t i = 0; i < kin->nReactions(); i++) {
        if (kin->reactantStoichCoeff(kCO2, i) || kin->productStoichCoeff(kCO2, i)) {
            auto rxn = kin->reaction(i);
            writelog("{:3d}  {:30s}  {: .8e}\n", i, rxn->equation(), wdot[i]);
        }
    }
    writelog("\n");

    // Access the Transport object. Based on the phase definition in the input file,
    // this will reference a MixTransport object.
    auto trans = sol->transport();
    writelog("T        viscosity     thermal conductivity\n");
    writelog("------   -----------   --------------------\n");
    for (size_t n = 0; n < 5; n++) {
        double T = 300 + 100 * n;
        gas->setState_TP(T, gas->pressure());
        writelog("{:.1f}    {:.4e}    {:.4e}\n",
            T, trans->viscosity(), trans->thermalConductivity());
    }
}

// the main program just calls function simple_demo2 within a 'try' block, and
// catches exceptions that might be thrown
int main()
{
    try {
        simple_demo2();
    } catch (std::exception& err) {
        std::cout << err.what() << std::endl;
        return 1;
    }
    return 0;
}

This program produces the output below:

Net reaction rates for reactions involving CO2
 11  CO + O (+M) <=> CO2 (+M)         3.54150687e-08
 13  HCO + O <=> CO2 + H              1.95679990e-11
 29  CH2CO + O <=> CH2 + CO2          3.45366954e-17
 30  CO + O2 <=> CO2 + O              2.70102741e-13
 98  CO + OH <=> CO2 + H              6.46935827e-03
119  CO + HO2 <=> CO2 + OH            1.86807592e-10
131  CH + CO2 <=> CO + HCO            9.41365868e-14
151  CH2(S) + CO2 <=> CH2 + CO2       3.11161343e-12
152  CH2(S) + CO2 <=> CH2O + CO       2.85339294e-11
225  NCO + O2 <=> CO2 + NO            3.74127381e-19
228  NCO + NO <=> CO2 + N2            6.25672710e-14
261  HNCO + O <=> CO2 + NH            6.84524918e-13
267  HNCO + OH <=> CO2 + NH2          7.78871222e-10
279  CO2 + NH <=> CO + HNO           -3.30333709e-09
281  NCO + NO2 <=> CO2 + N2O          2.14286657e-20
282  CO2 + N <=> CO + NO              6.42658345e-10
289  CH2 + O2 => CO2 + 2 H            1.51032305e-18
304  CH2CHO + O => CH2 + CO2 + H      1.00331721e-19

T        viscosity     thermal conductivity
------   -----------   --------------------
300.0    1.6701e-05    4.2143e-02
400.0    2.0896e-05    5.2797e-02
500.0    2.4704e-05    6.2827e-02
600.0    2.8230e-05    7.2625e-02
700.0    3.1536e-05    8.2311e-02