
: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/python/kinetics/mechanism_reduction.py"
.. LINE NUMBERS ARE GIVEN BELOW.

.. only:: html

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

        :ref:`Go to the end <sphx_glr_download_examples_python_kinetics_mechanism_reduction.py>`
        to download the full example code.

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

.. _sphx_glr_examples_python_kinetics_mechanism_reduction.py:


Mechanism reduction
===================

A simplistic approach to mechanism reduction which demonstrates Cantera's
features for dynamically manipulating chemical mechanisms.

Here, we use the detailed NUI Galway mechanism to simulate adiabatic, constant pressure
ignition of a lean n-hexane/air mixture. We track the maximum reaction rates for
each reaction to determine which reactions are the most important, according to
a simple metric based on the relative net reaction rate.

We then create a sequence of reduced mechanisms including only the top reactions
and the associated species, and run the simulations again with these mechanisms
to see whether the reduced mechanisms with a certain number of species are able
to adequately simulate the ignition delay problem.

Requires: cantera >= 3.2.0, matplotlib >= 2.0

.. tags:: Python, kinetics, combustion, reactor network, editing mechanisms,
          ignition delay, plotting

.. GENERATED FROM PYTHON SOURCE LINES 23-45

.. code-block:: Python


    import cantera as ct
    import numpy as np
    import matplotlib.pyplot as plt
    import dataclasses

    @dataclasses.dataclass
    class Result:
        species: int
        reactions: int
        time: list
        temperature: list


    ct.suppress_thermo_warnings()
    T0 = 975
    P0 = 5 * ct.one_atm
    gas = ct.Solution('example_data/n-hexane-NUIG-2015.yaml')
    gas.set_equivalence_ratio(0.8, 'NC6H14:1.0', 'O2:1.0, N2:3.76')
    X0 = gas.mole_fraction_dict()
    gas.TP = T0, P0








.. GENERATED FROM PYTHON SOURCE LINES 46-47

Run a simulation with the full mechanism

.. GENERATED FROM PYTHON SOURCE LINES 47-65

.. code-block:: Python

    r = ct.IdealGasConstPressureMoleReactor(gas, clone=False)
    sim = ct.ReactorNet([r])
    sim.preconditioner = ct.AdaptivePreconditioner()

    tt = []
    TT = []
    # Rmax is the maximum relative reaction rate at any timestep.
    Rmax = np.zeros(gas.n_reactions)
    while sim.time < 0.04:
        sim.step()
        tt.append(1000 * sim.time)
        TT.append(r.T)
        rnet = abs(r.phase.net_rates_of_progress)
        rnet /= max(rnet)
        Rmax = np.maximum(Rmax, rnet)

    baseline = Result(gas.n_species, gas.n_reactions, tt, TT)








.. GENERATED FROM PYTHON SOURCE LINES 66-67

Get the reaction objects, and sort them so the most active reactions are first.

.. GENERATED FROM PYTHON SOURCE LINES 67-69

.. code-block:: Python

    R = sorted(zip(Rmax, gas.reactions()), key=lambda x: -x[0])








.. GENERATED FROM PYTHON SOURCE LINES 70-72

Test reduced mechanisms with different numbers of reactions and collect results for
plotting.

.. GENERATED FROM PYTHON SOURCE LINES 72-109

.. code-block:: Python

    results = []
    for i, N in enumerate([100, 200, 300, 400, 600, 800]):
        # Get the N most active reactions
        reactions = [r[1] for r in R[:N]]

        # find the species involved in these reactions. At a minimum, include all
        # species in the reactant mixture
        species_names = {'N2', 'NC6H14', 'O2'}
        for reaction in reactions:
            species_names.update(reaction.reactants)
            species_names.update(reaction.products)

        # Get the species objects
        species = [gas.species(name) for name in species_names]

        # create the new reduced mechanism
        gas2 = ct.Solution(thermo='ideal-gas', kinetics='gas',
                           species=species, reactions=reactions)

        # save the reduced mechanism for later use
        gas2.write_yaml(f"hexane-reduced-{N}-reaction.yaml")

        # Re-run the ignition problem with the reduced mechanism
        gas2.TPX = T0, P0, X0
        r = ct.IdealGasConstPressureMoleReactor(gas2, clone=False)
        sim = ct.ReactorNet([r])
        sim.preconditioner = ct.AdaptivePreconditioner()

        tt = []
        TT = []
        while sim.time < 0.04:
            sim.step()
            tt.append(1000 * sim.time)
            TT.append(r.T)

        results.append(Result(gas2.n_species, gas2.n_reactions, tt, TT))








.. GENERATED FROM PYTHON SOURCE LINES 110-111

Plot the results.

.. GENERATED FROM PYTHON SOURCE LINES 111-126

.. code-block:: Python

    C = plt.cm.winter(np.linspace(0, 1, 6))

    fig, ax = plt.subplots()
    ax.plot(baseline.time, baseline.temperature, color='k', lw=3, zorder=100,
            label=f'K={baseline.species}, R={baseline.reactions}')

    for i, result in enumerate(results):
        ax.plot(result.time, result.temperature, lw=2, color=C[i],
                label=f'K={result.species}, R={result.reactions}')
    ax.set(xlabel='Time (ms)', ylabel='Temperature (K)')
    ax.legend(loc='upper left')
    ax.set_title('Reduced mechanism ignition delay times\n'
                 'K: number of species; R: number of reactions')
    ax.set_xlim(0, 40)
    plt.show()



.. image-sg:: /examples/python/kinetics/images/sphx_glr_mechanism_reduction_001.png
   :alt: Reduced mechanism ignition delay times K: number of species; R: number of reactions
   :srcset: /examples/python/kinetics/images/sphx_glr_mechanism_reduction_001.png, /examples/python/kinetics/images/sphx_glr_mechanism_reduction_001_2_00x.png 2.00x
   :class: sphx-glr-single-img






.. rst-class:: sphx-glr-timing

   **Total running time of the script:** (0 minutes 24.424 seconds)


.. _sphx_glr_download_examples_python_kinetics_mechanism_reduction.py:

.. only:: html

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

    .. container:: sphx-glr-download sphx-glr-download-jupyter

      :download:`Download Jupyter notebook: mechanism_reduction.ipynb <mechanism_reduction.ipynb>`

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

      :download:`Download Python source code: mechanism_reduction.py <mechanism_reduction.py>`

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

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


.. only:: html

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

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