
: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/reactors/preconditioned_integration.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_reactors_preconditioned_integration.py>`
        to download the full example code.

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

.. _sphx_glr_examples_python_reactors_preconditioned_integration.py:


Acceleration of reactor integration using a sparse preconditioned solver
========================================================================

Ideal gas, constant-pressure, adiabatic kinetics simulation that compares preconditioned
and non-preconditioned integration of n-hexane.

Requires: cantera >= 3.2.0, matplotlib >= 2.0

.. tags:: Python, combustion, reactor network, preconditioner

.. GENERATED FROM PYTHON SOURCE LINES 12-17

.. code-block:: Python

    import cantera as ct
    import matplotlib.pyplot as plt
    plt.rcParams['figure.constrained_layout.use'] = True
    from timeit import default_timer








.. GENERATED FROM PYTHON SOURCE LINES 18-24

Simulation setup
----------------

Create a reactor network for simulating the constant pressure ignition of a
stoichiometric n-hexane/air mixture, with or without the use of the preconditioned
solver.

.. GENERATED FROM PYTHON SOURCE LINES 24-60

.. code-block:: Python

    def integrate_reactor(preconditioner=True):
        # Use a detailed n-hexane mechanism with 1268 species
        gas = ct.Solution('example_data/n-hexane-NUIG-2015.yaml')
        gas.TP = 1000, ct.one_atm
        gas.set_equivalence_ratio(1, 'NC6H14', 'N2:3.76, O2:1.0')
        reactor = ct.IdealGasConstPressureMoleReactor(gas, clone=False)
        # set volume for reactors
        reactor.volume = 0.1
        # Create reactor network
        sim = ct.ReactorNet([reactor])
        # Add preconditioner
        if preconditioner:
            sim.derivative_settings = {"skip-third-bodies":True, "skip-falloff":True}
            sim.preconditioner = ct.AdaptivePreconditioner()
        sim.initialize()
        # Advance to steady state
        integ_time = default_timer()
        # solution array for state data
        states = ct.SolutionArray(reactor.phase, extra=['time'])
        # advance to steady state manually
        while (sim.time < 0.1):
            states.append(reactor.phase.state, time=sim.time)
            sim.step()
        integ_time = default_timer() - integ_time
        # Return time to integrate
        if preconditioner:
            print(f"Preconditioned Integration Time: {integ_time:f}")
        else:
            print(f"Non-preconditioned Integration Time: {integ_time:f}")
        # Get and output solver stats
        for key, value in sim.solver_stats.items():
            print(f"{key:>24s}: {value}")
        print("\n")
        # return some variables for plotting
        return states.time, states.T, states('CO2').Y, states('NC6H14').Y








.. GENERATED FROM PYTHON SOURCE LINES 61-63

Integrate with sparse, preconditioned solver
--------------------------------------------

.. GENERATED FROM PYTHON SOURCE LINES 63-65

.. code-block:: Python

    timep, Tp, CO2p, NC6H14p = integrate_reactor(preconditioner=True)





.. rst-class:: sphx-glr-script-out

 .. code-block:: none

    Preconditioned Integration Time: 6.648215
            step_solve_fails: 10
                       steps: 991
                   rhs_evals: 1199
             nonlinear_iters: 1196
        nonlinear_conv_fails: 27
              err_test_fails: 30
                  last_order: 2
       stab_order_reductions: 0
                   jac_evals: 0
            lin_solve_setups: 138
               lin_rhs_evals: 3860
                   lin_iters: 3860
              lin_conv_fails: 264
                  prec_evals: 37
                 prec_solves: 4963
          jt_vec_setup_evals: 0
           jt_vec_prod_evals: 3860






.. GENERATED FROM PYTHON SOURCE LINES 66-68

Integrate with direct linear solver
-----------------------------------

.. GENERATED FROM PYTHON SOURCE LINES 68-70

.. code-block:: Python

    timenp, Tnp, CO2np, NC6H14np  = integrate_reactor(preconditioner=False)





.. rst-class:: sphx-glr-script-out

 .. code-block:: none

    Non-preconditioned Integration Time: 421.477391
            step_solve_fails: 0
                       steps: 2795
                   rhs_evals: 4456
             nonlinear_iters: 4453
        nonlinear_conv_fails: 40
              err_test_fails: 184
                  last_order: 4
       stab_order_reductions: 0
                   jac_evals: 62
            lin_solve_setups: 416
               lin_rhs_evals: 78678
                   lin_iters: 0
              lin_conv_fails: 0
                  prec_evals: 0
                 prec_solves: 0
          jt_vec_setup_evals: 0
           jt_vec_prod_evals: 0






.. GENERATED FROM PYTHON SOURCE LINES 71-73

Plot selected state variables
-----------------------------

.. GENERATED FROM PYTHON SOURCE LINES 73-93

.. code-block:: Python

    fig, (ax1, ax2, ax3) = plt.subplots(3, 1, figsize=(5, 8))
    # temperature plot
    ax1.set_xlabel("Time")
    ax1.set_ylabel("Temperature")
    ax1.plot(timenp, Tnp, linewidth=2)
    ax1.plot(timep, Tp, linewidth=2, linestyle=":")
    ax1.legend(["Normal", "Preconditioned"])
    # CO2 plot
    ax2.set_xlabel("Time")
    ax2.set_ylabel("CO2")
    ax2.plot(timenp, CO2np, linewidth=2)
    ax2.plot(timep, CO2p, linewidth=2, linestyle=":")
    ax2.legend(["Normal", "Preconditioned"])
    # C12H26 plot
    ax3.set_xlabel("Time")
    ax3.set_ylabel("NC6H14")
    ax3.plot(timenp, NC6H14np, linewidth=2)
    ax3.plot(timep, NC6H14p, linewidth=2, linestyle=":")
    ax3.legend(["Normal", "Preconditioned"])
    plt.show()



.. image-sg:: /examples/python/reactors/images/sphx_glr_preconditioned_integration_001.png
   :alt: preconditioned integration
   :srcset: /examples/python/reactors/images/sphx_glr_preconditioned_integration_001.png, /examples/python/reactors/images/sphx_glr_preconditioned_integration_001_2_00x.png 2.00x
   :class: sphx-glr-single-img






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

   **Total running time of the script:** (7 minutes 9.906 seconds)


.. _sphx_glr_download_examples_python_reactors_preconditioned_integration.py:

.. only:: html

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

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

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

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

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

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

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


.. only:: html

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

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