From 1ab939ee566988612700d4306764b7c297b2c61a Mon Sep 17 00:00:00 2001
From: ssun30 <sun.su@northeastern.edu>
Date: Wed, 18 Jun 2025 14:18:16 -0400
Subject: [PATCH 01/14] [MATLAB] Fixed Sim1D.display Changed refine criteria
 for Diff_flame

---
 interfaces/matlab_experimental/OneDim/Sim1D.m | 9 ++-------
 samples/matlab_experimental/diff_flame.m      | 2 +-
 2 files changed, 3 insertions(+), 8 deletions(-)

diff --git a/interfaces/matlab_experimental/OneDim/Sim1D.m b/interfaces/matlab_experimental/OneDim/Sim1D.m
index b506f736b8b..83207ad0c27 100644
--- a/interfaces/matlab_experimental/OneDim/Sim1D.m
+++ b/interfaces/matlab_experimental/OneDim/Sim1D.m
@@ -54,14 +54,9 @@ function delete(s)
 
         %% Sim1D Utility Methods
 
-        function display(s, fname)
+        function display(s)
             % Show all domains.
-
-            if nargin == 1
-                fname = '-';
-            end
-
-            ctFunc('sim1D_show', s.stID, fname);
+            ctFunc('sim1D_show', s.stID);
         end
 
         function restore(s, fname, id)
diff --git a/samples/matlab_experimental/diff_flame.m b/samples/matlab_experimental/diff_flame.m
index 021e54f0ef9..769fa42b3c0 100644
--- a/samples/matlab_experimental/diff_flame.m
+++ b/samples/matlab_experimental/diff_flame.m
@@ -104,7 +104,7 @@
 % ``help setRefineCriteria``.
 
 f.energyEnabled = true;
-fl.setRefineCriteria(2, 200.0, 0.1, 0.2);
+fl.setRefineCriteria(2, 4, 0.2, 0.3, 0.04);
 fl.solve(loglevel, refine_grid);
 
 %% Show statistics of solution and elapsed time

From 3256fca38b58c341e684122344d242f62a2903fe Mon Sep 17 00:00:00 2001
From: Colin Willits <willitscb@MacBook-Air.local>
Date: Thu, 19 Jun 2025 16:01:32 -0400
Subject: [PATCH 02/14] [MATLAB] add correct units for axis for periodic_cstr.m
 sample

---
 samples/matlab_experimental/periodic_cstr.m | 2 ++
 1 file changed, 2 insertions(+)

diff --git a/samples/matlab_experimental/periodic_cstr.m b/samples/matlab_experimental/periodic_cstr.m
index 380bc513006..e35f3f77433 100644
--- a/samples/matlab_experimental/periodic_cstr.m
+++ b/samples/matlab_experimental/periodic_cstr.m
@@ -109,6 +109,8 @@
     plot(tm, y)
     legend('H2', 'O2', 'H2O')
     title('Mass Fractions')
+    ylabel('Mass Fractions')
+    xlabel('Time (s)')
 
     toc
 end

From 326043cec4136944fbe1bc88accb11fa6382c29e Mon Sep 17 00:00:00 2001
From: ssun30 <sun.su@northeastern.edu>
Date: Fri, 20 Jun 2025 10:38:19 -0400
Subject: [PATCH 03/14] [MATLAB] Added massDensity property to ThermoPhase
 Modified several samples to include this change

---
 interfaces/matlab_experimental/Base/ThermoPhase.m | 6 ++++++
 samples/matlab_experimental/conhp.m               | 2 +-
 samples/matlab_experimental/conuv.m               | 2 +-
 samples/matlab_experimental/diamond_cvd.m         | 2 +-
 samples/matlab_experimental/flame.m               | 2 +-
 samples/matlab_experimental/ignite.m              | 2 +-
 samples/matlab_experimental/periodic_cstr.m       | 2 +-
 samples/matlab_experimental/plug_flow_reactor.m   | 2 +-
 8 files changed, 13 insertions(+), 7 deletions(-)

diff --git a/interfaces/matlab_experimental/Base/ThermoPhase.m b/interfaces/matlab_experimental/Base/ThermoPhase.m
index 8f42f53db8d..78b7cfc1667 100644
--- a/interfaces/matlab_experimental/Base/ThermoPhase.m
+++ b/interfaces/matlab_experimental/Base/ThermoPhase.m
@@ -66,6 +66,8 @@
         %     Scalar double mean molecular weight. Units: kg/kmol.
         meanMolecularWeight
 
+        massDensity % Mass basis density. Units: kg/m^3.
+
         molarDensity % Molar basis density. Units: kmol/m^3.
 
         molecularWeights % Molecular weights of the species. Units: kg/kmol.
@@ -903,6 +905,10 @@ function display(tp)
             mmw = ctFunc('thermo_meanMolecularWeight', tp.tpID);
         end
 
+        function density = get.massDensity(tp)
+            density = ctFunc('thermo_density', tp.tpID);
+        end
+
         function density = get.molarDensity(tp)
             density = ctFunc('thermo_molarDensity', tp.tpID);
         end
diff --git a/samples/matlab_experimental/conhp.m b/samples/matlab_experimental/conhp.m
index 8483a40d0a0..d4776825063 100644
--- a/samples/matlab_experimental/conhp.m
+++ b/samples/matlab_experimental/conhp.m
@@ -6,6 +6,7 @@
     % It assumes that the ``gas`` object represents a reacting ideal gas mixture.
 
     % Set the state of the gas, based on the current solution vector.
+    gas.basis = 'mass';
     gas.Y = y(2:end);
     gas.TP = {y(1), gas.P};
     nsp = gas.nSpecies;
@@ -13,7 +14,6 @@
     % energy equation
     wdot = gas.netProdRates;
     H = gas.partialMolarEnthalpies';
-    gas.basis = 'mass';
     tdot =- 1 / (gas.D * gas.cp) .* wdot * H;
 
     % set up column vector for dydt
diff --git a/samples/matlab_experimental/conuv.m b/samples/matlab_experimental/conuv.m
index 84f47687777..b183f25e5bf 100644
--- a/samples/matlab_experimental/conuv.m
+++ b/samples/matlab_experimental/conuv.m
@@ -6,6 +6,7 @@
     % It assumes that the ``gas`` object represents a reacting ideal gas mixture.
 
     % Set the state of the gas, based on the current solution vector.
+    gas.basis = 'mass';
     gas.Y = y(2:end);
     gas.TD = {y(1), gas.D};
     nsp = gas.nSpecies;
@@ -13,7 +14,6 @@
     % energy equation
     wdot = gas.netProdRates;
     U = gas.partialMolarIntEnergies';
-    gas.basis = 'mass';
     tdot =- 1 / (gas.D * gas.cv) .* wdot * U;
 
     % set up column vector for dydt
diff --git a/samples/matlab_experimental/diamond_cvd.m b/samples/matlab_experimental/diamond_cvd.m
index 651d4ec6f62..dd9499e5af2 100644
--- a/samples/matlab_experimental/diamond_cvd.m
+++ b/samples/matlab_experimental/diamond_cvd.m
@@ -69,7 +69,7 @@
     r = surf_phase.netProdRates;
     carbon_dot = r(iC);
     mdot = mw * carbon_dot;
-    rate = mdot / dbulk.D;
+    rate = mdot / dbulk.massDensity;
     xx = [xx; x(ih)];
     rr = [rr; rate * 1.0e6 * 3600.0];
     cov = [cov; surf_phase.coverages];
diff --git a/samples/matlab_experimental/flame.m b/samples/matlab_experimental/flame.m
index 12ba9beeb56..ded43474c3b 100644
--- a/samples/matlab_experimental/flame.m
+++ b/samples/matlab_experimental/flame.m
@@ -32,7 +32,7 @@
     f = Sim1D({left flow right});
 
     % set default initial profiles.
-    rho0 = gas.D;
+    rho0 = gas.massDensity;
 
     % find the adiabatic flame temperature and corresponding
     % equilibrium composition
diff --git a/samples/matlab_experimental/ignite.m b/samples/matlab_experimental/ignite.m
index 239b6c0aa3d..13f7e00625a 100644
--- a/samples/matlab_experimental/ignite.m
+++ b/samples/matlab_experimental/ignite.m
@@ -25,7 +25,7 @@
     gas.TPX = {1001.0, OneAtm, 'H2:2,O2:1,N2:4'};
     gas.basis = 'mass';
     y0 = [gas.U
-          1.0 / gas.D
+          1.0 / gas.massDensity
           gas.Y'];
 
     time_interval = [0 0.001];
diff --git a/samples/matlab_experimental/periodic_cstr.m b/samples/matlab_experimental/periodic_cstr.m
index e35f3f77433..18414fc598e 100644
--- a/samples/matlab_experimental/periodic_cstr.m
+++ b/samples/matlab_experimental/periodic_cstr.m
@@ -72,7 +72,7 @@
     % controller (constant mdot). Set the mass flow rate to 1.25 sccm.
     sccm = 1.25;
     vdot = sccm * 1.0e-6/60.0 * ((OneAtm / gas.P) * (gas.T / 273.15)); % m^3/s
-    mdot = gas.D * vdot; % kg/s
+    mdot = gas.massDensity * vdot; % kg/s
     mfc = MassFlowController(upstream, cstr);
     mfc.massFlowRate = mdot;
 
diff --git a/samples/matlab_experimental/plug_flow_reactor.m b/samples/matlab_experimental/plug_flow_reactor.m
index d026fcf0315..94e16b03519 100644
--- a/samples/matlab_experimental/plug_flow_reactor.m
+++ b/samples/matlab_experimental/plug_flow_reactor.m
@@ -99,7 +99,7 @@
 
 T_calc(1) = gas_calc.T;
 Y_calc(1, :) = gas_calc.Y;
-rho_calc(1) = gas_calc.D;
+rho_calc(1) = gas_calc.massDensity;
 
 for i = 2:length(x_calc)
 

From 9d3180be0c11ab5d5b5b71de8f3daef270ca88cf Mon Sep 17 00:00:00 2001
From: Colin Willits <willitscb@MacBook-Air.local>
Date: Fri, 11 Jul 2025 17:31:39 -0400
Subject: [PATCH 04/14] [MATLAB] Fixed some samples and added try catch to
 test_examples.m

---
 samples/matlab_experimental/equil.m         |  4 +-
 samples/matlab_experimental/isentropic.m    |  2 +-
 samples/matlab_experimental/periodic_cstr.m | 31 ++++---
 samples/matlab_experimental/prandtl1.m      |  2 +-
 samples/matlab_experimental/prandtl2.m      |  2 +-
 samples/matlab_experimental/rankine.m       | 10 ++-
 samples/matlab_experimental/reactor1.m      | 11 ++-
 samples/matlab_experimental/reactor2.m      |  6 +-
 samples/matlab_experimental/test_examples.m | 89 ++++++++++++++-------
 9 files changed, 101 insertions(+), 56 deletions(-)

diff --git a/samples/matlab_experimental/equil.m b/samples/matlab_experimental/equil.m
index e06aaae1dd5..21e87baaa1b 100644
--- a/samples/matlab_experimental/equil.m
+++ b/samples/matlab_experimental/equil.m
@@ -20,7 +20,7 @@ function equil(g)
 
     nsp = gas.nSpecies;
 
-    % find methane, nitrogen, and oxygen indices
+    %% find methane, nitrogen, and oxygen indices
     ich4 = gas.speciesIndex('CH4');
     io2 = gas.speciesIndex('O2');
     in2 = gas.speciesIndex('N2');
@@ -41,7 +41,7 @@ function equil(g)
         xeq(:, i) = gas.X;
     end
 
-    % make plots
+    %% make plots
     clf;
     subplot(1, 2, 1);
     plot(phi, tad);
diff --git a/samples/matlab_experimental/isentropic.m b/samples/matlab_experimental/isentropic.m
index d0ac15ddf13..fee820e3acb 100644
--- a/samples/matlab_experimental/isentropic.m
+++ b/samples/matlab_experimental/isentropic.m
@@ -1,5 +1,5 @@
 function isentropic(g)
-    %% Isentropic, adiabatic flow
+    %% Adiabatic, isentropic flow
     %
     % In this example, the area ratio vs. Mach number curve is computed for a
     % hydrogen/nitrogen gas mixture.
diff --git a/samples/matlab_experimental/periodic_cstr.m b/samples/matlab_experimental/periodic_cstr.m
index 18414fc598e..ccc5471b5e2 100644
--- a/samples/matlab_experimental/periodic_cstr.m
+++ b/samples/matlab_experimental/periodic_cstr.m
@@ -30,7 +30,7 @@
     tic
     help periodic_cstr
 
-    % create the gas mixture
+    %% create the gas mixture
     gas = Solution('h2o2.yaml', 'ohmech');
 
     % pressure = 60 Torr, T = 770 K
@@ -39,26 +39,27 @@
 
     gas.TPX = {300, p, 'H2:2, O2:1'};
 
-    % create an upstream reservoir that will supply the reactor.  The
-    % temperature, pressure, and composition of the upstream reservoir are
+    %%
+    % Create an upstream reservoir that will supply the reactor.  
+    % The temperature, pressure, and composition of the upstream reservoir are
     % set to those of the 'gas' object at the time the reservoir is
     % created.
     upstream = Reservoir(gas);
-
+    %%
     % Now set the gas to the initial temperature of the reactor, and create
     % the reactor object.
     gas.TP = {t, p};
     cstr = IdealGasReactor(gas);
-
+    %%
     % Set its volume to 10 cm^3. In this problem, the reactor volume is
     % fixed, so the initial volume is the volume at all later times.
     cstr.V = 10.0 * 1.0e-6;
-
+    %%
     % We need to have heat loss to see the oscillations. Create a
     % reservoir to represent the environment, and initialize its
     % temperature to the reactor temperature.
     env = Reservoir(gas);
-
+    %%
     % Create a heat-conducting wall between the reactor and the
     % environment. Set its area, and its overall heat transfer
     % coefficient. Larger U causes the reactor to be closer to isothermal.
@@ -67,7 +68,7 @@
     w = Wall(cstr, env);
     w.area = 1.0;
     w.heatTransferCoeff = 0.02;
-
+    %%
     % Connect the upstream reservoir to the reactor with a mass flow
     % controller (constant mdot). Set the mass flow rate to 1.25 sccm.
     sccm = 1.25;
@@ -75,19 +76,19 @@
     mdot = gas.massDensity * vdot; % kg/s
     mfc = MassFlowController(upstream, cstr);
     mfc.massFlowRate = mdot;
-
+    %%
     % now create a downstream reservoir to exhaust into.
     downstream = Reservoir(gas);
-
+    %%
     % connect the reactor to the downstream reservoir with a valve, and
     % set the coefficient sufficiently large to keep the reactor pressure
     % close to the downstream pressure of 60 Torr.
     v = Valve(cstr, downstream);
     v.valveCoeff = 1.0e-9;
-
+    %%
     % create the network
     network = ReactorNet({cstr});
-
+    %%
     % now integrate in time
     tme = 0.0;
     dt = 0.1;
@@ -104,6 +105,12 @@
         y(3, n) = cstr.massFraction('H2O');
     end
 
+    for i = 2:length(tm)
+        if abs(y(3, i-1) - y(3, i)) > (y(3, i-1) * 5)
+            disp(tm(i))
+        end
+    end
+    %% plot
     clf
     figure(1)
     plot(tm, y)
diff --git a/samples/matlab_experimental/prandtl1.m b/samples/matlab_experimental/prandtl1.m
index a6c038cdc59..2b9b031f374 100644
--- a/samples/matlab_experimental/prandtl1.m
+++ b/samples/matlab_experimental/prandtl1.m
@@ -56,7 +56,7 @@ function prandtl1(g)
 
     disp(['CPU time = ' num2str(cputime - t0)]);
 
-    % plot results
+    %% plot results
 
     clf;
     %figure(1);
diff --git a/samples/matlab_experimental/prandtl2.m b/samples/matlab_experimental/prandtl2.m
index bc93a6992bc..1b56b63d987 100644
--- a/samples/matlab_experimental/prandtl2.m
+++ b/samples/matlab_experimental/prandtl2.m
@@ -55,7 +55,7 @@ function prandtl2(g)
 
     disp(['CPU time = ' num2str(cputime - t0)]);
 
-    % plot results
+    %% plot results
 
     clf;
     subplot(2, 2, 1);
diff --git a/samples/matlab_experimental/rankine.m b/samples/matlab_experimental/rankine.m
index 44de1aedd45..311c597e952 100644
--- a/samples/matlab_experimental/rankine.m
+++ b/samples/matlab_experimental/rankine.m
@@ -10,35 +10,39 @@
 
 help rankine
 
-% Initialize parameters
+%% Initialize parameters
 eta_pump = 0.6;
 eta_turbine = 0.8;
 p_max = 8.0 * OneAtm;
 t1 = 300.0;
 
-% create an object representing water
+%% create an object representing water
 w = Water;
 
+%%
 % start with saturated liquid water at t1
 basis = 'mass';
 w.TQ = {t1, 0};
 h1 = w.H;
 p1 = w.P;
 
+%%
 % pump it to p2
 pump_work = pump(w, p_max, eta_pump);
 h2 = w.H;
 
+%%
 % heat to saturated vapor
 w.PQ = {p_max, 1.0};
 h3 = w.H;
 
 heat_added = h3 - h2;
 
+%%
 % expand adiabatically back to the initial pressure
 turbine_work = expand(w, p1, eta_turbine);
 
-% compute the efficiency
+%% compute the efficiency
 efficiency = (turbine_work - pump_work) / heat_added;
 disp(sprintf('efficiency = %.2f%%', efficiency*100));
 
diff --git a/samples/matlab_experimental/reactor1.m b/samples/matlab_experimental/reactor1.m
index f6d2232a6dc..cd52e6ed232 100644
--- a/samples/matlab_experimental/reactor1.m
+++ b/samples/matlab_experimental/reactor1.m
@@ -20,7 +20,7 @@ function reactor1(g)
     end
 
     P = 101325.0;
-    % set the initial conditions
+    %% set the initial conditions
     gas.TP = {1001.0, P};
     nsp = gas.nSpecies;
     xx = zeros(nsp, 1);
@@ -29,26 +29,29 @@ function reactor1(g)
     xx(48) = 0.573;
     gas.X = xx;
 
-    % create a reactor, and insert the gas
+    %% create a reactor, and insert the gas
     r = IdealGasReactor(gas);
 
-    % create a reservoir to represent the environment
+    %% create a reservoir to represent the environment
     a = Solution('air.yaml', 'air', 'none');
     a.TP = {a.T, P};
     env = Reservoir(a);
 
+    %%
     % Define a wall between the reactor and the environment and
     % make it flexible, so that the pressure in the reactor is held
     % at the environment pressure.
     w = Wall(r, env);
 
+    %%
     % set expansion parameter. dV/dt = KA(P_1 - P_2)
     w.expansionRateCoeff = 1.0e6;
 
+    %%
     % set wall area
     w.area = 1.0;
 
-    % create a reactor network and insert the reactor:
+    %% create a reactor network and insert the reactor:
     network = ReactorNet({r});
 
     nSteps = 100;
diff --git a/samples/matlab_experimental/reactor2.m b/samples/matlab_experimental/reactor2.m
index 06cc112b26c..bc26abd2391 100644
--- a/samples/matlab_experimental/reactor2.m
+++ b/samples/matlab_experimental/reactor2.m
@@ -19,13 +19,13 @@ function reactor2(g)
         gas = Solution('gri30.yaml', 'gri30', 'none');
     end
 
-    % set the initial conditions
+    %% set the initial conditions
     gas.TPX = {1001.0, OneAtm, 'H2:2,O2:1,N2:4'};
 
-    % create a reactor, and insert the gas
+    %% create a reactor, and insert the gas
     r = IdealGasReactor(gas);
 
-    % create a reactor network and insert the reactor
+    %% create a reactor network and insert the reactor
     network = ReactorNet({r});
 
     nSteps = 100;
diff --git a/samples/matlab_experimental/test_examples.m b/samples/matlab_experimental/test_examples.m
index 00f389af144..6cb424a2015 100644
--- a/samples/matlab_experimental/test_examples.m
+++ b/samples/matlab_experimental/test_examples.m
@@ -1,31 +1,62 @@
 % runs selected examples without pausing
+run_test_examples();
 
-ctLoad
-clear all
-close all
-
-equil();
-isentropic();
-reactor1();
-reactor2();
-surf_reactor;
-periodic_cstr;
-% plug_flow_reactor;
-lithium_ion_battery
-rankine;
-prandtl1();
-prandtl2();
-% flame1;
-flame2;
-catcomb;
-diff_flame;
-ignite;
-ignite_hp;
-ignite_uv;
-diamond_cvd;
-
-clear all
-close all
-ctUnload
-
-disp('Test example run successful.')
+function run_test_examples()
+    clear all
+    close all
+    ctLoad
+
+    examples = {
+        'equil', 'isentropic', 'reactor1', 'reactor2', 'surf_reactor', ...
+        'periodic_cstr', 'plug_flow_reactor', 'lithium_ion_battery', ...
+        'rankine', 'prandtl1', 'prandtl2', 'flame1', 'flame2', ...
+        'catcomb', 'diff_flame', 'ignite', 'ignite_hp', 'ignite_uv', 'diamond_cvd'
+    };
+
+    passed = {};
+    failed = {};
+
+    for idx = 1:length(examples)
+        scriptName = examples{idx};
+        try
+            % Run script in base workshop to protect local variables like
+            % examples
+            evalin('base', sprintf('run(''%s.m'')', scriptName));
+            passed{end+1} = scriptName;
+        catch ME
+            fprintf('An error occurred while running %s: %s\n', scriptName, ME.message);
+            fprintf('Identifier: %s\n', ME.identifier);
+            if strcmp(ME.identifier, 'Cantera:ctError')
+                fprintf('Caught a CanteraError. Continuing execution...\n');
+            end
+            failed{end+1} = scriptName;
+        end
+    end
+    
+
+        % Summary report
+    disp(' ');
+    disp('============================');
+    disp('Summary of example runs');
+    disp('============================');
+
+    if ~isempty(passed)
+        fprintf('✅ Passed: %s\n', strjoin(passed, ', '));
+    else
+        disp('✅ Passed: (none)');
+    end
+
+    if ~isempty(failed)
+        fprintf('❌ Failed: %s\n', strjoin(failed, ', '));
+    else
+        disp('❌ Failed: (none)');
+    end
+
+    disp(' ');
+
+    clear all
+    close all
+    ctUnload
+    
+    disp('Test example run successfully.');
+end

From 9e3cc28f51611adbc39fd5bd3fe5026871d45418 Mon Sep 17 00:00:00 2001
From: Colin Willits <willitscb@MacBook-Air.local>
Date: Fri, 11 Jul 2025 17:43:55 -0400
Subject: [PATCH 05/14] [MATLAB] Fixed test_examples.m to use the try catch
 format

---
 samples/matlab_experimental/periodic_cstr.m | 2 +-
 samples/matlab_experimental/test_examples.m | 6 +++---
 2 files changed, 4 insertions(+), 4 deletions(-)

diff --git a/samples/matlab_experimental/periodic_cstr.m b/samples/matlab_experimental/periodic_cstr.m
index ccc5471b5e2..6d5ddd34cfb 100644
--- a/samples/matlab_experimental/periodic_cstr.m
+++ b/samples/matlab_experimental/periodic_cstr.m
@@ -40,7 +40,7 @@
     gas.TPX = {300, p, 'H2:2, O2:1'};
 
     %%
-    % Create an upstream reservoir that will supply the reactor.  
+    % Create an upstream reservoir that will supply the reactor.
     % The temperature, pressure, and composition of the upstream reservoir are
     % set to those of the 'gas' object at the time the reservoir is
     % created.
diff --git a/samples/matlab_experimental/test_examples.m b/samples/matlab_experimental/test_examples.m
index 6cb424a2015..6acc689627e 100644
--- a/samples/matlab_experimental/test_examples.m
+++ b/samples/matlab_experimental/test_examples.m
@@ -27,12 +27,12 @@ function run_test_examples()
             fprintf('An error occurred while running %s: %s\n', scriptName, ME.message);
             fprintf('Identifier: %s\n', ME.identifier);
             if strcmp(ME.identifier, 'Cantera:ctError')
-                fprintf('Caught a CanteraError. Continuing execution...\n');
+                disp('Caught a CanteraError. Continuing execution...\n');
             end
             failed{end+1} = scriptName;
         end
     end
-    
+
 
         % Summary report
     disp(' ');
@@ -57,6 +57,6 @@ function run_test_examples()
     clear all
     close all
     ctUnload
-    
+
     disp('Test example run successfully.');
 end

From f1868d605cedb8997ffba2f200bfa2ce59cd0bf7 Mon Sep 17 00:00:00 2001
From: ssun30 <sun.su@northeastern.edu>
Date: Wed, 11 Jun 2025 21:39:17 -0400
Subject: [PATCH 06/14] [MATLAB] Changed how search paths are set Replaced
 ctRoot with ctPaths, this will allow the user to set paths to the shared
 libraries, header files, and toolboxes more freely, instead of relying on
 relative paths to a fixed path for a certain conda environment.

---
 .../matlab_experimental/Utility/ctLoad.m      | 41 ++++++-----
 .../matlab_experimental/Utility/ctPaths.m     | 73 +++++++++++++++++++
 test/matlab_experimental/ctTestPath.m         | 23 ------
 test/matlab_experimental/ctTestSetUp.m        | 25 +------
 test/matlab_experimental/ctTestTearDown.m     | 11 +--
 5 files changed, 101 insertions(+), 72 deletions(-)
 create mode 100644 interfaces/matlab_experimental/Utility/ctPaths.m
 delete mode 100644 test/matlab_experimental/ctTestPath.m

diff --git a/interfaces/matlab_experimental/Utility/ctLoad.m b/interfaces/matlab_experimental/Utility/ctLoad.m
index cdedd769c0c..783910292b5 100644
--- a/interfaces/matlab_experimental/Utility/ctLoad.m
+++ b/interfaces/matlab_experimental/Utility/ctLoad.m
@@ -1,32 +1,39 @@
 function ctLoad()
+    % ctLoad
     % Load the Cantera C Library into Memory
 
+    paths = ctPaths();
+
+    if any(cellfun(@isempty, {paths.libPath, paths.includePath, paths.toolboxPath}))
+        error('ctLoad:MissingPath', ...
+              'Library, header, and toolbox paths must be configured with ctPaths(libPath,includePath, toolboxPath).');
+    end
+
     if ispc
-        ctName = '/bin/cantera_shared.dll';
+        libName = 'cantera_shared.dll';
     elseif ismac
-        ctName = '/Lib/libcantera_shared.dylib';
+        libName = 'libcantera_shared.dylib';
     elseif isunix
-        ctName = '/lib/libcantera_shared.so';
+        libName = 'libcantera_shared.so';
     else
-        error('Operating System Not Supported!');
-        return;
+        error('ctLoad:UnsupportedPlatform', 'Operating system not supported.');
     end
 
-    root = ctRoot;
-    if ispc
-        root = [ctRoot, '/Library'];
-    end
+    fullLibPath = fullfile(paths.libPath, libName);
+    fullHeaderPath = fullfile(paths.includePath, 'cantera', 'clib', 'ctmatlab.h');
 
     if ~libisloaded(ctLib)
-        [~, warnings] = loadlibrary([root, ctName], ...
-                                    [root, '/include/cantera/clib/ctmatlab.h'], ...
-                                    'includepath', [root, '/include'], ...
-                                    'addheader', 'ct', 'addheader', 'ctfunc', ...
-                                    'addheader', 'ctmultiphase', 'addheader', ...
-                                    'ctonedim', 'addheader', 'ctreactor', ...
-                                    'addheader', 'ctrpath', 'addheader', 'ctsurf');
+        [~, warnings] = loadlibrary(fullLibPath, fullHeaderPath, ...
+                                    'includepath', paths.includePath, ...
+                                    'addheader', 'ct', ...
+                                    'addheader', 'ctfunc', ...
+                                    'addheader', 'ctmultiphase', ...
+                                    'addheader', 'ctonedim', ...
+                                    'addheader', 'ctreactor', ...
+                                    'addheader', 'ctrpath', ...
+                                    'addheader', 'ctsurf');
     end
 
-    disp(sprintf('Cantera %s is ready for use.', ctVersion))
+    fprintf('Cantera %s is ready for use.\n', ctVersion);
 
 end
diff --git a/interfaces/matlab_experimental/Utility/ctPaths.m b/interfaces/matlab_experimental/Utility/ctPaths.m
new file mode 100644
index 00000000000..d5d95c795fc
--- /dev/null
+++ b/interfaces/matlab_experimental/Utility/ctPaths.m
@@ -0,0 +1,73 @@
+function paths = ctPaths(varargin)
+    % ctPaths ::
+    % Configure or retrieve the library/header/toolbox paths for Cantera.
+    % The paths are stored as MATLAB preferences.
+    %
+    %   >> ctPaths()                       % Get current config as struct
+    %   >> ctPaths(libPath, includePath)   % Set library/header/toolbox paths
+    %   >> ctPaths('clear')                % Clear saved paths
+    %
+    % :return:
+    %   paths: struct with fields 'libPath', 'includePath', and
+    %   'toolboxPath'
+
+    if nargin == 1 && strcmp(varargin{1}, 'clear')
+        if ispref('Cantera', 'Paths')
+            paths = getpref('Cantera', 'Paths');
+            subDirs = strsplit(genpath(paths.toolboxPath), pathsep);
+            currentPaths = strsplit(path, pathsep);
+            subdirsToRemove = intersect(subDirs, currentPaths);
+            if ~isempty(subdirsToRemove)
+                rmpath(subdirsToRemove{:});
+            end
+            rmpref('Cantera', 'Paths');
+        end
+        paths = struct('libPath', '', 'includePath', '', 'toolboxPath', '');
+        return
+    elseif nargin == 3 && all(cellfun(@ischar, varargin))
+        paths = struct('libPath', varargin{1}, ...
+                       'includePath', varargin{2}, ...
+                       'toolboxPath', varargin{3});
+        setpref('Cantera', 'Paths', paths);
+    else
+        % Load from saved preferences if available
+        if ispref('Cantera', 'Paths')
+            paths = getpref('Cantera', 'Paths');
+            return
+        else
+            paths = struct('libPath', '', 'includePath', '', 'toolboxPath', '');
+        end
+    end
+
+    if any(cellfun(@isempty, {paths.libPath, paths.includePath, paths.toolboxPath}))
+        error('ctPaths:MissingPath', ...
+              'Library, header, and toolbox paths must be configured with ctPaths(libPath,includePath, toolboxPath).');
+    end
+
+    mapping = {
+        'interfaces/matlab_experimental', 'toolbox';
+        'samples/matlab_experimental', 'samples';
+        'test/matlab_experimental', 'test/matlab_toolbox';
+        'test/data', 'test/data';
+        'data', 'data'
+        };
+
+    % Check whether user is using Cantera source code or MLTBX based on folder structure
+    if isfolder(fullfile(paths.toolboxPath, 'interfaces'))
+        col = 1;
+    else
+        col = 2;
+    end
+
+    for i = 1:size(mapping, 1)
+        subdir = fullfile(paths.toolboxPath, mapping{i, col});
+        if isfolder(subdir)
+            addpath(genpath(subdir));
+        else
+            warning('ctPaths:MissingDirectory', ...
+                    'Directory not found: %s', subdir);
+        end
+    end
+
+    savepath();
+end
diff --git a/test/matlab_experimental/ctTestPath.m b/test/matlab_experimental/ctTestPath.m
deleted file mode 100644
index f7fa212efc0..00000000000
--- a/test/matlab_experimental/ctTestPath.m
+++ /dev/null
@@ -1,23 +0,0 @@
-% ctTestPath.m
-% Set up environment for testing the Cantera Matlab interface
-% from within the Cantera source tree. Run this file from the
-% root of the Cantera source tree, for example:
-%
-%    cd ~/src/cantera
-%    run interfaces/matlab_experimental/testpath.m
-
-% get the list of directories on the Matlab path
-dirs = split(path, pathsep);
-
-% if 'cantera' is already in the path, remove it
-for i = 1:length(dirs)
-    if contains(dirs{i}, 'CANTERA', 'IgnoreCase', true)
-        rmpath(dirs{i});
-    end
-end
-
-cantera_root = getenv('CANTERA_ROOT');
-
-% Add the Cantera toolbox to the Matlab path
-addpath(genpath([cantera_root, '/interfaces/matlab_experimental']));
-addpath(genpath([cantera_root, '/test/matlab_experimental']));
diff --git a/test/matlab_experimental/ctTestSetUp.m b/test/matlab_experimental/ctTestSetUp.m
index 9590f149da7..c58f82d71e1 100644
--- a/test/matlab_experimental/ctTestSetUp.m
+++ b/test/matlab_experimental/ctTestSetUp.m
@@ -1,24 +1,3 @@
-clear all
-
-% Copy library to test folder
-ctTestPath;
-
-if ispc
-    ctName = '/build/lib/cantera_shared.dll';
-elseif ismac
-    ctName = '/build/lib/libcantera_shared.dylib';
-elseif isunix
-    ctName = '/build/lib/libcantera_shared.so';
+function ctTestSetUp()
+    ctLoad();
 end
-% Load Cantera
-if ~libisloaded('libcantera_shared')
-    [~, warnings] = loadlibrary([cantera_root, ctName], ...
-                                [cantera_root, '/include/cantera/clib/ctmatlab.h'], ...
-                                'includepath', [cantera_root, '/include'], ...
-                                'addheader', 'ct', 'addheader', 'ctfunc', ...
-                                'addheader', 'ctmultiphase', 'addheader', ...
-                                'ctonedim', 'addheader', 'ctreactor', ...
-                                'addheader', 'ctrpath', 'addheader', 'ctsurf');
-end
-
-disp('Cantera is loaded for test');
diff --git a/test/matlab_experimental/ctTestTearDown.m b/test/matlab_experimental/ctTestTearDown.m
index bab612e0866..2174435ec89 100644
--- a/test/matlab_experimental/ctTestTearDown.m
+++ b/test/matlab_experimental/ctTestTearDown.m
@@ -1,10 +1,3 @@
-clear all
-% Unload Cantera
-if ispc
-    lib = 'cantera_shared';
-else
-    lib = 'libcantera_shared';
+function ctTestTearDown()
+    ctUnload();
 end
-
-unloadlibrary(lib);
-disp('Cantera has been unloaded');

From 60d34e63d565f0b1a6ceff6fcaf0f879f9e64a3d Mon Sep 17 00:00:00 2001
From: ssun30 <sun.su@northeastern.edu>
Date: Wed, 11 Jun 2025 21:41:35 -0400
Subject: [PATCH 07/14] [MATLAB] Changed CI workflow to account for change to
 how paths are set

---
 .github/workflows/main.yml | 9 +++++++++
 1 file changed, 9 insertions(+)

diff --git a/.github/workflows/main.yml b/.github/workflows/main.yml
index 009cfcf3e39..f9ecb43ad3a 100644
--- a/.github/workflows/main.yml
+++ b/.github/workflows/main.yml
@@ -66,6 +66,15 @@ jobs:
           echo "LD_PRELOAD=$LIB_STDCXX:$LIB_OPENBLAS:$LIB_LAPACK" >> $GITHUB_ENV
       - name: Set up MATLAB
         uses: matlab-actions/setup-matlab@718d4320188c73c86eb94ce76b553cbf89778487 # v2.5.0
+      - name: Set MATLAB search paths for tests and build MLTBX toolbox # v2.2.1
+        uses: matlab-actions/run-command@67f012f1ee4bc1627a7a95a0ccdeec745c9f6f36
+        with:
+          command: |
+            toolboxPath = fullfile(getenv('CANTERA_ROOT'));
+            libPath = fullfile(toolboxPath, 'build', 'lib');
+            includePath = fullfile(toolboxPath, 'include');
+            addpath(fullfile(toolboxPath, 'interfaces', 'matlab_experimental', 'Utility'));
+            ctPaths(libPath, includePath, toolboxPath);
       - name: Run tests
         uses: matlab-actions/run-tests@f9fc3d8ca29fadef6227fa52884b144b9011fa2f # v2.1.1
         with:

From c2abf1eb8e026c4199b325f69f409d5662a44b0a Mon Sep 17 00:00:00 2001
From: Colin Willits <willitscb@MacBook-Air.local>
Date: Fri, 11 Jul 2025 18:03:45 -0400
Subject: [PATCH 08/14] [MATLAB] Added test examples job to CI

---
 .github/workflows/main.yml                  | 6 ++++++
 samples/matlab_experimental/test_examples.m | 2 +-
 2 files changed, 7 insertions(+), 1 deletion(-)

diff --git a/.github/workflows/main.yml b/.github/workflows/main.yml
index f9ecb43ad3a..a2b14c60f3a 100644
--- a/.github/workflows/main.yml
+++ b/.github/workflows/main.yml
@@ -79,6 +79,12 @@ jobs:
         uses: matlab-actions/run-tests@f9fc3d8ca29fadef6227fa52884b144b9011fa2f # v2.1.1
         with:
           select-by-folder: /home/runner/work/cantera/cantera/test/matlab_experimental
+      - name: Run MATLAB test_examples
+        uses: matlab-actions/run-command@67f012f1ee4bc1627a7a95a0ccdeec745c9f6f36 # v2.2.1
+        with:
+          command: |
+            addpath('/home/runner/work/cantera/cantera/samples/matlab_experimental');
+            test_examples;
 
   ubuntu-multiple-pythons:
     name: ${{ matrix.os }} with Python ${{ matrix.python-version }}, Numpy ${{ matrix.numpy || 'latest' }}, Cython ${{ matrix.cython || 'latest' }}
diff --git a/samples/matlab_experimental/test_examples.m b/samples/matlab_experimental/test_examples.m
index 6acc689627e..b8d3a3d8ec1 100644
--- a/samples/matlab_experimental/test_examples.m
+++ b/samples/matlab_experimental/test_examples.m
@@ -52,7 +52,7 @@ function run_test_examples()
         disp('❌ Failed: (none)');
     end
 
-    disp(' ');
+    disp('  ');
 
     clear all
     close all

From 898bded4f8fcdfc9690296dc4189edcedfa05971 Mon Sep 17 00:00:00 2001
From: Colin Willits <willitscb@MacBook-Air.local>
Date: Fri, 18 Jul 2025 18:00:26 -0400
Subject: [PATCH 09/14] [MATLAB] Modified the code for test_examples.m

---
 samples/matlab_experimental/test_examples.m | 28 ++++++++++++---------
 1 file changed, 16 insertions(+), 12 deletions(-)

diff --git a/samples/matlab_experimental/test_examples.m b/samples/matlab_experimental/test_examples.m
index b8d3a3d8ec1..4f593055885 100644
--- a/samples/matlab_experimental/test_examples.m
+++ b/samples/matlab_experimental/test_examples.m
@@ -1,10 +1,17 @@
 % runs selected examples without pausing
-run_test_examples();
+clear all
+close all
+ctLoad
+if run_test_examples()
+    disp('test_examples run successfully');
+end
+ctUnload
+clear all
+close all
+
+function success = run_test_examples()
 
-function run_test_examples()
-    clear all
-    close all
-    ctLoad
+    success = false;
 
     examples = {
         'equil', 'isentropic', 'reactor1', 'reactor2', 'surf_reactor', ...
@@ -50,13 +57,10 @@ function run_test_examples()
         fprintf('❌ Failed: %s\n', strjoin(failed, ', '));
     else
         disp('❌ Failed: (none)');
+        success = true;
     end
 
-    disp('  ');
-
-    clear all
-    close all
-    ctUnload
-
-    disp('Test example run successfully.');
+    disp(' ');
+    
+    
 end

From f793750ae2cca81cc3fd965aec38e480df56d5a0 Mon Sep 17 00:00:00 2001
From: Colin Willits <willitscb@MacBook-Air.local>
Date: Mon, 21 Jul 2025 09:04:35 -0400
Subject: [PATCH 10/14] [MATLAB] Added experimental samples crit_properties and
 shock_tube based on python examples

---
 samples/matlab_experimental/crit_properites.m |  61 +++++++++
 samples/matlab_experimental/shock_tube.m      | 121 ++++++++++++++++++
 2 files changed, 182 insertions(+)
 create mode 100644 samples/matlab_experimental/crit_properites.m
 create mode 100644 samples/matlab_experimental/shock_tube.m

diff --git a/samples/matlab_experimental/crit_properites.m b/samples/matlab_experimental/crit_properites.m
new file mode 100644
index 00000000000..27dd3bc8cd1
--- /dev/null
+++ b/samples/matlab_experimental/crit_properites.m
@@ -0,0 +1,61 @@
+%% Critical state properties
+% Print the critical state properties for the fluids for which Cantera has 
+% built-in liquid/vapor equations of state.
+%
+% .. tags:: Matlab, thermodynamics, multiphase, non-ideal fluid
+
+clear all;
+close all;
+
+tic % total runtime of script
+help crit_properites
+
+%% Create a pure fluid object
+fluids = struct('water', Water(), ...
+    'nitrogen', Nitrogen(), ...
+    'methane', Methane(), ...
+    'hydrogen', Hydrogen(), ...
+    'oxygen', Oxygen(), ...
+    'carbon_dioxide', CarbonDioxide(), ...
+    'heptane', Heptane(), ...
+    'HFC_134a', HFC134a());
+
+names = fieldnames(fluids);
+
+%% Plot critical properties and print tabulated values
+
+figure;
+hold on;
+title('Critical Properties of Pure Fluids');
+xlabel('Critical Temperature [K]');
+ylabel('Critical Pressure [Pa]');
+xlim([0 750]);
+grid on;
+
+% Print header
+fprintf('Critical State Properties\n');
+fprintf('%-16s   %-7s   %-10s   %-7s\n', 'Fluid', 'Tc [K]', 'Pc [Pa]', 'Zc');
+fprintf('%s   %s   %s   %s\n', repmat('-',1,16), repmat('-',1,7), repmat('-',1,10), repmat('-',1,7));
+
+% Loop through fluids
+for i = 1:length(names)
+    name = names{i};
+    f = fluids.(name);
+
+    Tc = f.critTemperature;
+    Pc = f.critPressure;
+    rhoc = f.critDensity;
+    mw = f.meanMolecularWeight;
+    R = GasConstant;
+
+    Zc = Pc * mw / (rhoc * R * Tc);
+
+    % Plot
+    plot(Tc, Pc, 'o');
+    text(Tc + 4, Pc + 2e5, strrep(name, '_', ' '), 'FontSize', 9);
+
+    % Print table row
+    fprintf('%-16s   %7.2f   %10.4g   %7.4f\n', strrep(name, '_', ' '), Tc, Pc, Zc);
+end
+
+toc
\ No newline at end of file
diff --git a/samples/matlab_experimental/shock_tube.m b/samples/matlab_experimental/shock_tube.m
new file mode 100644
index 00000000000..7e63b4c2c90
--- /dev/null
+++ b/samples/matlab_experimental/shock_tube.m
@@ -0,0 +1,121 @@
+%% Shock-tube species profiles as a function of time
+% Simulate species profiles for a shock tube as a function of time, and 
+% observe the impact of incorporating the reduced-pressure linear mixture 
+% rule (LMR-R) in such calculations.
+%
+% Here we predict the H2O mole fraction time profiles for a mixture of 1163
+% ppm H2O2/1330 ppm H2O/665 ppm O2/20% CO2/Ar following reflected shock 
+% waves (1196 K, 2.127 atm) and compare results against the experimental 
+% measurements of Shao et al. [1] Two models are compared in this example:
+%
+% 1. A 2023 model of H2 and NH3 chemistry published by Alzueta et al. [2]
+% 2. An adapted version of this model that has applied the reduced-pressure
+% linear mixture rule (LMR-R) and ab initio third-body efficiencies. [3]
+% References:
+%
+% [1] J. Shao, R. Choudhary, D. F. Davidson, R. K. Hanson, Shock tube/laser
+% absorption measurement of the rate constant of the reaction: H2O2+CO2 = 
+% 2OH+CO2, Proc. Combust. Inst. 39 (2023) 735 – 743.
+%
+% [2] M. U. Alzueta, I. Salas, H. Hashemi, P. Glarborg, CO-assisted NH3 
+% oxidation, Combust. Flame 257 (2023) 112438.
+%
+% [3] P. J. Singal, J. Lee, L. Lei, R. L. Speth, M. P. Burke, 
+% Implementation of New Mixture Rules Has a Substantial Impact on 
+% Combustion Predictions for H2 and NH3, Proc. Combust. Inst. 40 (2024) 
+% 105779.
+%
+% .. tags:: Matlab, shock tube, kinetics, combustion
+
+clear all;
+close all;
+
+tic % total running time of the script
+help shock_tube
+
+file = 'example_data/ammonia-CO-H2-Alzueta-2023.yaml';
+%% Models and colors
+models = struct('Original', 'baseline', 'LMR_R', 'linear-Burke');
+colors = struct('Original', [0.6, 0.6, 0.6], 'LMR_R', [0.5, 0, 0.5]); % grey and purple
+
+results = struct();
+
+%% Experimental data from Shao et al.
+expData.t = [12.3, 20.3, 26.4, 39.6, 58.5, 79.2, 96.1, 113.8, 131.6, ...
+             145.7, 161.2, 181.6, 195.3, 219.9, 237.2, 248.6, 262.4, ...
+             272.2, 280.9];  % microseconds
+expData.X_H2O = [1.47E-03, 1.59E-03, 1.66E-03, 1.78E-03, 1.98E-03, ...
+                 2.06E-03, 2.15E-03, 2.22E-03, 2.26E-03, 2.30E-03, ...
+                 2.39E-03, 2.38E-03, 2.40E-03, 2.42E-03, 2.47E-03, ...
+                 2.53E-03, 2.51E-03, 2.50E-03, 2.47E-03];
+
+%% Loop over mechanisms
+model_names = fieldnames(models);
+for k = 1:length(model_names)
+    model_key = model_names{k};
+    mech_name = models.(model_key);
+    
+    % Define mixture composition
+    X_H2O2 = 1163e-6;
+    X_H2O = 1330e-6;
+    X_O2   = 665e-6;
+    X_CO2  = 0.2 * (1 - X_H2O2 - X_H2O - X_O2);
+    X_Ar   = 1 - X_CO2;
+
+    gasComp = 'H2O2: 1163e-6, H2O: 1330e-6, O2: 665e-6, CO2: 0.1994, AR: 0.8006';
+
+    gas = Solution(file, mech_name);
+    gas.TPX = {1196, 2.127 * OneAtm, gasComp};
+
+    r = Reactor(gas,'Reactor');
+    r.energy = 'on';
+    sim = ReactorNet({r});
+
+    time = 0.0;
+    estIgnDelay = 1.0; % seconds
+    counter = 0;
+    dt = 5.0e-7;
+    
+    time_data = [];
+    H2O_X_data = [];
+
+    while time < estIgnDelay
+        time = time + dt;
+        sim.advance(time);
+        time = sim.time;
+        if mod(counter, 10) == 0
+            time_data(end+1) = time * 1e6;  % convert to µs
+            H2O_X_data(end+1) = gas.moleFraction({'H2O'});
+        end
+        counter = counter + 1;
+    end
+
+    % Store results
+    results.(model_key).t = time_data;
+    results.(model_key).X_H2O = H2O_X_data;
+end
+
+%% Plotting
+figure;
+hold on;
+
+for k = 1:length(model_names)
+    model_key = model_names{k};
+    plot(results.(model_key).t, 100 * results.(model_key).X_H2O, ...
+         'LineWidth', 2, 'Color', colors.(model_key), 'DisplayName', strrep(model_key, '_', '-'));
+end
+
+plot(expData.t, 100 * expData.X_H2O, 'ko', 'MarkerSize', 6, ...
+     'MarkerFaceColor', 'white', 'DisplayName', 'Shao et al.');
+
+xlabel('Time [\mus]');
+ylabel('H_2O Mole Fraction [%]');
+legend('Location', 'southwest');
+xlim([0 300]);
+title('Ignition Simulation vs. Experiment');
+grid on;
+
+
+toc
+
+

From 18b14f3ffb27101a36267e102012804dace0ed94 Mon Sep 17 00:00:00 2001
From: Colin Willits <willitscb@MacBook-Air.local>
Date: Tue, 29 Jul 2025 11:24:02 -0400
Subject: [PATCH 11/14] [MATLAB] Fixed shock tube example

---
 samples/matlab_experimental/shock_tube.m | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/samples/matlab_experimental/shock_tube.m b/samples/matlab_experimental/shock_tube.m
index 7e63b4c2c90..f19d94b994f 100644
--- a/samples/matlab_experimental/shock_tube.m
+++ b/samples/matlab_experimental/shock_tube.m
@@ -33,7 +33,7 @@
 tic % total running time of the script
 help shock_tube
 
-file = 'example_data/ammonia-CO-H2-Alzueta-2023.yaml';
+file = 'example_data/ammonia-CO-H2-Alzueta-2023.yaml'; % Had to add example/data here to find the file
 %% Models and colors
 models = struct('Original', 'baseline', 'LMR_R', 'linear-Burke');
 colors = struct('Original', [0.6, 0.6, 0.6], 'LMR_R', [0.5, 0, 0.5]); % grey and purple

From 8632f7833acaa4aec9a2b57769d731f3b0fa0e39 Mon Sep 17 00:00:00 2001
From: Colin Willits <willitscb@MacBook-Air.local>
Date: Wed, 30 Jul 2025 11:40:27 -0400
Subject: [PATCH 12/14] [MATLAB] Add vapor dome sample

[MATLAB] Fix vapor_dome and shock_tube and crit_properties samples
---
 samples/matlab_experimental/crit_properites.m |   2 +-
 samples/matlab_experimental/shock_tube.m      |   6 +-
 samples/matlab_experimental/vapor_dome.m      | 103 ++++++++++++++++++
 3 files changed, 106 insertions(+), 5 deletions(-)
 create mode 100644 samples/matlab_experimental/vapor_dome.m

diff --git a/samples/matlab_experimental/crit_properites.m b/samples/matlab_experimental/crit_properites.m
index 27dd3bc8cd1..fe55ed40882 100644
--- a/samples/matlab_experimental/crit_properites.m
+++ b/samples/matlab_experimental/crit_properites.m
@@ -58,4 +58,4 @@
     fprintf('%-16s   %7.2f   %10.4g   %7.4f\n', strrep(name, '_', ' '), Tc, Pc, Zc);
 end
 
-toc
\ No newline at end of file
+toc
diff --git a/samples/matlab_experimental/shock_tube.m b/samples/matlab_experimental/shock_tube.m
index f19d94b994f..f03cdc2a0df 100644
--- a/samples/matlab_experimental/shock_tube.m
+++ b/samples/matlab_experimental/shock_tube.m
@@ -74,15 +74,13 @@
     time = 0.0;
     estIgnDelay = 1.0; % seconds
     counter = 0;
-    dt = 5.0e-7;
+    dt = 1e-5;
     
     time_data = [];
     H2O_X_data = [];
 
     while time < estIgnDelay
-        time = time + dt;
-        sim.advance(time);
-        time = sim.time;
+        time = sim.step();
         if mod(counter, 10) == 0
             time_data(end+1) = time * 1e6;  % convert to µs
             H2O_X_data(end+1) = gas.moleFraction({'H2O'});
diff --git a/samples/matlab_experimental/vapor_dome.m b/samples/matlab_experimental/vapor_dome.m
new file mode 100644
index 00000000000..2adcd24e39b
--- /dev/null
+++ b/samples/matlab_experimental/vapor_dome.m
@@ -0,0 +1,103 @@
+%% Vapor Dome
+% This example generates a saturated steam table and plots the vapor dome.
+% The steam table corresponds to data typically found in thermodynamic text
+% books and uses the same customary units
+%
+% .. tags:: Matlab, thermodynamics, non-ideal fluid, plotting
+
+clear all
+close all
+
+tic
+help vapor_dome
+
+
+w = Water;
+w.basis = {'mass'};
+
+%% Temperature vector in degC
+degC = [w.minTemp - 273.15, 4, 5, 6, 8, ...
+        10:36, 38, ...
+        40:5:95, 100:10:290, ...
+        300:20:360, 370:1:373, ...
+        w.critTemperature - 273.15];
+
+T_vec = degC + 273.15;
+n = length(T_vec);
+
+%% Preallocate arrays
+T = zeros(n,1); P = zeros(n,1);
+vf = zeros(n,1); vg = zeros(n,1); vfg = zeros(n,1);
+uf = zeros(n,1); ug = zeros(n,1); ufg = zeros(n,1);
+hf = zeros(n,1); hg = zeros(n,1); hfg = zeros(n,1);
+sf = zeros(n,1); sg = zeros(n,1); sfg = zeros(n,1);
+
+%% Loop through temperatures and get saturated liquid/vapor properties
+for i = 1:n
+    T_i = T_vec(i);
+    T(i) = T_i - 273.15; % store in degC
+
+    % Saturated vapor
+    w.TQ = {T_i, 1};
+    vg(i) = w.V;
+    ug(i) = w.U/1e3;
+    hg(i) = w.H/1e3;
+    sg(i) = w.S/1e3;
+
+    % Saturated liquid
+    w.TQ = {T_i, 0};
+    vf(i) = w.V;
+    uf(i) = w.U/1e3;
+    hf(i) = w.H/1e3;
+    sf(i) = w.S/1e3;
+    
+    % Pressure (same for both states)
+    P(i) = w.P/1e5;
+end
+
+%% Delta values
+vfg = vg - vf;
+ufg = ug - uf;
+hfg = hg - hf;
+sfg = sg - sf;
+
+%% Reference state: triple point liquid
+w.TQ = {w.minTemp, 0};
+uf0 = w.U/1e3;
+hf0 = w.H/1e3;
+sf0 = w.S/1e3;
+pv0 = w.P * w.V/1e3;
+
+%% Adjust reference state
+uf = uf - uf0;
+ug = ug - uf0;
+hf = hf - hf0 + pv0;
+hg = hg - hf0 + pv0;
+sf = sf - sf0;
+sg = sg - sf0;
+
+%% Print and write saturated steam table to csv file
+data = table(T, P, vf, vfg, vg, uf, ufg, ug, hf, hfg, hg, sf, sfg, sg);
+disp(data);
+writetable(data, 'saturated_steam_T.csv');
+%% Plot P-v Diagram
+figure;
+semilogx(vf, P, 'b', 'LineWidth', 1.5); hold on;
+semilogx(vg, P, 'r', 'LineWidth', 1.5);
+plot(vg(end), P(end), 'ko', 'MarkerFaceColor', 'w');
+xlabel('Specific Volume v [m^3/kg]');
+ylabel('Pressure P [bar]');
+legend('Saturated Liquid', 'Saturated Vapor', 'Critical Point');
+title('Vapor Dome - P-v Diagram');
+grid on;
+
+%% Plot T-s Diagram
+figure;
+plot(sf, T, 'b', 'LineWidth', 1.5); hold on;
+plot(sg, T, 'r', 'LineWidth', 1.5);
+plot(sg(end), T(end), 'ko', 'MarkerFaceColor', 'w');
+xlabel('Specific Entropy s [kJ/kg-K]');
+ylabel('Temperature T [°C]');
+legend('Saturated Liquid', 'Saturated Vapor', 'Critical Point');
+title('Vapor Dome - T-s Diagram');
+grid on;

From 86bcfdc069e9a488f6ced02155a49df61939576e Mon Sep 17 00:00:00 2001
From: Colin Willits <willitscb@MacBook-Air.local>
Date: Mon, 4 Aug 2025 11:51:46 -0400
Subject: [PATCH 13/14] [MATLAB] Fix samples and add equations_of_state.m

[MATLAB] Fix equations_of_state.m
---
 samples/matlab_experimental/crit_properites.m |   4 +-
 .../matlab_experimental/equations_of_state.m  | 228 ++++++++++++++++++
 samples/matlab_experimental/shock_tube.m      |   6 +-
 samples/matlab_experimental/vapor_dome.m      |   2 +-
 4 files changed, 234 insertions(+), 6 deletions(-)
 create mode 100644 samples/matlab_experimental/equations_of_state.m

diff --git a/samples/matlab_experimental/crit_properites.m b/samples/matlab_experimental/crit_properites.m
index fe55ed40882..2c367439018 100644
--- a/samples/matlab_experimental/crit_properites.m
+++ b/samples/matlab_experimental/crit_properites.m
@@ -7,7 +7,7 @@
 clear all;
 close all;
 
-tic % total runtime of script
+tic
 help crit_properites
 
 %% Create a pure fluid object
@@ -37,7 +37,7 @@
 fprintf('%-16s   %-7s   %-10s   %-7s\n', 'Fluid', 'Tc [K]', 'Pc [Pa]', 'Zc');
 fprintf('%s   %s   %s   %s\n', repmat('-',1,16), repmat('-',1,7), repmat('-',1,10), repmat('-',1,7));
 
-% Loop through fluids
+%% Loop through fluids
 for i = 1:length(names)
     name = names{i};
     f = fluids.(name);
diff --git a/samples/matlab_experimental/equations_of_state.m b/samples/matlab_experimental/equations_of_state.m
new file mode 100644
index 00000000000..fbfb45b2122
--- /dev/null
+++ b/samples/matlab_experimental/equations_of_state.m
@@ -0,0 +1,228 @@
+%% Non-ideal equations of state
+%
+% This example demonstrates a comparison between ideal and non-ideal equations of
+% state (EoS) using Cantera.
+%
+% The following equations of state are used to evaluate thermodynamic properties
+% in this example:
+%
+% 1. Ideal-gas EoS from Cantera
+% 2. Non-ideal Redlich-Kwong EoS (R-K EoS) from Cantera
+%
+% .. tags:: Matlab, thermodynamics, non-ideal fluid
+
+clear all
+close all
+
+tic
+help equations_of_state
+
+%% Helper Functions
+% This examples uses \mathrm{ CO_2 } as the only species. The function
+% "get_thermo_Cantera" calculates thermodynamic properties based on the
+% thermodynamic state (T, p) of the species using Cantera. Applicable phases
+% are "Ideal-gas" and "Redlich-Kwong". The ideal-gas equation can be stated
+% as
+%
+% .. math:: 
+%     \mathrm{ pv = RT },
+%
+% where p, v, and T represent thermodynamic pressure, molar volume, and the
+% temperature of the gas-phase. R is the universal gas constant. The
+% Redlich-Kwong equation is a cubic, non-ideal equation of state,
+% represented as
+%
+% .. math::
+%     \mathrm{ p = \frac{RT}{v-b^\ast}-\frac{a^\ast}{v\sqrt{T}(v+b^\ast)} }.
+%
+% In this expression, R is the universal gas constant and v is the molar
+% volume. The temperature-dependent van der Waals attraction parameter
+% \mathrm{ a^\ast } and volume correction parameter (repulsive parameter)
+% \mathrm{ b^\ast } represent molecular interactions.
+%
+% To plot the comparision of thermodynamic properties among the tree EoS,
+% the "plotEoS" function is used.
+
+function [h, u, s, cp, cv] = get_thermo_Cantera(phase, T, p)
+    phase.basis = {'mass'};
+    n = length(p);
+    u = zeros(1, n);
+    h = zeros(1, n);
+    s = zeros(1, n);
+    cp = zeros(1, n);
+    cv = zeros(1, n);
+
+    for i = 1:n
+        phase.TPX = {T, p(i), 'CO2:1.0'};
+        u(i) = phase.U / 1000;
+        h(i) = phase.H / 1000;
+        s(i) = phase.S / 1000;
+        cp(i) = phase.cp / 1000;
+        cv(i) = phase.cv / 1000;
+    end
+
+    % Reference to first point
+    u = u - u(1);
+    h = h - h(1);
+    s = s - s(1);
+end
+
+function plotEoS(p, ideal, rk, y_label)
+    figure;
+    plot(p / 1e5, ideal, 'b-', 'LineWidth', 2); 
+    hold on;
+    plot(p / 1e5, rk, 'r-', 'LineWidth', 2);
+    xlabel('Pressure [bar]');
+    ylabel(y_label);
+    legend('Ideal EoS', 'R-K EoS');
+    grid on;
+end
+
+%% EoS Comparison based on thermodynamic properties
+%
+% This is th emain subroutine that compares the plots and thermodynamic
+% values obtained using three equations of state.
+
+% Input parameters
+T = 300; % K
+p = 1e5 * linspace(1, 100, 1000); % Pa
+
+% Read the ideal-gas phase
+idealGasPhase = Solution('example_data/co2-thermo.yaml', 'CO2-Ideal');
+[h_ideal, u_ideal, s_ideal, cp_ideal, cv_ideal] = get_thermo_Cantera(idealGasPhase, T, p);
+
+% Read the Redlich-Kwong phase
+redlichKwongPhase = Solution('example_data/co2-thermo.yaml', 'CO2-RK');
+[h_RK, u_RK, s_RK, cp_RK, cv_RK] = get_thermo_Cantera(redlichKwongPhase, T, p);
+
+% Plot the results
+plotEoS(p, u_ideal, u_RK, "Relative Internal Energy [kJ/kg]");
+plotEoS(p, h_ideal, h_RK, "Relative Enthalpy [kJ/kg]"); 
+plotEoS(p, s_ideal, s_RK, "Relative Entropy [kJ/kg-K]");
+
+%%
+% The thermodynamic properties such as internal energy, enthalpy, and 
+% entropy are plotted against the operating pressure at a constant 
+% temperature T = 300 K. The three equations follow each other closely at low 
+% pressures (P < 10bar). However, the ideal gas EoS departs significantly 
+% from the observed behavior of gases near the critical regime 
+% (Pcrit = 73.77 bar).
+%
+% The ideal gas EoS does not consider inter-molecular interactions and the 
+% volume occupied by individual gas particles. At low temperatures and high
+% pressures, inter-molecular forces become particularly significant due to 
+% a reduction in inter-molecular distances. Additionally, at high density, 
+% the volume of individual molecules becomes significant. Both of these 
+% factors contribute to the deviation from ideal behavior at high pressures. 
+% The cubic Redlich-Kwong EoS, on the other hand, predicts thermodynamic 
+% properties accurately near the critical regime.
+
+% Specific heat at constant pressure 
+plotEoS(p, cp_ideal, cp_RK, "C_p [kJ/kg-K]");
+% Specific heat at constant volume
+plotEoS(p, cv_ideal, cv_RK, "C_v [kJ/kg-K]");
+%%
+% In the case of Ideal gas EoS, the specific heats at constant pressure
+% (Cp) and constant volume (Cv) are independent of the pressure. Hence, Cp
+% and Cv for ideal EoS do not change as the pressure is varied from 1 bar
+% to 100 bar in this study.
+%
+% Cp for the R-K EoS follows the trend closely with the Helmholtz EoS from
+% CoolProp up to the critical regime. Alhtough Cp shows reasonable
+% agreement with the Helmholtz EoS in sub-critical and supercritical
+% reginmes, it inaccurately predicts a very high value near the critical
+% point. However, Cp at the critical point is finite for the real fluid.
+% The sudden rise in Cp in the case of the R-K EoS is just a numerical
+% artifact, due to the EoS yielding infinite values in the limiting case,
+% and not a real singularity.
+%
+% Cv, on the other hand, predicts smaller values in the subcritical and
+% critical regime. However, it shows completely incorrect values in the
+% super-critical region, making it invalid at very high pressures. It is
+% well known that the cubic equations typically fail to predics accurate
+% constant-volume heat capacity in the transcritical region [2]_. Certain
+% cubic EoS models have been extended to resolve the discrepancy using
+% crossover models. For further information sese the work of Span [2]_ and
+% Saeed et al. [3]_.
+
+%% Temperature-Density plots
+%
+% The following function plots the T-'rho' diagram over a wide pressure and
+% temperature range. The temperature is varied from 250 K to 400 K. The
+% pressure is changed from 1 bar to 600 bar.
+
+% Input parameters
+% Set up arrays for pressure and temperature
+p_arr = logspace(1, log10(600), 10);
+T_arr = linspace(250, 401, 20);
+p_arr = 1e5 * p_arr(:);
+
+% Initialize matrices to hold densities
+density_ideal = zeros(length(p_arr), length(T_arr));
+density_RK = zeros(length(p_arr), length(T_arr));
+%density_CP = zeros(length(p_arr), length(T_arr));
+
+% Loop over pressure and temperature to compute densities
+for i = 1:length(p_arr)
+    for j = 1:length(T_arr)
+        T = T_arr(j);
+        P = p_arr(i);
+
+        % Ideal gas
+        idealGasPhase.TP = {T, P};
+        density_ideal(i, j) = idealGasPhase.D;
+
+        % Redlich-Kwong gas
+        redlichKwongPhase.TP = {T, P};
+        density_RK(i, j) = redlichKwongPhase.D;
+
+    end
+end
+
+% Plotting
+figure;
+hold on;
+colors = parula(length(p_arr));
+
+for i = 1:length(p_arr)
+    % Ideal gas lines (dashed)
+    ideal_line = plot(density_ideal(i,:), T_arr, '--', 'Color', colors(i,:), 'HandleVisibility', 'on');
+    % RK EoS lines (circles)
+    RK_line = plot(density_RK(i,:), T_arr, 'o', 'Color', colors(i,:), 'HandleVisibility', 'on');
+end
+xlabel('Density [kg/m^3]');
+
+ylabel('Temperature [K]');
+legend([ideal_line(1), RK_line(1)], {"Ideal EoS", "R-K EoS"}, 'Location', 'northeast'); 
+title('T vs Density for CO2 at Various Pressures');
+grid on;
+
+% Text annotations
+text(30, 320, 'p = 1 bar', 'Color', colors(1,:), 'Rotation', 90);
+text(430, 318, 'p = 97 bar', 'Color', colors(6,:), 'Rotation', -12);
+text(960, 320, 'p = 600 bar', 'Color', colors(10,:), 'Rotation', -68);
+
+toc
+%%
+% The figure compares T-\rho plots for ideal, R-K, and Helmholtz EoS at
+% different operating pressures. All three EoS yield the same plots at low pressures (0
+% bar and 10 bar). However, the Ideal gas EoS departs significantly at high pressures
+% (P > 10 bar), where non-ideal effects are prominent. The R-K EoS closely
+% matches the Helmholtz EoS at supercritical pressures (P > 70 bar). However,
+% it does depart in the liquid-vapor region that exists at P < P_crit
+% and low temperatures (T_crit).
+%
+% .. [1] I.H. Bell, J. Wronski, S. Quoilin, V. Lemort, "Pure and Pseudo-pure Fluid
+%    Thermophysical Property Evaluation and the Open-Source Thermophysical Property
+%    Library CoolProp," Industrial & Engineering Chemistry Research 53 (2014),
+%    https://pubs.acs.org/doi/10.1021/ie4033999
+%
+% .. [2] R. Span, "Multiparameter Equations of State - An Accurate Source of
+%    Thermodynamic Property Data," Springer Berlin Heidelberg (2000),
+%    https://dx.doi.org/10.1007/978-3-662-04092-8
+%
+% .. [3] A. Saeed, S. Ghader, "Calculation of density, vapor pressure and heat capacity
+%    near the critical point by incorporating cubic SRK EoS and crossover translation,"
+%    Fluid Phase Equilibria (2019) 493, https://doi.org/10.1016/j.fluid.2019.03.027
+%
+% sphinx_gallery_thumbnail_number = -1
diff --git a/samples/matlab_experimental/shock_tube.m b/samples/matlab_experimental/shock_tube.m
index f03cdc2a0df..54fffe50044 100644
--- a/samples/matlab_experimental/shock_tube.m
+++ b/samples/matlab_experimental/shock_tube.m
@@ -30,13 +30,13 @@
 clear all;
 close all;
 
-tic % total running time of the script
+tic
 help shock_tube
 
-file = 'example_data/ammonia-CO-H2-Alzueta-2023.yaml'; % Had to add example/data here to find the file
+file = 'example_data/ammonia-CO-H2-Alzueta-2023.yaml';
 %% Models and colors
 models = struct('Original', 'baseline', 'LMR_R', 'linear-Burke');
-colors = struct('Original', [0.6, 0.6, 0.6], 'LMR_R', [0.5, 0, 0.5]); % grey and purple
+colors = struct('Original', [0.6, 0.6, 0.6], 'LMR_R', [0.5, 0, 0.5]);
 
 results = struct();
 
diff --git a/samples/matlab_experimental/vapor_dome.m b/samples/matlab_experimental/vapor_dome.m
index 2adcd24e39b..c7a4b6e9fb4 100644
--- a/samples/matlab_experimental/vapor_dome.m
+++ b/samples/matlab_experimental/vapor_dome.m
@@ -35,7 +35,7 @@
 %% Loop through temperatures and get saturated liquid/vapor properties
 for i = 1:n
     T_i = T_vec(i);
-    T(i) = T_i - 273.15; % store in degC
+    T(i) = T_i - 273.15;
 
     % Saturated vapor
     w.TQ = {T_i, 1};

From 8d3c894a94c579b2d6e8460a79e0e66e235cd424 Mon Sep 17 00:00:00 2001
From: Colin Willits <willitscb@MacBook-Air.local>
Date: Fri, 8 Aug 2025 18:55:40 -0400
Subject: [PATCH 14/14] [MATLAB] Trimming whitespaces and trailing comments

---
 samples/matlab_experimental/crit_properites.m |  2 +-
 .../matlab_experimental/equations_of_state.m  |  1 -
 samples/matlab_experimental/shock_tube.m      | 25 ++++++++-----------
 samples/matlab_experimental/test_examples.m   |  4 +--
 samples/matlab_experimental/vapor_dome.m      |  2 +-
 5 files changed, 14 insertions(+), 20 deletions(-)

diff --git a/samples/matlab_experimental/crit_properites.m b/samples/matlab_experimental/crit_properites.m
index 2c367439018..0b8e3d20044 100644
--- a/samples/matlab_experimental/crit_properites.m
+++ b/samples/matlab_experimental/crit_properites.m
@@ -1,5 +1,5 @@
 %% Critical state properties
-% Print the critical state properties for the fluids for which Cantera has 
+% Print the critical state properties for the fluids for which Cantera has
 % built-in liquid/vapor equations of state.
 %
 % .. tags:: Matlab, thermodynamics, multiphase, non-ideal fluid
diff --git a/samples/matlab_experimental/equations_of_state.m b/samples/matlab_experimental/equations_of_state.m
index fbfb45b2122..c29d1670ac2 100644
--- a/samples/matlab_experimental/equations_of_state.m
+++ b/samples/matlab_experimental/equations_of_state.m
@@ -160,7 +160,6 @@ function plotEoS(p, ideal, rk, y_label)
 % Initialize matrices to hold densities
 density_ideal = zeros(length(p_arr), length(T_arr));
 density_RK = zeros(length(p_arr), length(T_arr));
-%density_CP = zeros(length(p_arr), length(T_arr));
 
 % Loop over pressure and temperature to compute densities
 for i = 1:length(p_arr)
diff --git a/samples/matlab_experimental/shock_tube.m b/samples/matlab_experimental/shock_tube.m
index 54fffe50044..aa546cbd4fa 100644
--- a/samples/matlab_experimental/shock_tube.m
+++ b/samples/matlab_experimental/shock_tube.m
@@ -1,11 +1,11 @@
 %% Shock-tube species profiles as a function of time
-% Simulate species profiles for a shock tube as a function of time, and 
-% observe the impact of incorporating the reduced-pressure linear mixture 
+% Simulate species profiles for a shock tube as a function of time, and
+% observe the impact of incorporating the reduced-pressure linear mixture
 % rule (LMR-R) in such calculations.
 %
 % Here we predict the H2O mole fraction time profiles for a mixture of 1163
-% ppm H2O2/1330 ppm H2O/665 ppm O2/20% CO2/Ar following reflected shock 
-% waves (1196 K, 2.127 atm) and compare results against the experimental 
+% ppm H2O2/1330 ppm H2O/665 ppm O2/20% CO2/Ar following reflected shock
+% waves (1196 K, 2.127 atm) and compare results against the experimental
 % measurements of Shao et al. [1] Two models are compared in this example:
 %
 % 1. A 2023 model of H2 and NH3 chemistry published by Alzueta et al. [2]
@@ -14,15 +14,15 @@
 % References:
 %
 % [1] J. Shao, R. Choudhary, D. F. Davidson, R. K. Hanson, Shock tube/laser
-% absorption measurement of the rate constant of the reaction: H2O2+CO2 = 
+% absorption measurement of the rate constant of the reaction: H2O2+CO2 =
 % 2OH+CO2, Proc. Combust. Inst. 39 (2023) 735 – 743.
 %
-% [2] M. U. Alzueta, I. Salas, H. Hashemi, P. Glarborg, CO-assisted NH3 
+% [2] M. U. Alzueta, I. Salas, H. Hashemi, P. Glarborg, CO-assisted NH3
 % oxidation, Combust. Flame 257 (2023) 112438.
 %
-% [3] P. J. Singal, J. Lee, L. Lei, R. L. Speth, M. P. Burke, 
-% Implementation of New Mixture Rules Has a Substantial Impact on 
-% Combustion Predictions for H2 and NH3, Proc. Combust. Inst. 40 (2024) 
+% [3] P. J. Singal, J. Lee, L. Lei, R. L. Speth, M. P. Burke,
+% Implementation of New Mixture Rules Has a Substantial Impact on
+% Combustion Predictions for H2 and NH3, Proc. Combust. Inst. 40 (2024)
 % 105779.
 %
 % .. tags:: Matlab, shock tube, kinetics, combustion
@@ -54,7 +54,7 @@
 for k = 1:length(model_names)
     model_key = model_names{k};
     mech_name = models.(model_key);
-    
+
     % Define mixture composition
     X_H2O2 = 1163e-6;
     X_H2O = 1330e-6;
@@ -74,8 +74,7 @@
     time = 0.0;
     estIgnDelay = 1.0; % seconds
     counter = 0;
-    dt = 1e-5;
-    
+
     time_data = [];
     H2O_X_data = [];
 
@@ -115,5 +114,3 @@
 
 
 toc
-
-
diff --git a/samples/matlab_experimental/test_examples.m b/samples/matlab_experimental/test_examples.m
index 4f593055885..1bcc1b10948 100644
--- a/samples/matlab_experimental/test_examples.m
+++ b/samples/matlab_experimental/test_examples.m
@@ -40,7 +40,6 @@
         end
     end
 
-
         % Summary report
     disp(' ');
     disp('============================');
@@ -61,6 +60,5 @@
     end
 
     disp(' ');
-    
-    
+
 end
diff --git a/samples/matlab_experimental/vapor_dome.m b/samples/matlab_experimental/vapor_dome.m
index c7a4b6e9fb4..f8eb21d0f05 100644
--- a/samples/matlab_experimental/vapor_dome.m
+++ b/samples/matlab_experimental/vapor_dome.m
@@ -50,7 +50,7 @@
     uf(i) = w.U/1e3;
     hf(i) = w.H/1e3;
     sf(i) = w.S/1e3;
-    
+
     % Pressure (same for both states)
     P(i) = w.P/1e5;
 end