1% -------------------------------------------------------------------------
2% EXAMPLE: ROCKET Propellants considering an Infinite-Area-Chamber (IAC)
3%
4% Compute rocket propellant performance considering an Infinite-Area-Chamber
5% for lean to rich LH2-LOX mixtures at 101.325 bar, a set of 11 species
6% considered, a set of equivalence ratios phi contained in (1, 5) [-], and
7% an exit area ratio A_exit/A_throat = 3
8%
9% HYDROGEN_L == {'H','H2O','OH','H2','O','O3','O2','HO2','H2O2',...
10% 'H2bLb','O2bLb'}
11%
12% See wiki or setListspecies method from ChemicalSystem class for more
13% predefined sets of species
14%
15% @author: Alberto Cuadra Lara
16% Postdoctoral researcher - Group Fluid Mechanics
17% Universidad Carlos III de Madrid
18%
19% Last update Jul 29 2024
20% -------------------------------------------------------------------------
21
22% Import packages
23import combustiontoolbox.databases.NasaDatabase
24import combustiontoolbox.core.*
25import combustiontoolbox.rocket.*
26import combustiontoolbox.utils.display.*
27
28% Get Nasa database
29DB = NasaDatabase();
30
31% Define chemical system
32system = ChemicalSystem(DB, 'hydrogen_L');
33
34% Initialize mixture
35mix = Mixture(system);
36
37% Define chemical state
38set(mix, {'H2bLb'}, 'fuel', 1);
39set(mix, {'O2bLb'}, 'oxidizer', 1);
40
41% Define properties
42mixArray1 = setProperties(mix, 'temperature', 90, 'pressure', 100 * 1.01325, 'equivalenceRatio', 1:0.05:5, 'areaRatio', 3);
43
44% Initialize solver
45solver = RocketSolver('problemType', 'ROCKET_IAC');
46
47% Solve problem
48[mixArray1, mixArray2, mixArray3, mixArray4] = solver.solveArray(mixArray1);
49
50% Plot properties
51ax1 = plotProperties(repmat({'equivalenceRatio'}, 1, 12), mixArray4, {'T', 'p', 'h', 'e', 'g', 'cp', 's', 'gamma_s', 'sound', 'u', 'I_sp', 'I_vac'}, mixArray4, 'basis', {[], [], 'mi', 'mi', 'mi', 'mi', 'mi', [], [], [], [], []});
52ax1 = plotProperties(repmat({'equivalenceRatio'}, 1, 12), mixArray3, {'T', 'p', 'h', 'e', 'g', 'cp', 's', 'gamma_s', 'sound', 'u', 'I_sp', 'I_vac'}, mixArray3, 'basis', {[], [], 'mi', 'mi', 'mi', 'mi', 'mi', [], [], [], [], []}, 'ax', ax1);
53ax1 = plotProperties(repmat({'equivalenceRatio'}, 1, 10), mixArray2, {'T', 'p', 'h', 'e', 'g', 'cp', 's', 'gamma_s', 'sound', 'u'}, mixArray2, 'basis', {[], [], 'mi', 'mi', 'mi', 'mi', 'mi', [], [], []}, 'ax', ax1);
54leg = legend(ax1.Children(end), {'Exit', 'Throat', 'Chamber'}, 'Interpreter', 'latex', 'FontSize', ax1.Children(end).FontSize);
55
56% Plot molar fractions
57plotComposition(mixArray2(1), mixArray1, 'equivalenceRatio', 'Xi', 'mintol', 1e-14, 'y_var', mixArray2, 'title', 'Chamber');
58plotComposition(mixArray3(1), mixArray1, 'equivalenceRatio', 'Xi', 'mintol', 1e-14, 'y_var', mixArray3, 'title', 'Throat');
59plotComposition(mixArray4(1), mixArray1, 'equivalenceRatio', 'Xi', 'mintol', 1e-14, 'y_var', mixArray4, 'title', 'Exit');