1% -------------------------------------------------------------------------
2% EXAMPLE: SHOCK_I_THERMO
3%
4% Influence of caloric models on jump conditions in normal shocks
5%
6% This script examines the effects of different caloric models on the jump
7% conditions (changes in temperature, density, and adiabatic index)
8% encountered in normal shock waves. The models under consideration are:
9%
10% 1. Thermochemical frozen: assumes a calorically perfect gas, where
11% specific heat values remain constant.
12%
13% 2. Frozen: assumes a thermally perfect gas, where specific heat values
14% vary with temperature.
15%
16% 3. Equilibrium: assumes a calorically imperfect gas, where both
17% thermal and caloric properties vary, accounting for chemical
18% reactions in equilibrium.
19%
20% @author: Alberto Cuadra Lara
21%
22% Last update Sep 05 2024
23% -------------------------------------------------------------------------
24
25% Import packages
26import combustiontoolbox.databases.NasaDatabase
27import combustiontoolbox.core.*
28import combustiontoolbox.shockdetonation.ShockSolver
29import combustiontoolbox.utils.display.*
30
31% Get Nasa database
32DB = NasaDatabase();
33
34% Define chemical system
35system = ChemicalSystem(DB);
36
37% Initialize mixture
38mix = Mixture(system);
39
40% Define chemical state
41set(mix, {'N2', 'O2'}, [79/21, 1]);
42
43% Initialize figure
44plotConfig = PlotConfig();
45plotConfig.innerposition = [0.05, 0.05, 0.45, 0.55];
46plotConfig.outerposition = [0.05, 0.05, 0.45, 0.55];
47ax1 = setFigure(plotConfig);
48ax2 = setFigure(plotConfig);
49ax3 = setFigure(plotConfig);
50
51% Calculations
52for i = 1:3
53
54 switch i
55 case 1
56 FLAG_TCHEM_FROZEN = true;
57 FLAG_FROZEN = false;
58 linestyle = ':';
59 case 2
60 FLAG_TCHEM_FROZEN = false;
61 FLAG_FROZEN = true;
62 linestyle = '--';
63 case 3
64 FLAG_TCHEM_FROZEN = false;
65 FLAG_FROZEN = false;
66 linestyle = '-';
67 end
68
69 % Define properties
70 mixArray1 = setProperties(mix, 'temperature', 300, 'pressure', 1, 'M1', 1:0.1:10);
71
72 % Initialize solver
73 solver = ShockSolver('problemType', 'SHOCK_I', 'FLAG_TCHEM_FROZEN', FLAG_TCHEM_FROZEN, 'FLAG_FROZEN', FLAG_FROZEN);
74
75 % Solve problem
76 [mixArray1, mixArray2] = solver.solveArray(mixArray1);
77
78 % Plots
79 ax1 = plotFigure('M1', [mixArray1.mach], 'T_2/T_1', [mixArray2.T] ./ [mixArray1.T], 'linestyle', linestyle, 'color', [0 0 0], 'ax', ax1);
80 ax2 = plotFigure('M1', [mixArray1.mach], '\rho_2/\rho_1', [mixArray2.rho] ./ [mixArray1.rho], 'linestyle', linestyle, 'color', [0 0 0], 'ax', ax2);
81 ax3 = plotFigure('M1', [mixArray1.mach], '\gamma_s', [mixArray2.gamma_s], 'linestyle', linestyle, 'color', [0 0 0], 'ax', ax3);
82
83 if i ~= 3
84 continue
85 end
86
87 ax1 = plotFigure('M1', [5 5], 'T_2/T_1', [ax1.YLim(1), ax1.YLim(2)], 'linestyle', '--', 'color', [0.5 0.5 0.5], 'ax', ax1);
88 ax2 = plotFigure('M1', [5 5], '\rho_2/\rho_1', [ax2.YLim(1), ax2.YLim(2)], 'linestyle', '--', 'color', [0.5 0.5 0.5], 'ax', ax2);
89 ax3 = plotFigure('M1', [5 5], '\gamma_s', [ax3.YLim(1), ax3.YLim(2)], 'linestyle', '--', 'color', [0.5 0.5 0.5], 'ax', ax3);
90end