! Fortran 90 Demo
!
! This program illustrates using Cantera in Fortran 90 to compute
! thermodynamic, kinetic, and transport properties of a gas mixture.
!
! Keywords: tutorial, equilibrium, thermodynamics, kinetics, transport
! 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.
program main
! use the Cantera module
use cantera
implicit none
! objects representing phases of matter have type 'phase_t'
type(phase_t) :: gas, surf
type(interface_t) iface
integer nsp, nrxns
double precision :: t, p
character*1000 err_buf
write(*,*)
write(*,*) '******** Fortran 90 Test Program ********'
! Read in a definition of the 'gas' phase.
! This will take the definition with name 'ohmech' from file
! 'h2o2.yaml', located in the Cantera data directory
gas = importPhase('h2o2.yaml','ohmech')
t = 1200.0 ! K
p = 101325.0 ! Pa
! set the temperature, pressure, and mole fractions.
call setState_TPX(gas, t, p, 'H2:1, O2:1, AR:2')
nsp = nSpecies(gas) ! number of species
nrxns = nReactions(gas) ! number of reactions
call demo(gas, nsp, nrxns)
! Demo of importing a phase with surface kinetics
gas = importPhase('ptcombust.yaml', 'gas')
iface = importInterface('ptcombust.yaml', 'Pt_surf', gas)
surf = iface%surf
call setState_TPX(gas, t, p, 'H2:0.1, O2:0.65, H2O:0.2, CO:0.05')
call setState_TPX(surf, t, p, 'O(S):0.01, PT(S): 0.8, CO(S): 0.19')
call demo_surf(surf, nSpecies(surf), nReactions(surf))
end program main
!--------------------------------------------------------
subroutine demo(gas, MAXSP, MAXRXNS)
! use the Cantera module
use cantera
implicit none
! declare the arguments
type(phase_t), intent(inout) :: gas
integer, intent(in) :: MAXSP
integer, intent(in) :: MAXRXNS
double precision q(MAXRXNS), qf(MAXRXNS), qr(MAXRXNS)
double precision diff(MAXSP)
double precision molar_cp(MAXSP), molar_h(MAXSP), g_rt(MAXSP)
character*80 eq
character*20 name
double precision :: dnu, dlam
integer :: i, irxns, nsp, k
write(*,*) 'Initial state properties:'
write(*,10) temperature(gas), pressure(gas), density(gas), &
enthalpy_mole(gas), entropy_mole(gas), cp_mole(gas)
! compute the equilibrium state holding the specific
! enthalpy and pressure constant
call equilibrate(gas, 'HP')
write(*,*) 'Equilibrium state properties:'
write(*,10) temperature(gas), pressure(gas), density(gas), &
enthalpy_mole(gas), entropy_mole(gas), cp_mole(gas)
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'//)
! Reaction information
irxns = nReactions(gas)
! forward and reverse rates of progress should be equal
! in equilibrium states
call getFwdRatesOfProgress(gas, qf)
call getRevRatesOfProgress(gas, qr)
! net rates of progress should be zero in equilibrium states
call getNetRatesOfProgress(gas, q)
! for each reaction, print the equation and the rates of progress
do i = 1,irxns
call getReactionString(gas, i,eq)
write(*,20) eq,qf(i),qr(i),q(i)
20 format(a27,3e14.5,' kmol/m3/s')
end do
! transport properties
dnu = viscosity(gas)
dlam = thermalConductivity(gas)
call getMixDiffCoeffs(gas, 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(gas)
do k = 1, nsp
call getSpeciesName(gas, k, name)
write(*,40) name, diff(k)
40 format(' ',a20,e14.5,' m2/s')
end do
! Thermodynamic properties
CALL getPartialMolarCp(gas, molar_cp)
CALL getPartialMolarEnthalpies(gas, molar_h)
CALL getGibbs_RT(gas, g_rt)
write(*,*) 'Species molar cp, molar enthalpy, normalized Gibbs free energy'
do k = 1, nsp
call getSpeciesName(gas, k, name)
write(*,50) name, molar_cp(k), molar_h(k), g_rt(k)
50 format(' ',a20,g14.5,' J/kmol-K',g14.5,' J/kmol',g14.5,'-')
end do
return
end subroutine demo
subroutine demo_surf(surf, nsp, nrxn)
use cantera
implicit none
type(phase_t) surf
integer :: nsp, nrxn, i
character*40 equation
double precision :: ropf(nrxn), ropnet(nrxn)
write(*,*)
write(*,*) '******** Interface Kinetics Test ********'
write(*,*)
call getFwdRatesOfProgress(surf, ropf)
call getNetRatesOfProgress(surf, ropnet)
write(*,*) 'Equation Fwd rate Net rate'
do i = 1, nrxn
call getReactionString(surf, i, equation)
write(*,60) equation, ropf(i), ropnet(i)
60 format(' ',a40, e14.5, e14.5)
end do
end subroutine demo_surf