
:html_theme.sidebar_secondary.remove:

.. py:currentmodule:: cantera


.. DO NOT EDIT.
.. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY.
.. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE:
.. "examples/cxx/demo.cpp"
.. LINE NUMBERS ARE GIVEN BELOW.

.. only:: html

    .. note::
        :class: sphx-glr-download-link-note

        :ref:`Go to the end <sphx_glr_download_examples_cxx_demo.cpp>`
        to download the full example code.

.. rst-class:: sphx-glr-example-title

.. _sphx_glr_examples_cxx_demo.cpp:

Property calculation demo
=========================

This demonstration program builds an object representing a reacting gas
mixture, and uses it to compute thermodynamic properties, chemical
equilibrium, and transport properties.

.. tags:: C++, tutorial, thermodynamics, equilibrium, transport

.. GENERATED FROM PYTHON SOURCE LINES 10-127

.. code-block:: C++

    // This file is part of Cantera. See License.txt in the top-level directory or
    // at https://cantera.org/license.txt for license and copyright information.

    // Include cantera header files. They should be included in the form
    // "cantera/*.h". These headers are designed for use in C++ programs and
    // provide a simplified interface to the Cantera header files. If you need
    // to include core headers directly, use the format "cantera/module/*.h".

    #include "cantera/core.h"
    #include "cantera/kinetics/Reaction.h"
    #include "cantera/base/global.h" // provides Cantera::writelog
    #include <iostream>

    // All Cantera kernel names are in namespace Cantera. You can either
    // reference everything as Cantera::<name>, or include the following
    // 'using namespace' line.
    using namespace Cantera;

    // The program is put into a function so that error handling code can
    // be conveniently put around the whole thing. See main() below.

    void demoprog()
    {
        writelog("\n**** C++ Test Program ****\n");

        auto sol = newSolution("h2o2.yaml", "ohmech");
        auto gas = sol->thermo();
        double temp = 1200.0;
        double pres = OneAtm;
        gas->setState_TPX(temp, pres, "H2:1, O2:1, AR:2");


        // Thermodynamic properties

        writelog("\n\nInitial state:\n\n");
        writelog(
            "Temperature:    {:14.5g} K\n"
            "Pressure:       {:14.5g} Pa\n"
            "Density:        {:14.5g} kg/m3\n"
            "Molar Enthalpy: {:14.5g} J/kmol\n"
            "Molar Entropy:  {:14.5g} J/kmol-K\n"
            "Molar cp:       {:14.5g} J/kmol-K\n",
            gas->temperature(), gas->pressure(), gas->density(),
            gas->enthalpy_mole(), gas->entropy_mole(), gas->cp_mole());

        // set the gas to the equilibrium state with the same specific
        // enthalpy and pressure
        gas->equilibrate("HP");

        writelog("\n\nEquilibrium state:\n\n");
        writelog(
            "Temperature:    {:14.5g} K\n"
            "Pressure:       {:14.5g} Pa\n"
            "Density:        {:14.5g} kg/m3\n"
            "Molar Enthalpy: {:14.5g} J/kmol\n"
            "Molar Entropy:  {:14.5g} J/kmol-K\n"
            "Molar cp:       {:14.5g} J/kmol-K\n",
            gas->temperature(), gas->pressure(), gas->density(),
            gas->enthalpy_mole(), gas->entropy_mole(), gas->cp_mole());


        // Reaction information

        auto kin = sol->kinetics();
        int irxns = kin->nReactions();
        vector<double> qf(irxns);
        vector<double> qr(irxns);
        vector<double> q(irxns);

        // since the gas has been set to an equilibrium state, the forward
        // and reverse rates of progress should be equal for all
        // reversible reactions, and the net rates should be zero.
        kin->getFwdRatesOfProgress(&qf[0]);
        kin->getRevRatesOfProgress(&qr[0]);
        kin->getNetRatesOfProgress(&q[0]);

        writelog("\n\n");
        for (int i = 0; i < irxns; i++) {
            const auto& rxn = kin->reaction(i);
            writelog("{:30s} {:14.5g} {:14.5g} {:14.5g}  kmol/m3/s\n",
                   rxn->equation(), qf[i], qr[i], q[i]);
        }


        // transport properties

        // create a transport manager for the gas that computes
        // mixture-averaged properties
        sol->setTransportModel("mixture-averaged");
        auto tr = sol->transport();

        // print the viscosity, thermal conductivity, and diffusion
        // coefficients
        writelog("\n\nViscosity:            {:14.5g} Pa-s\n", tr->viscosity());
        writelog("Thermal conductivity: {:14.5g} W/m/K\n", tr->thermalConductivity());

        int nsp = gas->nSpecies();
        vector<double> diff(nsp);
        tr->getMixDiffCoeffs(&diff[0]);
        int k;
        writelog("\n\n{:20s}  {:21s}\n", "Species", "Diffusion Coefficient");
        for (k = 0; k < nsp; k++) {
            writelog("{:20s}  {:14.5g} m2/s\n", gas->speciesName(k), diff[k]);
        }
    }



    int main()
    {
        try {
            demoprog();
            return 0;
        } catch (CanteraError& err) {
            std::cout << err.what() << std::endl;
            return -1;
        }
    }

.. _sphx_glr_download_examples_cxx_demo.cpp:

.. only:: html

  .. container:: sphx-glr-footer sphx-glr-footer-example

    .. container:: sphx-glr-download sphx-glr-download-python

      :download:`Download C++ source code: demo.cpp <demo.cpp>`

    .. container:: sphx-glr-download sphx-glr-download-zip

      :download:`Download zipped: demo.zip <demo.zip>`


.. only:: html

 .. rst-class:: sphx-glr-signature

    `Gallery generated by Sphinx-Gallery <https://sphinx-gallery.github.io>`_
