💻 Fortran Source Code Library
We currently offer 172 open-source, production-grade Fortran codes for offline testing. Run calculations locally on your own machine, view code structure, read technical explanations, and download compilation packages including sample input files.
Adiabatic Flame Temperature
Core Numerical Engine in Fortran 90 • 32 total downloads
combustion_aft.f90
! =========================================================================
! Source File: combustion_aft.f90
! =========================================================================
program combustion_aft
implicit none
integer :: fuel_type, iostat_val, i, n_sweep, iter
double precision :: excess_air, T_inlet, P_kPa
double precision :: C_at, H_at, O_at, M_fuel, LHV
double precision :: a_stoich, air_stoich, AFR_stoich, AFR_actual
double precision :: actual_O2, n_CO2, n_H2O, n_N2, n_O2_ex, n_total
double precision :: y_CO2, y_H2O, y_N2, y_O2
double precision :: T_ad, T_low, T_high, T_mid
double precision :: H_react, H_prod, hf_fuel
double precision :: cp_CO2, cp_H2O, cp_N2, cp_O2
double precision :: hf_CO2, hf_H2O
double precision :: P_H2O, T_dew
double precision :: exc_sw, T_ad_sw, AFR_sw
! Average cp in J/(mol.K)
cp_CO2 = 54.0d0
cp_H2O = 41.0d0
cp_N2 = 32.0d0
cp_O2 = 34.0d0
! Enthalpies of formation in kJ/mol
hf_CO2 = -393.5d0
hf_H2O = -241.8d0
! ── Read inputs ─────────────────────────────────────────────
read(*,*,iostat=iostat_val) fuel_type
if (iostat_val /= 0) then
write(*,*) 'ERROR: Invalid fuel type input.'
stop
end if
read(*,*,iostat=iostat_val) excess_air
read(*,*,iostat=iostat_val) T_inlet
read(*,*,iostat=iostat_val) P_kPa
read(*,*,iostat=iostat_val) C_at
read(*,*,iostat=iostat_val) H_at
read(*,*,iostat=iostat_val) O_at
if (iostat_val /= 0) then
write(*,*) 'ERROR: Failed to read all inputs.'
stop
end if
! ── Input validation ────────────────────────────────────────
if (excess_air < 0.0d0) excess_air = 0.0d0
if (T_inlet <= 0.0d0) T_inlet = 298.15d0
if (P_kPa <= 0.0d0) P_kPa = 101.325d0
! ── Fuel properties ─────────────────────────────────────────
if (fuel_type == 1) then
C_at = 1.0d0; H_at = 4.0d0; O_at = 0.0d0
M_fuel = 16.04d0; LHV = 50050.0d0
else if (fuel_type == 2) then
C_at = 3.0d0; H_at = 8.0d0; O_at = 0.0d0
M_fuel = 44.10d0; LHV = 46357.0d0
else if (fuel_type == 3) then
C_at = 8.0d0; H_at = 18.0d0; O_at = 0.0d0
M_fuel = 114.23d0; LHV = 44427.0d0
else if (fuel_type == 4) then
C_at = 0.0d0; H_at = 2.0d0; O_at = 0.0d0
M_fuel = 2.016d0; LHV = 119960.0d0
else
! Custom fuel
if (C_at < 0.0d0) C_at = 1.0d0
if (H_at < 0.0d0) H_at = 4.0d0
if (O_at < 0.0d0) O_at = 0.0d0
M_fuel = C_at*12.011d0 + H_at*1.008d0 + O_at*16.0d0
if (M_fuel < 1.0d0) M_fuel = 16.04d0
! Estimate LHV from Dulong-like formula (kJ/kg)
LHV = 33900.0d0*C_at*12.011d0/M_fuel + &
142300.0d0*H_at*1.008d0/M_fuel - &
9090.0d0*O_at*16.0d0/M_fuel
if (LHV < 1000.0d0) LHV = 50000.0d0
end if
! ── Stoichiometry ───────────────────────────────────────────
! CxHyOz + a(O2 + 3.76 N2) -> C CO2 + H/2 H2O + 3.76a N2
! a_stoich = C + H/4 - O_fuel/2
a_stoich = C_at + H_at / 4.0d0 - O_at / 2.0d0
if (a_stoich < 0.0d0) a_stoich = 0.001d0
air_stoich = a_stoich / 0.21d0 ! moles air per mole fuel
AFR_stoich = air_stoich * 28.97d0 / M_fuel ! mass basis
AFR_actual = AFR_stoich * (1.0d0 + excess_air / 100.0d0)
! Actual O2 supplied
actual_O2 = a_stoich * (1.0d0 + excess_air / 100.0d0)
! Product moles per mole fuel
n_CO2 = C_at
n_H2O = H_at / 2.0d0
n_N2 = 3.76d0 * actual_O2
n_O2_ex = a_stoich * excess_air / 100.0d0
n_total = n_CO2 + n_H2O + n_N2 + n_O2_ex
! Mole fractions
y_CO2 = n_CO2 / max(n_total, 1.0d-10)
y_H2O = n_H2O / max(n_total, 1.0d-10)
y_N2 = n_N2 / max(n_total, 1.0d-10)
y_O2 = n_O2_ex / max(n_total, 1.0d-10)
! ── Adiabatic flame temperature (bisection) ────────────────
! Enthalpy of formation of fuel from LHV:
! LHV = -hf_fuel + C*(-hf_CO2) + H/2*(-hf_H2O_g) [per mole fuel, roughly]
! Actually: hf_fuel = C*hf_CO2 + H/2*hf_H2O - (-LHV*M_fuel/1000)
hf_fuel = n_CO2 * hf_CO2 + n_H2O * hf_H2O + LHV * M_fuel / 1000.0d0
! Energy balance: H_react = H_prod(T_ad)
! H_react = hf_fuel + actual_O2*(0) + n_N2_react*(0) + cp contributions at T_inlet
! H_react preheat = [actual_O2*cp_O2 + 3.76*actual_O2*cp_N2] * (T_inlet-298.15) / 1000
H_react = hf_fuel + (actual_O2 * cp_O2 + 3.76d0*actual_O2 * cp_N2) * &
(T_inlet - 298.15d0) / 1000.0d0
! H_prod(T) = n_CO2*hf_CO2 + n_H2O*hf_H2O + [n_CO2*cp_CO2 + n_H2O*cp_H2O + n_N2*cp_N2 + n_O2_ex*cp_O2]*(T-298.15)/1000
! Solve: H_react = H_prod(T_ad)
! => hf_fuel + preheat = n_CO2*hf_CO2 + n_H2O*hf_H2O + cp_prod_total*(T_ad-298.15)/1000
! => T_ad = 298.15 + 1000 * (hf_fuel + preheat - n_CO2*hf_CO2 - n_H2O*hf_H2O) / cp_prod_total
T_low = 300.0d0
T_high = 5000.0d0
do iter = 1, 100
T_mid = 0.5d0 * (T_low + T_high)
H_prod = n_CO2*hf_CO2 + n_H2O*hf_H2O + &
(n_CO2*cp_CO2 + n_H2O*cp_H2O + n_N2*cp_N2 + n_O2_ex*cp_O2) * &
(T_mid - 298.15d0) / 1000.0d0
if (H_prod < H_react) then
T_low = T_mid
else
T_high = T_mid
end if
if (abs(T_high - T_low) < 0.1d0) exit
end do
T_ad = 0.5d0 * (T_low + T_high)
! ── Dew point ───────────────────────────────────────────────
P_H2O = y_H2O * P_kPa ! partial pressure of H2O in kPa
! Simplified Antoine for water: T_dew(K) ~ 273.15 + 100 * (P_H2O/101.325)^0.25
! More accurate simplified: T_dew = 273.15 + 45.0*log(P_H2O) + 1.0 (rough fit for low P)
if (P_H2O > 0.0d0) then
T_dew = 273.15d0 + 45.0d0 * log(max(P_H2O, 0.001d0)) + 7.0d0
else
T_dew = 273.15d0
end if
if (T_dew < 273.15d0) T_dew = 273.15d0
! ── Output ──────────────────────────────────────────────────
write(*,'(A)') '============================================================'
write(*,'(A)') ' COMBUSTION & ADIABATIC FLAME TEMPERATURE'
write(*,'(A)') '============================================================'
write(*,*)
write(*,'(A)') '--- INPUTS --------------------------------------------------'
write(*,'(A,I8)') ' Fuel Type = ', fuel_type
write(*,'(A,F10.2,A)') ' Fuel Formula : C', C_at, ''
write(*,'(A,F10.2)') ' H', H_at
if (O_at > 0.0d0) &
write(*,'(A,F10.2)') ' O', O_at
write(*,'(A,F10.2,A)') ' Molar Mass M_fuel = ', M_fuel, ' g/mol'
write(*,'(A,F12.2,A)') ' LHV = ', LHV, ' kJ/kg'
write(*,'(A,F10.2,A)') ' Excess Air = ', excess_air, ' percent'
write(*,'(A,F12.2,A)') ' Inlet Temperature = ', T_inlet, ' K'
write(*,'(A,F12.4,A)') ' Pressure = ', P_kPa, ' kPa'
write(*,*)
write(*,'(A)') '--- STOICHIOMETRY -------------------------------------------'
write(*,'(A,F10.4)') ' O2 stoichiometric (mol) = ', a_stoich
write(*,'(A,F10.4)') ' Air stoich (mol/mol fuel) = ', air_stoich
write(*,'(A,F10.4)') ' AFR stoichiometric (mass) = ', AFR_stoich
write(*,'(A,F10.4)') ' AFR actual (mass) = ', AFR_actual
write(*,'(A,F10.4)') ' Equivalence Ratio phi = ', 1.0d0/(1.0d0+excess_air/100.0d0)
write(*,*)
write(*,'(A)') '--- FLUE GAS COMPOSITION (per mole fuel) -------------------'
write(*,'(A,F10.4,A)') ' CO2 = ', n_CO2, ' mol'
write(*,'(A,F10.4,A)') ' H2O = ', n_H2O, ' mol'
write(*,'(A,F10.4,A)') ' N2 = ', n_N2, ' mol'
write(*,'(A,F10.4,A)') ' O2 (excess) = ', n_O2_ex, ' mol'
write(*,'(A,F10.4,A)') ' Total products = ', n_total, ' mol'
write(*,*)
write(*,'(A)') '--- FLUE GAS MOLE FRACTIONS ---------------------------------'
write(*,'(A,F10.4,A)') ' y_CO2 = ', y_CO2*100.0d0, ' percent'
write(*,'(A,F10.4,A)') ' y_H2O = ', y_H2O*100.0d0, ' percent'
write(*,'(A,F10.4,A)') ' y_N2 = ', y_N2*100.0d0, ' percent'
write(*,'(A,F10.4,A)') ' y_O2 = ', y_O2*100.0d0, ' percent'
write(*,*)
write(*,'(A)') '--- ADIABATIC FLAME TEMPERATURE -----------------------------'
write(*,'(A,F12.2,A)') ' T_ad = ', T_ad, ' K'
write(*,'(A,F12.2,A)') ' T_ad = ', T_ad - 273.15d0, ' C'
write(*,*)
write(*,'(A)') '--- DEW POINT -----------------------------------------------'
write(*,'(A,F10.4,A)') ' P_H2O (partial) = ', P_H2O, ' kPa'
write(*,'(A,F12.2,A)') ' T_dew = ', T_dew, ' K'
write(*,'(A,F12.2,A)') ' T_dew = ', T_dew - 273.15d0, ' C'
write(*,*)
! ── Sensitivity sweep: excess air 0% to 200% ──────────────
n_sweep = 40
write(*,'(A)') '--- SENSITIVITY: T_AD VS EXCESS AIR -------------------------'
write(*,'(A)') ' ExcessAir[%] T_ad[K] AFR_actual'
write(*,'(A)') ' -----------------------------------------------------------'
do i = 1, n_sweep
exc_sw = dble(i-1) * 200.0d0 / dble(n_sweep - 1)
! Recompute for this excess air
actual_O2 = a_stoich * (1.0d0 + exc_sw / 100.0d0)
n_CO2 = C_at
n_H2O = H_at / 2.0d0
n_N2 = 3.76d0 * actual_O2
n_O2_ex = a_stoich * exc_sw / 100.0d0
n_total = n_CO2 + n_H2O + n_N2 + n_O2_ex
H_react = hf_fuel + (actual_O2 * cp_O2 + 3.76d0*actual_O2 * cp_N2) * &
(T_inlet - 298.15d0) / 1000.0d0
T_low = 300.0d0; T_high = 5000.0d0
do iter = 1, 80
T_mid = 0.5d0*(T_low+T_high)
H_prod = n_CO2*hf_CO2 + n_H2O*hf_H2O + &
(n_CO2*cp_CO2+n_H2O*cp_H2O+n_N2*cp_N2+n_O2_ex*cp_O2) * &
(T_mid-298.15d0)/1000.0d0
if (H_prod < H_react) then; T_low=T_mid; else; T_high=T_mid; end if
if (abs(T_high-T_low)<0.1d0) exit
end do
T_ad_sw = 0.5d0*(T_low+T_high)
AFR_sw = AFR_stoich * (1.0d0 + exc_sw/100.0d0)
write(*,'(F10.2,4X,F12.2,4X,F10.4)') exc_sw, T_ad_sw, AFR_sw
end do
write(*,*)
write(*,'(A)') '--- CORRELATIONS USED ---------------------------------------'
write(*,'(A)') ' Stoichiometric: CxHyOz + a(O2+3.76N2) -> CO2+H2O+N2'
write(*,'(A)') ' a = C + H/4 - O_fuel/2'
write(*,'(A)') ' AFR = a/0.21 * 28.97 / M_fuel (mass basis)'
write(*,'(A)') ' T_ad by bisection on energy balance (constant cp)'
write(*,'(A)') ' Dew point from partial pressure of H2O (simplified)'
end program combustion_aft
Solver Description
Solves combustion stoichiometry and calculates the Adiabatic Flame Temperature (AFT) under constant-pressure or constant-volume conditions. Uses temperature-dependent polynomial specific heat equations for reactants and products. Computes flue gas molar composition, air-fuel ratios (gravimetric and molar), equivalence ratio, and combustion products' dew point temperature.
Key Numerical Methods & Architecture
- Input Redirection: Reads parameters sequentially from standard input (`stdin`) using Fortran sequential read (`read(*,*)`), ensuring modular integration.
- Modular Design: Formulated using pure mathematical routines, separation of equations from output formatting, and precise numerical solvers (e.g. bisection, Newton-Raphson).
- Standard Compliant: Written in clean, standards-compliant Fortran 90 to ensure cross-compiler compatibility.
🛠️ Local Compilation
To test this code on your machine, compile the source code file(s) using a standard Fortran compiler (e.g., `gfortran`).
Compilation Command:
gfortran -O3 combustion_aft.f90 -o combustion_aft
Execution Command:
Execute the program by feeding the sample input file into the program using stdin redirection:
combustion_aft < input.txt
📥 Downloads & Local Files
Preview of the required input file (input.txt):
! Fuel type (1=Methane, 2=Propane, 3=Octane, 4=Hydrogen, 5=Custom CxHyOz)
1
! Excess Air [%]
0.0
! Inlet Temperature [K]
298.15
! Pressure [kPa]
101.325
! Carbon atoms (x) for custom fuel
1.0
! Hydrogen atoms (y) for custom fuel
4.0
! Oxygen atoms (z) for custom fuel
0.0
1
! Excess Air [%]
0.0
! Inlet Temperature [K]
298.15
! Pressure [kPa]
101.325
! Carbon atoms (x) for custom fuel
1.0
! Hydrogen atoms (y) for custom fuel
4.0
! Oxygen atoms (z) for custom fuel
0.0