From 55ccdc63b902c9bc51adf6d70eafa4de5913b78c Mon Sep 17 00:00:00 2001 From: ssun30 Date: Mon, 16 May 2022 13:53:20 -0400 Subject: [PATCH] [Matlab] Update 1D documentation Replace magic numbers for 1D domain types with string constants --- .../matlab_experimental/1D/AxiStagnFlow.m | 2 +- .../matlab_experimental/1D/AxisymmetricFlow.m | 2 +- interfaces/matlab_experimental/1D/Domain1D.m | 439 ++++++---- interfaces/matlab_experimental/1D/FreeFlame.m | 2 +- interfaces/matlab_experimental/1D/Inlet.m | 2 +- interfaces/matlab_experimental/1D/Outlet.m | 2 +- interfaces/matlab_experimental/1D/OutletRes.m | 2 +- interfaces/matlab_experimental/1D/Sim1D.m | 2 +- interfaces/matlab_experimental/1D/Surface.m | 4 +- interfaces/matlab_experimental/1D/SymmPlane.m | 2 +- .../matlab_experimental/Base/Interface.m | 77 +- .../matlab_experimental/Base/Kinetics.m | 3 + interfaces/matlab_experimental/Base/Mixture.m | 342 +++++--- .../matlab_experimental/Base/Solution.m | 55 +- .../matlab_experimental/Base/ThermoPhase.m | 814 +++++++++++++----- interfaces/matlab_experimental/Func/Func.m | 161 ++-- interfaces/matlab_experimental/LoadCantera.m | 3 +- .../matlab_experimental/Reactor/FlowDevice.m | 130 +-- .../matlab_experimental/Reactor/ReactorNet.m | 3 - samples/matlab_experimental/PFR.m | 2 +- samples/matlab_experimental/catcomb.m | 2 +- samples/matlab_experimental/isentropic.m | 2 +- .../matlab_experimental/lithium_ion_battery.m | 4 +- 23 files changed, 1343 insertions(+), 714 deletions(-) diff --git a/interfaces/matlab_experimental/1D/AxiStagnFlow.m b/interfaces/matlab_experimental/1D/AxiStagnFlow.m index 4a02ed74bf..8a342bb15e 100644 --- a/interfaces/matlab_experimental/1D/AxiStagnFlow.m +++ b/interfaces/matlab_experimental/1D/AxiStagnFlow.m @@ -1,4 +1,4 @@ function m = AxiStagnFlow(gas) % Get an axisymmetric stagnation flow domain. - m = Domain1D(1, gas); + m = Domain1D('StagnationFlow', gas); end diff --git a/interfaces/matlab_experimental/1D/AxisymmetricFlow.m b/interfaces/matlab_experimental/1D/AxisymmetricFlow.m index a1a25c7e11..aebce2adfa 100644 --- a/interfaces/matlab_experimental/1D/AxisymmetricFlow.m +++ b/interfaces/matlab_experimental/1D/AxisymmetricFlow.m @@ -2,7 +2,7 @@ % 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 diff --git a/interfaces/matlab_experimental/1D/Domain1D.m b/interfaces/matlab_experimental/1D/Domain1D.m index d83f44fc76..525e3da01b 100644 --- a/interfaces/matlab_experimental/1D/Domain1D.m +++ b/interfaces/matlab_experimental/1D/Domain1D.m @@ -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' - % - % :parameter b: - % Instance of class 'Solution' (for a == 'StagnationFlow') - % or 'Interface' (for a == 'ReactingSurface'). - % Not used for all other valid values of 'a'. - % - % :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'. - + % `StagnationFlow` + % `Inlet1D` + % `Surf1D` + % `Symm1D` + % `Outlet1D` + % `ReactingSurface` + % `Sim1D` + % `OutletRes` + % + % :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``. + % checklib; d.domainID = -1; @@ -96,37 +96,54 @@ %% 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 @@ function clear(d) %% 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 @@ function clear(d) 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 @@ function clear(d) 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 @@ function clear(d) 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 @@ function clear(d) 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 @@ function clear(d) 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 @@ function clear(d) 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 @@ function clear(d) 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 @@ function clear(d) 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 @@ function clear(d) 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 @@ function clear(d) 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 @@ function clear(d) 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 @@ function setBounds(d, component, lower, upper) 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 @@ function setCoverageEqs(d, onoff) 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 @@ function setFixedTempProfile(d, profile) 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 @@ function setProfileD(d, n, p) 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 @@ function setSteadyTolerances(d, component, rtol, atol) 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 @@ function setTransientTolerances(d, component, rtol, atol) 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 diff --git a/interfaces/matlab_experimental/1D/FreeFlame.m b/interfaces/matlab_experimental/1D/FreeFlame.m index a069866b3e..63355df161 100644 --- a/interfaces/matlab_experimental/1D/FreeFlame.m +++ b/interfaces/matlab_experimental/1D/FreeFlame.m @@ -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 diff --git a/interfaces/matlab_experimental/1D/Inlet.m b/interfaces/matlab_experimental/1D/Inlet.m index 8e59b2f337..a2037fe2fc 100644 --- a/interfaces/matlab_experimental/1D/Inlet.m +++ b/interfaces/matlab_experimental/1D/Inlet.m @@ -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 diff --git a/interfaces/matlab_experimental/1D/Outlet.m b/interfaces/matlab_experimental/1D/Outlet.m index 30afdc7839..a2e602a0a8 100644 --- a/interfaces/matlab_experimental/1D/Outlet.m +++ b/interfaces/matlab_experimental/1D/Outlet.m @@ -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 diff --git a/interfaces/matlab_experimental/1D/OutletRes.m b/interfaces/matlab_experimental/1D/OutletRes.m index 813b0bf394..b57b879fc6 100644 --- a/interfaces/matlab_experimental/1D/OutletRes.m +++ b/interfaces/matlab_experimental/1D/OutletRes.m @@ -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 diff --git a/interfaces/matlab_experimental/1D/Sim1D.m b/interfaces/matlab_experimental/1D/Sim1D.m index 89ed1ac6cf..96f0ee09ef 100644 --- a/interfaces/matlab_experimental/1D/Sim1D.m +++ b/interfaces/matlab_experimental/1D/Sim1D.m @@ -26,7 +26,7 @@ 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 diff --git a/interfaces/matlab_experimental/1D/Surface.m b/interfaces/matlab_experimental/1D/Surface.m index 03c7544c0f..2e88fdfb98 100644 --- a/interfaces/matlab_experimental/1D/Surface.m +++ b/interfaces/matlab_experimental/1D/Surface.m @@ -4,14 +4,14 @@ % 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 diff --git a/interfaces/matlab_experimental/1D/SymmPlane.m b/interfaces/matlab_experimental/1D/SymmPlane.m index 32388544bb..0ffa2e3eb6 100644 --- a/interfaces/matlab_experimental/1D/SymmPlane.m +++ b/interfaces/matlab_experimental/1D/SymmPlane.m @@ -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 diff --git a/interfaces/matlab_experimental/Base/Interface.m b/interfaces/matlab_experimental/Base/Interface.m index d3dd77a4f8..2233fd7e0e 100644 --- a/interfaces/matlab_experimental/Base/Interface.m +++ b/interfaces/matlab_experimental/Base/Interface.m @@ -9,15 +9,27 @@ %% 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 `__. + % + % 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 @@ %% 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 @@ 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 @@ 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 @@ 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. diff --git a/interfaces/matlab_experimental/Base/Kinetics.m b/interfaces/matlab_experimental/Base/Kinetics.m index bad326163c..993048eb05 100644 --- a/interfaces/matlab_experimental/Base/Kinetics.m +++ b/interfaces/matlab_experimental/Base/Kinetics.m @@ -6,8 +6,11 @@ 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 diff --git a/interfaces/matlab_experimental/Base/Mixture.m b/interfaces/matlab_experimental/Base/Mixture.m index 67312c37b3..f6f77d9642 100644 --- a/interfaces/matlab_experimental/Base/Mixture.m +++ b/interfaces/matlab_experimental/Base/Mixture.m @@ -2,42 +2,48 @@ 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}); - % - % Phases may also be added later using the addPhase method: - % >> water = Solution('water.cti'); - % >> mix.addPhase(water, 3.0); - % - % 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. - % - % 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. - % - % :parameter phases: - % Cell array of phases and mole numbers. + % 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:: + % + % >> gas = Solution('gas.yaml'); + % >> graphite = Solution('graphite.yaml'); + % >> mix = Mixture({gas, 1.0; graphite, 0.1}); + % + % Phases may also be added later using the addPhase method:: + % + % >> water = Solution('water.yaml'); + % >> mix.addPhase(water, 3.0); + % + % 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 @@ %% 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 @@ function display(m) 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 @@ function addPhase(m, phase, moles) %% 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 @@ function addPhase(m, phase, moles) 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 @@ function addPhase(m, phase, moles) 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 @@ function addPhase(m, phase, moles) 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 @@ function addPhase(m, phase, moles) %% 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 @@ function setSpeciesMoles(m, moles) 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); - % - % :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. + % 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. :: + % + % >> 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 diff --git a/interfaces/matlab_experimental/Base/Solution.m b/interfaces/matlab_experimental/Base/Solution.m index 5d4d2a87ba..88f3672d23 100644 --- a/interfaces/matlab_experimental/Base/Solution.m +++ b/interfaces/matlab_experimental/Base/Solution.m @@ -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 @@ % 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; diff --git a/interfaces/matlab_experimental/Base/ThermoPhase.m b/interfaces/matlab_experimental/Base/ThermoPhase.m index c53c8c2a66..640789994b 100644 --- a/interfaces/matlab_experimental/Base/ThermoPhase.m +++ b/interfaces/matlab_experimental/Base/ThermoPhase.m @@ -11,11 +11,11 @@ S % entropy U % internal energy G % Gibbs free energy - V % specific volume basis end properties (Dependent) + V % specific volume DP DPX DPY @@ -64,6 +64,15 @@ %% ThermoPhase class constructor function tp = ThermoPhase(src, id) + % THERMOPHASE ThermoPhase class constructor. + % t = ThermoPhase(src, id) + % :param src: + % Input string of YAML, CTI, or XML file name. + % :param id: + % ID of the phase to import as specified in the input file. (optional) + % :return: + % Instance of class :mat:func:`ThermoPhase` + % checklib; if nargin > 2 error('ThermoPhase expects 1 or 2 input arguments.'); @@ -87,14 +96,26 @@ function display(tp, threshold) end function tpClear(tp) - % Delete the kernel object. - + % CLEAR Delete the kernel object. + % tp.tpClear + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) + % calllib(ct, 'thermo_del', tp.tpID); end function tp = set.basis(tp, b) - % Set basis of thermodynamic properties - % Default is molar + % BASIS Determines whether intensive thermodynamic properties are + % treated on a mass (per kg) or molar (per kmol) basis. This + % affects the values returned by the properties H, U, S, G, V, + % Density, Cv, and Cp, as well as the values used with the + % state-setting properties such as HPX and UV. + % tp.Basis + % :param tp: + % Instance of class :mat:func:`ThermoPhase`. + % :param b: + % String. Can be 'mole'/'molar'/'Molar'/Mole' or 'mass'/'Mass'. if strcmp(b, 'mole') || strcmp(b, 'molar') ... || strcmp(b, 'Mole') || strcmp(b, 'Molar') @@ -105,14 +126,67 @@ function tpClear(tp) end end + function tp = equilibrate(tp, xy, solver, rtol, maxsteps, ... + maxiter, loglevel) + % EQUILIBRATE Set the phase to a state of chemical equilibrium. + % tp.equilibrate(xy, solver, rtol, maxsteps, maxiter, loglevel) + % :param XY: + % A two-letter string, which must be one of the set + % ``['TP','TV','HP','SP','SV','UV','UP']``, + % indicating which pair of properties should be held constant. + % Not all of the properties to be held constant are available with + % all of the solvers. + % :param solver: + % Specifies the equilibrium solver to use. If solver = 0, a fast + % solver using the element potential method will be used. If + % solver = 1, a slower but more robust Gibbs minimization solver + % will be used. If solver >= 2, a version of the VCS algorithm will + % be used. If solver < 0 or is unspecified, the fast solver + % will be tried first, then if it fails the Gibbs minimization solver + % will be tried. + % :param rtol: + % The relative error tolerance. + % :param maxsteps: + % Maximum number of steps in composition to take to find a + % converged solution. + % :param maxiter: + % For the Gibbs minimization solver only, this specifies the number + % of 'outer' iterations on T or P when some property pair other than + % TP is specified. + % :param loglevel: + % Set to a value > 0 to write diagnostic output. Larger values + % generate more detailed information. + % + if nargin < 3 + solver = -1; + end + if nargin < 4 + rtol = 1.0e-9; + end + if nargin < 5 + maxsteps = 1000; + end + if nargin < 6 + maxiter = 100; + end + if nargin < 7 + loglevel = 0; + end + calllib(ct, 'thermo_equilibrate', tp.tpID, xy, solver, ... + rtol, maxsteps, maxiter, loglevel); + end + %% PhaseGet single methods - function amu = atomicWeights(tp) - % Get the atomic masses of the elements. - % + function amu = atomicMasses(tp) + % ATOMICMASSES Get the atomic masses of the elements. + % x = tp.atomicMasses + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase). % :return: - % Vector of element atomic masses. Unit: kg/kmol - + % Vector of element atomic masses. Units: kg/kmol + % nel = tp.nElements; aa = zeros(1, nel); pt = libpointer('doublePtr', aa); @@ -122,11 +196,15 @@ function tpClear(tp) end function e = charges(tp) - % Get the array of species charges. + % CHARGES Get the array of species charges + % x = tp.charges % + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) % :return: - % Vector of species charges. Unit: elem. charge - + % Vector of species charges. Units: elem. charge + % nsp = tp.nSpecies; yy = zeros(1, nsp); pt = libpointer('doublePtr', yy); @@ -136,19 +214,32 @@ function tpClear(tp) end function k = elementIndex(tp, name) - % Get the index of an element given the name. - % The index is an integer assigned to each element in sequence - % as it is read in from the input file. + % ELEMENTINDEX Get the index of an element given its name. + % k = tp.elementIndex(name) + % The index is an integer assigned to each element in sequence as it + % is read in from the input file. + % + % If ``name`` is a single string, the return value will be a integer + % containing the corresponding index. If it is an cell array of + % strings, the output will be an array of the same shape + % containing the indices. % - % Note: In keeping with the conventions used by Matlab, the - % indices start from 1 instead of 0 as in Cantera C++ and - % Python interfaces. + % NOTE: In keeping with the conventions used by Matlab, this method + % returns 1 for the first element. In contrast, the corresponding + % method elementIndex in the Cantera C++ and Python interfaces + % returns 0 for the first element, 1 for the second one, etc. :: % - % :parameter name: - % String or cell array of elements whose index is requested + % >> ic = gas.elementIndex('C'); + % >> ih = gas.elementIndex('H'); + % + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) + % :param name: + % String or cell array of strings of elements to look up % :return: - % Integer number of elements in the phase. - + % Integer or vector of integers of element indices + % if iscell(name) [m, n] = size(name); k = zeros(m, n); @@ -177,8 +268,17 @@ function tpClear(tp) end function elMassFrac = elementalMassFraction(tp, element) - % Determine the elemental mass fraction in gas object. - % Check input :parameters. + % ELEMENTALMASSFRACTION Determine the elemental mass fraction in gas object. + % elMassFrac = tp.elementalMassFraction(element) + % :param tp: + % Object representing the gas, instance of class :mat:func:`Solution`, + % and an ideal gas. The state of this object should be set to an + % estimate of the gas state before calling elementalMassFraction. + % :param element: + % String representing the element name. + % :return: + % Elemental mass fraction within a gas object. + % if nargin ~= 2 error('elementalMassFraction expects two input arguments.'); end @@ -216,11 +316,17 @@ function tpClear(tp) end function mmw = meanMolecularWeight(tp) - % Get the mean molecular weight. + % MEANMOLECULARWEIGHT Get the mean molecular weight. + % mmw = tp.meanMolecularWeight + % The mean molecular weight is the mole-fraction-weighted sum of the + % molar masses of the individual species in the phase. % + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) % :return: - % Double mean molecular weight. Unit: kg/kmol - + % Scalar double mean molecular weight. Units: kg/kmol + % mmw = calllib(ct, 'thermo_meanMolecularWeight', tp.tpID); end @@ -231,11 +337,15 @@ function tpClear(tp) end function mw = MolecularWeights(tp) - % Get the array of molecular weights of all species. + % MOLECULARWEIGHTS Get the molecular weights of the species. + % mw = tp.molecularWeights % + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) % :return: - % Vector of species molecular weights. Unit: kg/kmol - + % Vector of species molecular weights. Units: kg/kmol + % nsp = tp.nSpecies; yy = zeros(1, nsp); pt = libpointer('doublePtr', yy); @@ -245,15 +355,18 @@ function tpClear(tp) end function n = nAtoms(tp, species, element) - % Get the number of atoms of an element in a species. - % - % :parameter k: - % String species name or integer species number. - % :parameter m: - % String element name or integer element number. + % NATOMS Get the number of atoms of an element in a species. + % n = tp.nAtoms(k,m) + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) + % :param k: + % String species name or integer species number + % :param m: + % String element name or integer element number % :return: - % Integer number of atoms of the element in the species. - + % Number of atoms of element ``m`` in species ``k``. + % if nargin == 3 if ischar(species) k = tp.speciesIndex(species); @@ -278,37 +391,55 @@ function tpClear(tp) end function nel = nElements(tp) - % Get the number of elements. - % + % NELEMENTS Get the number of elements. + % n = tp.nElements + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) % :return: - % Integer number of elements in the phase. - + % Number of elements in the phase. + % nel = calllib(ct, 'thermo_nElements', tp.tpID); end function nsp = nSpecies(tp) - % Get the number of species. - % + % NSPECIES Get the number of species. + % n = tp.nSpecies + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) % :return: - % Integer number of species in the phase. - + % Number of species in the phase. + % nsp = calllib(ct, 'thermo_nSpecies', tp.tpID); end function k = speciesIndex(tp, name) - % Get the index of a species given the name. - % The index is an integer assigned to each species in sequence - % as it is read in from the input file. + % SPECIESINDEX Get the index of a species given the name. + % k = tp.speciesIndex(name) + % The index is an integer assigned to each species in sequence as it + % is read in from the input file. + % + % NOTE: In keeping with the conventions used by Matlab, this method + % returns 1 for the first species, 2 for the second, etc. In + % contrast, the corresponding method speciesIndex in the Cantera C++ + % and Python interfaces returns 0 for the first species, 1 for the + % second one, etc. :: % - % Note: In keeping with the conventions used by Matlab, the - % indices start from 1 instead of 0 as in Cantera C++ and - % Python interfaces. + % >> ich4 = gas.speciesIndex('CH4'); + % >> iho2 = gas.speciesIndex('HO2'); % - % :parameter name: - % String or cell array of species whose index is requested. + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % class derived from ThermoPhase) + % :param name: + % If name is a single string, the return value will be a integer + % containing the corresponding index. If it is an cell array of + % strings, the output will be an array of the same shape + % containing the indices. % :return: - % Integer number of species in the phase. - + % Scalar or array of integers + % if iscell(name) [m, n] = size(name); k = zeros(m, n); @@ -337,13 +468,16 @@ function tpClear(tp) end function nm = speciesName(tp, k) - % Get the name of a species given the index. - % - % :parameter k: - % Scalar or array of integer species index. + % SPECIESNAME Get the name of a species given the index. + % nm = tp.speciesName(k) + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % class derived from ThermoPhase) + % :param k: + % Scalar or array of integer species numbers % :return: - % Cell array of strings species name. - + % Cell array of strings + % [m, n] = size(k); nm = cell(m, n); for i = 1:m @@ -362,51 +496,75 @@ function tpClear(tp) end function n = speciesNames(tp) - % Get all species names. - + % SPECIESNAMES Get the species names. + % n = tp.speciesNames + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % class derived from ThermoPhase) + % :return: + % Cell array of strings of all of the species names + % n = tp.speciesName(1:tp.nSpecies); end function temperature = get.T(tp) - % Get the temperature. - % + % GET.T Get the temperature. + % temperature = tp.T + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % class derived from ThermoPhase) % :return: - % Double temperature. Unit: K - + % Temperature. Units: K + % temperature = calllib(ct, 'thermo_temperature', tp.tpID); end function pressure = get.P(tp) - % Get the pressure. - % + % GET.P Get the pressure. + % pressure = tp.P + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) % :return: - % Double pressure. Unit: Pa - + % Pressure. Units: Pa + % pressure = calllib(ct, 'thermo_pressure', tp.tpID); end function density = get.D(tp) - % Get the mass basis density in kg/m^3. - + % GET.D Get the density. + % density = tp.D + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) + % :return: + % Mass density. Units: kg/m**3 + % density = calllib(ct, 'thermo_density', tp.tpID); end function volume = get.V(tp) - % Get the specific volume depending on the basis. - % + % GET.V Get the specific volume. + % volume = tp.V + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) % :return: - % Density depending on the basis. Units: + % Specific volume depending on the basis. Units: % m^3/kmol (molar) m^3/kg (mass). - volume = 1/tp.D; end function moleFractions = get.X(tp) - % Get the mole fractions of all species. - % + % GET.X Get the mole fractions of all species. + % moleFractions = tp.X + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) % :return: - % Vector of species mole fractions. - + % Vector of species mole fractions for input phase. If + % no output argument is specified, a bar plot is produced. + % nsp = tp.nSpecies; xx = zeros(1, nsp); pt = libpointer('doublePtr', xx); @@ -426,14 +584,17 @@ function tpClear(tp) end function x = moleFraction(tp, species) - % Get the mole fraction of one or a list of species. - % - % :parameter species: - % String or cell array of species whose mole fraction is - % requested. + % MOLEFRACTION Get the mole fraction of one or a list of species. + % x = tp.moleFraction(species) + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) + % :param species: + % String or cell array of strings of species whose mole + % fraction is desired % :return: - % Scalar or vector of species mole fractions. - + % Scalar or vector double mole fractions + % xarray = tp.X; if isa(species, 'char') k = tp.speciesIndex(species); @@ -455,11 +616,15 @@ function tpClear(tp) end function massFractions = get.Y(tp) - % Get the mass fractions of all species. - % + % GET.Y Get the mass fractions of all species. + % massFractions = tp.Y + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) % :return: - % Vector of species mass fractions. - + % Vector of species mass fractions for input phase. If + % no output argument is specified, a bar plot is produced. + % nsp = tp.nSpecies; yy = zeros(1, nsp); pt = libpointer('doublePtr', yy); @@ -479,14 +644,17 @@ function tpClear(tp) end function y = massFraction(tp, species) - % Get the mass fraction of one or a list of species. - % - % :parameter species: - % String or cell array of species whose mass fraction is - % requested. + % MASSFRACTION Get the mass fraction of one or a list of species. + % y = tp.massFraction(species) + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) + % :param species: + % String or cell array of strings of species whose mass + % fraction is desired % :return: - % Scalar or vector of species mass fractions. - + % Scalar or vector double mass fractions + % yy = tp.Y; if isa(species, 'char') k = tp.speciesIndex(species); @@ -509,12 +677,24 @@ function tpClear(tp) %% ThermoGet single methods - function mu = chemical_potentials(tp) - % Get the chemical potentials of the species. + function mu = chemicalPotentials(tp) + % CHEMPOTENTIALS Get the chemical potentials of the species. + % mu = tp.chemPotentials + % The expressions used to compute the chemical potential + % depend on the model implemented by the underlying kernel + % thermo manager. % + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase). % :return: - % Vector of species chemical potentials. Unit: J/kmol. - + % Vector of species chemical potentials. Units: J/kmol + % + % This method returns an array containing the species + % chemical potentials [J/kmol]. The expressions used to + % compute these depend on the model implemented by the + % underlying kernel thermo manager. + % nsp = tp.nSpecies; xx = zeros(1, nsp); pt = libpointer('doublePtr', xx); @@ -524,12 +704,15 @@ function tpClear(tp) end function c = cv(tp) - % Get the specific heat at constant volume. - % + % CV Get the basis-dependent specific heat at constant volume. + % c = tp.cv + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) % :return: - % Specific heat of the mixture at constant volume depending - % on the basis. Units: J/kmol-K (molar) J/kg-K (mass). - + % specific heat of the mixture at constant volume. + % Units: J/kg-K (mass basis) or J/kmol-K (molar basis). + % if strcmp(tp.basis, 'molar') c = calllib(ct, 'thermo_cv_mole', tp.tpID); elseif strcmp(tp.basis, 'mass') @@ -539,12 +722,15 @@ function tpClear(tp) end function c = cp(tp) - % Get the specific heat at constant pressure. - % + % CP Get the basis-dependent specific heat at constant pressure. + % v = tp.cp + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) % :return: - % Specific heat of the mixture at constant pressure depending - % on the basis. Units: J/kmol-K (molar) J/kg-K (mass). - + % Specific heat of the mixture at constant pressure. + % Units: J/kg-K (mass-basis) or J/mol-K (molar-basis). + % if strcmp(tp.basis, 'molar') c = calllib(ct, 'thermo_cp_mole', tp.tpID); elseif strcmp(tp.basis, 'mass') @@ -554,44 +740,62 @@ function tpClear(tp) end function d = critDensity(tp) - % Get the critical density. - % + % CRITDENSITY Get the critical density. + % v = tp.critDensity + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) % :return: - % Critical density. Unit: K. - + % Critical density. Units: kg/m**3 + % d = calllib(ct, 'thermo_critDensity', tp.tpID); end function p = critPressure(tp) - % Get the critical pressure. - % + % CRITPRESSURE Get the critical pressure. + % v = tp.critPressure + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) % :return: - % Critical temperature. Unit: Pa. - + % Critical pressure. Units: Pa + % p = calllib(ct, 'thermo_critPressure', tp.tpID); end function t = critTemperature(tp) - % Get the critical temperature. - % + % CRITTEMPERATURE Get the critical temperature. + % v = tp.critTemperature + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) % :return: - % Critical temperature. Unit: K. - + % Critical temperature. Units: K + % t = calllib(ct, 'thermo_critTemperature', tp.tpID); end function v = electricPotential(tp) - % Get the electric potential - % + % ELECTRICPOTENTIAL Get the electric potential. + % v = tp.electricPotential + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) % :return: - % Electric potential of the phase. Unit: V. - + % The electric potential of the phase. Units: V + % v = calllib(ct, 'thermo_electricPotential', tp.tpID); end function e = eosType(tp) - % Get the type of the equation of state - + % EOSTYPE Get the type of the equation of state. + % e = tp.eosType + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) + % :return: + % An string identifying the equation of state. + % buflen = calllib(ct, 'thermo_getEosType', tp.tpID, 0, ''); if buflen > 0 aa = char(zeros(1, buflen)); @@ -602,8 +806,15 @@ function tpClear(tp) end function v = isIdealGas(tp) - % Get a flag indicating whether the phase is an ideal gas. - + % ISIDEALGAS Get a flag indicating whether the phase is an ideal gas. + % v = tp.isIdealGas + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) + % :return: + % True (1) if the phase is an ideal gas or ideal gas + % mixture, and false (0) otherwise. + % if strcmp(tp.eosType, 'IdealGas') v = 1; else @@ -613,68 +824,123 @@ function tpClear(tp) end function b = isothermalCompressibility(tp) - % Get the isothermal compressibility - % + % ISOTHERMALCOMPRESSIBILITY Get the isothermal compressibility. + % b = tp.isothermalCompressibility + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) % :return: - % Isothermal compressibility. Unit: 1/Pa. - + % Isothermal Compressibility. Units: 1/Pa + % b = calllib(ct, 'thermo_isothermalCompressibility', tp.tpID); end function t = maxTemp(tp) - % Get the maximum temperature of the :parameter fits. + % MAXTEMP Get the maximum temperature of the parameter fits. + % v = tp.maxTemp + % The parameterizations used to represent the temperature-dependent + % species thermodynamic properties are generally only valid in some + % finite temperature range, which may be different for each species + % in the phase. % + % See also: :mat:func:`minTemp` + % + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) % :return: - % Vector of maximum temperatures of all species. - + % Vector of maximum temperatures of all species + % t = calllib(ct, 'thermo_maxTemp', tp.tpID, -1); end function t = minTemp(tp) - % Get the minimum temperature of the :parameter fits. + % MINTEMP Get the minimum temperature of the parameter fits. + % v = tp.minTemp + % The parameterizations used to represent the temperature-dependent + % species thermodynamic properties are generally only valid in some + % finite temperature range, which may be different for each species + % in the phase. % + % See also: :mat:func:`maxTemp` + % + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) % :return: - % Vector of minimum temperatures of all species. - + % Vector of minimum temperatures of all species + % t = calllib(ct, 'thermo_minTemp', tp.tpID, -1); end - function p = P_sat(tp, t) - % Get the saturation pressure for a given temperature. - % - % :parameter t: - % Temperature. Unit: K. + function p = refPressure(tp) + % REFPRESSURE Get the reference pressure. + % v = tp.refPressure + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) % :return: - % Saturation pressure for temperature t. Unit: Pa. + % Reference pressure in Pa. All standard-state + % thermodynamic properties are for this pressure. + % + p = calllib(ct, 'thermo_refPressure', tp.tpID, -1); + end + function p = satPressure(tp, t) + % SATPRESSURE Get the saturation pressure for a given temperature. + % v = tp.satPressure(T) + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) + % :param T: + % Temperature Units: K + % :return: + % Saturation pressure for temperature T. Units: Pa + % p = calllib(ct, 'thermo_satPressure', tp.tpID, t); end - function p = refPressure(tp) - % Get the reference pressure. - % + function t = satTemperature(tp, p) + % SATTEMPERATURE Get the saturation temperature for a given pressure. + % v = tp.satTemperature(p) + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) + % :param p: + % Pressure. Units: Pa % :return: - % Reference pressure. Unit: Pa. - - p = calllib(ct, 'thermo_refPressure', tp.tpID, -1); + % Saturation temperature for pressure p. Units: K + % + t = calllib(ct, 'thermo_satTemperature', tp.tpID, p); end function c = soundspeed(tp) - % Get the speed of sound - % If the phase is an ideal gas, the speed of sound is - % calculated by: + % SOUNDSPEED Get the speed of sound. + % c = tp.soundspeed + % If the phase is an ideal gas, the speed of sound + % is calculated by: + % + % .. math:: c = \sqrt{\gamma * R * T} + % + % where :math:`\gamma` is the ratio of specific heats, :math:`R` is + % the specific gas constant, and :math:`T` is the temperature. If the + % phase is not an ideal gas, the speed of sound is calculated by % - % c = sqrt(gamma * R * T); + % .. math:: c = \sqrt{\left(\frac{\partial p}{\partial \rho}\right)_s} % - % where gamma is the ratio of specific heat, R is the specific - % gas constant, and T is the temperature. If the phase is not - % an ideal gas, the speed of sound is calculated by: + % where :math:`p` is the pressure and :math:`\rho` is the density, + % and the subscript :math:`s` indicates constant entropy. This is + % approximated by slightly increasing the density at constant entropy + % and computing the change in pressure. % - % c = sqrt( + % .. math:: c = \sqrt{\frac{p_1 - p_0}{\rho_1-\rho_0}} % + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % class derived from ThermoPhase) % :return: - % The speed of sound. Unit: m/s - + % The speed of sound. Units: m/s + % if tp.isIdealGas tp.basis = 'mass'; gamma = tp.cp/tp.cv; @@ -694,37 +960,35 @@ function tpClear(tp) end function a = thermalExpansionCoeff(tp) - % Get the thermal expansion coefficient. - % + % THERMALEXPANSIONCOEFF Get the thermal expansion coefficient. + % a = tp.thermalExpansionCoeff + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % class derived from ThermoPhase) % :return: - % Thermal expansion coefficient. Unit: 1/K. - - a = calllib(ct, 'thermo_thermalExpansionCoeff', tp.tpID); - end - - function t = T_sat(tp, p) - % Get the saturation temperature for a given pressure. + % Thermal Expansion Coefficient. Units: 1/K % - % :parameter p: - % Pressure. Unit: Pa. - % :return: - % Saturation temperature for pressure p. Unit: K. - - t = calllib(ct, 'thermo_satTemperature', tp.tpID, p); + a = calllib(ct, 'thermo_thermalExpansionCoeff', tp.tpID); end function v = vaporFraction(tp) - % Get the vapor fractions. - % + % VAPORFRACTION Get the vapor fraction. + % v = tp.vaporFraction + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % class derived from ThermoPhase) % :return: - % Vapor fraction. - + % Vapor fraction. + % v = calllib(ct, 'thermo_vaporFraction', tp.tpID); end function enthalpy = get.H(tp) - % Get the enthalpy. - % + % GET.H Get the mass specific enthalpy. + % enthalpy = tp.H + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) % :return: % Enthalpy of the mixture depending on the basis. % Units: J/kmol (molar) J/kg (mass). @@ -738,11 +1002,17 @@ function tpClear(tp) end function enthalpy = enthalpies_RT(tp) - % Get the non-dimensional enthalpy. - % + % ENTHALPIES_RT Get the non-dimensional enthalpies. + % v = tp.enthalpies_RT + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) % :return: - % Vector of standard-state species enthalpies divided by RT. - + % Vector of standard-state species enthalpies + % divided by RT, where R is the universal gas + % constant and T is the temperature. For gaseous species, these + % values are ideal gas enthalpies. + % nsp = tp.nSpecies; xx = zeros(1, nsp); pt = libpointer('doublePtr', xx); @@ -751,12 +1021,15 @@ function tpClear(tp) end function entropy = get.S(tp) - % Get the entropy. - % + % GET.S Get the mass specific entropy. + % entropy = tp.S + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) % :return: - % Entropy of the mixture depending on the basis. - % Units: J/kmol-K (molar) J/kg-K (mass). - + % Entropy of the mixture depending on the basis. + % Units: J/kmol-K (molar) J/kg-K (mass) + % if strcmp(tp.basis, 'molar') entropy = calllib(ct, 'thermo_entropy_mole', tp.tpID); elseif strcmp(tp.basis, 'mass') @@ -766,8 +1039,11 @@ function tpClear(tp) end function intEnergy = get.U(tp) - % Get the internal energy. - % + % GET.U Get the mass specific internal energy. + % intEnergy = tp.U + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) % :return: % Internal energy of the mixture depending on the basis. % Units: J/kmol (molar) J/kg (mass). @@ -781,8 +1057,11 @@ function tpClear(tp) end function gibbs = get.G(tp) - % Get the Gibss free energy. - % + % GET.G Get the mass specific Gibbs function. + % gibbs = tp.G + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % object that derives from ThermoPhase) % :return: % Gibbs free energy of the mixture depending on the basis. % Units: J/kmol (molar) J/kg (mass). @@ -968,24 +1247,62 @@ function tpClear(tp) %% PhaseSet single methods function tp = setElectricPotential(tp, phi) - % Set the electric potential in V. - + % SETELECTRICPOTENTIAL Set the electric potential. + % tp.setElectricPotential(phi) + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % class derived from ThermoPhase) + % :param phi: + % Electric potential. Units: V + % calllib(ct, 'thermo_setElectricPotential', tp.tpID, phi); end function tp = setState_Psat(tp, p, q) - % Set saturated vapor - + % SETSTATEPSAT Set the fluid state using the given pressure and quality. + % tp.setState_Psat(p, q) + % The fluid state will be set to a saturated liquid-vapor state using the + % input pressure and vapor fraction (quality) as the independent, + % intensive variables. + % + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % class derived from ThermoPhase) + % :param p: + % Pressure (Pa) + % :param q: + % Vapor fraction + % calllib(ct, 'thermo_setState_Psat', tp.tpID, p, q); end function tp = setState_Tsat(tp, t, q) - % Set saturated liquid - + % SETSTATETSAT Set the fluid state using the given temperature and quality. + % tp.setState_Tsat(t, q) + % The fluid state will be set to a saturated liquid-vapor state using the + % input temperature and vapor fraction (quality) as the independent, + % intensive variables. + % + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % class derived from ThermoPhase) + % :param t: + % Temperature (K) + % :param q: + % Vapor fraction + % calllib(ct, 'thermo_setState_Tsat', tp.tpID, t, 1 - q); end function set.T(tp, temperature) + % SET.T Set the temperature. + % tp.T = temperature + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % class derived from ThermoPhase) + % :param temperature: + % Temperature. Units: K + % if temperature <= 0 error('The temperature must be positive'); end @@ -993,6 +1310,17 @@ function tpClear(tp) end function set.P(tp, pressure) + % SET.P Set the pressure. + % tp.P = pressure + % The pressure is set by changing the density holding the + % temperature and chemical composition fixed. + % + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % class derived from ThermoPhase) + % :param pressure: + % Pressure. Units: Pa + % if pressure <= 0 error('The pressure must be positive'); end @@ -1000,6 +1328,14 @@ function tpClear(tp) end function set.D(tp, density) + % SET.D Set the density. + % tp.D = density + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % class derived from ThermoPhase) + % :param density: + % Density. Units: kg/m**3 + % if density <= 0 error('The density must be positive'); end @@ -1007,10 +1343,29 @@ function tpClear(tp) end function set.X(tp, xx) - lim = 1e-9; + % SET.X Set the species mole fractions. + % tp.X = xx + % Note that calling :mat:func:`setMoleFractions` leaves the temperature and + % density unchanged, and therefore the pressure changes if the new + % composition has a molar mass that is different than the old + % composition. If it is desired to change the composition and hold + % the pressure fixed, use method :mat:func:`set` and specify the mole + % fractions and the pressure, or call :mat:func:`setPressure` + % after calling :mat:func:`setmoleFractions`. + % + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % class derived from ThermoPhase) + % :param xx: + % Vector of mole fractions whose length must be the same as + % the number of species. Alternatively, a string in the format + % ``'SPEC:Y,SPEC2:Y2'`` that specifies the mole fraction of + % specific species. + % + tol = 1e-9; if isa(xx, 'double') nsp = tp.nSpecies; - if sum(xx) - 1 <= 1e-9 + if sum(xx) - 1 <= tol norm = 0; else norm = 1; end @@ -1022,6 +1377,25 @@ function tpClear(tp) end function set.Y(tp, yy) + % SET.Y Set the species mass fractions. + % tp.Y = yy + % Note that calling :mat:func:`setMassFractions` leaves the temperature and + % density unchanged, and therefore the pressure changes if the new + % composition has a molar mass that is different than the old + % composition. If it is desired to change the composition and hold + % the pressure fixed, use method :mat:func:`set` and specify the mass + % fractions and the pressure, or call :mat:func:`setPressure` + % after calling :mat:func:`setMassFractions`. + % + % :param tp: + % Instance of class :mat:func:`ThermoPhase` (or another + % class derived from ThermoPhase) + % :param yy: + % Vector of mass fractions whose length must be the same as + % the number of species. Alternatively, a string in the format + % ``'SPEC:Y,SPEC2:Y2'`` that specifies the mass fraction of + % specific species. + % if isa(yy, 'double') nsp = tp.nSpecies; if sum(yy) -1 <= 1e-9 @@ -1317,27 +1691,5 @@ function tpClear(tp) tp.VH = input(1:2); end - function tp = equilibrate(tp, xy, solver, rtol, maxsteps, ... - maxiter, loglevel) - % use the ChemEquil solver by default - if nargin < 3 - solver = -1; - end - if nargin < 4 - rtol = 1.0e-9; - end - if nargin < 5 - maxsteps = 1000; - end - if nargin < 6 - maxiter = 100; - end - if nargin < 7 - loglevel = 0; - end - calllib(ct, 'thermo_equilibrate', tp.tpID, xy, solver, ... - rtol, maxsteps, maxiter, loglevel); - end - end end diff --git a/interfaces/matlab_experimental/Func/Func.m b/interfaces/matlab_experimental/Func/Func.m index e3d15446dd..fc23554704 100644 --- a/interfaces/matlab_experimental/Func/Func.m +++ b/interfaces/matlab_experimental/Func/Func.m @@ -12,54 +12,58 @@ %% 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. + % + % See also: :mat:func:`polynom`, :mat:func:`gaussian`, :mat:func:`plus`, + % :mat:func:`rdivide`, :mat:func:`times` % - % Note: this MATLAB class shadows the underlying C++ Cantera - % class 'Func1'. See the Cantera C++ documentation for more - % details. + % :param typ: + % String indicating type of functor to create. Possible values are: % - % :parameter typ: - % String indicating type of functor to create. Possible - % values are: - % * 'polynomial' - % * 'fourier' - % * 'gaussian' - % * 'arrhenius' - % * 'sum' - % * 'diff' - % * 'ratio' - % * 'composite' - % * 'periodic' + % * ``'polynomial'`` + % * ``'fourier'`` + % * ``'gaussian'`` + % * ``'arrhenius'`` + % * ``'sum'`` + % * ``'diff'`` + % * ``'ratio'`` + % * ``'composite'`` + % * ``'periodic'`` % - % :parameter n: - % Number of :parameters required for the functor - % :parameter p: - % Vector of :parameters + % :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 @@ %% 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 @@ function display(a) 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 @@ function display(a) 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 @@ function display(a) end end else - s = ['*** char not yet implemented for' p.typ ' ***']; + s = ['*** char not yet implemented for' f.typ ' ***']; end end end diff --git a/interfaces/matlab_experimental/LoadCantera.m b/interfaces/matlab_experimental/LoadCantera.m index b89cf40275..e764568715 100644 --- a/interfaces/matlab_experimental/LoadCantera.m +++ b/interfaces/matlab_experimental/LoadCantera.m @@ -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 @@ 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'], ... diff --git a/interfaces/matlab_experimental/Reactor/FlowDevice.m b/interfaces/matlab_experimental/Reactor/FlowDevice.m index 9fb5087e7b..c6f3a1c82d 100644 --- a/interfaces/matlab_experimental/Reactor/FlowDevice.m +++ b/interfaces/matlab_experimental/Reactor/FlowDevice.m @@ -11,12 +11,27 @@ %% 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 @@ 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 @@ %% 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 @@ function install(f, upstream, downstream) 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 @@ function setMaster(f, d) 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 diff --git a/interfaces/matlab_experimental/Reactor/ReactorNet.m b/interfaces/matlab_experimental/Reactor/ReactorNet.m index 527022cf46..67c6185f53 100644 --- a/interfaces/matlab_experimental/Reactor/ReactorNet.m +++ b/interfaces/matlab_experimental/Reactor/ReactorNet.m @@ -34,9 +34,6 @@ end r.id = calllib(ct, 'reactornet_new'); -% if r.id < 0 -% error(geterr); -% end % add reactors nr = length(reactors); diff --git a/samples/matlab_experimental/PFR.m b/samples/matlab_experimental/PFR.m index df24cfb73d..c3165188e2 100644 --- a/samples/matlab_experimental/PFR.m +++ b/samples/matlab_experimental/PFR.m @@ -37,7 +37,7 @@ 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'); diff --git a/samples/matlab_experimental/catcomb.m b/samples/matlab_experimental/catcomb.m index 2b786ed583..6b2ec88fdb 100644 --- a/samples/matlab_experimental/catcomb.m +++ b/samples/matlab_experimental/catcomb.m @@ -95,7 +95,7 @@ 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{:}); diff --git a/samples/matlab_experimental/isentropic.m b/samples/matlab_experimental/isentropic.m index ffc8791d4e..c96e5f2610 100644 --- a/samples/matlab_experimental/isentropic.m +++ b/samples/matlab_experimental/isentropic.m @@ -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 diff --git a/samples/matlab_experimental/lithium_ion_battery.m b/samples/matlab_experimental/lithium_ion_battery.m index 8140afcde1..515d99fe84 100644 --- a/samples/matlab_experimental/lithium_ion_battery.m +++ b/samples/matlab_experimental/lithium_ion_battery.m @@ -118,7 +118,7 @@ % 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 @@ % 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); %