[Matlab] Update 1D documentation

Replace magic numbers for 1D domain types with string constants
This commit is contained in:
ssun30 2022-05-16 13:53:20 -04:00 committed by Ray Speth
parent 2318dd9c21
commit 55ccdc63b9
23 changed files with 1343 additions and 714 deletions

View File

@ -1,4 +1,4 @@
function m = AxiStagnFlow(gas)
% Get an axisymmetric stagnation flow domain.
m = Domain1D(1, gas);
m = Domain1D('StagnationFlow', gas);
end

View File

@ -2,7 +2,7 @@ function m = AxisymmetricFlow(gas, id)
% Create an axisymmetric flow domain.
% :param id:
% String ID of the flow.
m = Domain1D(1, gas);
m = Domain1D('StagnationFlow', gas);
if nargin == 1
m.setID('flow');
else

View File

@ -1,39 +1,39 @@
classdef Domain1D < handle
properties
domainID
type
T
P
domainID % ID of the domain
type % Type of the domain
T % Temperature
P % Pressure
end
methods
%% Domain1D class constructor.
function d = Domain1D(a, b, c)
% DOMAIN1D Domain1D class constructor.
% d = Domain1D(a, b, c)
%
% :parameter a:
% String type of domain. Possible values are:
% 'StagnationFlow'
% 'Inlet1D'
% 'Surf1D'
% 'Symm1D'
% 'Outlet1D'
% 'ReactingSurface'
% 'Sim1D'
% 'OutletRes'
% `StagnationFlow`
% `Inlet1D`
% `Surf1D`
% `Symm1D`
% `Outlet1D`
% `ReactingSurface`
% `Sim1D`
% `OutletRes`
%
% :parameter b:
% Instance of class 'Solution' (for a == 'StagnationFlow')
% or 'Interface' (for a == 'ReactingSurface').
% Not used for all other valid values of 'a'.
% :param b:
% Instance of class :mat:func:`Solution` (for ``a == 1``)
% or :mat:func:`Interface` (for ``a == 6``). Not used for
% all other valid values of ``a``.
% :param c:
% Integer, either 1 or 2, indicating whether an axisymmetric
% stagnation flow or a free flame should be created. If not
% specified, defaults to 1. Ignored if ``a != 1``.
%
% :parameter c:
% A string indicating whether an axisymmetric stagnation
% flow (for c == 'AxisymmetricFlow') or a free flame
% (for c == 'FreeFlame') should be created. If not specified,
% defaults to 'Axisymmetric'. Ignored if parameter "a" is not
% type 'StagnationFlow'.
checklib;
d.domainID = -1;
@ -96,37 +96,54 @@ classdef Domain1D < handle
%% Utility Methods
function clear(d)
% Delete the Domain1D object
% CLEAR Delete the C++ Domain1D object.
% d.clear
% :param d:
% Instance of class :mat:func:`Domain1D` (or another
% object that derives from Domain1D)
%
calllib(ct, 'domain_del', d.domainID);
end
function d = disableEnergy(d)
% Disable the energy equation.
% DISABLEENERGY Disable the energy equation.
% d.disableEnergy
% :param d:
% Instance of class :mat:func:`Domain1D`
%
disp(' ');
disp('Disabling the energy equation...');
calllib(ct, 'stflow_solveEnergyEqn', d.domainID, 0);
end
function d = enableEnergy(d)
% Enable the energy equation.
% ENABLEENERGY Enable the energy equation.
% d.enableEnergy
% :param d:
% Instance of class :mat:func:`Domain1D`
%
disp(' ');
disp('Enabling the energy equation...');
calllib(ct, 'stflow_solveEnergyEqn', d.domainID, 1);
end
function d = disableSoret(d)
% Disable the diffusive mass fluxes due to the Soret effect.
% DISABLESORET Disable the diffusive mass fluxes due to the Soret effect.
% d.disableSoret
% :param d:
% Instance of class :mat:func:`Domain1D`
%
disp(' ');
disp('Disabling the Soret effect...');
calllib(ct, 'stflow_enableSoret', d.domainID, 0);
end
function d = enableSoret(d)
% Enable the diffusive mass fluxes due to the Soret effect.
% ENABLESORET Enable the diffusive mass fluxes due to the Soret effect.
% d.enableSoret
% :param d:
% Instance of class :mat:func:`Domain1D`
%
disp(' ');
disp('Disabling the Soret effect...');
calllib(ct, 'stflow_enableSoret', d.domainID, 1);
@ -135,14 +152,14 @@ classdef Domain1D < handle
%% Domain Get Methods
function b = bounds(d, component)
% Return the (lower, upper) bounds for a solution component.
%
% BOUNDS Return the (lower, upper) bounds for a solution component.
% d.bounds(componoent)
% :param component:
% String name of the component for which the bounds are
% returned.
% :return:
% 1x2 Vector of the lower and upper bounds.
%
n = d.componentIndex(component);
lower = calllib(ct, 'domain_lowerBound', d.domainID, n);
upper = calllib(ct, 'domain_upperBound', d.domainID, n);
@ -150,16 +167,16 @@ classdef Domain1D < handle
end
function n = componentIndex(d, name)
% Get the index of a component given its name
%
% :parameter d:
% Instance of class 'Domain1D'.
% :parameter name:
% String name of the component to look up. If a numeric
% value is passed, it will be :returned directly.
% COMPONENTINDEX Get the index of a component given its name.
% n = d.componentIndex(name)
% :param d:
% Instance of class :mat:func:`Domain1D`
% :param name:
% String name of the component to look up. If a numeric value
% is passed, it will be returned.
% :return:
% Index of the component.
% Index of the component, or input numeric value.
%
if isa(name, 'double')
n = name;
else
@ -175,15 +192,16 @@ classdef Domain1D < handle
end
function n = componentName(d, index)
% Get the name of a component given its index.
%
% :parameter d:
% Instance of class 'Domain1D'.
% :parameter n:
% Integer or vector of integers of components' names.
% COMPONENTNAME Get the name of a component given its index.
% n = d.componentName(index)
% :param d:
% Instance of class :mat:func:`Domain1D`
% :param index:
% Integer or vector of integers of components' names
% to get.
% :return:
% Cell array of component names.
% Cell array of component names.
%
n = length(index);
s = cell(m);
for i = 1:n
@ -200,11 +218,14 @@ classdef Domain1D < handle
end
function i = domainIndex(d)
% Get the domain index.
% DOMAININDEX Get the domain index.
% i = d.domainIndex
% :param d:
% Instance of class :mat:func:`Domain1D`
% :return:
% This function :returns an integer flag denoting the
% location of the domain, beginning with 1 at the left.
% This function returns an integer flag denoting the location
% of the domain, beginning with 1 at the left.
%
i = calllib(ct, 'domain_index', d.domainID);
if i >= 0
i = i + 1;
@ -215,11 +236,14 @@ classdef Domain1D < handle
end
function i = domainType(d)
% Get the type of domain.
% DOMAINTYPE Get the type of domain.
% i = d.domainType
% :param d:
% Instance of class :mat:func:`Domain1D`
% :return:
% This function :returns an integer flag denoting the domain
% type.
% This function returns an integer flag denoting the domain
% type.
%
i = calllib(ct, 'domain_type', d.domainID);
end
@ -248,12 +272,13 @@ classdef Domain1D < handle
end
function a = isFlow(d)
% Determine whether a domain is a flow.
% :parameter d:
% Instance of class 'Domain1D'.
% ISFLOW Determine whether a domain is a flow.
% a = d.isFlow
% :param d:
% Instance of class :mat:func:`Domain1D`
% :return:
% 1 if the domain is a flow domain, and 0 otherwise.
% 1 if the domain is a flow domain, and 0 otherwise.
%
t = d.domainType;
% See Domain1D.h for definitions of constants.
if t < 100
@ -263,12 +288,13 @@ classdef Domain1D < handle
end
function a = isInlet(d)
% Determine whether a domain is an inlet.
% :parameter d:
% Instance of class 'Domain1D'.
% ISINLET Determine whether a domain is an inlet.
% a = d.isInlet
% :param d:
% Instance of class :mat:func:`Domain1D`
% :return:
% 1 if the domain is an inlet, and 0 otherwise.
% 1 if the domain is an inlet, and 0 otherwise.
%
t = d.domainType;
% See Domain1D.h for definitions of constants.
if t == 104
@ -278,12 +304,13 @@ classdef Domain1D < handle
end
function a = isSurface(d)
% Determine whether a domain is a surface.
% :parameter d:
% Instance of class 'Domain1D'.
% ISSURFACE Determine if a domain is a surface.
% a = d.isSurface
% :param d:
% Instance of class :mat:func:`Domain1D`
% :return:
% 1 if the domain is a surface, and 0 otherwise.
% 1 if the domain is a surface, and 0 otherwise.
%
t = d.domainType;
% See Domain1D.h for definitions of constants.
if t == 102
@ -293,26 +320,30 @@ classdef Domain1D < handle
end
function mdot = massFlux(d)
% Get the mass flux.
% MASSFLUX Get the mass flux.
% mdot = d.massFlux
% :param d:
% Instance of class :mat:func:`Domain1D`
% :return:
% The mass flux in the domain.
% The mass flux in the domain.
%
mdot = calllib(ct, 'bdry_mdot', d.domainID);
end
function y = massFraction(d, k)
% Get the mass fraction of a species given its integer index.
% This method :returns the mass fraction of species 'k', where k
% is the integer index of the species in the flow domain to
% which the boundary domain is attached.
% MASSFRACTION Get the mass fraction of a species given its integer index.
% y = d.massFraction(k)
% This method returns the mass fraction of species ``k``, where
% k is the integer index of the species in the flow domain
% to which the boundary domain is attached.
%
% :parameter d:
% Instance of class 'Domain1D'
% :parameter k:
% Integer species index
% :param d:
% Instance of class :mat:func:`Domain1D`
% :param k:
% Integer species index
% :return:
% Mass fraction of species
% Mass fraction of species
%
if d.domainIndex == 0
error('No flow domain attached!')
end
@ -324,25 +355,31 @@ classdef Domain1D < handle
end
function n = nComponents(d)
% Get the number of components.
% NCOMPONENTS Get the number of components.
% n = d.nComponents
% :param d:
% Instance of class :mat:func:`Domain1D`
% :return:
% Number of variables at each grid point.
% Number of variables at each grid point
%
n = calllib(ct, 'domain_nComponents', d.domainID);
end
function n = nPoints(d)
% Get the number of grid points.
% NPOINTS Get the number of grid points.
% n = d.nPoints
% :param d:
% Instance of class :mat:func:`Domain1D`
% :return:
% Integer number of grid points.
% Integer number of grid points.
%
n = calllib(ct, 'domain_nPoints', d.domainID);
end
function tol = tolerances(d, component)
% Return the (relative, absolute) error tolerances for a
% TOLERANCES Return the (relative, absolute) error tolerances for a
% solution component.
%
% tol = d.tolerances(component)
% :param component:
% String name of the component for which the bounds are
% returned.
@ -356,22 +393,37 @@ classdef Domain1D < handle
end
function temperature = get.T(d)
% Get the boundary temperature (K).
% GET.T Get the boundary temperature.
% temperature = d.T
% :param d:
% Instance of class :mat:func:`Domain1D`
% :return:
% Temperature. Units: K
%
temperature = calllib(ct, 'bdry_temperature', d.domainID);
end
function pressure = get.P(d)
% Get the pressure (Pa).
% GET.P Get the boundary pressure.
% pressure = d.P
% :param d:
% Instance of class :mat:func:`Domain1D`
% :return:
% Pressure. Units: Pa
%
pressure = calllibt(ct, 'stflow_pressure', d.domainID);
end
%% Domain Set Methods
function set.T(d, t)
% Set the temperature (K).
% SET.T Set the temperature.
% d.T = t
% :param d:
% Instance of class :mat:func:`Domain1D`
% :param t:
% Temperature to be set. Units: K
%
if t <= 0
error('The temperature must be positive');
end
@ -379,8 +431,13 @@ classdef Domain1D < handle
end
function set.P(d, p)
% Set the pressure (K).
% SET.P Set the pressure.
% d.P = p
% :param d:
% Instance of class :mat:func:`Domain1D`
% :param p:
% Pressure to be set. Units: Pa
%
if p <= 0
error('The pressure must be positive');
end
@ -403,15 +460,18 @@ classdef Domain1D < handle
end
function setCoverageEqs(d, onoff)
% Enable or disable solving the coverage equations.
% :parameter d:
% Instance of class 'Domain1D'
% :parameter onoff:
% string, one of 'on' or 'yes' to turn solving the coverage
% equations on. One of 'off' or 'no' to turn off the
% coverage equation.
if d.domainType ~= 6
% SETBOUNDS Set bounds on the solution components.
% d.setBounds(component, lower, upper)
% :param d:
% Instance of class :mat:func:`Domain1D`
% :param component:
% String, component to set the bounds on
% :param lower:
% Lower bound
% :param upper:
% Upper bound
%
if ~strcmp(d.type, 'ReactingSurface')
error('Wrong domain type. Expected a reacting surface domain.');
end
@ -431,17 +491,19 @@ classdef Domain1D < handle
end
function setFixedTempProfile(d, profile)
% Set a fixed temperature profile to use when the energy
% equation is not being solved. The profile must be entered as
% an array of positions/temperatures, which may be in rows or
% SETFIXEDTEMPPROFILE Set a fixed temperature profile.
% d.setFixedTempProfile(profile)
% Set the temperature profile to use when the
% energy equation is not being solved. The profile must be entered
% as an array of positions / temperatures, which may be in rows or
% columns.
%
% :parameter d:
% Instance of class 'Domain1D'.
% :parameter profile:
% n x 2 or 2 x n array of 'n' points at which the
% temperature is specified.
% :param d:
% Instance of class :mat:func:`Domain1D`
% :param profile:
% n x 2 or 2 x n array of ``n`` points at which the temperature
% is specified.
%
sz = size(profile);
if sz(1) == 2
l = length(profile(1, :));
@ -456,47 +518,55 @@ classdef Domain1D < handle
end
function setID(d, id)
% Set the ID tag for a domain.
% :parameter id:
% String ID to assign.
% SETID Set the ID tag for a domain.
% d.setID(id)
% :param d:
% Instance of class :mat:func:`Domain1D`
% :param id:
% String ID to assign
%
calllib(ct, 'domain_setID', d.domainID, id);
end
function setMdot(d, mdot)
% Set the mass flow rate.
% :parameter mdot:
% Mass flow rate.
% SETMDOT Set the mass flow rate.
% d.setMdot(mdot)
% :param d:
% Instance of class :mat:func:`Domain1D`
% :param mdot:
% Mass flow rate
%
calllib(ct, 'bdry_setMdot', d.domainID, mdot);
end
function setMoleFractions(d, x)
% Set the mole fractions.
% :parameter x:
% String specifying the species and mole fractions in the
% format "'Spec:X,Spec2:X2'".
% SETMOLEFRACTIONS Set the mole fractions.
% d.setMoleFractions(x)
% :param d:
% Instance of class :mat:func:`Domain1D`
% :param x:
% String specifying the species and mole fractions in
% the format ``'SPEC:X,SPEC2:X2'``.
%
calllib(ct, 'bdry_setMoleFractions', d.domainID, x);
end
function setProfileD(d, n, p)
% Set the profile of a component.
% Convenience function to allow an instance of 'Domain1D' to
% have a profile of its components set when it is part of a
% 'Stack'.
% SETPROFILED Set the profile of a component.
% d.setProfileD(n, p)
% Convenience function to allow an instance of :mat:func:`Domain1D` to
% have a profile of its components set when it is part of a :mat:func:`Stack`.
%
% :param d:
% Instance of class :mat:func:`Domain1D`
% :param n:
% Integer index of component, vector of component indices, string
% of component name, or cell array of strings of component names.
% :param p:
% n x 2 array, whose columns are the relative (normalized) positions
% and the component values at those points. The number of positions
% ``n`` is arbitrary.
%
% :parameter d:
% Instance of class 'Domain1d'.
% :parameter n:
% Integer index of component, vector of component indices,
% string of component name, or cell array of strings of
% component names.
% :parameter p:
% n x 2 array, whose columns are the relative (normalized)
% positions and the component values at those points. The
% number of positions 'n' is arbitrary.
if d.stack == 0
error('Install domain in stack before calling setProfile.');
end
@ -504,18 +574,19 @@ classdef Domain1D < handle
end
function setSteadyTolerances(d, component, rtol, atol)
% Set the steady-state tolerances.
% :parameter d:
% Instance of class 'Domain1D'.
% :parameter component:
% String or cell array of strings of component values whose
% tolerances should be set. If 'default' is specified, the
% tolerance of all components will be set.
% :parameter rtol:
% Relative tolerance.
% :parameter atol:
% Absolute tolerance.
% SETSTEADYTOLERANCES Set the steady-state tolerances.
% d.setSteadyTolerances(component, rtol, atol)
% :param d:
% Instance of class :mat:func:`Domain1D`
% :param component:
% String or cell array of strings of component values
% whose tolerances should be set. If ``'default'`` is
% specified, the tolerance of all components will be set.
% :param rtol:
% Relative tolerance
% :param atol:
% Absolute tolerance
%
if strcmp(component, 'default')
nc = d.nComponents;
for ii = 1:nc
@ -537,18 +608,19 @@ classdef Domain1D < handle
end
function setTransientTolerances(d, component, rtol, atol)
% Set the transient tolerances.
% :parameter d:
% Instance of class 'Domain1D'.
% :parameter component:
% String or cell array of strings of component values whose
% tolerances should be set. If 'default' is specified, the
% tolerance of all components will be set.
% :parameter rtol:
% Relative tolerance.
% :parameter atol:
% Absolute tolerance.
% SETTRANSIENTTOLERANCES Set the transient tolerances.
% d.setTransientTolerances(component, rtol, atol)
% :param d:
% Instance of class :mat:func:`Domain1D`
% :param component:
% String or cell array of strings of component values
% whose tolerances should be set. If ``'default'`` is
% specified, the tolerance of all components will be set.
% :param rtol:
% Relative tolerance
% :param atol:
% Absolute tolerance
%
if strcmp(component, 'default')
nc = d.nComponents;
for ii = 1:nc
@ -570,19 +642,22 @@ classdef Domain1D < handle
end
function setTransport(d, itr)
% Set the solution object used for calculating transport
% properties.
%
% SETTRANSPORT Set the solution object used for calculating transport properties.
% d.setTransport(itr)
% :param itr:
% ID of the solution object for which transport properties
% are calculated.
%
calllib(ct, 'stflow_setTransport', d.domainID, itr);
end
function setupGrid(d, grid)
% Set up the solution grid.
% SETUPGRID Set up the solution grid.
% d.setupGrid(grid)
% :param d:
% Instance of class :mat:func:`Domain1D`
% :param grid:
%
calllib(ct, 'domain_setupGrid', d.domainID, numel(grid), grid);
end

View File

@ -1,6 +1,6 @@
function m = FreeFlame(gas, id)
% Create a freely-propagating flat flame.
m = Domain1D(1, gas, 2);
m = Domain1D('StagnationFlow', gas, 2);
if nargin == 1
m.setID('flame');
else

View File

@ -1,6 +1,6 @@
function m = Inlet(id)
% Create an inlet domain.
m = Domain1D(2);
m = Domain1D('Inlet1D');
if nargin == 0
m.setID('inlet');
else

View File

@ -1,6 +1,6 @@
function m = Outlet(id)
% Create an outlet domain.
m = Domain1D(5);
m = Domain1D('Outlet1D');
if nargin == 0
m.setID('outlet');
else

View File

@ -1,6 +1,6 @@
function m = OutletRes(id)
% Create an outlet reservoir domain.
m = Domain1D(-2);
m = Domain1D('OutletRes');
if nargin == 0
m.setID('outletres');
else

View File

@ -26,7 +26,7 @@ classdef Stack < handle
nd = length(domains);
ids = zeros(1, nd);
for n=1:nd
ids(n) = domains(n).domID;
ids(n) = domains(n).domainID;
end
s.stID = calllib(ct, 'sim1D_new', nd, ids);
else

View File

@ -4,14 +4,14 @@ function m = Surface(id, surface_mech)
% Instance of class 'Interface' defining the surface reaction
% mechanism to be used. Optional.
if nargin < 2
m = Domain1D(3);
m = Domain1D('Surf1D');
if nargin == 0
m.setID('surface');
elseif nargin == 1
m.setID(id);
end
else
m = Domain1D(6, surface_mech);
m = Domain1D('ReactingSurface', surface_mech);
m.setID(id);
end
end

View File

@ -1,6 +1,6 @@
function m = SymmPlane(id)
% Create an symmetry plane domain.
m = Domain1D(4);
m = Domain1D(Symm1D);
if nargin == 0
m.setID('symmetry_plane');
else

View File

@ -9,15 +9,27 @@ classdef Interface < handle & ThermoPhase & Kinetics
%% Interface class constructor
function s = Interface(src, id, p1, p2, p3, p4)
% :parameter src:
% CTI or CTML file containing the interface or edge phase.
% :parameter id:
% Name of the interface or edge phase in the source file.
% :parameter p1/p2/p3/p4:
% Adjoining phase to the interface;
% INTERFACE Interface class constructor.
% s = Interface(src, id, p1, p2, p3, p4)
% See `Interfaces <https://cantera.org/tutorials/cti/phases.html#interfaces>`__.
%
% See also: :mat:func:`importEdge`, :mat:func:`importInterface`
%
% :param src:
% CTI or CTML file containing the interface or edge phase.
% :param id:
% Name of the interface or edge phase in the CTI or CTML file.
% :param p1:
% Adjoining phase to the interface.
% :param p2:
% Adjoining phase to the interface.
% :param p3:
% Adjoining phase to the interface.
% :param p4:
% Adjoining phase to the interface.
% :return:
% Instance of class 'Interface'.
% Instance of class :mat:func:`Interface`
%
checklib;
t = ThermoPhase(src, id);
s@ThermoPhase(src, id);
@ -39,11 +51,15 @@ classdef Interface < handle & ThermoPhase & Kinetics
%% Interface methods
function c = get.coverages(s)
% Get the surface coverages of the species on an interface.
%
% GET.COVERAGES Get the surface coverages of the species on an interface.
% c = s.coverages
% :param s:
% Instance of class :mat:func:`Interface` with surface species
% :return:
% Vector of length "n_surf_species" for coverage.
% If no output value is assigned, a bar graph will be plotted.
% Otherwise, a vector of length ``n_surf_species`` will be
% returned.
%
surfID = s.tpID;
nsp = s.nSpecies;
xx = zeros(1, nsp);
@ -53,12 +69,14 @@ classdef Interface < handle & ThermoPhase & Kinetics
end
function d = get.siteDensity(s)
% Get the site density.
%
% GET.SITEDENSITY Get the surface coverages of the species on an interface.
% c = s.siteDensity
% :param s:
% Instance of class :mat:func:`Interface` with surface species
% :return:
% Double site density. Unit: kmol/m^2 for surface phases,
% kmol/m for edge phases.
%
surfID = s.tpID;
d = calllibt(ct, 'surf_siteDensity', surfID);
end
@ -78,13 +96,22 @@ classdef Interface < handle & ThermoPhase & Kinetics
end
function set.coverages(s, cov, norm)
% Set surface coverages of the species on an interface.
% SET.COVERAGES Set surface coverages of the species on an interface.
% s.coverages = (cov, norm)
% :param s:
% Instance of class :mat:func:`Interface`
% :param cov:
% Coverage of the species. ``cov`` can be either a vector of
% length ``n_surf_species``, or a string in the format
% ``'Species:Coverage, Species:Coverage'``
% :param norm:
% Optional flag that denotes whether or not to normalize the species
% coverages. ``norm`` is either of the two strings ``'nonorm'``` or
% ``'norm'``. If left unset, the default is `norm`. This only works if
% ``s`` is a vector, not a string. Since unnormalized coverages can lead to
% unphysical results, ``'nonorm'`` should be used only in rare cases, such
% as computing partial derivatives with respect to a species coverage.
%
% :parameter cov:
% Coverage of the species. "Cov" can be either a vector of
% length "n_surf_species", or a string in the format
% "Species:Coverage, Species:Coverage".
if nargin == 3 && strcmp(norm, 'nonorm')
norm_flag = 0;
else
@ -111,9 +138,11 @@ classdef Interface < handle & ThermoPhase & Kinetics
end
function set.siteDensity(s, d)
% Set surface site densities.
%
% :parameter density:
% SET.SITEDENSITY Set the site density of a phase on an interface.
% s.siteDensity = d
% :param s:
% Instance of class :mat:func:`Interface`
% :parameter d
% Double site density. Unit: kmol/m^2 for surface phases,
% kmol/m for edge phases.

View File

@ -6,8 +6,11 @@ classdef Kinetics < handle
Kf % forward reaction rate
Kr % reverse reaction rate
dH % enthalpy of reaction
dHss % standard state enthalpy of reaction
dS % entropy of reaction
dSss % standard state entropy of reaction
dG % gibbs free energy of reaction
dGss % standard state gibbs free energy of reaction
end
methods

View File

@ -2,42 +2,48 @@ classdef Mixture < handle
properties
mixID
phases
T
P
phases % Phases in the mixture
T % Temperature
P % Pressure
end
methods
%% Mixture class constructor
function m = Mixture(phases)
% To construct a mixture, supply a cell array of phases and mole
% numbers:
% >> air = Solution('air.yaml');
% >> graphite = Solution('graphite.yaml');
% >> mix = Mixture({air, 1.0; graphite, 0.1});
% MIXTURE Multiphase mixture class constructor.
% m = Mixture(phases)
% Class :mat:func:`Mixture` represents mixtures of one or more phases of matter.
% To construct a mixture, supply a cell array of phases and
% mole numbers::
%
% Phases may also be added later using the addPhase method:
% >> water = Solution('water.cti');
% >> mix.addPhase(water, 3.0);
% >> gas = Solution('gas.yaml');
% >> graphite = Solution('graphite.yaml');
% >> mix = Mixture({gas, 1.0; graphite, 0.1});
%
% Note that the objects representing each phase compute only
% the intensive state of the phase - they do not store any
% information on the amount of this phase. Mixture objects, on
% the other hand, represent full extensive state.
% Phases may also be added later using the addPhase method::
%
% Mixture objects are lightweight in the sense that they do
% not store :parameters needed to compute thermodynamic or
% kinetic properties of the phases. These are contained in the
% ('heavyweight') phase objects. Multiple mixture objects are
% constructed using the same set of phase objects. Each one
% stores its own state information locally, and synchronizes the
% phase objects whenever it requires phase properties.
% >> water = Solution('water.yaml');
% >> mix.addPhase(water, 3.0);
%
% :parameter phases:
% Cell array of phases and mole numbers.
% Note that the objects representing each phase compute only the
% intensive state of the phase - they do not store any information
% on the amount of this phase. Mixture objects, on the other hand,
% represent the full extensive state.
%
% Mixture objects are 'lightweight' in the sense that they do not
% store parameters needed to compute thermodynamic or kinetic
% properties of the phases. These are contained in the
% ('heavyweight') phase objects. Multiple mixture objects may be
% constructed using the same set of phase objects. Each one stores
% its own state information locally, and synchronizes the phase
% objects whenever it requires phase properties.
%
% :param phases:
% Cell array of phases and mole numbers
% :return:
% Instance of class 'Mixture'.
% Instance of class :mat:func:`Mixture`
%
checklib;
@ -72,8 +78,11 @@ classdef Mixture < handle
%% Utility methods
function display(m)
% Display the state of the mixture on the terminal.
% DISPLAY Display the state of the mixture on the terminal.
% m.display
% :param self:
% Instance of class :mat:func:`Mixture`
%
calllib(ct, 'mix_updatePhases', m.mixID);
[np, nc] = size(m.phases);
for n = 1:np
@ -86,21 +95,26 @@ classdef Mixture < handle
end
function clear(m)
% Delete the MultiPhase object.
% CLEAR Delete the C++ MultiPhase object.
% m.clear
% :param m:
% Instance of class :mat:func:`Mixture`
%
calllib(ct, 'mix_del', m.mixID);
end
function addPhase(m, phase, moles)
% Add a phase to the mixture
% ADDPHASE Add a phase to a mixture.
% addPhase(self, phase, moles)
% :param self:
% Instance of class :mat:func:`Mixture` to which phases should be
% added
% :param phase:
% Instance of class :mat:func:`ThermoPhase` which should be added
% :param moles:
% Number of moles of the ``phase`` to be added to this mixture.
% Units: kmol
%
% :parameter m:
% Instance of class 'Mixture' to which phases is added.
% :parameter phase:
% Instance of class 'ThermoPhase' which should be added.
% :parameter moles:
% Number of moles of the phase to be added. Unit: kmol.
if ~isa(phase, 'ThermoPhase')
error('Phase object of wrong type.');
end
@ -121,84 +135,123 @@ classdef Mixture < handle
%% Mixture Get methods
function temperature = get.T(m)
% Get the temperature of the mixture.
%
% GET.T Get the temperature of a mixture.
% temperature = m.T
% :param m:
% Instance of class :mat:func:`Mixture`
% :return:
% Temperature in K.
% Temperature (K)
%
temperature = calllib(ct, 'mix_temperature', m.mixID);
end
function pressure = get.P(m)
% Get the pressure of themixture.
%
% GET.P Get the pressure of the mixture.
% pressure = m.P
% :param m:
% Instance of class :mat:func:`Mixture`
% :return:
% Pressure in Pa.
% Pressure. Units: Pa
%
pressure = calllib(ct, 'mix_pressure', m.mixID);
end
function n = nAtoms(m, k, e)
% Get the number of atoms of element e in species k.
function n = nAtoms(m, e)
% NATOMS Get the number of atoms of an element in a mixture.
% n = m.nAtoms(e)
% :param m:
% Instance of class :mat:func:`Mixture`
% :param e:
% Index of element
% :return:
% Number of atoms for element e
%
% Note: In keeping with the conventions used by Matlab, the
% indices start from 1 instead of 0 as in Cantera C++ and
% Python interfaces.
%
n = calllib(ct, 'mix_nPhases', m.mixID, k-1, e-1);
% Check back on this one!
end
function n = nElements(m)
% Get the number of elements in the mixture.
% NELEMENTS Get the number of elements in a mixture.
% n = m.nElements
% :param m:
% Instance of class :mat:func:`Mixture`
% :return:
% Number of elements in the input
%
n = calllib(ct, 'mix_nElements', m.mixID);
end
function n = nPhases(m)
% Get the number of phases in the mixture.
% NPHASES Get the number of phases in a mixture.
% n = m.nPhases
% :param m:
% Instance of class :mat:func:`Mixture`
% :return:
% Number of phases in the input
%
n = calllib(ct, 'mix_nPhases', m.mixID);
end
function n = nSpecies(m)
% Get the number of species in the mixture.
% NSPECIES Get the number of species in a mixture.
% n = m.nSpecies
% :param m:
% Instance of class :mat:func:`Mixture`
% :return:
% Number of species in the input
%
n = calllib(ct, 'mix_nSpecies', m.mixID);
end
function n = elementIndex(m, name)
% Get the index of element 'name'.
% ELEMENTINDEX Get the index of an element.
% n = m.elementIndex(name)
% :param m:
% Instance of class :mat:func:`Mixture`
% :param name:
% Name of the element whose index is desired
% :return:
% Index of element with name ``name``
%
% Note: In keeping with the conventions used by Matlab, the
% indices start from 1 instead of 0 as in Cantera C++ and
% Python interfaces.
%
n = calllib(ct, 'mix_elementIndex', m.mixID, name) + 1;
end
function n = speciesIndex(m, k, p)
% Get the index of species k in phase p.
% SPECIESINDEX Get the index of a species in a mixture.
% n = m.speciesIndex(k, p)
% :param m:
% Instance of class :mat:func:`Mixture`
% :param name:
% Name of the speces whose index is desired
% :return:
% Index of species with name ``name``
%
% :param k:
% Double
% Note: In keeping with the conventions used by Matlab, the
% indices start from 1 instead of 0 as in Cantera C++ and
% Python interfaces.
%
n = calllib(ct, 'mix_speciesIndex', m.mixID, k-1, p-1) + 1;
% Check back on this one!
end
function moles = elementMoles(m, e)
% Get the number of moles of an element in a mixture.
%
% :parameter e:
% ELEMENTMOLES Get the number of moles of an element in a mixture.
% moles = m.elementMoles(e)
% :param m:
% Instance of class :mat:func:`Mixture`
% :param e:
% Integer element number.
% :return:
% Moles of element number 'e'. If input 'e' is empty, return
% moles of every element in the mixture. Unit: kmol.
%
if nargin == 2
moles = calllib(ct, 'mix_elementMoles', m.mixID, e)
elseif nargin == 1
@ -212,14 +265,16 @@ classdef Mixture < handle
end
function moles = phaseMoles(m, n)
% Get the number of moles of a phase in a mixture.
%
% :parameter n:
% PHASEMOLES Get the number of moles of a phase in a mixture.
% moles = m.phaseMoles(n)
% :param m:
% Instance of class :mat:func:`Mixture`
% :param n:
% Integer phase number.
% :return:
% Moles of phase number 'n'. If input 'n' is empty, return
% moles of phase element in the mixture. Unit: kmol.
% moles of every phase in the mixture. Unit: kmol.
%
if nargin == 2
moles = calllib(ct, 'mix_phaseMoles', m.mixID, n);
elseif nargin == 1
@ -234,9 +289,11 @@ classdef Mixture < handle
end
function moles = speciesMoles(m, k)
% Get the number of moles of a species in a mixture.
%
% :parameter k:
% SPECIESMOLES Get the number of moles of a species in a mixture.
% moles = m.speciesMoles(n)
% :param m:
% Instance of class :mat:func:`Mixture`
% :param k:
% Integer species number.
% :return:
% Moles of species number 'k'. If input 'k' is empty, return
@ -256,11 +313,13 @@ classdef Mixture < handle
end
function mu = chemPotentials(m)
% Get the chemical potentials of species in the mixture.
%
% CHEMPOTENTIALS Get the chemical potentials of species in a mixture.
% mu = m.chemPotentials
% :param m:
% Instance of class :mat:func:`Mixture`
% :return:
% Vector of chemical potentials. Unit: J/kmol.
% Vector of chemical potentials. Units: J/kmol
%
nsp = m.nSpecies;
xx = zeros(1, nsp);
ptr = libpointer('doublePtr', xx);
@ -271,44 +330,56 @@ classdef Mixture < handle
%% Mixture Set methods
function m = set.T(m, temp)
% Set the mixture temperature.
% SET.T Set the mixture temperature.
% m.T = temp
% :param m:
% Instance of class :mat:func:`Mixture`
% :param temp:
% Temperature. Units: K
%
% :parameter temp:
% Temperature to set. Unit: K.
calllib(ct, 'mix_setTemperature', m.mixID, temp);
end
function m = set.P(m, pressure)
% Set the mixture pressure.
% SET.P Set the pressure of the mixture.
% m.P = pressure
% :param m:
% Instance of class :mat:func:`Mixture`
% :param pressure:
% Pressure. Units: Pa
%
% :parameter pressure:
% Pressure to set. Unit: Pa.
calllib(ct, 'mix_setPressure', m.mixID, pressure);
end
function setPhaseMoles(m, n, moles)
% Set the number of moles of phase n in the mixture.
% SETPHASEMOLES Set the number of moles of a phase in a mixture.
% m.setPhaseMoles(n, moles)
% :param m:
% Instance of class :mat:func:`Mixture`
% :param n:
% Phase number in the input
% :param moles:
% Number of moles to add. Units: kmol
%
% :parameter n:
% Phase number.
% :parameter moles:
% Number of moles to set. Unit: kmol.
calllib(ct, 'mix_setPhaseMoles', m.mixID, n-1, moles);
end
function setSpeciesMoles(m, moles)
% Set the moles of the species. The moles may be specified as a
% string or as a vector. If a vector is used, it must be
% dimensioned at least as large as the total number of species
% in the mixture. Note that the species may belong to any
% phase, and unspecified species are set to zero.
% SETSPECIESMOLES Set the moles of the species.
% m.setSpeciesMoles(moles)
% Set the moles of the species in kmol. The moles may
% be specified either as a string, or as an vector. If a vector is
% used, it must be dimensioned at least as large as the total number
% of species in the mixture. Note that the species may belong to any
% phase, and unspecified species are set to zero. ::
%
% >> mix.setSpeciesMoles('C(s):1.0, CH4:2.0, O2:0.2');
%
% :param m:
% Instance of class :mat:func:`Mixture`
% :param moles:
% Vector or string specifying the moles of species
%
% :parameter moles:
% Vector or string specifying the moles of species.
if isa(moles, 'double')
l = length(moles);
calllib(ct, 'mix_setMoles', m.mixID, l, moles);
@ -320,44 +391,49 @@ classdef Mixture < handle
end
function r = equilibrate(m, XY, err, maxsteps, maxiter, loglevel)
% Set the mixture to a state of chemical equilibrium.
%
% EQUILIBRATE Set the mixture to a state of chemical equilibrium.
% m.equilibrate(XY, err, maxsteps, maxiter, loglevel)
% This method uses a version of the VCS algorithm to find the
% composition that minimizes the total Gibbs free energy of the
% mixture, subject to element conservation constraints. For a
% description of the theory, see Smith and Missen, "Chemical
% Reaction Equilibrium". The VCS algorithm is implemented in
% Reaction Equilibrium." The VCS algorithm is implemented in
% Cantera kernel class MultiPhaseEquil.
%
% The VCS algorithm solves for the equilibrium composition for
% specified temperature and pressure. If any other property
% pair other than 'TP' is specified, then an outer iteration
% loop is used to adjust T and/or P so that the specified
% property values are obtained:
% >> equilibrate(mix, 'TP);
% >> equilibrate('TP', 1.0e-6, 500);
% specified temperature and pressure. If any other property pair
% other than "TP" is specified, then an outer iteration loop is
% used to adjust T and/or P so that the specified property
% values are obtained. ::
%
% :parameter XY:
% Two-letter string specifying the two properties to hold
% fixed. Currently 'TP', 'HP', 'TV', and 'SP' have been
% implemented. Default: 'TP'.
% :parameter err:
% Error tolerance. Iteration will continue until delta_Mu/RT
% is less than this value for each reaction. Default:
% 1.0e-9.
% :parameter maxsteps:
% Maximum number of steps to take while solving the
% equilibrium problem for specified T and P. Default: 1000.
% :parameter maxiter:
% Maximum number of temperature and/or pressure iterations.
% This is only relevant if a property pair other than (T,
% P)is specified. Default: 200.
% :parameter loglevel:
% Set to a value > 0 to write diagnostic output. Larger
% values generate more detailed information.
% >> mix.equilibrate('TP')
% >> mix.equilibrate('TP', 1.0e-6, 500)
%
% :param m:
% Instance of class :mat:func:`Mixture`
% :param XY:
% Two-letter string specifying the two properties to hold
% fixed. Currently, ``'TP'``, ``'HP'``, ``'TV'``, and ``'SP'`` are
% implemented. Default: ``'TP'``.
% :param err:
% Error tolerance. Iteration will continue until :math:`\Delta\mu)/RT`
% is less than this value for each reaction. Default:
% 1.0e-9. Note that this default is very conservative, and good
% equilibrium solutions may be obtained with larger error
% tolerances.
% :param maxsteps:
% Maximum number of steps to take while solving the
% equilibrium problem for specified T and P. Default: 1000.
% :param maxiter:
% Maximum number of temperature and/or pressure
% iterations. This is only relevant if a property pair other
% than (T,P) is specified. Default: 200.
% :param loglevel:
% Set to a value > 0 to write diagnostic output.
% Larger values generate more detailed information.
% :return:
% The error in the solution.
% The error in the solution
%
if nargin < 6
loglevel = 0;
end

View File

@ -1,10 +1,58 @@
classdef Solution < handle & ThermoPhase & Kinetics & Transport
properties
tp
end
methods
% Solution class constructor
%% Solution class constructor
function s = Solution(src, id, trans)
% SOLUTION Solution class constructor.
% s = Solution(src, id, trans)
% Class :mat:func:`Solution` represents solutions of multiple species. A
% solution is defined as a mixture of two or more constituents
% (species) that are completely mixed on molecular length
% scales. The macroscopic intensive thermodynamic state of a
% solution is specified by two thermodynamic properties (for
% example, the temperature and pressure), and the relative amounts
% of each species, which may be given as mole fractions or mass
% fractions. ::
%
% >> s = Solution('input.yaml'[, phase_name[, transport_model]])
%
% constructs a :mat:func:`Solution` object from a specification contained in
% file ``input.yaml``. Optionally, the name of the phase to be imported
% can be specified with ``phase_name``. If a :mat:func:`Transport` model is
% included in ``input.yaml``, it will be included in the :mat:func:`Solution`
% instance with the default transport modeling as set
% in the input file. To specify the transport modeling, set the input
% argument ``trans`` to one of ``'default'``, ``'None'``, ``'Mix'``, or ``'Multi'``.
% In this case, the phase name must be specified as well. Alternatively,
% change the ``transport`` node in the CTML file, or ``transport``
% property in the CTI file before loading the phase. The transport
% modeling cannot be changed once the phase is loaded.
%
% Class :mat:func:`Solution` derives from three more basic classes, and most of
% its methods are inherited from these classes. These are:
%
% * class :mat:func:`ThermoPhase` - composition information and thermodynamic properties
% * class :mat:func:`Kinetics` - homogeneous kinetics
% * class :mat:func:`Transport` - transport properties
%
% See also: :mat:func:`ThermoPhase`, :mat:func:`Kinetics`, :mat:func:`Transport`
%
% :param src:
% Input string of CTI or CTML file name.
% :param id:
% Optional unless ``trans`` is specified. ID of the phase to
% import as specified in the CTML or CTI file.
% :param trans:
% String, transport modeling. Possible values are ``'default'``, ``'None'``,
% ``'Mix'``, or ``'Multi'``. If not specified, ``'default'`` is used.
% :return:
% Instance of class :mat:func:`Solution`
if nargin == 1
id = '-';
end
@ -25,6 +73,11 @@ classdef Solution < handle & ThermoPhase & Kinetics & Transport
% Delete the kernel objects associated with a solution
function clear(s)
% CLEAR Delete the kernel objects associated with a Solution.
% s.clear
% :param s:
% Instance of class :mat:func:`Solution`
%
s.tpClear;
s.kinClear;
s.trClear;

File diff suppressed because it is too large Load Diff

View File

@ -12,54 +12,58 @@ classdef Func < handle
%% Functor class constructor
function x = Func(typ, n, p)
% Functor class constructor.
%
% FUNC Func class constructor.
% 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.
% 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 be 3.
% 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;
%
% 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.
% 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.
% Note: this MATLAB class shadows the underlying C++ Cantera class
% "Func1". See the Cantera C++ documentation for more details.
%
% :parameter typ:
% String indicating type of functor to create. Possible
% values are:
% * 'polynomial'
% * 'fourier'
% * 'gaussian'
% * 'arrhenius'
% * 'sum'
% * 'diff'
% * 'ratio'
% * 'composite'
% * 'periodic'
% See also: :mat:func:`polynom`, :mat:func:`gaussian`, :mat:func:`plus`,
% :mat:func:`rdivide`, :mat:func:`times`
%
% :parameter n:
% Number of :parameters required for the functor
% :parameter p:
% Vector of :parameters
% :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:func:`Func`
% Instance of class :mat:func:`Func`
%
if ~isa(typ, 'char')
error('Function type must be a string');
end
@ -125,34 +129,37 @@ classdef Func < handle
%% Utility methods
function clear(f)
% Clear the functor from memory.
% CLEAR Delete the C++ Func1 object.
% f.clear
% :param f:
% Instance of class :mat:func:`Func`
%
calllib(ct, 'func_del', f.id);
end
function display(a)
% Display the equation of the input function on the terminal.
function display(f)
% DISPLAY Display the equation of the input function on the terminal.
% f.display
% :param f:
% Instance of class :mat:func:`Func`
%
% :parameter a:
% Instance of class 'Func'
disp(' ');
disp([inputname(1),' = '])
disp(' ');
disp([' ' char(a)])
disp([' ' f.char])
disp(' ');
end
function b = subsref(a, s)
% Redefine subscripted references for functors.
%
% :parameter a:
% Instance of class 'Func'
% :parameter s:
% Value at which the function should be evaluated.
% SUBSREF Redefine subscripted references for functors.
% b = a.subsref(s)
% :param a:
% Instance of class :mat:func:`Func`
% :param s:
% Value at which the function should be evaluated.
% :return:
% The value of the function evaluated at 's'.
% Returns the value of the function evaluated at ``s``
%
if strcmp(s.type, '()')
ind= s.subs{:};
b = zeros(1, length(ind));
@ -163,25 +170,31 @@ classdef Func < handle
end
end
function s = char(p)
% Get the formatted string to display the function.
if strcmp(p.typ,'sum')
s = ['(' char(p.f1) ') + (' char(p.f2) ')'];
elseif strcmp(p.typ,'diff')
s = ['(' char(p.f1) ') - (' char(p.f2) ')'];
elseif strcmp(p.typ,'prod')
s = ['(' char(p.f1) ') * (' char(p.f2) ')'];
elseif strcmp(p.typ,'ratio')
s = ['(' char(p.f1) ') / (' char(p.f2) ')'];
function s = char(f)
% CHAR Get the formatted string to display the function.
% s = f.char
% :param f:
% Instance of class :mat:func:`Func`
% :return:
% Formatted string displaying the function
%
if strcmp(f.typ,'sum')
s = ['(' char(f.f1) ') + (' char(f.f2) ')'];
elseif strcmp(f.typ,'diff')
s = ['(' char(f.f1) ') - (' char(f.f2) ')'];
elseif strcmp(f.typ,'prod')
s = ['(' char(f.f1) ') * (' char(f.f2) ')'];
elseif strcmp(f.typ,'ratio')
s = ['(' char(f.f1) ') / (' char(f.f2) ')'];
elseif all(p.coeffs == 0)
s = '0';
else
if strcmp(p.typ,'polynomial')
if strcmp(f.typ,'polynomial')
d = length(p.coeffs) - 1;
s = [];
nn = 0;
for b = p.coeffs;
cc(d+1-nn) = b;
for b = f.coeffs;
cc(d + 1 - nn) = b;
nn = nn + 1;
end
for a = cc;
@ -208,12 +221,12 @@ classdef Func < handle
end
d = d - 1;
end
elseif strcmp(p.typ, 'gaussian')
s = ['Gaussian(' num2str(p.coeffs(1)) ',' ...
num2str(p.coeffs(2)) ',' ...
num2str(p.coeffs(3)) ')'];
elseif strcmp(p.typ, 'fourier')
c = reshape(p.coeffs, [], 2);
elseif strcmp(f.typ, 'gaussian')
s = ['Gaussian(' num2str(f.coeffs(1)) ',' ...
num2str(f.coeffs(2)) ',' ...
num2str(f.coeffs(3)) ')'];
elseif strcmp(f.typ, 'fourier')
c = reshape(f.coeffs, [], 2);
Ao = c(1, 1);
w = c(1, 2);
A = c(2:end, 1);
@ -249,7 +262,7 @@ classdef Func < handle
end
end
else
s = ['*** char not yet implemented for' p.typ ' ***'];
s = ['*** char not yet implemented for' f.typ ' ***'];
end
end
end

View File

@ -1,6 +1,6 @@
function LoadCantera
addpath('Class', 'Utility', 'PresetFlowDevices', 'PresetFunctors', ...
'PresetMixtures', 'PresetReactors', 'PresetObjects', '1D');
'PresetMixtures', 'PresetReactors', '1D');
if ispc
ctname = 'cantera_shared.dll';
elseif ismac
@ -14,6 +14,7 @@ function LoadCantera
end
if ~libisloaded(ct)
load('Utility/cantera_root.mat');
addpath([cantera_root, '/samples/matlab_new']);
[~,warnings] = loadlibrary([cantera_root '/lib/' ctname], ...
[cantera_root '/include/cantera/clib/ctmatlab.h'], ...
'includepath', [cantera_root '/include'], ...

View File

@ -11,12 +11,27 @@ classdef FlowDevice < handle
%% FlowDevice class constructor
function x = FlowDevice(typ)
% Flow Device class constructor.
% FLOWDEVICE FlowDevice class constructor.
% x = FlowDevice(typ)
% Base class for devices that allow flow between reactors.
% :mat:func:`FlowDevice` objects are assumed to be adiabatic,
% non-reactive, and have negligible internal volume, so that they are
% internally always in steady-state even if the upstream and downstream
% reactors are not. The fluid enthalpy, chemical composition, and mass
% flow rate are constant across a :mat:func:`FlowDevice`, and the
% pressure difference equals the difference in pressure between the
% upstream and downstream reactors.
%
% See also: :mat:func:`MassFlowController`, :mat:func:`Valve`
%
% :param typ:
% Type of :mat:func:`FlowDevice` to be created. ``typ='MassFlowController'``
% for :mat:func:`MassFlowController`, ``typ='PressureController'`` for
% :mat:func:`PressureController` and ``typ='Valve'`` for
% :mat:func:`Valve`
% :return:
% Instance of class :mat:func:`FlowDevice`
%
% :parameter typ:
% Type of flow device to be created. Type =
% 'MassFlowController', 'PressureController' or 'Valve'.
checklib;
if nargin == 0
@ -25,9 +40,6 @@ classdef FlowDevice < handle
x.type = typ;
x.id = calllib(ct, 'flowdev_new', typ);
% if x.id < 0
% error(geterr);
% end
x.upstream = -1;
x.downstream = -1;
end
@ -35,21 +47,28 @@ classdef FlowDevice < handle
%% Utility methods
function clear(f)
% Clear the specified flow device from memory.
% CLEAR Clear the specified flow device from memory.
% f.clear
% :param f:
% Instance of :mat:func:`FlowDevice` to be cleared.
%
calllib(ct, 'flowdev_del', f.id);
end
%% FlowDevice methods
function install(f, upstream, downstream)
% Install a flow device between reactors or reservoirs.
% INSTALL Install a flow device between reactors or reservoirs.
% f.install(upstream, downstream)
% :param f:
% Instance of class :mat:func:`FlowDevice` to install
% :param upstream:
% Upstream :mat:func:`Reactor` or :mat:func:`Reservoir`
% :param downstream:
% Downstream :mat:func:`Reactor` or :mat:func:`Reservoir`
% :return:
% Instance of class :mat:func:`FlowDevice`
%
% :parameter upstream:
% Upsteram 'Reactor' or 'Reservoir'.
% :parameter downstream:
% Downstream 'Reactor' or 'Reservoir'.
if nargin == 3
if ~isa(upstream, 'Reactor') || ~isa(downstream, 'Reactor')
error(['Flow devices can only be installed between',...
@ -58,59 +77,67 @@ classdef FlowDevice < handle
i = upstream.id;
j = downstream.id;
calllib(ct, 'flowdev_install', f.id, i, j);
% if ok < 0
% error(geterr)
% end
else error('install requires 3 arguments');
end
end
function mdot = massFlowRate(f)
% Get the mass flow rate.
%
% MASSFLOWRATE Get the mass flow rate.
% mdot = f.massFlowRate
% :param f:
% Instance of class :mat:func:`MassFlowController`
% :return:
% The mass flow rate through the flow device
% The mass flow rate through the :mat:func:`FlowDevice` at the current time
%
mdot = calllib(ct, 'flowdev_massFlowRate2', f.id);
end
function setFunction(f, mf)
% Set the time function with class 'func'.
% SETFUNCTION Set the mass flow rate with class :mat:func:`Func`.
% f.setFunction(mf)
%
% See also: :mat:func:`MassFlowController`, :mat:func:`Func`
%
% :param f:
% Instance of class :mat:func:`MassFlowController`
% :param mf:
% Instance of class :mat:func:`Func`
%
% :parameter mf:
% Instance of class 'func'.
if strcmp(f.type, 'MassFlowController')
k = calllib(ct, 'flowdev_setTimeFunction', f.id, ...
mf.id);
% if k < 0
% error(geterr);
% end
else
error('Time function can only be set for mass flow controllers.');
end
end
function setMassFlowRate(f, mdot)
% Set the mass flow rate to a constant value.
% SETMASSFLOWRATE Set the mass flow rate to a constant value.
% f.setMassFlowRate(mdot)
%
% See also: :mat:func:`MassFlowController`
%
% :param f:
% Instance of class :mat:func:`MassFlowController`
% :param mdot:
% Mass flow rate
%
% :parameter mdot:
% Mass flow rate
if strcmp(f.type, 'MassFlowController')
k = calllib(ct, 'flowdev_setMassFlowCoeff', f.id, mdot);
% if k < 0
% error(geterr);
% end
else
error('Mass flow rate can only be set for mass flow controllers.');
end
end
function setMaster(f, d)
% Set the Master flow device used to compute this device's mass
% SETMASTER Set the Master flow device used to compute this device's mass
% flow rate.
% f.setMaster(d)
% :param f:
% Instance of class :mat:func:`MassFlowController`
% :param mf:
% Instance of class :mat:func:`Func`
%
if strcmp(f.type, 'PressureController')
k = calllib(ct, 'flowdev_setMaster', f.id, d);
else
@ -119,23 +146,26 @@ classdef FlowDevice < handle
end
function setValveCoeff(f, k)
% Set the valve coefficient 'K'.
%
% SETVALVECOEFF Set the valve coefficient :math:`K`.
% f.setValveCoeff(k)
% The mass flow rate [kg/s] is computed from the expression
% mdot = K(P_upstream - P_downstream)
% as long as this produces a positive value. If thsi expression
% is negative, zero is :returned.
%
% :parameter k:
% Value fo the valve coefficient. Unit: kg/Pa-s.
% .. math:: \dot{m} = K(P_{upstream} - P_{downstream})
%
% as long as this produces a positive value. If this expression is
% negative, zero is returned.
%
% See also: :mat:func:`Valve`
%
% :param f:
% Instance of class :mat:func:`Valve`
% :param k:
% Value of the valve coefficient. Units: kg/Pa-s
%
if ~strcmp(f.type, 'Valve')
error('Valve coefficient can only be set for valves.');
end
ok = calllib(ct, 'flowdev_setValveCoeff', f.id, k);
% if k < 0
% error(geterr);
% end
end
end

View File

@ -34,9 +34,6 @@ classdef ReactorNet < handle
end
r.id = calllib(ct, 'reactornet_new');
% if r.id < 0
% error(geterr);
% end
% add reactors
nr = length(reactors);

View File

@ -37,7 +37,7 @@ P0 = 4.47*101325;
Phi = 0.2899;
% Import the gas phase, read out key species indices:
gas_calc = Solution('gri30.yaml');
gas_calc = Solution('gri30.yaml', 'gri30');
ich4 = gas_calc.speciesIndex('CH4');
io2 = gas_calc.speciesIndex('O2');
in2 = gas_calc.speciesIndex('N2');

View File

@ -95,7 +95,7 @@ surf_phase.advanceCoverages(1.0);
flow = AxisymmetricFlow(gas, 'flow');
% set some parameters for the flow
flow.setPressure(p);
flow.P = p;
flow.setupGrid(initial_grid);
flow.setSteadyTolerances('default', tol_ss{:});
flow.setTransientTolerances('default', tol_ts{:});

View File

@ -9,7 +9,7 @@ function isentropic(g)
if nargin == 1
gas = g;
else
gas = Solution('gri30.yaml');
gas = Solution('gri30.yaml', 'gri30');
end
% set the stagnation state

View File

@ -118,7 +118,7 @@ function anCurr = anode_curr(phi_s, phi_l, X_Li_an, anode, elde, elyt, anode_int
% Get the net reaction rate at the anode-side interface
% Reaction according to cti file: Li+[elyt] + V[anode] + electron <=> Li[anode]
r = anode_interface.rop_net; % [kmol/m2/s]
r = anode_interface.ropNet; % [kmol/m2/s]
% Calculate the current. Should be negative for cell discharge.
anCurr = r*faradayconstant*S_an; %
@ -135,7 +135,7 @@ function caCurr = cathode_curr(phi_s, phi_l, X_Li_ca, cathode, elde, elyt, catho
% Get the net reaction rate at the cathode-side interface
% Reaction according to cti file: Li+[elyt] + V[cathode] + electron <=> Li[cathode]
r = cathode_interface.rop_net; % [kmol/m2/s]
r = cathode_interface.ropNet; % [kmol/m2/s]
% Calculate the current. Should be negative for cell discharge.
caCurr = r*faradayconstant*S_ca*(-1); %