1% -------------------------------------------------------------------------
2% EXAMPLE: EXOPLANET WASP-43b - METALLICITY 1
3%
4% Compute equilibrium vertical composition with a metallicity 1 of WASP-43b
5%
6% URL RESULTS TEA:
7% https://github.com/dzesmin/RRC-BlecicEtal-2015a-ApJS-TEA/tree/master/Fig6/WASP43b-solar
8%
9% @author: Alberto Cuadra Lara
10% Postdoctoral researcher - Group Fluid Mechanics
11% Universidad Carlos III de Madrid
12%
13% Last update April 02 2024
14% -------------------------------------------------------------------------
15
16% Import packages
17import combustiontoolbox.databases.SolarAbundances
18import combustiontoolbox.databases.NasaDatabase
19import combustiontoolbox.core.*
20import combustiontoolbox.equilibrium.*
21import combustiontoolbox.utils.display.*
22
23% Definitions
24listSpecies = {'C2H2_acetylene', 'C2H4', 'C', 'CH4', 'CO2', 'CO', 'H2', 'H2O', 'H2S', 'H', 'HCN', 'He', 'HS_M', 'N2', 'N', 'NH3', 'O', 'S'};
25species = {'H', 'He', 'C', 'N', 'O', 'S'};
26metallicity = 1;
27
28% Get initial composition from solar abundances
29DB_solar = SolarAbundances();
30moles = DB_solar.abundances2moles(species, metallicity);
31
32% Get Nasa database
33DB = NasaDatabase();
34
35% Define chemical system
36system = ChemicalSystem(DB, listSpecies);
37
38% Initialize mixture
39mix = Mixture(system);
40
41% Define chemical state
42set(mix, species, moles);
43
44% Define properties
45mixArray = setProperties(mix, 'temperature', linspace(100, 4000, 300), 'pressure', logspace(-5, 2, 300));
46
47% Initialize solver
48solver = EquilibriumSolver('problemType', 'TP');
49
50% Solve problem
51solver.solveArray(mixArray);
52
53% Plot molar fractions
54plotComposition(mixArray(1), mixArray, 'Xi', 'p', 'mintol', 1e-14, 'ydir', 'reverse', 'xscale', 'log');