1% -------------------------------------------------------------------------
2% EXAMPLE: HP - Humid air
3%
4% Compute adiabatic temperature and equilibrium composition at constant
5% pressure (e.g., 1.01325 bar) for lean to rich CH4-ideal_air mixtures at
6% standard conditions, a set of 10 species considered, and excess of air
7% of 15%, and a set of relative humidity (0, 20, 40, 60, 80, 100) [%]
8%
9% LS == {'CO2', 'CO', 'H2O', 'H2', 'O2', 'N2', 'NO', 'OH', 'O', 'H'}
10%
11% See wiki or list_species() for more predefined sets of species
12%
13% Results compared with [1]
14%
15% [1] Sánchez, Y. A. C., Piedrahita, C. A. R., & Serrano, E. G. F. (2014).
16% Influencia del aire húmedo en la combustión del metano. Scientia et
17% technica, 19(4), 364-370.
18%
19% @author: Alberto Cuadra Lara
20% Postdoctoral researcher - Group Fluid Mechanics
21% Universidad Carlos III de Madrid
22%
23% Last update Feb 21 2024
24% -------------------------------------------------------------------------
25
26% User inputs
27air_excess = 15; % [%]
28humidity_relative = 0:20:100; % [%]
29T = [10:2:60] + 273.15; % [K]
30p = 1 * 1.01325; % [bar]
31
32% Initialization
33self = App();
34DB_master = self.DB_master;
35DB = self.DB;
36LS = {'CO2', 'CO', 'H2O', 'H2', 'O2', 'N2', 'NO', 'OH', 'O', 'H'};
37
38% Definitions
39phi_t = 2; % stoichiometric number of moles of air per mole of fuel
40phi = 1 / (1 + air_excess / 100); % equivalence ratio [-]
41mm_air = self.DB.Air.mm; % [g/mol]
42mm_H2O = self.DB.H2O.mm; % [g/mol]
43n_air = 4.76; % total number of moles of stoichiometric ideal dry air (79% N2 and 21% O2 in volume)
44
45% Miscellaneous
46config = self.Misc.config;
47config.innerposition = [0.2 0.2 0.35 0.5]; % Set figure inner position [normalized]
48config.outerposition = [0.2 0.2 0.35 0.5]; % Set figure outer position [normalized]
49colors = brewermap(length(humidity_relative), 'Greys');
50FLAG_FIRST = true;
51
52for j = length(humidity_relative):-1:1
53 % Get specific humidity
54 W = humidity_specific(T, p, humidity_relative(j));
55 % Get moles H2O
56 n_H20 = W * n_air * (mm_air / mm_H2O);
57 % Solve problem
58 for i = length(T):-1:1
59 % Initialization
60 self = App('fast', DB_master, DB, LS);
61 % Initial conditions
62 self = set_prop(self, 'TR', T(i), 'pR', p, 'phi', phi);
63 self.PD.S_Fuel = {'CH4'};
64 self.PD.S_Oxidizer = {'N2', 'O2', 'H2O'};
65 self.PD.N_Oxidizer = phi_t ./ phi .* [79/21, 1, n_H20(i)];
66 % self.PD.ratio_oxidizers_O2 = [79, 21, nH20(i) * 21] / 21;
67 % Set additional inputs (depends of the problem selected)
68 self = set_prop(self, 'pP', self.PD.pR.value);
69 % Solve problem
70 self = solve_problem(self, 'HP');
71 % Get data
72 mix1{i} = self.PS.strR{1};
73 mix2{i} = self.PS.strP{1};
74 end
75
76 if FLAG_FIRST
77 index_CO2 = find_ind(self.S.LS, 'CO2');
78 index_H2O = find_ind(self.S.LS, 'H2O');
79 index_CH4 = find_ind(self.S.LS, 'CH4');
80 FLAG_FIRST = false;
81 end
82
83 % Extract data
84 mR = cell2vector(mix1, 'mi');
85 TR = cell2vector(mix1, 'T');
86 Yi_R = cell2vector(mix1, 'Yi');
87 Yi_R_CH4 = Yi_R(index_CH4, :);
88
89 mP = cell2vector(mix2, 'mi');
90 TP = cell2vector(mix2, 'T');
91 Yi_P = cell2vector(mix2, 'Yi');
92 Xi_P = cell2vector(mix2, 'Xi');
93 Yi_P_H2O = Yi_P(index_H2O, :);
94 Xi_P_CO2 = Xi_P(index_CO2, :);
95
96 m_CH4 = mR .* Yi_R_CH4;
97 m_H20 = mP .* Yi_P_H2O;
98
99 % Plot results
100 if ~exist('ax1', 'var')
101 ax1 = set_figure(config);
102 ax2 = set_figure(config);
103 ax3 = set_figure(config);
104 end
105
106 % Adiabatic flame temperature vs. reactant's temperature
107 ax1 = plot_figure('T_R\ [\rm K]', TR, 'T_P\ [\rm K]', TP, 'ax', ax1, 'color', colors(j, :));
108 % Molar fraction of CO2 vs. reactant's temperature
109 ax2 = plot_figure('T_R\ [\rm K]', TR, 'X_{\rm CO_2}', Xi_P_CO2, 'ax', ax2, 'color', colors(j, :));
110 % Water content after the combustion process (kg H2O / kg fuel) vs. reactant's temperature
111 ax3 = plot_figure('T_R\ [\rm K]', TR, '\rm Water\ content\ [kg_{H_2O}/kg_{CH_4}]', m_H20./m_CH4, 'ax', ax3, 'color', colors(j, :));
112end