
: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/matlab/ignite_hp.m"
.. LINE NUMBERS ARE GIVEN BELOW.

.. only:: html

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

        :ref:`Go to the end <sphx_glr_download_examples_matlab_ignite_hp.m>`
        to download the full example code.

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

.. _sphx_glr_examples_matlab_ignite_hp.m:

Constant pressure ignition with user-specified equations
========================================================


Solves the same ignition problem as :doc:`reactor1.m <reactor1>`, but uses
the local function ``REACTOR_ODE`` to implement the governing equations for
an adiabatic, constant-pressure, zero-dimensional reactor.

Requires: cantera >= 3.2.0

.. tags:: Matlab, combustion, user-defined model, ignition delay, plotting

.. GENERATED FROM PYTHON SOURCE LINES 0-39

.. code-block:: Matlab

    function ignite_hp(gas)
        tic
        help ignite_hp

        if nargin == 0
            gas = ct.Solution('gri30.yaml', 'gri30');
        end

        mw = gas.molecularWeights;
        gas.TPX = {1001.0, ct.OneAtm, 'H2:2,O2:1,N2:4'};
        gas.basis = 'mass';
        nsp = gas.nSpecies;

        y0 = [gas.T
              gas.Y'];
        tel = [0, 0.001];
        options = odeset('RelTol', 1.e-5, 'AbsTol', 1.e-12, 'Stats', 'on');
        t0 = cputime;
        out = ode15s(@reactor_ode, tel, y0, options, gas, mw, nsp);
        disp(['CPU time = ' num2str(cputime - t0)]);

        if nargout == 0
            % plot the temperature and OH mole fractions.
            figure(1);
            plot(out.x, out.y(1, :));
            xlabel('time');
            ylabel('Temperature');
            title(['Final T = ' num2str(out.y(1, end)), ' K']);

            figure(2);
            ioh = gas.speciesIndex('OH');
            plot(out.x, out.y(1 + ioh, :));
            xlabel('time');
            ylabel('Mass Fraction');
            title('OH Mass Fraction');
        end

        toc
    end

.. GENERATED FROM PYTHON SOURCE LINES 50-52

Local Functions
---------------

.. GENERATED FROM PYTHON SOURCE LINES 53-79

.. code-block:: Matlab

    function dydt = reactor_ode(t, y, gas, mw, nsp)
        % ODE system for a constant-pressure, adiabatic reactor
        %
        % Function ``REACTOR_ODE`` evaluates the system of ordinary differential equations
        % for an adiabatic, constant-pressure, zero-dimensional reactor.
        % It assumes that the ``gas`` object represents a reacting ideal gas mixture.

        % Set the state of the gas, based on the current solution vector.
        gas.TPY = {y(1), gas.P, y(2:end)};

        % energy equation
        wdot = gas.netProdRates;
        H = gas.partialMolarEnthalpies';
        tdot =- 1 / (gas.D * gas.cp) .* wdot * H;

        % set up column vector for dydt
        dydt = [tdot
                zeros(nsp, 1)];

        % species equations
        rrho = 1.0 / gas.D;

        for i = 1:nsp
            dydt(i + 1) = rrho * mw(i) * wdot(i);
        end

    end

.. _sphx_glr_download_examples_matlab_ignite_hp.m:

.. only:: html

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

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

      :download:`Download Matlab source code: ignite_hp.m <ignite_hp.m>`

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

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


.. only:: html

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

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