diff --git a/Cantera/clib/src/ctfunc.cpp b/Cantera/clib/src/ctfunc.cpp index 2bdb85926..7f1680540 100755 --- a/Cantera/clib/src/ctfunc.cpp +++ b/Cantera/clib/src/ctfunc.cpp @@ -41,6 +41,12 @@ extern "C" { r = new Fourier1(n, params[n+1], params[0], params + 1, params + n + 2); } + else if (type == GaussianFuncType) { + if (lenp < 3) + throw CanteraError("func_new", + "not enough Gaussian coefficients"); + r = new Gaussian(params[0], params[1], params[2]); + } else if (type == PolyFuncType) { if (lenp < n + 1) throw CanteraError("func_new", diff --git a/Cantera/python/Cantera/Func.py b/Cantera/python/Cantera/Func.py index caf5c82d2..efceaf973 100644 --- a/Cantera/python/Cantera/Func.py +++ b/Cantera/python/Cantera/Func.py @@ -83,7 +83,13 @@ class Polynomial(Func1): Func1.__init__(self, 2, len(coeffs)-1, coeffs) - +class Gaussian(Func1): + """A Gaussian pulse. + """ + def __init__(self, A = 0.0, t0 = 0.0, FWHM = 0.0): + coeffs = array([A, t0, 0.5*FWHM], 'd') + Func1.__init__(self, 4, 0, coeffs) + class Fourier(Func1): """ Fourier series. diff --git a/Cantera/python/Cantera/Reactor.py b/Cantera/python/Cantera/Reactor.py index 3a274cf16..66a554758 100644 --- a/Cantera/python/Cantera/Reactor.py +++ b/Cantera/python/Cantera/Reactor.py @@ -26,8 +26,8 @@ class ReactorBase: self._name = name self._verbose = verbose self.insert(contents) - self.setInitialVolume(volume) - self.setEnergy(energy) + self._setInitialVolume(volume) + self._setEnergy(energy) def __del__(self): @@ -49,7 +49,7 @@ class ReactorBase: return s def name(self): - """Reactor name.""" + """The name of the reactor specified when it was constructed.""" return self._name def reactor_id(self): @@ -68,22 +68,23 @@ class ReactorBase: _cantera.reactor_setKineticsMgr(self.__reactor_id, contents.ckin) - def setInitialTime(self, t0): - """Set the initial time. Restarts integration from this time + def setInitialTime(self, T0): + """Deprecated. + Set the initial time. Restarts integration from this time using the current state as the initial condition. Default: 0.0 s""" - _cantera.reactor_setInitialTime(self.__reactor_id, t0) + _cantera.reactor_setInitialTime(self.__reactor_id, T0) - def setInitialVolume(self, t0): - """Set the initial reactor volume. Default: 1.0 m^3.""" - _cantera.reactor_setInitialVolume(self.__reactor_id, t0) + def _setInitialVolume(self, V0): + """Set the initial reactor volume. """ + _cantera.reactor_setInitialVolume(self.__reactor_id, V0) - def setEnergy(self, e): + def _setEnergy(self, eflag): """Turn the energy equation on or off. If the argument is the string 'off' or the number 0, the energy equation is disabled, and the reactor temperature is held constant at its initial value.""" ie = 1 - if e == 'off' or e == 0: + if eflag == 'off' or eflag == 0: ie = 0 if self._verbose: if ie: @@ -93,58 +94,66 @@ class ReactorBase: _cantera.reactor_setEnergy(self.__reactor_id, ie) def temperature(self): - """Temperature [K].""" + """The temperature in the reactor [K].""" return _cantera.reactor_temperature(self.__reactor_id) def density(self): - """Density [kg/m^3].""" + """The density of the fluid in the reactor [kg/m^3].""" return _cantera.reactor_density(self.__reactor_id) def volume(self): - """Volume [m^3].""" + """The total reactor volume [m^3]. The volume may change with time + if non-rigid walls are installed on the reactor.""" return _cantera.reactor_volume(self.__reactor_id) def time(self): - """Time [s]. The reactor time is set by method advance.""" + """Deprecated. The current time [s].""" return _cantera.reactor_time(self.__reactor_id) def mass(self): - """The total mass of the reactor contents [kg].""" + """The total mass of fluid in the reactor [kg].""" return _cantera.reactor_mass(self.__reactor_id) def enthalpy_mass(self): - """The specific enthalpy [J/kg].""" + """The specific enthalpy of the fluid in the reactor [J/kg].""" return _cantera.reactor_enthalpy_mass(self.__reactor_id) def intEnergy_mass(self): - """The specific interhal energy [J/kg].""" + """The specific internal energy of the fluid in the reactor [J/kg].""" return _cantera.reactor_intEnergy_mass(self.__reactor_id) def pressure(self): - """The pressure [Pa].""" + """The pressure in the reactor [Pa].""" return _cantera.reactor_pressure(self.__reactor_id) def advance(self, time): - """Advance the state of the reactor in time from the current + """Deprecated. + Advance the state of the reactor in time from the current time to time 'time'. Note: this method is deprecated. See class ReactorNet.""" return _cantera.reactor_advance(self.__reactor_id, time) def step(self, time): - """Take one internal time step from the current time toward + """Deprecated. + Take one internal time step from the current time toward time 'time'. Note: this method is deprecated. See class ReactorNet.""" return _cantera.reactor_step(self.__reactor_id, time) - def massFraction(self, k): - """Mass fraction of species k.""" - if type(k) == types.StringType: - kk = self._contents.speciesIndex(k) + def massFraction(self, s): + """The mass fraction of species s, specified either by name or + index number. + >>> y1 = r.massFraction(7) + >>> y2 = r.massFraction('CH3O') + """ + if type(s) == types.StringType: + kk = self._contents.speciesIndex(s) else: - kk = k + kk = s return _cantera.reactor_massFraction(self.__reactor_id, kk) def massFractions(self): + """Return an array of the species mass fractions.""" nsp = self._contents.nSpecies() y = zeros(nsp,'d') for k in range(nsp): @@ -152,30 +161,49 @@ class ReactorBase: return y def moleFractions(self): + """Return an array of the species mole fractions.""" y = self.massFractions() self._contents.setMassFractions(y) return self._contents.moleFractions() - def moleFraction(self, k): - """Mole fraction of species k.""" - if type(k) == types.StringType: - kk = self._contents.speciesIndex(k) + def moleFraction(self, s): + """The mole fraction of species s, specified either by name or + index number. + >>> x1 = r.moleFraction(7) + >>> x2 = r.moleFraction('CH3O') + """ + if type(s) == types.StringType: + kk = self._contents.speciesIndex(s) else: - kk = k + kk = s x = self.moleFractions() return x[kk] def inlets(self): - """Return the list of flow devices installed on inlets to this reactor.""" + """Return the list of flow devices installed on inlets to this reactor. + This method can be used to access information about the flows entering + the reactor: + >>> for n in r.inlets(): + ... print n.name(), n.massFlowRate() + See MassFlowController, Valve. + """ return self._inlets def outlets(self): """Return the list of flow devices installed on outlets - on this reactor.""" + on this reactor. + >>> for o in r.outlets(): + ... print o.name(), o.massFlowRate() + See MassFlowController, Valve. + """ return self._outlets def walls(self): - """Return the list of walls installed on this reactor.""" + """Return the list of walls installed on this reactor. + >>> for w in r.walls(): + ... print w.name() + See Wall. + """ return self._walls def _addInlet(self, inlet): @@ -193,15 +221,40 @@ class ReactorBase: so that it will not be deleted before this object.""" self._walls.append(wall) - def updateContents(self): + def syncContents(self): """Set the state of the object representing the reactor contents - to the current reactor state.""" + to the current reactor state. + >>> r = Reactor(gas) + >>> (statements that change the state of object 'gas') + >>> r.syncContents(self) + After this statement, the state of object 'gas' is synchronized + with the reactor state. + See 'contents'. + """ self._contents.setState_TRY(self.temperature(), self.density(), self.massFractions()) def contents(self): - updateContents() + """Return an object representing the reactor contents, after first + synchronizing its state with the current reactor state. This method + is useful when some property of the fluid in the reactor is + needed that is not provided by a method of class Reactor. + >>> r = Reactor(gas) + >>> (statements that change the state of object 'gas') + >>> c = r.contents() + >>> print c.gibbs_mole(), c.chemPotentials() + + Note that after calling method 'contents', object 'c' + references the same underlying kernel object as object 'gas' + does. Therefore, all properties of 'c' and 'gas' are + identical. (Remember that Python objects are really C + pointers; at the C level, both point to the same data + structure.) + It is also allowed to write + >>> gas = r.contents() + """ + syncContents() return self._contents @@ -210,14 +263,46 @@ _reservoircount = 0 class Reactor(ReactorBase): """ - A reactor. + Zero-dimensional reactors. Instances of class Reactor represent + zero-dimensional reactors. By default, they are closed (no inlets + or outlets), have fixed volume, and have adiabatic, chemically-intert + walls. These properties may all be changed by adding appropriate + components. + See classes 'Wall', 'MassFlowController', and 'Valve'. """ def __init__(self, contents = None, name = '', volume = 1.0, energy = 'on', verbose = 0): """ - Create a Reactor instance, and if 'contents' is specified, - insert it. + contents - Reactor contents. If not specified, the reactor is + initially empty. In this case, call method insert to specify + the contents. + + name - Used only to identify this reactor in output. If not + specified, defaults to 'Reactor_n', where n is an integer + assigned in the order Reactor objects are created. + + volume - Initial reactor volume. Defaults to 1 m^3. + + energy - Set to 'on' or 'off'. If set to 'off', the energy + equation is not solved, and the temperature is held at its + initial value. The default in 'on'. + + verbose - if set to a non-zero value, additional diagnostic + information will be printed. + + Some examples showing how to create Reactor objects are shown below. + >>> gas = GRI30() + >>> r1 = Reactor(gas) + This is equivalent to: + >>> r1 = Reactor() + >>> r1.insert(gas) + Arguments may be specified using keywords in any order: + >>> r2 = Reactor(contents = gas, energy = 'off', + ... name = 'isothermal_reactor') + >>> r3 = Reactor(contents = gas, name = 'adiabatic_reactor') + Here's an array of reactors: + >>> reactor_array = [Reactor(), Reactor(gas), Reactor(Air())] """ global _reactorcount if name == '': @@ -230,11 +315,34 @@ class Reactor(ReactorBase): class Reservoir(ReactorBase): """ - A reservoir is a reactor with a constant state. Class Reservoir - derives from class ReactorBase, and overloads method advance to do - nothing. + A reservoir is a reactor with a constant state. The temperature, + pressure, and chemical composition in a reservoir never change from + their initial values. """ - def __init__(self, contents = None, name = '', verbose = 0): + def __init__(self, contents = None, name = '', verbose = 0): + """ + contents - Reservoir contents. If not specified, the reservoir is + initially empty. In this case, call method insert to specify + the contents. + + name - Used only to identify this reservoir in output. If not + specified, defaults to 'Reservoir_n', where n is an integer + assigned in the order Reservoir objects are created. + + verbose - if set to a non-zero value, additional diagnostic + information will be printed. + + Some examples showing how to create Reservoir objects are shown below. + >>> gas = GRI30() + >>> res1 = Reservoir(gas) + This is equivalent to: + >>> res1 = Reactor() + >>> res1.insert(gas) + Arguments may be specified using keywords in any order: + >>> res2 = Reservoir(contents = Air(), + ... name = 'environment') + >>> res3 = Reservoir(contents = gas, name = 'upstream_state') + """ global _reservoircount if name == '': name = 'Reservoir_'+`_reservoircount` @@ -243,7 +351,7 @@ class Reservoir(ReactorBase): name = name, verbose = verbose, type = 2) def advance(self, time): - """Do nothing.""" + """Deprecated. Do nothing.""" pass @@ -271,36 +379,36 @@ class FlowDevice: _cantera.flowdev_del(self.__fdev_id) def name(self): + """The name specified when initially constructed.""" return self._name def ready(self): """ - Returns true if the device is ready to use. + Deprecated. Returns true if the device is ready to use. """ return _cantera.flowdev_ready(self.__fdev_id) def massFlowRate(self): - """ - Mass flow rate (kg/s). - """ + """Mass flow rate (kg/s). """ return _cantera.flowdev_massFlowRate(self.__fdev_id) def setSetpoint(self, v): """ - Set the set point. + Deprecated. Set the set point. """ _cantera.flowdev_setSetpoint(self.__fdev_id, v) def setpoint(self): """ - The setpoint value. + Deprecated. The setpoint value. """ return _cantera.flowdev_setpoint(self.__fdev_id) def install(self, upstream, downstream): """ Install the device between the upstream and downstream - reactors. + reactors or reservoirs. + >>> f.install(upstream = reactor1, downstream = reservoir2) """ if self._verbose: print @@ -309,15 +417,61 @@ class FlowDevice: downstream._addInlet(self) _cantera.flowdev_install(self.__fdev_id, upstream.reactor_id(), downstream.reactor_id()) - def setParameters(self, c): + def _setParameters(self, c): params = array(c,'d') n = len(params) return _cantera.flowdev_setParameters(self.__fdev_id, n, params) + _mfccount = 0 class MassFlowController(FlowDevice): - def __init__(self, upstream=None, downstream=None, name='', verbose=0): + + """Mass flow controllers. A mass flow controller maintains a + constant mass flow rate independent of upstream and downstream + conditions. The equation used to compute the mass flow rate is + \f[ \dot m = \dot m_0, \f] where \f$ \dot m_0 \f$ is a + non-negative value specified when the object is constructed or set + by calling method setMassFlowRate. + + Unlike a real mass flow controller, a MassFlowController object + will maintain the flow even if the downstream pressure is greater + than the upstream pressure. This allows simple implementation of + loops, in which exhaust gas from a reactor is fed back into it + through an inlet. But note that this capability should be used + with caution, since no account is taken of the work required to do + this. + + A mass flow controller is assumed to be adiabatic, non-reactive, + and have negligible volume, so that it is 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 mass flow controller, and the pressure + difference equals the difference in pressure between the upstream + and downstream reactors. + """ + def __init__(self, upstream=None, + downstream=None, + name='', + verbose=0, mdot = 0.0): + """ + upstream - upstream reactor or reservoir. + + downstream - downstream reactor or reservoir. + + name - name used to identify the mass flow controller in output. + If no name is specified, it defaults to 'MFC_n', where n is an + integer assigned in the order the MassFlowController object + was created. + + mdot - Mass flow rate [kg/s]. This mass flow rate will be + maintained, independent of unstream and downstream conditions, + unless reset by calling method 'setMassFlowRate'. + + verbose - if set to a positive integer, additional diagnostic + information will be printed. + + """ global _mfccount if name == '': name = 'MFC_'+`_mfccount` @@ -325,18 +479,79 @@ class MassFlowController(FlowDevice): FlowDevice.__init__(self,1,name,verbose) if upstream and downstream: self.install(upstream, downstream) + if mdot: + self.set(mdot = mdot) - def setMassFlowRate(self, mdot): + def _setMassFlowRate(self, mdot): + """Set or reset the mass flow rate to 'mdot' [kg/s]. + """ if self._verbose: print self._name+': setting mdot to '+`mdot`+' kg/s' self.setSetpoint(mdot) + def set(self, mdot = 0.0): + """Set the mass flow rate [kg/s]. + + >>> mfc.set(mdot = 0.2) + """ + self.setSetpoint(mdot) + _valvecount = 0 class Valve(FlowDevice): + """Valves. In Cantera, a Valve object is a flow devices with mass + flow rate proportional to the pressure drop across it. The equation + used to compute the mass flow rate is + \f[ \dot m = K_v (P_1 - P_2) \f] + if \f$ P_1 > P_2. \f$ + Otherwise, + \f$ \dot m = 0 \f$. It is never possible for the flow to reverse + and go from the downstream to the upstream reactor/reservoir through + a line containing a Valve object. + + 'Valve' objects are often used between an upstream reactor and a + downstream reactor or reservoir to maintain them both at nearly the + same pressure. By setting the constant \f$ K_v \f$ to a + sufficiently large value, very small pressure differences will + result in flow between the reactors that counteracts the pressure + difference. + + Since the mass flow rate is assumed to be linear in \f$ \Delta P \f$, + these objects do not model real, physical valves, in which the flow rate + is proportional to \f$ \sqrt(\Delta P) \f$ for small pressure + differences, and becomes independent of \f$ \Delta P \f$ when + it becomes large (choked flow). Perhaps the name of this class should + be changed to avoid confusion with real valves -- if you have suggestions, + post a comment at the Cantera User's Group site. + + A Valve is assumed to be adiabatic, non-reactive, and have + negligible internal volume, so that it is 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 Valve, and the pressure difference equals + the difference in pressure between the upstream and downstream + reactors. + + """ def __init__(self, upstream=None, downstream=None, - name='', K = 0.0, verbose=0): + name='', Kv = 0.0, verbose=0): + """ + upstream - upstream reactor or reservoir. + + downstream - downstream reactor or reservoir. + + name - name used to identify the valve in output. + If no name is specified, it defaults to 'Valve_n', where n is an + integer assigned in the order the Valve object + was created. + + Kv - the constant in the mass flow rate equation. + + verbose - if set to a positive integer, additional diagnostic + information will be printed. + + """ global _valvecount if name == '': name = 'Valve_'+`_valvecount` @@ -344,15 +559,17 @@ class Valve(FlowDevice): FlowDevice.__init__(self,3,name,verbose) if upstream and downstream: self.install(upstream, downstream) - self.setValveCoeff(K) + self.setValveCoeff(Kv) + def setValveCoeff(self, v): + """Set or reset the valve coefficient \f$ K_v \f$.""" vv = zeros(1,'d') vv[0] = v if self._verbose: print print self._name+': setting valve coefficient to '+`v`+' kg/Pa-s' - self.setParameters(vv) + self._setParameters(vv) #------------- Wall --------------------------- @@ -361,8 +578,8 @@ _wallcount = 0 class Wall: """ - A Wall separates two reactors. Any number of walls may be created - between any pair of reactors. + Reactor walls. + A Wall separates two reactors, or a reactor and a reservoir. """ def __init__(self, left=None, right=None, name = '', A = 1.0, K = 0.0, U = 0.0, @@ -434,6 +651,8 @@ class Wall: return _cantera.wall_setHeatFlux(self.__wall_id, n) def setExpansionRateCoeff(self, k): + """Set the coefficient K that determines the expansion rate + resulting from a unit pressure drop.""" _cantera.wall_setExpansionRateCoeff(self.__wall_id, k) def setExpansionRate(self, vfunc=None): @@ -443,7 +662,19 @@ class Wall: n = 0 if vfunc: n = vfunc.func_id() _cantera.wall_setExpansionRate(self.__wall_id, n) - + + def vdot(self): + """Rate of volume change [m^3]. A positive value corresponds + to the left-hand reactor volume increasing, and the right-hand + reactor volume decreasing.""" + return _cantera.wall_vdot(self.__wall_id) + + def heatFlowRate(self): + """Rate of heat flow through the wall. A positive value + corresponds to heat flowing from the left-hand reactor to the + right-hand one.""" + return _cantera.wall_Q(self.__wall_id) + def install(self, left, right): left._addWall(self) right._addWall(self) diff --git a/Cantera/src/Func1.h b/Cantera/src/Func1.h index 4bb03ab39..93b32e496 100644 --- a/Cantera/src/Func1.h +++ b/Cantera/src/Func1.h @@ -23,6 +23,7 @@ namespace Cantera { const int FourierFuncType = 1; const int PolyFuncType = 2; const int ArrheniusFuncType = 3; + const int GaussianFuncType = 4; const int SumFuncType = 20; const int DiffFuncType = 25; const int ProdFuncType = 30; @@ -44,6 +45,24 @@ namespace Cantera { }; + class Gaussian : public Func1 { + public: + Gaussian(double A, double t0, double tau) { + m_A = A; + m_t0 = t0; + m_tau = tau; + } + virtual ~Gaussian() {} + virtual doublereal eval(doublereal t) { + doublereal x = (t - m_t0)/m_tau; + return m_A*exp(-x*x); + } + protected: + doublereal m_A, m_t0, m_tau; + private: + }; + + /** * Polynomial of degree n. */ diff --git a/Cantera/src/zeroD/Wall.cpp b/Cantera/src/zeroD/Wall.cpp index ca83472dd..5b1c61fe6 100644 --- a/Cantera/src/zeroD/Wall.cpp +++ b/Cantera/src/zeroD/Wall.cpp @@ -71,7 +71,7 @@ namespace Cantera { /** * The heat flux is given by - * \f[ Q = h A (T_{left} - T_{right}) + G(t) \f] + * \f[ Q = h A (T_{left} - T_{right}) + A G(t) \f] * where h is the heat transfer coefficient, and * \f$ G(t) \f$ is a specified function of time. */