
: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/fortran/f77_demo.f"
.. LINE NUMBERS ARE GIVEN BELOW.

.. only:: html

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

        :ref:`Go to the end <sphx_glr_download_examples_fortran_f77_demo.f>`
        to download the full example code.

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

.. _sphx_glr_examples_fortran_f77_demo.f:

Fortran 77 demo
===============

This program uses functions defined in :doc:`demo_ftnlib.cpp
<demo_ftnlib>` to create an ideal gas mixture and print some of
its properties.

For a C++ version of this program, see :doc:`demo.cpp
<../cxx/demo>`.

Replace this sample main program with your program

.. tags:: Fortran 77, tutorial, thermodynamics, kinetics,
         transport

.. GENERATED FROM PYTHON SOURCE LINES 15-39

.. code-block:: FortranFixed

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

          program demo
          implicit double precision (a-h,o-z)
          parameter (MAXSP = 20, MAXRXNS = 100)
          double precision q(MAXRXNS), qf(MAXRXNS), qr(MAXRXNS)
          double precision diff(MAXSP)
          character*80 eq
          character*20 name
    c
          write(*,*)
          write(*,*) '********   Fortran 77 Test Program   ********'

          call newIdealGasMix('h2o2.yaml','ohmech','mixture-averaged')
          t = 1200.0
          p = 101325.0
          call setState_TPX_String(t, p,
         $     'H2:2, O2:1, OH:0.01, H:0.01, O:0.01')

    c
          write(*,*) 'Initial state properties:'
          write(*,10) temperature(), pressure(), density(),
         $     enthalpy_mole(), entropy_mole(), cp_mole()

.. GENERATED FROM PYTHON SOURCE LINES 40-44

Compute the equilibrium state
-----------------------------

Hold the specific enthalpy and pressure constant.

.. GENERATED FROM PYTHON SOURCE LINES 45-58

.. code-block:: FortranFixed
   :dedent: 1

          call equilibrate('HP')

          write(*,*) 'Equilibrium state properties:'
          write(*,10) temperature(), pressure(), density(),
         $     enthalpy_mole(), entropy_mole(), cp_mole()

     10   format(//'Temperature:   ',g14.5,' K'/
         $         'Pressure:      ',g14.5,' Pa'/
         $         'Density:       ',g14.5,' kg/m3'/
         $         'Molar Enthalpy:',g14.5,' J/kmol'/
         $         'Molar Entropy: ',g14.5,' J/kmol-K'/
         $         'Molar cp:      ',g14.5,' J/kmol-K'//)


.. GENERATED FROM PYTHON SOURCE LINES 59-61

Reaction information
--------------------

.. GENERATED FROM PYTHON SOURCE LINES 63-79

.. code-block:: FortranFixed

          irxns = nReactions()

    c     forward and reverse rates of progress should be equal
    c     in equilibrium states
          call getFwdRatesOfProgress(qf)
          call getRevRatesOfProgress(qr)

    c     net rates of progress should be zero in equilibrium states
          call getNetRatesOfProgress(q)

    c     for each reaction, print the equation and the rates of progress
          do i = 1,irxns
             call getReactionEqn(i,eq)
             write(*,20) eq,qf(i),qr(i),q(i)
     20      format(a27,3e14.5,' kmol/m3/s')
          end do

.. GENERATED FROM PYTHON SOURCE LINES 80-82

Transport properties
--------------------

.. GENERATED FROM PYTHON SOURCE LINES 82-99

.. code-block:: FortranFixed
   :dedent: 1

          dnu = viscosity()
          dlam = thermalConductivity()
          call getMixDiffCoeffs(diff)

          write(*,30) dnu, dlam
     30   format(//'Viscosity:             ',g14.5,'  Pa-s'/
         $         'Thermal conductivity:  ',g14.5,'  W/m/K'/)
          write(*,*) 'Species            ',
         $     '    Diffusion Coefficient'
          nsp = nSpecies()
          do k = 1, nsp
             call getSpeciesName(k, name)
             write(*,40) name, diff(k)
     40      format(' ',a20,e14.5,' m2/s')
          end do

          stop
          end

.. _sphx_glr_download_examples_fortran_f77_demo.f:

.. only:: html

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

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

      :download:`Download FortranFixed source code: f77_demo.f <f77_demo.f>`

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

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


.. only:: html

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

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