diff --git a/interfaces/matlab_experimental/1D/AxisymmetricFlow.m b/interfaces/matlab_experimental/1D/AxisymmetricFlow.m
index c57986fc70..8391a7e032 100644
--- a/interfaces/matlab_experimental/1D/AxisymmetricFlow.m
+++ b/interfaces/matlab_experimental/1D/AxisymmetricFlow.m
@@ -1,7 +1,7 @@
 classdef AxisymmetricFlow < Domain1D
-    % Create an axisymmetric flow domain.
+    % Create an axisymmetric flow domain. ::
     %
-    % m = AxisymmetricFlow(gas, id)
+    %     >> m = AxisymmetricFlow(gas, id)
     %
     % :param gas:
     %     Instance of class :mat:class:`Solution`
@@ -13,10 +13,10 @@
 
     methods
 
-        % Constructor
         function m = AxisymmetricFlow(gas, id)
+            % Constructor
 
-            m = m@Domain1D('StagnationFlow', gas);
+            m@Domain1D('StagnationFlow', gas);
 
             if nargin == 1
                 m.setID('flow');
diff --git a/interfaces/matlab_experimental/1D/CounterFlowDiffusionFlame.m b/interfaces/matlab_experimental/1D/CounterFlowDiffusionFlame.m
index 50af94e498..c5b6153c52 100644
--- a/interfaces/matlab_experimental/1D/CounterFlowDiffusionFlame.m
+++ b/interfaces/matlab_experimental/1D/CounterFlowDiffusionFlame.m
@@ -1,6 +1,7 @@
-function flame = CounterFlowDiffusionFlame(left, flow, right, tp_f, tp_o, oxidizer)
-    % Create a counter flow diffusion flame stack.
-    % flame = CounterFlowDiffusionFlame(left, flow, right, tp_f, tp_o, oxidizer)
+classdef CounterFlowDiffusionFlame < Sim1D
+    % Create a counter flow diffusion flame stack. ::
+    %
+    %     >> flame = CounterFlowDiffusionFlame(left, flow, right, tp_f, tp_o, oxidizer)
     %
     % :param left:
     %     Object representing the left inlet, which must be
@@ -20,207 +21,222 @@
     % :param oxidizer:
     %     String representing the oxidizer species. Most commonly O2.
     % :return:
-    %     Instance of :mat:class:`Sim1D` object representing the left
-    %     inlet, flow, and right inlet.
+    %     Instance of :mat:class:`CounterFlowDiffusionFlame` object
+    %     representing the left inlet, flow, and right inlet.
     %
 
-    %% Check input parameters
+    methods
 
-    if nargin ~= 6
-        error('CounterFlowDiffusionFlame expects six input arguments.');
-    end
+        function flame = CounterFlowDiffusionFlame(left, flow, right, ...
+                                                   tp_f, tp_o, oxidizer)
 
-    if ~tp_f.isIdealGas
-        error('Fuel gas object must represent an ideal gas mixture.');
-    end
+            % Constructor
 
-    if ~tp_o.isIdealGas
-        error('Oxidizer gas object must represent an ideal gas mixture.');
-    end
 
-    if ~left.isInlet
-        error('Left inlet object of wrong type.');
-    end
+            %% Check input parameters
+            if nargin ~= 6
+                error('CounterFlowDiffusionFlame expects six input arguments.');
+            end
 
-    if ~flow.isFlow
-        error('Flow object of wrong type.');
-    end
+            if ~tp_f.isIdealGas
+                error('Fuel gas object must represent an ideal gas mixture.');
+            end
 
-    if ~right.isInlet
-        error('Right inlet object of wrong type.');
-    end
+            if ~tp_o.isIdealGas
+                error('Oxidizer gas object must represent an ideal gas mixture.');
+            end
 
-    if ~ischar(oxidizer)
-        error('Oxidizer name must be of format character.');
-    end
+            if ~left.isInlet
+                error('Left inlet object of wrong type.');
+            end
 
-    %% Get the density of both fuel and oxidizer streams.
-    % To be used in determining velocity of each stream. Also get the
-    % temperature of both inlet streams.
+            if ~flow.isFlow
+                error('Flow object of wrong type.');
+            end
 
-    rhof = tp_f.D;
-    rho0 = tp_o.D;
-    tf = left.T;
-    tox = right.T;
+            if ~right.isInlet
+                error('Right inlet object of wrong type.');
+            end
 
-    %% Find the species index of the oxidizer.
-    % To be used in determining initial strain rate.
+            if ~ischar(oxidizer)
+                error('Oxidizer name must be of format character.');
+            end
 
-    ioxidizer = tp_o.speciesIndex(oxidizer);
+            %% Define elMoles function
 
-    %% Calculate the stoichiometric mixture fraction.
-    % Needed for determining location of flame edges and composition.
-    % elMoles function used to calculate the number of moles of C, H, and O
-    % atoms in the fuel and oxidizer streams:
-    %
-    %    elMoles = elementalMassFraction/element atomic weight.
-    %
-    % From this, the stoichiometric Air/Fuel ratio can be determined.
-    % 1 Mole of O needs 2 Moles of C and 0.5 Moles of H for stoichiometric
-    % conditions. The stoichiometric mixture fraction, Zst, is then
-    % calculated.
-
-    sFuel = elMoles(tp_f, 'O') - 2 * elMoles(tp_f, 'C') - 0.5 * elMoles(tp_f, 'H');
-    sOx = elMoles(tp_o, 'O') - 2 * elMoles(tp_o, 'C') - 0.5 * elMoles(tp_o, 'H');
-    phi = sFuel / sOx;
-    zst = 1.0 / (1.0 - phi);
-
-    %% Compute the stoichiometric mass fractions of each species.
-    % Use this to set the fuel gas object and calculate adiabatic flame
-    % temperature and equilibrium composition.
-
-    spec = tp_f.speciesNames; % Get all of the species names in gas object.
-    nsp = tp_f.nSpecies; % Get total number of species in gas object.
-    % Get the current mass fractions of both fuel and inlet streams.
-    yox = tp_o.Y;
-    yf = tp_f.Y;
-    ystoich_double = zeros(1, nsp); % Create empty vector for stoich mass frac.
-
-    for n = 1:nsp
-        % Calculate stoichiometric mass fractions.
-        ystoich_double(n) = zst * yf(n) + (1.0 - zst) * yox(n);
-        % Convert mass fraction vector to string vector.
-        ystoich_str{n} = num2str(ystoich_double(n));
-        % Convert string vector to cell with SPECIES:MASS FRACTION format.
-        y_stoich{n} = [spec{n}, ':', ystoich_str{n}];
-    end
+            function moles = elMoles(tp, element)
+                % Determine the elemental moles in a gas object per unit mass.
 
-    % Initialize stoichiometric mass fraction cell with first SP:Y value.
-    ystoich = [y_stoich{1}];
+                % Check input parameters
+                if nargin ~= 2
+                    error('elMoles expects two input arguments.');
+                end
 
-    for i = 2:nsp
-        % Update cell to have format similar to N2:Yst,O2:Yst,...
-        ystoich = [ystoich ',', y_stoich{i}];
-    end
+                if ~tp.isIdealGas
+                    error('Gas object must represent an ideal gas mixture.');
+                end
 
-    % Set the fuel gas object as stoichiometric values and use equilibrate
-    % function to determine stoichiometric equilibrium temperature and mass
-    % fractions.
-    tp_f.TPY = {tf, tp_f.P, ystoich};
-    tp_f.equilibrate('HP');
-    teq = tp_f.T;
-    yeq = tp_f.Y;
-
-    %% Estimate the strain rate.
-    % Based on the inlet stream velocities and determine initial 'guess'
-    % for mixture fraction based on mass flux ratio.
-
-    zz = flow.gridPoints;
-    dz = zz(end) - zz(1);
-    mdotl = left.massFlux;
-    mdotr = right.massFlux;
-    uleft = mdotl / rhof;
-    uright = mdotr / rho0;
-    a = (abs(uleft) + abs(uright)) / dz;
-    diff = tp_f.mixDiffCoeffs;
-    f = sqrt(a / (2.0 * diff(ioxidizer)));
-    x0num = sqrt(uleft * mdotl) * dz;
-    x0den = sqrt(uleft * mdotr) + sqrt(uright * mdotr);
-    x0 = x0num / x0den;
-
-    %% Calculate initial values of temperature and mass fractions.
-    % These values to be used for energy equation solution. Method is based
-    % on the Burke-Schumann model.
-
-    nz = flow.nPoints;
-    zm = zeros(1, nz);
-    u = zeros(1, nz);
-    v = zeros(1, nz);
-    y = zeros(nz, nsp);
-    t = zeros(1, nz);
-
-    for j = 1:nz
-        x = zz(j);
-        zeta = f * (x - x0);
-        zmix = 0.5 * (1.0 - erf(zeta)); % Mixture fraction in flame.
-        zm(j) = zmix;
-        u(j) = a * (x0 - zz(j)); % Axial velocity.
-        v(j) = a; % Radial velocity.
-
-        if zmix > zst
+                if ~ischar(element)
+                    error('Element name must be of format character.');
+                end
 
-            for n = 1:nsp
-                y(j, n) = yeq(n) + (zmix - zst) * (yf(n) - yeq(n)) / (1.0 - zst);
-            end
+                % Calculate the moles per mass of mixture of an element within a gas
+                % object. The equation used is:
+                %
+                %    elMoles = elMassFrac/Mel
+                %
+                % where elMassFrac is the elemental mass fraction within the gas object
+                % using the elementalMassFraction function; Mel is the atomic mass of
+                % the element.
 
-            t(j) = teq + (tf - teq) * (zmix - zst) / (1.0 - zst);
-        else
+                elMassFrac = tp.elementalMassFraction(element);
+                eli = tp.elementIndex(element);
+                M = tp.atomicMasses;
+                Mel = M(eli);
+                moles = elMassFrac / Mel;
 
-            for n = 1:nsp
-                y(j, n) = yox(n) + zmix * (yeq(n) - yox(n)) / zst;
             end
 
-            t(j) = tox + zmix * (teq - tox) / zst;
-        end
+            %% Get the density of both fuel and oxidizer streams.
+            % To be used in determining velocity of each stream. Also get the
+            % temperature of both inlet streams.
+
+            rhof = tp_f.D;
+            rho0 = tp_o.D;
+            tf = left.T;
+            tox = right.T;
+
+            %% Find the species index of the oxidizer.
+            % To be used in determining initial strain rate.
+
+            ioxidizer = tp_o.speciesIndex(oxidizer);
+
+            %% Calculate the stoichiometric mixture fraction.
+            % Needed for determining location of flame edges and composition.
+            % elMoles function used to calculate the number of moles of C, H, and O
+            % atoms in the fuel and oxidizer streams:
+            %
+            %    elMoles = elementalMassFraction/element atomic weight.
+            %
+            % From this, the stoichiometric Air/Fuel ratio can be determined.
+            % 1 Mole of O needs 2 Moles of C and 0.5 Moles of H for stoichiometric
+            % conditions. The stoichiometric mixture fraction, Zst, is then
+            % calculated.
+
+            sFuel = elMoles(tp_f, 'O') - 2 * elMoles(tp_f, 'C') - ...
+                    0.5 * elMoles(tp_f, 'H');
+            sOx = elMoles(tp_o, 'O') - 2 * elMoles(tp_o, 'C') - ...
+                    0.5 * elMoles(tp_o, 'H');
+            phi = sFuel / sOx;
+            zst = 1.0 / (1.0 - phi);
+
+            %% Compute the stoichiometric mass fractions of each species.
+            % Use this to set the fuel gas object and calculate adiabatic flame
+            % temperature and equilibrium composition.
+
+            spec = tp_f.speciesNames; % Get all of the species names in gas object.
+            nsp = tp_f.nSpecies; % Get total number of species in gas object.
+            % Get the current mass fractions of both fuel and inlet streams.
+            yox = tp_o.Y;
+            yf = tp_f.Y;
+            ystoich_double = zeros(1, nsp); % Create empty vector for stoich mass frac.
 
-    end
+            for n = 1:nsp
+                % Calculate stoichiometric mass fractions.
+                ystoich_double(n) = zst * yf(n) + (1.0 - zst) * yox(n);
+                % Convert mass fraction vector to string vector.
+                ystoich_str{n} = num2str(ystoich_double(n));
+                % Convert string vector to cell with SPECIES:MASS FRACTION format.
+                y_stoich{n} = [spec{n}, ':', ystoich_str{n}];
+            end
 
-    zrel = zz / dz;
+            % Initialize stoichiometric mass fraction cell with first SP:Y value.
+            ystoich = [y_stoich{1}];
 
-    %% Create the flame stack.
-    % Set the profile of the flame with the estimated axial velocities,
-    % radial velocities, temperature, and mass fractions calculated above.
-    flame = Sim1D([left flow right]);
-    flame.setProfile(2, {'velocity', 'spread_rate'}, [zrel; u; v]);
-    flame.setProfile(2, 'T', [zrel; t]);
+            for i = 2:nsp
+                % Update cell to have format similar to N2:Yst,O2:Yst,...
+                ystoich = [ystoich ',', y_stoich{i}];
+            end
 
-    for n = 1:nsp
-        nm = tp_f.speciesName(n);
-        flame.setProfile(2, nm, [zrel; transpose(y(:, n))])
-    end
+            % Set the fuel gas object as stoichiometric values and use equilibrate
+            % function to determine stoichiometric equilibrium temperature and mass
+            % fractions.
+            tp_f.TPY = {tf, tp_f.P, ystoich};
+            tp_f.equilibrate('HP');
+            teq = tp_f.T;
+            yeq = tp_f.Y;
+
+            %% Estimate the strain rate.
+            % Based on the inlet stream velocities and determine initial 'guess'
+            % for mixture fraction based on mass flux ratio.
+
+            zz = flow.gridPoints;
+            dz = zz(end) - zz(1);
+            mdotl = left.massFlux;
+            mdotr = right.massFlux;
+            uleft = mdotl / rhof;
+            uright = mdotr / rho0;
+            a = (abs(uleft) + abs(uright)) / dz;
+            diff = tp_f.mixDiffCoeffs;
+            f = sqrt(a / (2.0 * diff(ioxidizer)));
+            x0num = sqrt(uleft * mdotl) * dz;
+            x0den = sqrt(uleft * mdotr) + sqrt(uright * mdotr);
+            x0 = x0num / x0den;
+
+            %% Calculate initial values of temperature and mass fractions.
+            % These values to be used for energy equation solution. Method is based
+            % on the Burke-Schumann model.
+
+            nz = flow.nPoints;
+            zm = zeros(1, nz);
+            u = zeros(1, nz);
+            v = zeros(1, nz);
+            y = zeros(nz, nsp);
+            t = zeros(1, nz);
+
+            for j = 1:nz
+                x = zz(j);
+                zeta = f * (x - x0);
+                zmix = 0.5 * (1.0 - erf(zeta)); % Mixture fraction in flame.
+                zm(j) = zmix;
+                u(j) = a * (x0 - zz(j)); % Axial velocity.
+                v(j) = a; % Radial velocity.
+
+                if zmix > zst
+
+                    for n = 1:nsp
+                        y(j, n) = yeq(n) + (zmix - zst) * (yf(n) - yeq(n)) /...
+                                 (1.0 - zst);
+                    end
+
+                    t(j) = teq + (tf - teq) * (zmix - zst) / (1.0 - zst);
+                else
+
+                    for n = 1:nsp
+                        y(j, n) = yox(n) + zmix * (yeq(n) - yox(n)) / zst;
+                    end
+
+                    t(j) = tox + zmix * (teq - tox) / zst;
+                end
 
-end
+            end
 
-%% Define elMoles function
+            zrel = zz / dz;
 
-function moles = elMoles(tp, element)
-    % Determine the elemental moles in a gas object per unit mass.
+            %% Create the flame stack.
+            % Set the profile of the flame with the estimated axial velocities,
+            % radial velocities, temperature, and mass fractions calculated above.
+            flame@Sim1D({left flow right});
+            flame.setProfile(2, {'velocity', 'spread_rate'}, [zrel; u; v]);
+            flame.setProfile(2, 'T', [zrel; t]);
 
-    % Check input parameters
-    if nargin ~= 2
-        error('elMoles expects two input arguments.');
-    end
+            for n = 1:nsp
+                nm = tp_f.speciesName(n);
+                flame.setProfile(2, nm, [zrel; transpose(y(:, n))]);
+            end
 
-    if ~tp.isIdealGas
-        error('Gas object must represent an ideal gas mixture.');
-    end
+        end
 
-    if ~ischar(element)
-        error('Element name must be of format character.');
     end
 
-    % Calculate the moles per mass of mixture of an element within a gas
-    % object. The equation used is:
-    %
-    %    elMoles = elMassFrac/Mel
-    %
-    % where elMassFrac is the elemental mass fraction within the gas object
-    % using the elementalMassFraction function; Mel is the atomic mass of
-    % the element.
-
-    elMassFrac = tp.elementalMassFraction(element);
-    eli = tp.elementIndex(element);
-    M = tp.atomicMasses;
-    Mel = M(eli);
-    moles = elMassFrac / Mel;
 end
diff --git a/interfaces/matlab_experimental/1D/Domain1D.m b/interfaces/matlab_experimental/1D/Domain1D.m
index 9f22522e5e..917daf7d7f 100644
--- a/interfaces/matlab_experimental/1D/Domain1D.m
+++ b/interfaces/matlab_experimental/1D/Domain1D.m
@@ -1,7 +1,7 @@
 classdef Domain1D < handle
-    % Domain1D Class
+    % Domain1D Class ::
     %
-    % d = Domain1D(a, b)
+    %     >> d = Domain1D(a, b)
     %
     % :param a:
     %    String type of domain. Possible values are:
@@ -122,10 +122,13 @@
     methods
         %% Domain1D Class Constructor.
 
-        function d = Domain1D(a, b)
+        function d = Domain1D(varargin)
+            % Create a :mat:class:`Domain1D` object.
+
             checklib;
 
             if nargin == 1
+                a = varargin{1};
 
                 if strcmp(a, 'Inlet1D')
                     d.domainID = callct('inlet_new');
@@ -142,7 +145,9 @@
                 end
 
             elseif nargin == 2
-                % a stagnation flow
+                a = varargin{1};
+                b = varargin{2};
+
                 if strcmp(a, 'StagnationFlow')
 
                     if isa(b, 'Solution')
@@ -183,7 +188,7 @@
         %% Domain1D Class Destructor
 
         function delete(d)
-            % Delete the C++ Domain1D object.
+            % Delete the :mat:class:`Domain1D` object.
 
             callct('domain_del', d.domainID);
         end
diff --git a/interfaces/matlab_experimental/1D/FreeFlame.m b/interfaces/matlab_experimental/1D/FreeFlame.m
index e9018b77c8..5aa26e4985 100644
--- a/interfaces/matlab_experimental/1D/FreeFlame.m
+++ b/interfaces/matlab_experimental/1D/FreeFlame.m
@@ -1,7 +1,7 @@
 classdef FreeFlame < Domain1D
-    % Create a freely-propagating flat flame.
+    % Create a freely-propagating flat flame. ::
     %
-    % m = FreeFlame(gas, id)
+    %     >> m = FreeFlame(gas, id)
     %
     % :param gas:
     %     Instance of class :mat:class:`Solution`
@@ -13,10 +13,10 @@
     %
     methods
 
-        % Constructor
         function m = FreeFlame(gas, id)
+            % Constructor
 
-            m = m@Domain1D('StagnationFlow', gas, 2);
+            m@Domain1D('StagnationFlow', gas, 2);
 
             if nargin == 1
                 m.setID('flame');
diff --git a/interfaces/matlab_experimental/1D/Inlet.m b/interfaces/matlab_experimental/1D/Inlet.m
index f8ae095e87..8269843461 100644
--- a/interfaces/matlab_experimental/1D/Inlet.m
+++ b/interfaces/matlab_experimental/1D/Inlet.m
@@ -1,7 +1,7 @@
 classdef Inlet < Domain1D
-    % Create an inlet domain.
+    % Create an inlet domain. ::
     %
-    % m = Inlet(id)
+    %     >> m = Inlet(id)
     %
     % Note that an inlet can only be a terminal domain - it must be
     % either the leftmost or rightmost domain in a stack.
@@ -14,10 +14,10 @@
 
     methods
 
-        % Constructor
         function m = Inlet(id)
+            % Constructor
 
-            m = m@Domain1D('Inlet1D');
+            m@Domain1D('Inlet1D');
 
             if nargin == 0
                 m.setID('inlet');
diff --git a/interfaces/matlab_experimental/1D/Outlet.m b/interfaces/matlab_experimental/1D/Outlet.m
index 45bb4fa292..4032cdde89 100644
--- a/interfaces/matlab_experimental/1D/Outlet.m
+++ b/interfaces/matlab_experimental/1D/Outlet.m
@@ -1,7 +1,7 @@
 classdef Outlet < Domain1D
-    % Create an outlet domain.
+    % Create an outlet domain. ::
     %
-    % m = Outlet(id)
+    %     >> m = Outlet(id)
     %
     % :param id:
     %     String ID of the outlet.
@@ -11,10 +11,10 @@
 
     methods
 
-        % Constructor
         function m = Outlet(id)
+            % Constructor
 
-            m = m@Domain1D('Outlet1D');
+            m@Domain1D('Outlet1D');
 
             if nargin == 0
                 m.setID('outlet');
diff --git a/interfaces/matlab_experimental/1D/OutletRes.m b/interfaces/matlab_experimental/1D/OutletRes.m
index 32aa19b8db..e84751236d 100644
--- a/interfaces/matlab_experimental/1D/OutletRes.m
+++ b/interfaces/matlab_experimental/1D/OutletRes.m
@@ -1,7 +1,7 @@
 classdef OutletRes < Domain1D
-    % Create an outlet reservoir domain.
+    % Create an outlet reservoir domain. ::
     %
-    % m = OutletRes(id)
+    %     >> m = OutletRes(id)
     %
     % :param id:
     %     String ID of the outlet reservoir.
@@ -11,10 +11,10 @@
 
     methods
 
-        % Constructor
         function m = OutletRes(id)
+            % Constructor
 
-            m = m@Domain1D('OutletRes');
+            m@Domain1D('OutletRes');
 
             if nargin == 0
                 m.setID('outletres');
diff --git a/interfaces/matlab_experimental/1D/Sim1D.m b/interfaces/matlab_experimental/1D/Sim1D.m
index c0bedba5cd..db95081568 100644
--- a/interfaces/matlab_experimental/1D/Sim1D.m
+++ b/interfaces/matlab_experimental/1D/Sim1D.m
@@ -1,7 +1,7 @@
 classdef Sim1D < handle
-    % Sim1D Class
+    % Sim1D Class. ::
     %
-    % s = Sim1D(domains)
+    %     >> s = Sim1D(domains)
     %
     % A Sim1D object is a container for one-dimensional domains,
     % which are instances of class Domain1D. The domains are of two
@@ -10,9 +10,9 @@
     % See also: :mat:class:`Domain1D`
     %
     % :param domains:
-    %     Vector of domain instances
+    %     Cell array of instances of :mat:class:`Domain1D` and its subclasses.
     % :return:
-    %     Instance of class :mat:class:`Sim1D`
+    %     Instance of class :mat:class:`Sim1D`.
     %
 
     properties (SetAccess = immutable)
@@ -24,9 +24,12 @@
     end
 
     methods
+
         %% Sim1D Class Constructor
 
         function s = Sim1D(domains)
+            % Create a :mat:class:`Sim1D` object.
+
             checklib;
 
             s.stID = -1;
@@ -41,7 +44,7 @@
             ids = zeros(1, nd);
 
             for n = 1:nd
-                ids(n) = domains(n).domainID;
+                ids(n) = domains{n}.domainID;
             end
 
             s.stID = callct('sim1D_new', nd, ids);
@@ -51,7 +54,7 @@
         %% Sim1D Class Destructor
 
         function delete(s)
-            % Delete the C++ Sim1D object.
+            % Delete the :mat:class:`Sim1D` object.
 
             callct('sim1D_del', s.stID);
         end
@@ -83,7 +86,7 @@ function plotSolution(s, domain, component)
             %
 
             n = s.stackIndex(domain);
-            d = s.domains(n);
+            d = s.domains{n};
             z = d.gridPoints;
             x = s.solution(domain, component);
             plot(z, x);
@@ -162,7 +165,7 @@ function saveSoln(s, fname, id, desc)
             %
 
             idom = s.stackIndex(domain);
-            d = s.domains(idom);
+            d = s.domains{idom};
             np = d.nPoints;
 
             if nargin == 3
@@ -275,7 +278,7 @@ function getInitialSoln(s)
             %     The grid in domain name
 
             n = s.stackIndex(name);
-            d = s.domains(n);
+            d = s.domains{n};
             z = d.gridPoints;
         end
 
@@ -298,7 +301,7 @@ function getInitialSoln(s)
             end
 
             idom = s.stackIndex(domain);
-            d = s.domains(idom);
+            d = s.domains{idom};
 
             nc = d.nComponents;
             np = d.nPoints;
@@ -428,7 +431,7 @@ function setProfile(s, name, comp, p)
                 n = s.domainIndex(name);
             end
 
-            d = s.domains(n);
+            d = s.domains{n};
 
             if isa(comp, 'double') || isa(comp, 'cell')
                 c = comp;
diff --git a/interfaces/matlab_experimental/1D/Surface.m b/interfaces/matlab_experimental/1D/Surface.m
index 74e6172fce..020b1dd37a 100644
--- a/interfaces/matlab_experimental/1D/Surface.m
+++ b/interfaces/matlab_experimental/1D/Surface.m
@@ -1,7 +1,7 @@
 classdef Surface < Domain1D
-    % Create a surface domain.
+    % Create a surface domain. ::
     %
-    % m = Surface(id, surface_mech)
+    %     >> m = Surface(id, surface_mech)
     %
     % :param id:
     %     String ID of surface
@@ -15,17 +15,16 @@
 
     methods
 
-        % Constructor
         function m = Surface(id, surface_mech)
+            % Constructor
 
-            if nargin > 1
-                m = m@Domain1D('ReactingSurface', surface_mech);
-                m.setID(id);
-                return
+            if nargin < 2
+                param = {'Surf1D'};
+            else
+                param = {'ReactingSurface', surface_mech};
             end
 
-            m = m@Domain1D('Surf1D');
-
+            m@Domain1D(param{:});
             if nargin == 0
                 m.setID('surface');
             elseif nargin == 1
diff --git a/interfaces/matlab_experimental/1D/SymmPlane.m b/interfaces/matlab_experimental/1D/SymmPlane.m
index 6005bd5cec..09ea3bc311 100644
--- a/interfaces/matlab_experimental/1D/SymmPlane.m
+++ b/interfaces/matlab_experimental/1D/SymmPlane.m
@@ -1,7 +1,7 @@
 classdef SymmPlane < Domain1D
-    % Create a symmetry plane domain.
+    % Create a symmetry plane domain. ::
     %
-    % m = SymmPlane(id)
+    %     >> m = SymmPlane(id)
     %
     % :param id:
     %     String ID of the symmetry plane.
@@ -11,10 +11,10 @@
 
     methods
 
-        % Constructor
         function m = SymmPlane(id)
+            % Constructor
 
-            m = Domain1D(Symm1D);
+            m@Domain1D('Symm1D');
 
             if nargin == 0
                 m.setID('symmetry_plane');
diff --git a/interfaces/matlab_experimental/Base/CarbonDioxide.m b/interfaces/matlab_experimental/Base/CarbonDioxide.m
index 7377d54e51..6b16653949 100644
--- a/interfaces/matlab_experimental/Base/CarbonDioxide.m
+++ b/interfaces/matlab_experimental/Base/CarbonDioxide.m
@@ -1,7 +1,7 @@
 function c = CarbonDioxide()
-    % Return an object representing carbon dioxide.
+    % Return an object representing carbon dioxide. ::
     %
-    % c = CarbonDioxide
+    %     >> c = CarbonDioxide
     %
     % The object returned by this method implements an accurate equation of
     % state for carbon dioxide that can be used in the liquid, vapor, saturated
diff --git a/interfaces/matlab_experimental/Base/HFC134a.m b/interfaces/matlab_experimental/Base/HFC134a.m
index 2beef6bdf9..c463a2bcd9 100644
--- a/interfaces/matlab_experimental/Base/HFC134a.m
+++ b/interfaces/matlab_experimental/Base/HFC134a.m
@@ -1,7 +1,7 @@
 function h = HFC134a()
-    % Return an object representing refrigerant HFC134a.
+    % Return an object representing refrigerant HFC134a. ::
     %
-    % h = HFC134a()
+    %     >> h = HFC134a()
     %
     % The object returned by this method implements an accurate equation of
     % state for refrigerant HFC134a (R134a) that can be used in the liquid,
diff --git a/interfaces/matlab_experimental/Base/Heptane.m b/interfaces/matlab_experimental/Base/Heptane.m
index 1efcd95194..200ffcb107 100644
--- a/interfaces/matlab_experimental/Base/Heptane.m
+++ b/interfaces/matlab_experimental/Base/Heptane.m
@@ -1,7 +1,7 @@
 function h = Heptane
-    % Return an object representing n-heptane.
+    % Return an object representing n-heptane. ::
     %
-    % h = Heptane()
+    %     >> h = Heptane()
     %
     % The object returned by this method implements an accurate equation of
     % state for n-heptane that can be used in the liquid, vapor, saturated
diff --git a/interfaces/matlab_experimental/Base/Hydrogen.m b/interfaces/matlab_experimental/Base/Hydrogen.m
index 6aa03ea51f..51f1b5fc09 100644
--- a/interfaces/matlab_experimental/Base/Hydrogen.m
+++ b/interfaces/matlab_experimental/Base/Hydrogen.m
@@ -1,7 +1,7 @@
 function h = Hydrogen()
-    % Return an object representing hydrogen.
+    % Return an object representing hydrogen. ::
     %
-    % h = Hydrogen()
+    %     >> h = Hydrogen()
     %
     % The object returned by this method implements an accurate equation of
     % state for hydrogen that can be used in the liquid, vapor, saturated
diff --git a/interfaces/matlab_experimental/Base/Interface.m b/interfaces/matlab_experimental/Base/Interface.m
index f5fe5aa766..7ec47f113c 100644
--- a/interfaces/matlab_experimental/Base/Interface.m
+++ b/interfaces/matlab_experimental/Base/Interface.m
@@ -1,7 +1,7 @@
 classdef Interface < handle & ThermoPhase & Kinetics
-    % Interface Class
+    % Interface Class ::
     %
-    % s = Interface(src, id, p1, p2, p3, p4)
+    %     >> s = Interface(src, id, p1, p2, p3, p4)
     %
     % See `ideal-surface <https://cantera.org/documentation/docs-2.6/sphinx/html/yaml/phases.html#sec-yaml-ideal-surface>`__
     % and `Declaring adjacent phases <https://cantera.org/tutorials/yaml/phases.html#declaring-adjacent-phases>`__.
@@ -38,7 +38,8 @@
         %% Interface Class Constructor
 
         function s = Interface(src, id, p1, p2, p3, p4)
-            % Interface class constructor
+            % Create an :mat:class:`Interface` object.
+
             checklib;
 
             t = ThermoPhase(src, id);
diff --git a/interfaces/matlab_experimental/Base/Kinetics.m b/interfaces/matlab_experimental/Base/Kinetics.m
index 9858238915..13f2ac95e2 100644
--- a/interfaces/matlab_experimental/Base/Kinetics.m
+++ b/interfaces/matlab_experimental/Base/Kinetics.m
@@ -73,16 +73,17 @@
         %% Kinetics Class Constructor
 
         function kin = Kinetics(ph, src, id, n1, n2, n3, n4)
-            % Kinetics Class
+            % Kinetics Class ::
             %
-            % k = Kinetics(ph, neighbor1, neighbor2, neighbor3, neighbor4)
+            %     >> k = Kinetics(ph, neighbor1, neighbor2, neighbor3, neighbor4)
             %
             % Class Kinetics represents kinetics managers, which are classes that manage
-            % reaction mechanisms.  The reaction mechanism attributes are specified in a YAML file.
+            % reaction mechanisms. The reaction mechanism attributes are specified in a
+            % YAML file.
             %
-            % Instances of class :mat:class:`Kinetics` are responsible for evaluating reaction rates
-            % of progress, species production rates, and other quantities pertaining to
-            % a reaction mechanism.
+            % Instances of class :mat:class:`Kinetics` are responsible for evaluating
+            % reaction rates of progress, species production rates, and other
+            % quantities pertaining to a reaction mechanism.
             %
             % :param ph:
             %     An instance of class :mat:class:`ThermoPhase` representing the phase
@@ -92,20 +93,21 @@
             % :param id:
             %     ID of the phase to import as specified in the input file. (optional)
             % :param neighbor1:
-            %     Instance of class :mat:class:`ThermoPhase` or :mat:class:`Solution` representing a
-            %     neighboring phase.
+            %     Instance of class :mat:class:`ThermoPhase` or :mat:class:`Solution`
+            %     representing a neighboring phase.
             % :param neighbor2:
-            %     Instance of class :mat:class:`ThermoPhase` or :mat:class:`Solution` representing a
-            %     neighboring phase.
+            %     Instance of class :mat:class:`ThermoPhase` or :mat:class:`Solution`
+            %     representing a neighboring phase.
             % :param neighbor3:
-            %     Instance of class :mat:class:`ThermoPhase` or :mat:class:`Solution` representing a
-            %     neighboring phase.
+            %     Instance of class :mat:class:`ThermoPhase` or :mat:class:`Solution`
+            %     representing a neighboring phase.
             % :param neighbor4:
-            %     Instance of class :mat:class:`ThermoPhase` or :mat:class:`Solution` representing a
-            %     neighboring phase.
+            %     Instance of class :mat:class:`ThermoPhase` or :mat:class:`Solution`
+            %     representing a neighboring phase.
             % :return:
             %      Instance of class :mat:class:`Kinetics`
             %
+
             checklib;
 
             % indices for bulk phases in a heterogeneous mechanism
@@ -145,7 +147,7 @@
         %% Kinetics Class Destructor
 
         function delete(kin)
-            % Delete the kernel object
+            % Delete the :mat:class:`Kinetics` object.
 
             callct('kin_del', kin.kinID);
         end
diff --git a/interfaces/matlab_experimental/Base/Methane.m b/interfaces/matlab_experimental/Base/Methane.m
index 88d7696a37..8bd13b261d 100644
--- a/interfaces/matlab_experimental/Base/Methane.m
+++ b/interfaces/matlab_experimental/Base/Methane.m
@@ -1,7 +1,7 @@
 function m = Methane()
-    % Return an object representing methane.
+    % Return an object representing methane. ::
     %
-    % m = Methane()
+    %     >> m = Methane()
     %
     % The object returned by this method implements an accurate equation of
     % state for methane that can be used in the liquid, vapor, saturated
diff --git a/interfaces/matlab_experimental/Base/Nitrogen.m b/interfaces/matlab_experimental/Base/Nitrogen.m
index c25f074e6a..600b988fa9 100644
--- a/interfaces/matlab_experimental/Base/Nitrogen.m
+++ b/interfaces/matlab_experimental/Base/Nitrogen.m
@@ -1,7 +1,7 @@
 function n = Nitrogen()
-    % Return an object representing nitrogen.
+    % Return an object representing nitrogen. ::
     %
-    % n = Nitrogen()
+    %     >> n = Nitrogen()
     %
     % The object returned by this method implements an accurate equation of
     % state for nitrogen that can be used in the liquid, vapor, saturated
diff --git a/interfaces/matlab_experimental/Base/Oxygen.m b/interfaces/matlab_experimental/Base/Oxygen.m
index b4c7e51ae0..8c48c90c9c 100644
--- a/interfaces/matlab_experimental/Base/Oxygen.m
+++ b/interfaces/matlab_experimental/Base/Oxygen.m
@@ -1,7 +1,7 @@
 function o = Oxygen()
-    % Return an object representing oxygen.
+    % Return an object representing oxygen. ::
     %
-    % o = Oxygen()
+    %     >> o = Oxygen()
     %
     % The object returned by this method implements an accurate equation of
     % state for oxygen that can be used in the liquid, vapor, saturated
diff --git a/interfaces/matlab_experimental/Base/Solution.m b/interfaces/matlab_experimental/Base/Solution.m
index c47034c085..036a1e447c 100644
--- a/interfaces/matlab_experimental/Base/Solution.m
+++ b/interfaces/matlab_experimental/Base/Solution.m
@@ -1,7 +1,7 @@
 classdef Solution < handle & ThermoPhase & Kinetics & Transport
-    % Solution Class
+    % Solution Class ::
     %
-    % s = Solution(src, id, trans)
+    %     >> s = Solution(src, id, trans)
     %
     % Class :mat:class:`Solution` represents solutions of multiple species. A
     % solution is defined as a mixture of two or more constituents
@@ -54,6 +54,8 @@
         %% Solution Class Constructor
 
         function s = Solution(src, id, trans)
+            % Create a :mat:class:`Solution` object.
+
             checklib;
 
             if nargin < 2 || nargin > 3
@@ -77,7 +79,7 @@
         %% Solution Class Destructor
 
         function delete(s)
-            % Delete the kernel objects associated with a Solution.
+            % Delete :mat:class:`Solution` object.
 
             disp('Solution class object has been deleted');
         end
diff --git a/interfaces/matlab_experimental/Base/ThermoPhase.m b/interfaces/matlab_experimental/Base/ThermoPhase.m
index 94185e6f5b..bac9e542ba 100644
--- a/interfaces/matlab_experimental/Base/ThermoPhase.m
+++ b/interfaces/matlab_experimental/Base/ThermoPhase.m
@@ -1,14 +1,14 @@
 classdef ThermoPhase < handle
-    % ThermoPhase Class
+    % ThermoPhase Class ::
     %
-    % t = ThermoPhase(src, id)
+    %     >> t = ThermoPhase(src, id)
     %
     % :param src:
     %     Input string of YAML file name.
     % :param id:
     %     ID of the phase to import as specified in the input file.
     % :return:
-    %     Instance of class :mat:class:`ThermoPhase`
+    %     Instance of class :mat:class:`ThermoPhase`.
     %
 
     properties (SetAccess = public)
@@ -251,6 +251,7 @@
         %% ThermoPhase Class Constructor
 
         function tp = ThermoPhase(src, id)
+            % Create a :mat:class:`ThermoPhase` object.
             checklib;
 
             if nargin ~= 2
@@ -264,7 +265,7 @@
         %% ThermoPhase Class Destructor
 
         function tpClear(tp)
-            % Delete the ThermoPhase object.
+            % Delete the :mat:class:`ThermoPhase` object.
 
             callct('thermo_del', tp.tpID);
         end
diff --git a/interfaces/matlab_experimental/Base/Transport.m b/interfaces/matlab_experimental/Base/Transport.m
index 72e0d6e760..5d4af80b63 100644
--- a/interfaces/matlab_experimental/Base/Transport.m
+++ b/interfaces/matlab_experimental/Base/Transport.m
@@ -30,12 +30,12 @@
         %% Transport Class Constructor
 
         function tr = Transport(tp, model, loglevel)
-            % Transport Class
+            % Transport Class ::
             %
-            % tr = Transport(r, th, model, loglevel)
+            %     >> tr = Transport(r, th, model, loglevel)
             %
-            % Create a new instance of class :mat:class:`Transport`. One to three arguments
-            % may be supplied. The first must be an instance of class
+            % Create a new instance of class :mat:class:`Transport`. One to three
+            % arguments may be supplied. The first must be an instance of class
             % :mat:class:`ThermoPhase`. The second (optional) argument is the type of
             % model desired, specified by the string ``'default'``, ``'Mix'`` or
             % ``'Multi'``. ``'default'`` uses the default transport specified in the
@@ -50,8 +50,9 @@
             % :param loglevel:
             %     Level of diagnostic logging. Default if not specified is 4.
             % :return:
-            %     Instance of class :mat:class:`Transport`
+            %     Instance of class :mat:class:`Transport`.
             %
+
             checklib;
             tr.trID = 0;
 
@@ -82,7 +83,7 @@
         %% Transport Class Destructor
 
         function delete(tr)
-            % Delete the kernel object.
+            % Delete the :mat:class:`Transport` object.
 
             callct('trans_del', tr.trID);
         end
diff --git a/interfaces/matlab_experimental/Base/Water.m b/interfaces/matlab_experimental/Base/Water.m
index f3ad5c1784..388170ffbd 100644
--- a/interfaces/matlab_experimental/Base/Water.m
+++ b/interfaces/matlab_experimental/Base/Water.m
@@ -1,7 +1,7 @@
 function w = Water()
-    % Return an object representing water.
+    % Return an object representing water. ::
     %
-    % w = Water()
+    %     >> w = Water()
     %
     % The object returned by this method implements an accurate equation of
     % state for water that can be used in the liquid, vapor, saturated
diff --git a/interfaces/matlab_experimental/Constants/faradayconstant.m b/interfaces/matlab_experimental/Constants/faradayconstant.m
index d79841d480..e64ff31342 100644
--- a/interfaces/matlab_experimental/Constants/faradayconstant.m
+++ b/interfaces/matlab_experimental/Constants/faradayconstant.m
@@ -1,6 +1,7 @@
 function r = faradayconstant
-    % Get the Faraday constant in C/kmol of electron.
-    % r = faradayconstant
+    % Get the Faraday constant in C/kmol of electron. ::
+    %
+    %     >> r = faradayconstant
     %
     % :return:
     %     The Faraday constant in C/kmol-e.
diff --git a/interfaces/matlab_experimental/Constants/gasconstant.m b/interfaces/matlab_experimental/Constants/gasconstant.m
index e610d2b102..b4235b5e1d 100644
--- a/interfaces/matlab_experimental/Constants/gasconstant.m
+++ b/interfaces/matlab_experimental/Constants/gasconstant.m
@@ -1,4 +1,10 @@
 function r = gasconstant
-    % Universal gas constant in J/kmol-K.
+    % Get the universal gas constant in J/kmol-K. ::
+    %
+    %     >> r = gasconstant
+    %
+    % :return:
+    %     The universal gas constant in J/kmol-K.
+    %
     r = 8314.4621;
 end
diff --git a/interfaces/matlab_experimental/Constants/oneatm.m b/interfaces/matlab_experimental/Constants/oneatm.m
index e8b9eec9b8..ab2104b949 100644
--- a/interfaces/matlab_experimental/Constants/oneatm.m
+++ b/interfaces/matlab_experimental/Constants/oneatm.m
@@ -1,4 +1,10 @@
 function p = oneatm
-    % Get 1 atm in Pa.
+    % Get one atmosphere in Pa. ::
+    %
+    %     >> p = oneatm
+    %
+    % :return:
+    %     One atmosphere in Pascals.
+    %
     p = 101325.0;
 end
diff --git a/interfaces/matlab_experimental/Func/Func.m b/interfaces/matlab_experimental/Func/Func.m
index 2cb73f5f62..94a15c63ca 100644
--- a/interfaces/matlab_experimental/Func/Func.m
+++ b/interfaces/matlab_experimental/Func/Func.m
@@ -1,57 +1,4 @@
 classdef Func < handle
-    % Func Class
-    %
-    % x = Func(typ, n, p)
-    %
-    % A class for functors.
-    % A functor is an object that behaves like a function. Cantera
-    % defines a set of functors to use to create arbitrary functions to
-    % specify things like heat fluxes, piston speeds, etc., in reactor
-    % network simulations. Of course, they can be used for other things
-    % too.
-    %
-    % The main feature of a functor class is that it overloads the ``()``
-    % operator to evaluate the function. For example, suppose object
-    % ``f`` is a functor that evaluates the polynomial :math:`2x^2 - 3x + 1`.
-    % Then writing ``f(2)`` would cause the method that evaluates the
-    % function to be invoked, and would pass it the argument ``2``. The
-    % return value would of course be 3.
-    %
-    % The types of functors you can create in Cantera are these:
-    %
-    % 1. A polynomial
-    % 2. A Fourier series
-    % 3. A sum of Arrhenius terms
-    % 4. A Gaussian.
-    %
-    % You can also create composite functors by adding, multiplying, or
-    % dividing these basic functors, or other composite functors.
-    %
-    % Note: this MATLAB class shadows the underlying C++ Cantera class
-    % "Func1". See the Cantera C++ documentation for more details.
-    %
-    % See also: :mat:class:`polynom`, :mat:class:`gaussian`, :mat:class:`fplus`,
-    % :mat:class:`frdivide`, :mat:class:`ftimes`
-    %
-    % :param typ:
-    %     String indicating type of functor to create. Possible values are:
-    %
-    %     * ``'polynomial'``
-    %     * ``'fourier'``
-    %     * ``'gaussian'``
-    %     * ``'arrhenius'``
-    %     * ``'sum'``
-    %     * ``'diff'``
-    %     * ``'ratio'``
-    %     * ``'composite'``
-    %     * ``'periodic'``
-    % :param n:
-    %     Number of parameters required for the functor
-    % :param p:
-    %     Vector of parameters
-    % :return:
-    %     Instance of class :mat:class:`Func`
-    %
 
     properties (SetAccess = immutable)
         f1
@@ -65,6 +12,59 @@
         %% Func Class Constructor
 
         function x = Func(typ, n, p)
+            % Func Class ::
+            %
+            %     >> x = Func(typ, n, p)
+            %
+            % A functor is an object that behaves like a function. Cantera
+            % defines a set of functors to use to create arbitrary functions to
+            % specify things like heat fluxes, piston speeds, etc., in reactor
+            % network simulations. Of course, they can be used for other things
+            % too.
+            %
+            % The main feature of a functor class is that it overloads the ``()``
+            % operator to evaluate the function. For example, suppose object
+            % ``f`` is a functor that evaluates the polynomial :math:`2x^2 - 3x + 1`.
+            % Then writing ``f(2)`` would cause the method that evaluates the
+            % function to be invoked, and would pass it the argument ``2``. The
+            % return value would of course be 3.
+            %
+            % The types of functors you can create in Cantera are these:
+            %
+            % 1. A polynomial
+            % 2. A Fourier series
+            % 3. A sum of Arrhenius terms
+            % 4. A Gaussian.
+            %
+            % You can also create composite functors by adding, multiplying, or
+            % dividing these basic functors, or other composite functors.
+            %
+            % Note: this MATLAB class shadows the underlying C++ Cantera class
+            % "Func1". See the Cantera C++ documentation for more details.
+            %
+            % See also: :mat:class:`polynom`, :mat:class:`gaussian`, :mat:class:`fplus`,
+            % :mat:class:`frdivide`, :mat:class:`ftimes`
+            %
+            % :param typ:
+            %     String indicating type of functor to create. Possible values are:
+            %
+            %     * ``'polynomial'``
+            %     * ``'fourier'``
+            %     * ``'gaussian'``
+            %     * ``'arrhenius'``
+            %     * ``'sum'``
+            %     * ``'diff'``
+            %     * ``'ratio'``
+            %     * ``'composite'``
+            %     * ``'periodic'``
+            % :param n:
+            %     Number of parameters required for the functor
+            % :param p:
+            %     Vector of parameters
+            % :return:
+            %     Instance of class :mat:class:`Func`
+            %
+
             checklib;
 
             if ~isa(typ, 'char')
@@ -77,8 +77,7 @@
             itype = -1;
 
             function nn = newFunc(itype, n, p)
-                % helper function to pass the correct parameters to the C
-                % library
+                % helper function to pass the correct parameters to the C library.
                 if itype < 20
                     [msize, nsize] = size(p);
                     lenp = msize * nsize;
@@ -135,7 +134,7 @@
         %% Func Class Destructor
 
         function delete(f)
-            % Delete the C++ Func1 object.
+            % Delete the :mat:class:`Func` object.
 
             callct('func_del', f.id);
         end
diff --git a/interfaces/matlab_experimental/Func/fplus.m b/interfaces/matlab_experimental/Func/fplus.m
index 778fbb6358..b974122b33 100644
--- a/interfaces/matlab_experimental/Func/fplus.m
+++ b/interfaces/matlab_experimental/Func/fplus.m
@@ -1,7 +1,7 @@
-function r = fplus(a, b)
-    % Get a functor representing the sum of two input functors.
+classdef fplus < Func
+    % Fplus - Get a functor representing the sum of two input functors. ::
     %
-    % r = fplus(a, b)
+    %     >> f = fplus(a, b)
     %
     % :param a:
     %     Instance of class :mat:class:`Func`
@@ -12,9 +12,10 @@
     %
     methods
 
-        % Constructor
         function f = fplus(a, b)
-            f = f@Func('sum', a, b);
+            % Constructor
+
+            f@Func('sum', a, b);
         end
 
     end
diff --git a/interfaces/matlab_experimental/Func/frdivide.m b/interfaces/matlab_experimental/Func/frdivide.m
index d4298ad79e..4de5b1c1a0 100644
--- a/interfaces/matlab_experimental/Func/frdivide.m
+++ b/interfaces/matlab_experimental/Func/frdivide.m
@@ -1,7 +1,7 @@
-function r = frdivide(a, b)
-    % Get a functor that is the ratio of the input functors.
+classdef frdivide < Func
+    % Frdivide - Get a functor representing the ratio of two input functors. ::
     %
-    % r = frdivide(a,b)
+    %     >> f = frdivide(a, b)
     %
     % :param a:
     %     Instance of class :mat:class:`Func`
@@ -12,9 +12,10 @@
     %
     methods
 
-        % Constructor
         function f = frdivide(a, b)
-            f = f@Func('ratio', a, b);
+            % Constructor
+
+            f@Func('ratio', a, b);
         end
 
     end
diff --git a/interfaces/matlab_experimental/Func/ftimes.m b/interfaces/matlab_experimental/Func/ftimes.m
index 7cb83ac4dd..3bb7b41de2 100644
--- a/interfaces/matlab_experimental/Func/ftimes.m
+++ b/interfaces/matlab_experimental/Func/ftimes.m
@@ -1,7 +1,7 @@
 classdef ftimes < Func
-    % Ftimes - Get a functor representing the product of two input functors.
+    % Ftimes - Get a functor representing the product of two input functors. ::
     %
-    % f = ftimes(a, b)
+    %     >> f = ftimes(a, b)
     %
     % :param a:
     %     Instance of class :mat:class:`Func`
@@ -12,9 +12,10 @@
     %
     methods
 
-        % Constructor
         function f = ftimes(a, b)
-            f = f@Func('prod', a, b);
+            % Constructor
+
+            f@Func('prod', a, b);
         end
 
     end
diff --git a/interfaces/matlab_experimental/Func/gaussian.m b/interfaces/matlab_experimental/Func/gaussian.m
index 502fdd7df6..d31851b41d 100644
--- a/interfaces/matlab_experimental/Func/gaussian.m
+++ b/interfaces/matlab_experimental/Func/gaussian.m
@@ -1,7 +1,7 @@
 classdef gaussian < Func
-    % Gaussian - Create a Gaussian :mat:class:`Func` instance.
+    % Gaussian - Create a Gaussian :mat:class:`Func` instance. ::
     %
-    % f = gaussian(peak, center, width)
+    %     >> f = gaussian(peak, center, width)
     %
     % :param peak:
     %     The peak value
@@ -16,9 +16,10 @@
     %
     methods
 
-        % Constructor
         function f = gaussian(peak, center, width)
-            f = f@Func('gaussian', 0, [peak, center, width]);
+            % Constructor
+
+            f@Func('gaussian', 0, [peak, center, width]);
         end
 
     end
diff --git a/interfaces/matlab_experimental/Func/polynom.m b/interfaces/matlab_experimental/Func/polynom.m
index dbb4e79484..308352f34c 100644
--- a/interfaces/matlab_experimental/Func/polynom.m
+++ b/interfaces/matlab_experimental/Func/polynom.m
@@ -1,7 +1,7 @@
 classdef polynom < Func
-    % Polynom - Create a polynomial :mat:class:`Func` instance.
+    % Polynom - Create a polynomial :mat:class:`Func` instance. ::
     %
-    % poly = polynom(coeffs)
+    %     >> f = polynom(coeffs)
     %
     % The polynomial coefficients are specified by a vector
     % ``[a0 a1 .... aN]``. Examples:
@@ -18,18 +18,22 @@
     %
     methods
 
-        % Constructor
         function f = polynom(coeffs)
+            % Constructor
+
             [n m] = size(coeffs);
 
-    [n m] = size(coeffs);
+            if n == 1
+                p = m - 1;
+            elseif m == 1
+                p = n - 1;
+            else
+                error('wrong shape for coefficient array');
+            end
+
+            f@Func('polynomial', p, coeffs);
+        end
 
-    if n == 1
-        poly = Func('polynomial', m - 1, coeffs);
-    elseif m == 1
-        poly = Func('polynomial', n - 1, coeffs);
-    else
-        error('wrong shape for coefficient array');
     end
 
 end
diff --git a/interfaces/matlab_experimental/Phases/Mixture.m b/interfaces/matlab_experimental/Phases/Mixture.m
index 181e2c52d9..d6701cef67 100644
--- a/interfaces/matlab_experimental/Phases/Mixture.m
+++ b/interfaces/matlab_experimental/Phases/Mixture.m
@@ -1,7 +1,7 @@
 classdef Mixture < handle
-    % Mixture Class.
+    % Mixture Class ::
     %
-    % m = Mixture(phases)
+    %     >> m = Mixture(phases)
     %
     % Class :mat:class:`Mixture` represents mixtures of one or more phases of matter.
     % To construct a mixture, supply a cell array of phases and
@@ -150,6 +150,8 @@
         %% Mixture Class Constructor
 
         function m = Mixture(phases)
+            % Create a :mat:class:`Mixture` object.
+
             checklib;
 
             if nargin > 1
@@ -189,7 +191,7 @@
         %% Mixture Class Destructor
 
         function delete(m)
-            % Delete the MultiPhase object.
+            % Delete the :mat:class:`Mixture` object.
 
             calllib(ct, 'mix_del', m.mixID);
         end
diff --git a/interfaces/matlab_experimental/Reactor/ConstPressureReactor.m b/interfaces/matlab_experimental/Reactor/ConstPressureReactor.m
index 22979809c8..dc4808bc99 100644
--- a/interfaces/matlab_experimental/Reactor/ConstPressureReactor.m
+++ b/interfaces/matlab_experimental/Reactor/ConstPressureReactor.m
@@ -1,7 +1,7 @@
 classdef ConstPressureReactor < Reactor
-    % Create a constant pressure reactor object.
+    % Create a constant pressure reactor object. ::
     %
-    % r = ConstPressureReactor(contents)
+    %     >> r = ConstPressureReactor(contents)
     %
     % A :mat:class:`ConstPressureReactor` is an instance of class
     % :mat:class:`Reactor` where the pressure is held constant. The volume
@@ -23,13 +23,14 @@
 
     methods
 
-        % Constructor
         function r = ConstPressureReactor(contents)
+            % Constructor
+
             if nargin == 0
                 contents = 0;
             end
 
-            r = r@Reactor(contents, 'ConstPressureReactor');
+            r@Reactor(contents, 'ConstPressureReactor');
         end
 
     end
diff --git a/interfaces/matlab_experimental/Reactor/FlowDevice.m b/interfaces/matlab_experimental/Reactor/FlowDevice.m
index bfb8c87b26..00273d1fcb 100644
--- a/interfaces/matlab_experimental/Reactor/FlowDevice.m
+++ b/interfaces/matlab_experimental/Reactor/FlowDevice.m
@@ -1,7 +1,7 @@
 classdef FlowDevice < handle
-    % FlowDevice Class.
+    % FlowDevice Class ::
     %
-    % x = FlowDevice(typ)
+    %     >> x = FlowDevice(typ)
     %
     % Base class for devices that allow flow between reactors.
     % :mat:class:`FlowDevice` objects are assumed to be adiabatic,
@@ -30,7 +30,7 @@
 
     end
 
-    properties (SetAccess = protected)
+    properties (SetAccess = public)
         %
         % Upstream object of type :mat:class:`Reactor` or :mat:class:`Reservoir`.
         upstream
@@ -59,6 +59,8 @@
         %% FlowDevice Class Constructor
 
         function x = FlowDevice(typ)
+            % Create a :mat:class:`FlowDevice` object.
+
             checklib;
 
             if nargin == 0
@@ -74,7 +76,7 @@
         %% FlowDevice Class Destructor
 
         function delete(f)
-            % Delete the specified flow device from memory.
+            % Delete the :mat:class:`FlowDevice` object.
 
             callct('flowdev_del', f.id);
         end
diff --git a/interfaces/matlab_experimental/Reactor/FlowReactor.m b/interfaces/matlab_experimental/Reactor/FlowReactor.m
index 6b46d1d3a6..5ad7819ccf 100644
--- a/interfaces/matlab_experimental/Reactor/FlowReactor.m
+++ b/interfaces/matlab_experimental/Reactor/FlowReactor.m
@@ -1,7 +1,7 @@
 classdef FlowReactor < Reactor
-    % Create a flow reactor object.
+    % Create a flow reactor object. ::
     %
-    % r = FlowReactor(contents)
+    %     >> r = FlowReactor(contents)
     %
     % A reactor representing adiabatic plug flow in a constant-area
     % duct. Examples:
@@ -21,13 +21,14 @@
 
     methods
 
-        % Constructor
         function r = FlowReactor(contents)
+            % Constructor
+
             if nargin == 0
                 contents = 0;
             end
 
-            r = r@Reactor(contents, 'FlowReactor');
+            r@Reactor(contents, 'FlowReactor');
         end
 
     end
diff --git a/interfaces/matlab_experimental/Reactor/IdealGasConstPressureReactor.m b/interfaces/matlab_experimental/Reactor/IdealGasConstPressureReactor.m
index 722719b61b..1b9c626992 100644
--- a/interfaces/matlab_experimental/Reactor/IdealGasConstPressureReactor.m
+++ b/interfaces/matlab_experimental/Reactor/IdealGasConstPressureReactor.m
@@ -1,7 +1,7 @@
 classdef IdealGasConstPressureReactor < Reactor
-    % Create a constant pressure reactor with an ideal gas.
+    % Create a constant pressure reactor with an ideal gas. ::
     %
-    % r = IdealGasConstPressureReactor(contents)
+    %     >> r = IdealGasConstPressureReactor(contents)
     %
     % An IdealGasConstPressureReactor is an instance of class Reactor where the
     % pressure is held constant. The volume is not a state variable, but
@@ -25,13 +25,14 @@
 
     methods
 
-        % Constructor
         function r = IdealGasConstPressureReactor(contents)
+            % Constructor
+
             if nargin == 0
                 contents = 0;
             end
 
-            r = r@Reactor(contents, 'IdealGasConstPressureReactor');
+            r@Reactor(contents, 'IdealGasConstPressureReactor');
         end
 
     end
diff --git a/interfaces/matlab_experimental/Reactor/IdealGasReactor.m b/interfaces/matlab_experimental/Reactor/IdealGasReactor.m
index 2f44e40c70..1747ec5ec5 100644
--- a/interfaces/matlab_experimental/Reactor/IdealGasReactor.m
+++ b/interfaces/matlab_experimental/Reactor/IdealGasReactor.m
@@ -1,7 +1,7 @@
 classdef IdealGasReactor < Reactor
-    % Create a reactor with an ideal gas.
+    % Create a reactor with an ideal gas. ::
     %
-    % r = IdealGasReactor(contents)
+    %     >> r = IdealGasReactor(contents)
     %
     % An IdealGasReactor is an instance of class Reactor where the governing
     % equations are specialized for the ideal gas equation of state (and do not
@@ -22,13 +22,14 @@
 
     methods
 
-        % Constructor
         function r = IdealGasReactor(contents)
+            % Constructor
+
             if nargin == 0
                 contents = 0;
             end
 
-            r = r@Reactor(contents, 'IdealGasReactor');
+            r@Reactor(contents, 'IdealGasReactor');
         end
 
     end
diff --git a/interfaces/matlab_experimental/Reactor/MassFlowController.m b/interfaces/matlab_experimental/Reactor/MassFlowController.m
index c5273efcf0..44ba2f2aee 100644
--- a/interfaces/matlab_experimental/Reactor/MassFlowController.m
+++ b/interfaces/matlab_experimental/Reactor/MassFlowController.m
@@ -1,7 +1,7 @@
 classdef MassFlowController < FlowDevice
-    % Create a mass flow controller.
+    % Create a mass flow controller. ::
     %
-    % m = MassFlowController(upstream, downstream)
+    %     >> m = MassFlowController(upstream, downstream)
     %
     % Creates an instance of class :mat:class:`FlowDevice` configured to
     % simulate a mass flow controller that maintains a constant mass flow
@@ -23,10 +23,10 @@
 
     methods
 
-        % Constructor
         function m = MassFlowController(upstream, downstream)
+            % Constructor
 
-            m = m@FlowDevice('MassFlowController');
+            m@FlowDevice('MassFlowController');
 
             if nargin == 2
                 m.install(upstream, downstream)
diff --git a/interfaces/matlab_experimental/Reactor/Reactor.m b/interfaces/matlab_experimental/Reactor/Reactor.m
index f963df29b3..055330a32e 100644
--- a/interfaces/matlab_experimental/Reactor/Reactor.m
+++ b/interfaces/matlab_experimental/Reactor/Reactor.m
@@ -1,7 +1,7 @@
 classdef Reactor < handle
-    % Reactor Class
+    % Reactor Class ::
     %
-    % r = Reactor(content, typ)
+    %     >> r = Reactor(content, typ)
     %
     % A 'Reactor' object simulates a perfectly-stirred reactor. It
     % has a time-dependent tstate, and may be coupled to other
@@ -30,7 +30,7 @@
 
     end
 
-    properties (SetAccess = protected)
+    properties (SetAccess = public)
 
         contents
         %
@@ -111,6 +111,8 @@
         %% Reactor Class Constructor
 
         function r = Reactor(content, typ)
+            % Create a :mat:class:`Reactor` object.
+
             checklib;
 
             if nargin == 0
@@ -126,8 +128,18 @@
             r.id = callct('reactor_new', typ);
 
             if isa(content, 'Solution')
-                r.insert(content);
-            elseif ~(isa(contents, 'double') && contents == 0)
+                r.contents = content;
+                if ~isa(content, 'ThermoPhase')
+                    error('Wrong object type');
+                end
+
+                callct('reactor_setThermoMgr', r.id, content.tpID);
+
+                if ~strcmp(r.type, 'Reservoir')
+                    callct('reactor_setKineticsMgr', r.id, content.kinID);
+                end
+
+            elseif ~(isa(content, 'double') && content == 0)
                 error('Reactor contents must be an object of type "Solution"');
             end
 
@@ -136,7 +148,7 @@
         %% Reactor Class Destructor
 
         function delete(r)
-            % Delete the Reactor object from memory.
+            % Delete the :mat:class:`Reactor` object.
 
             callct('reactor_del', r.id);
         end
@@ -157,26 +169,6 @@ function addSensitivityReaction(r, m)
             callct('reactor_addSensitivityReaction', r.id, m);
         end
 
-        function insert(r, gas)
-            % Insert a solution or mixture into a reactor.
-            %
-            % r.insert(gas)
-            %
-            % :param r:
-            %    Instance of class 'Reactor'.
-            % :param gas:
-            %    Instance of class 'Solution'.
-            %
-
-            r.contents = gas;
-            setThermoMgr(r, gas);
-
-            if ~strcmp(r.type, 'Reservoir')
-                setKineticsMgr(r, gas);
-            end
-
-        end
-
         %% Reactor Get Methods
 
         function temperature = get.T(r)
@@ -272,48 +264,6 @@ function insert(r, gas)
 
         end
 
-        function setThermoMgr(r, t)
-            % Set the thermodynamics manager.
-            %
-            % r.setThermoMgr(t)
-            %
-            % This method is used internally during Reactor initialization,
-            % but is usually not called by users.
-            %
-            % :param r:
-            %    Instance of class 'Reactor'.
-            % :param t:
-            %    Instance of class 'ThermoPhase' or another object
-            %    containing an instance of that class.
-
-            if ~isa(t, 'ThermoPhase')
-                error('Wrong object type');
-            end
-
-            callct('reactor_setThermoMgr', r.id, t.tpID);
-        end
-
-        function setKineticsMgr(r, k)
-            % Set the kinetics manager.
-            %
-            % r.setKineticsMgr(k)
-            %
-            % This method is used internally during Reactor initialization,
-            % but is usually not called by users.
-            %
-            % :param r:
-            %    Instance of class 'Reactor'.
-            % :param t:
-            %    Instance of class 'Kinetics' or another object
-            %    containing an instance of that class.
-
-            if ~isa(k, 'Kinetics')
-                error('Wrong object type');
-            end
-
-            callct('reactor_setKineticsMgr', r.id, k.kinID);
-        end
-
     end
 
 end
diff --git a/interfaces/matlab_experimental/Reactor/ReactorNet.m b/interfaces/matlab_experimental/Reactor/ReactorNet.m
index d107e20c27..7ffa5c6b55 100644
--- a/interfaces/matlab_experimental/Reactor/ReactorNet.m
+++ b/interfaces/matlab_experimental/Reactor/ReactorNet.m
@@ -1,7 +1,7 @@
 classdef ReactorNet < handle
-    % ReactorNet class
+    % ReactorNet class ::
     %
-    % r = ReactorNet(reacgtors)
+    %     >> r = ReactorNet(reactors)
     %
     % A 'ReactorNet' object is a container that holds one or more
     % 'Reactor' objects in a network. 'ReactorNet' objects are used
@@ -59,6 +59,8 @@
         %% ReactorNet Class Constructor
 
         function r = ReactorNet(reactors)
+            % Create a :mat:class:`ReactorNet` object.
+
             checklib;
 
             if nargin ~= 1
@@ -77,7 +79,7 @@
             nr = length(reactors);
 
             for i = 1:nr
-                r.addReactor(reactors{i});
+                callct('reactornet_addreactor', r.id, reactors{i}.id);
             end
 
         end
@@ -85,27 +87,13 @@
         %% ReactorNet Class Destructor
 
         function delete(r)
-            % Delete the ReactorNet object from the memory.
+            % Delete the :mat:class:`ReactorNet` object object.
 
             callct('reactornet_del', r.id);
         end
 
         %% ReactorNet Utility Methods
 
-        function addReactor(net, reactor)
-            % Add a reactor to a network.
-            %
-            % net.addReactor(reactor)
-            %
-            % :param net:
-            %    Instance of class 'ReactorNet'.
-            % :param reactor:
-            %    Instance of class 'Solution'.
-            %
-
-            callct('reactornet_addreactor', net.id, reactor.id);
-        end
-
         function advance(r, tout)
             % Advance the state of the reactor network in time.
             %
diff --git a/interfaces/matlab_experimental/Reactor/ReactorSurface.m b/interfaces/matlab_experimental/Reactor/ReactorSurface.m
index 59208f465b..72d767e716 100644
--- a/interfaces/matlab_experimental/Reactor/ReactorSurface.m
+++ b/interfaces/matlab_experimental/Reactor/ReactorSurface.m
@@ -1,7 +1,7 @@
 classdef ReactorSurface < handle
-    % ReactorSurface Class
+    % ReactorSurface Class ::
     %
-    % s = ReactorSurface(kleft, reactor, area)
+    %     >> s = ReactorSurface(surf, reactor, area)
     %
     % A surface on which heterogeneous reactions take place. The
     % mechanism object (typically an instance of class 'Interface')
@@ -14,7 +14,7 @@
     % after initial construction by using the various methods of
     % the 'ReactorSurface' class.
     %
-    % :param kleft:
+    % :param surf:
     %    Surface reaction mechanisms for the left-facing surface.
     %    This must bean instance of class 'Kinetics', or of a class
     %    derived from 'Kinetics', such as 'Interface'.
@@ -43,20 +43,29 @@
     methods
         %% ReactorSurface Class Constructor
 
-        function s = ReactorSurface(kleft, reactor, area)
+        function s = ReactorSurface(surf, reactor, area)
+            % Create a :mat:class:`ReactorSurface` object.
+
             checklib;
 
             s.surfID = callct('reactorsurface_new', 0);
             s.reactor = -1;
 
             if nargin >= 1
-                s.setKinetics(kleft);
+                ikin = 0;
+
+                if isa(surf, 'Kinetics')
+                    ikin = surf.kinID;
+                end
+
+                callct('reactorsurface_setkinetics', s.surfID, ikin);
             end
 
             if nargin >= 2
 
                 if isa(reactor, 'Reactor')
-                    s.install(reactor);
+                    s.reactor = reactor;
+                    callct('reactorsurface_install', s.surfID, reactor.id);
                 else
                     warning('Reactor was not installed due to incorrect type');
                 end
@@ -78,26 +87,13 @@
         %% ReactorSurface Class Destructor
 
         function delete(s)
-            % Delete the ReactorSurface object from the memory.
+            % Delete the :mat:class:`ReactorSurface` object.
 
             callct('reactorsurface_del', s.surfID);
         end
 
         %% ReactorSurface Utility Methods
 
-        function install(s, r)
-            % Install a ReactorSurface in a Reactor.
-            %
-            % s.install(r)
-            %
-            % :param r:
-            %    Instance of class 'Reactor'.
-            %
-
-            s.reactor = r;
-            callct('reactorsurface_install', s.surfID, r.id);
-        end
-
         function addSensitivityReaction(s, m)
             % Specifies that the sensitivity of the state variables with
             % respect to reaction m should be computed. The surface must be
@@ -123,26 +119,6 @@ function addSensitivityReaction(s, m)
             callct('reactorsurface_setArea', s.surfID, a);
         end
 
-        function setKinetics(s, kin)
-            % Setthe surface reaction mechanism on a reactor surface.
-            %
-            % s.setKinetics(kin)
-            %
-            % :param kin:
-            %    Instance of class 'Kinetics' (or another object derived
-            %    from kin) to be used as the kinetic mechanism for this
-            %    surface. Typically an instance of class 'Interface'.
-            %
-
-            ikin = 0;
-
-            if isa(kin, 'Kinetics')
-                ikin = kin.kinID;
-            end
-
-            callct('reactorsurface_setkinetics', s.surfID, ikin);
-        end
-
     end
 
 end
diff --git a/interfaces/matlab_experimental/Reactor/Reservoir.m b/interfaces/matlab_experimental/Reactor/Reservoir.m
index 6046d46180..bcab24d72b 100644
--- a/interfaces/matlab_experimental/Reactor/Reservoir.m
+++ b/interfaces/matlab_experimental/Reactor/Reservoir.m
@@ -1,7 +1,7 @@
 classdef Reservoir < Reactor
-    % Create a Reservoir object.
+    % Create a Reservoir object. ::
     %
-    % r = Reservoir(contents)
+    %     >> r = Reservoir(contents)
     %
     % A :mat:class:`Reservoir` is an instance of class :mat:class:`Reactor`
     % configured so that its intensive state is constant in time. A
@@ -28,13 +28,14 @@
 
     methods
 
-        % Constructor
         function r = Reservoir(contents)
+            % Constructor
+
             if nargin == 0
                 contents = 0;
             end
 
-            r = r@Reactor(contents, 'Reservoir');
+            r@Reactor(contents, 'Reservoir');
         end
 
     end
diff --git a/interfaces/matlab_experimental/Reactor/Valve.m b/interfaces/matlab_experimental/Reactor/Valve.m
index 52a35ee626..86682eef2a 100644
--- a/interfaces/matlab_experimental/Reactor/Valve.m
+++ b/interfaces/matlab_experimental/Reactor/Valve.m
@@ -1,6 +1,8 @@
 classdef Valve < FlowDevice
-    % Create a valve.
-    % v = Valve(upstream, downstream)
+    % Create a valve. ::
+    %
+    %     >> v = Valve(upstream, downstream)
+    %
     % Create an instance of class :mat:class:`FlowDevice` configured to
     % simulate a valve that produces a flow rate proportional to the
     % pressure difference between the upstream and downstream reactors.
@@ -30,10 +32,10 @@
 
     methods
 
-        % Constructor
         function v = Valve(upstream, downstream)
+            % Constructor
 
-            v = v@FlowDevice('Valve');
+            v@FlowDevice('Valve');
 
             if nargin == 2
                 v.install(upstream, downstream)
diff --git a/interfaces/matlab_experimental/Reactor/Wall.m b/interfaces/matlab_experimental/Reactor/Wall.m
index 884f3e019d..4b1b6e8570 100644
--- a/interfaces/matlab_experimental/Reactor/Wall.m
+++ b/interfaces/matlab_experimental/Reactor/Wall.m
@@ -1,7 +1,7 @@
 classdef Wall < handle
-    % Wall Class
+    % Wall Class ::
     %
-    % x = Wall()
+    %     >> x = Wall(l, r)
     %
     % A Wall separates two reactors, or a reactor and a reservoir.
     % A Wall has a finite area, may conduct heat between the two
@@ -34,6 +34,15 @@
     % after initial construction by using the various methods of
     % the 'Wall' class.
     %
+    % :param l:
+    %    Instance of class 'Reactor' to be used as the bulk phase
+    %    on the left side of the wall.
+    % :param r:
+    %    Instance of class 'Reactor' to be used as the bulk phase
+    %    on the right side of the wall.
+    % :return:
+    %    Instance of class 'Wall'.
+    %
 
     properties (SetAccess = immutable)
 
@@ -77,60 +86,44 @@
     methods
         %% Wall Class Constructor
 
-        function x = Wall()
+        function w = Wall(l, r)
+            % Create a :mat:class:`Wall` object.
             checklib;
 
             % At the moment, only one wall type is implemented
             typ = 'Wall';
 
-            x.type = char(typ);
-            x.id = callct('wall_new', x.type);
+            w.type = char(typ);
+            w.id = callct('wall_new', w.type);
+
+            % Install the wall between left and right reactors
+            w.left = l;
+            w.right = r;
+            callct('wall_install', w.id, l.id, r.id);
 
             % Set default values.
-            x.left = -1;
-            x.right = -1;
-            x.area = 1.0;
-            x.setExpansionRateCoeff(0.0);
-            x.setHeatTransferCoeff(0.0);
+            w.area = 1.0;
+            w.expansionRateCoeff = 0.0;
+            w.heatTransferCoeff = 0.0;
+
+            % Check whether the wall is ready.
+            ok = callct('wall_ready', w.id);
+            if ok
+                disp('The wall object is ready.');
+            else
+                error('The wall object is not ready.');
+            end
 
         end
 
         %% Wall Class Destructor
 
         function delete(w)
-            % Clear the Wall object from the memory.
+            % Clear the :mat:class:`Wall` object.
 
             callct('wall_del', w.id);
         end
 
-        %% Wall Class Utility Methods
-
-        function install(w, l, r)
-            % Install a wall between two reactors.
-            %
-            % w.install(l, r)
-            %
-            % :param l:
-            %    Instance of class 'Reactor' to be used as the bulk phase
-            %    on the left side of the wall.
-            % :param r:
-            %    Instance of class 'Reactor' to be used as the bulk phase
-            %    on the right side of the wall.
-
-            w.left = l;
-            w.right = r;
-            callct('wall_install', w.id, l.id, r.id);
-        end
-
-        function ok = ready(w)
-            % Check whether a wall is ready.
-            %
-            % ok = w.ready
-            %
-
-            ok = callct('wall_ready', w.id);
-        end
-
         %% ReactorNet get methods
 
         function a = get.area(w)
@@ -138,7 +131,7 @@ function install(w, l, r)
         end
 
         function q = qdot(w, t)
-            % Total heat transfer through a wall at a given time t.
+            % Total heat transfer rate through a wall at a given time t.
 
             q = callct('wall_Q', w.id, t);
         end
diff --git a/interfaces/matlab_experimental/Utility/canteraGitCommit.m b/interfaces/matlab_experimental/Utility/canteraGitCommit.m
index f9e964d66c..446a0dec57 100644
--- a/interfaces/matlab_experimental/Utility/canteraGitCommit.m
+++ b/interfaces/matlab_experimental/Utility/canteraGitCommit.m
@@ -1,7 +1,7 @@
 function v = canteraGitCommit()
-    % Get Cantera Git commit hash.
+    % Get Cantera Git commit hash. ::
     %
-    % canteraGitCommit()
+    %     >> canteraGitCommit()
     %
     % :return:
     %     A string containing the Git commit hash for the current version of Cantera.
diff --git a/interfaces/matlab_experimental/Utility/canteraVersion.m b/interfaces/matlab_experimental/Utility/canteraVersion.m
index edba03d55c..6dbe16c965 100644
--- a/interfaces/matlab_experimental/Utility/canteraVersion.m
+++ b/interfaces/matlab_experimental/Utility/canteraVersion.m
@@ -1,7 +1,7 @@
 function v = canteraVersion()
-    % Get Cantera version information.
+    % Get Cantera version information. ::
     %
-    % canteraVersion()
+    %     >> canteraVersion()
     %
     % :return:
     %     A string containing the Cantera version.
diff --git a/interfaces/matlab_experimental/Utility/getDataDirectories.m b/interfaces/matlab_experimental/Utility/getDataDirectories.m
index 735bcef603..ca7c3840f7 100644
--- a/interfaces/matlab_experimental/Utility/getDataDirectories.m
+++ b/interfaces/matlab_experimental/Utility/getDataDirectories.m
@@ -1,7 +1,7 @@
 function d = getDataDirectories()
-    % Get a cell array of the directories searched for data files.
+    % Get a cell array of the directories searched for data files. ::
     %
-    % getdatadirectories()
+    %     >> getdatadirectories()
     %
     % Get a cell array of the directories Cantera searches for data files
     %
diff --git a/samples/matlab_experimental/catcomb.m b/samples/matlab_experimental/catcomb.m
index 9ca8c3c797..3032745e47 100644
--- a/samples/matlab_experimental/catcomb.m
+++ b/samples/matlab_experimental/catcomb.m
@@ -18,7 +18,6 @@
 
 clear all
 close all
-clc
 
 t0 = cputime; % record the starting time
 
@@ -130,7 +129,7 @@
 % Once the component parts have been created, they can be assembled
 % to create the 1D simulation.
 
-stack = Sim1D([inlt, flow, surf]);
+stack = Sim1D({inlt, flow, surf});
 
 % set the initial profiles.
 stack.setProfile(2, {'velocity', 'spread_rate', 'T'}, ...
diff --git a/samples/matlab_experimental/diamond_cvd.m b/samples/matlab_experimental/diamond_cvd.m
index 1060728bb2..822ba3d1ba 100644
--- a/samples/matlab_experimental/diamond_cvd.m
+++ b/samples/matlab_experimental/diamond_cvd.m
@@ -14,12 +14,10 @@
 
 %% Initialization
 
-help diamond_cvd
-
 clear all
 close all
-clc
 
+help diamond_cvd
 t0 = cputime; % record the starting time
 
 %% Operation Parameters
diff --git a/samples/matlab_experimental/diff_flame.m b/samples/matlab_experimental/diff_flame.m
index 1f3de8a8d2..98f016716c 100644
--- a/samples/matlab_experimental/diff_flame.m
+++ b/samples/matlab_experimental/diff_flame.m
@@ -10,7 +10,6 @@
 
 clear all
 close all
-clc
 
 tic % total running time of the script
 help diff_flame
diff --git a/samples/matlab_experimental/equil.m b/samples/matlab_experimental/equil.m
index a8e09aea5f..69ac3ba52c 100644
--- a/samples/matlab_experimental/equil.m
+++ b/samples/matlab_experimental/equil.m
@@ -8,7 +8,6 @@ function equil(g)
 
     clear all
     close all
-    clc
 
     tic
     help equil
diff --git a/samples/matlab_experimental/flame.m b/samples/matlab_experimental/flame.m
index 0a0320b3f2..cd567bfd07 100644
--- a/samples/matlab_experimental/flame.m
+++ b/samples/matlab_experimental/flame.m
@@ -27,7 +27,7 @@
     end
 
     % create the container object
-    f = Sim1D([left flow right]);
+    f = Sim1D({left flow right});
 
     % set default initial profiles.
     rho0 = gas.D;
diff --git a/samples/matlab_experimental/flame1.m b/samples/matlab_experimental/flame1.m
index d04f174de5..5b13d048ea 100644
--- a/samples/matlab_experimental/flame1.m
+++ b/samples/matlab_experimental/flame1.m
@@ -9,7 +9,6 @@
 
 clear all
 close all
-clc
 
 tic
 help flame1
diff --git a/samples/matlab_experimental/flame2.m b/samples/matlab_experimental/flame2.m
index 7d6a2ad685..cf2957bd01 100644
--- a/samples/matlab_experimental/flame2.m
+++ b/samples/matlab_experimental/flame2.m
@@ -8,7 +8,6 @@
 
 clear all
 close all
-clc
 
 tic
 help flame2
diff --git a/samples/matlab_experimental/ignite.m b/samples/matlab_experimental/ignite.m
index 1d21b1eb11..6c8cdf09a2 100644
--- a/samples/matlab_experimental/ignite.m
+++ b/samples/matlab_experimental/ignite.m
@@ -9,7 +9,6 @@
 
     clear all
     close all
-    clc
 
     tic
     help ignite
diff --git a/samples/matlab_experimental/ignite_hp.m b/samples/matlab_experimental/ignite_hp.m
index 3f467c45e6..2141ce862e 100644
--- a/samples/matlab_experimental/ignite_hp.m
+++ b/samples/matlab_experimental/ignite_hp.m
@@ -6,7 +6,6 @@ function ignite_hp(gas)
 
     clear all
     close all
-    clc
 
     tic
     help ignite_hp
diff --git a/samples/matlab_experimental/ignite_uv.m b/samples/matlab_experimental/ignite_uv.m
index 5d521a3dd9..e544905e52 100644
--- a/samples/matlab_experimental/ignite_uv.m
+++ b/samples/matlab_experimental/ignite_uv.m
@@ -6,7 +6,6 @@ function ignite_uv(gas)
 
     clear all
     close all
-    clc
 
     tic
     help ignite_uv
diff --git a/samples/matlab_experimental/isentropic.m b/samples/matlab_experimental/isentropic.m
index 61fc686fa9..bba146849f 100644
--- a/samples/matlab_experimental/isentropic.m
+++ b/samples/matlab_experimental/isentropic.m
@@ -8,7 +8,6 @@ function isentropic(g)
 
     clear all
     close all
-    clc
 
     tic
     help isentropic
diff --git a/samples/matlab_experimental/lithium_ion_battery.m b/samples/matlab_experimental/lithium_ion_battery.m
index 29185a109d..9869fca14a 100644
--- a/samples/matlab_experimental/lithium_ion_battery.m
+++ b/samples/matlab_experimental/lithium_ion_battery.m
@@ -25,7 +25,6 @@
 
 clear all
 close all
-clc
 
 tic
 help lithium_ion_battery
diff --git a/samples/matlab_experimental/periodic_cstr.m b/samples/matlab_experimental/periodic_cstr.m
index 631b086ea0..e58c57e0fc 100644
--- a/samples/matlab_experimental/periodic_cstr.m
+++ b/samples/matlab_experimental/periodic_cstr.m
@@ -21,7 +21,6 @@
 
     clear all
     close all
-    clc
 
     tic
     help periodic_cstr
@@ -48,7 +47,7 @@
 
     % 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.setInitialVolume(10.0 * 1.0e-6);
+    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
@@ -60,10 +59,9 @@
     % coefficient. Larger U causes the reactor to be closer to isothermal.
     % If U is too small, the gas ignites, and the temperature spikes and
     % stays high.
-    w = Wall;
-    w.install(cstr, env);
+    w = Wall(cstr, env);
     w.area = 1.0;
-    w.setHeatTransferCoeff(0.02);
+    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.
@@ -72,7 +70,7 @@
     mdot = gas.D * vdot; % kg/s
     mfc = MassFlowController;
     mfc.install(upstream, cstr);
-    mfc.setMassFlowRate(mdot);
+    mfc.massFlowRate = mdot;
 
     % now create a downstream reservoir to exhaust into.
     downstream = Reservoir(gas);
@@ -82,7 +80,7 @@
     % close to the downstream pressure of 60 Torr.
     v = Valve;
     v.install(cstr, downstream);
-    v.setValveCoeff(1.0e-9);
+    v.valveCoeff = 1.0e-9;
 
     % create the network
     network = ReactorNet({cstr});
diff --git a/samples/matlab_experimental/plug_flow_reactor.m b/samples/matlab_experimental/plug_flow_reactor.m
index 6a58912b94..e467c279e4 100644
--- a/samples/matlab_experimental/plug_flow_reactor.m
+++ b/samples/matlab_experimental/plug_flow_reactor.m
@@ -26,7 +26,6 @@
 
 clear all
 close all
-clc
 
 tic
 help plug_flow_reactor
diff --git a/samples/matlab_experimental/prandtl1.m b/samples/matlab_experimental/prandtl1.m
index eda9995b4b..bd14236e2c 100644
--- a/samples/matlab_experimental/prandtl1.m
+++ b/samples/matlab_experimental/prandtl1.m
@@ -9,7 +9,6 @@ function prandtl1(g)
 
     clear all
     close all
-    clc
 
     tic
     help prandtl1
diff --git a/samples/matlab_experimental/prandtl2.m b/samples/matlab_experimental/prandtl2.m
index c0247a98b0..8f16c545cc 100644
--- a/samples/matlab_experimental/prandtl2.m
+++ b/samples/matlab_experimental/prandtl2.m
@@ -8,7 +8,6 @@ function prandtl2(g)
 
     clear all
     close all
-    clc
 
     tic
     help prandtl2
diff --git a/samples/matlab_experimental/rankine.m b/samples/matlab_experimental/rankine.m
index 28b7ecff40..e13ba1dd13 100644
--- a/samples/matlab_experimental/rankine.m
+++ b/samples/matlab_experimental/rankine.m
@@ -21,23 +21,19 @@
 w.setState_Tsat(t1, 1.0);
 h1 = w.H;
 p1 = w.P;
-w
 
 % pump it to p2
 pump_work = pump(w, p_max, eta_pump);
 h2 = w.H;
-w
 
 % heat to saturated vapor
 w.setState_Psat(p_max, 1.0);
 h3 = w.H;
-w
 
 heat_added = h3 - h2;
 
 % expand adiabatically back to the initial pressure
 turbine_work = expand(w, p1, eta_turbine);
-w
 
 % compute the efficiency
 efficiency = (turbine_work - pump_work) / heat_added;
diff --git a/samples/matlab_experimental/reactor1.m b/samples/matlab_experimental/reactor1.m
index 9b632123d9..0c0d790f41 100644
--- a/samples/matlab_experimental/reactor1.m
+++ b/samples/matlab_experimental/reactor1.m
@@ -9,7 +9,6 @@ function reactor1(g)
 
     clear all
     close all
-    clc
 
     tic
     help reactor1
@@ -42,11 +41,10 @@ function reactor1(g)
     % 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;
-    w.install(r, env);
+    w = Wall(r, env);
 
     % set expansion parameter. dV/dt = KA(P_1 - P_2)
-    w.setExpansionRateCoeff(1.0e6);
+    w.expansionRateCoeff = 1.0e6;
 
     % set wall area
     w.area = 1.0;
diff --git a/samples/matlab_experimental/reactor2.m b/samples/matlab_experimental/reactor2.m
index c3ccc83216..5e3599960a 100644
--- a/samples/matlab_experimental/reactor2.m
+++ b/samples/matlab_experimental/reactor2.m
@@ -9,7 +9,6 @@ function reactor2(g)
 
     clear all
     close all
-    clc
 
     tic
     help reactor2
diff --git a/samples/matlab_experimental/surf_reactor.m b/samples/matlab_experimental/surf_reactor.m
index ff8fef7b46..45bc23aac1 100644
--- a/samples/matlab_experimental/surf_reactor.m
+++ b/samples/matlab_experimental/surf_reactor.m
@@ -9,7 +9,6 @@
 
 clear all
 close all
-clc
 
 tic
 help surfreactor
@@ -32,7 +31,7 @@
 
 % create a reactor, and insert the gas
 r = IdealGasReactor(gas);
-r.setInitialVolume(1.0e-6)
+r.V = 1.0e-6;
 
 % create a reservoir to represent the environment
 a = Solution('air.yaml', 'air', 'None');
@@ -42,8 +41,7 @@
 % 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;
-w.install(r, env);
+w = Wall(r, env);
 
 A = 1e-4; % Wall area
 
@@ -52,10 +50,10 @@
 
 % set the wall area and heat transfer coefficient.
 w.area = A;
-w.setHeatTransferCoeff(1.0e1); % W/m2/K
+w.heatTransferCoeff = 1.0e1; % W/m2/K
 
 % set expansion rate parameter. dV/dt = KA(P_1 - P_2)
-w.setExpansionRateCoeff(1.0);
+w.expansionRateCoeff = 1.0;
 
 network = ReactorNet({r});
 % setTolerances(network, 1.0e-8, 1.0e-12);
diff --git a/samples/matlab_experimental/test_examples.m b/samples/matlab_experimental/test_examples.m
index be79bcecbf..e97f964b45 100644
--- a/samples/matlab_experimental/test_examples.m
+++ b/samples/matlab_experimental/test_examples.m
@@ -19,9 +19,13 @@
 flame2;
 catcomb;
 diff_flame;
+ignite;
 ignite_hp;
 ignite_uv;
+diamond_cvd;
 
 clear all
 close all
 UnloadCantera
+
+disp('Test example run successful.')