1% -------------------------------------------------------------------------
2% EXAMPLE: DET REFLECTED
3%
4% Compute pre-shock and post-shock state for a reflected planar detonation
5% considering Chapman-Jouguet (CJ) theory for lean to rich CH4-air mixtures
6% at standard conditions, a set of 24 species considered and a set of
7% equivalence ratios (phi) contained in (0.5, 5) [-]
8%
9%
10% Soot formation == {'CO2','CO','H2O','H2','O2','N2','Ar','Cbgrb',...
11% 'C2','C2H4','CH','CH3','CH4','CN','H',...
12% 'HCN','HCO','N','NH','NH2','NH3','NO','O','OH'}
13%
14% See wiki or setListspecies method from ChemicalSystem class for more
15% predefined sets of species
16%
17% @author: Alberto Cuadra Lara
18% Postdoctoral researcher - Group Fluid Mechanics
19% Universidad Carlos III de Madrid
20%
21% Last update July 28 2024
22% -------------------------------------------------------------------------
23
24% Import packages
25import combustiontoolbox.databases.NasaDatabase
26import combustiontoolbox.core.*
27import combustiontoolbox.shockdetonation.*
28import combustiontoolbox.utils.display.*
29
30% Get Nasa database
31DB = NasaDatabase();
32
33% Define chemical system
34system = ChemicalSystem(DB, 'soot formation');
35
36% Initialize mixture
37mix = Mixture(system);
38
39% Define chemical state
40set(mix, {'CH4'}, 'fuel', 1);
41set(mix, {'N2', 'O2', 'Ar', 'CO2'}, 'oxidizer', [78.084, 20.9476, 0.9365, 0.0319] / 20.9476);
42
43% Define properties
44mixArray1 = setProperties(mix, 'temperature', 300, 'pressure', 1.01325, 'equivalenceRatio', 0.5:0.01:5);
45
46% Initialize solver
47solver = DetonationSolver('problemType', 'DET_R');
48
49% Solve problem
50[mixArray1, mixArray2, mixArray3] = solver.solveArray(mixArray1);
51
52% Plot Hugoniot curves
53ax1 = plotFigure('\rho_1 / \rho_i', [mixArray1.rho] ./ [mixArray2.rho], 'p_i / p_1', [mixArray2.p] ./ [mixArray1.p], 'xScale', 'log', 'yScale', 'log', 'linestyle', '-');
54ax1 = plotFigure('\rho_1 / \rho_i', [mixArray1.rho] ./ [mixArray3.rho], 'p_i / p_1', [mixArray3.p] ./ [mixArray1.p], 'xScale', 'log', 'yScale', 'log', 'linestyle', '--', 'ax', ax1);
55legend(ax1, {'Incident', 'Reflected'}, 'Interpreter', 'latex', 'FontSize', ax1.FontSize);
56
57% Plot post-shock temperature
58ax2 = plotFigure('equivalenceRatio', mixArray1, 'T', mixArray2, 'linestyle', '-');
59ax2 = plotFigure('equivalenceRatio', mixArray1, 'T', mixArray3, 'linestyle', '--', 'ax', ax2);
60legend(ax2, {'Incident', 'Reflected'}, 'Interpreter', 'latex', 'FontSize', ax2.FontSize);
61
62% Plot molar fractions
63plotComposition(mixArray2(1), mixArray1, 'equivalenceRatio', 'Xi', 'mintol', 1e-14, 'y_var', mixArray2);
64plotComposition(mixArray3(1), mixArray1, 'equivalenceRatio', 'Xi', 'mintol', 1e-14, 'y_var', mixArray3);