From 4f0c424fb7e5e722513a35da219f089aa277b2dc Mon Sep 17 00:00:00 2001 From: Dave Goodwin Date: Sun, 22 Feb 2004 17:10:05 +0000 Subject: [PATCH] *** empty log message *** --- Cantera/clib/src/ct.cpp | 8 ++- Cantera/clib/src/ct.h | 2 +- Cantera/python/Cantera/Kinetics.py | 11 +++- Cantera/python/examples/catcomb.py | 33 ++++++++--- Cantera/python/examples/dustygas.py | 30 ++++++++++ Cantera/python/examples/flame1.py | 5 +- Cantera/python/examples/flame2.py | 22 ++++---- Cantera/python/src/ctkinetics_methods.cpp | 4 +- Cantera/python/src/methods.h | 4 ++ Cantera/python/src/pycantera.cpp | 8 ++- Cantera/src/InterfaceKinetics.cpp | 16 +++++- Cantera/src/InterfaceKinetics.h | 26 ++++++--- Cantera/src/Kinetics.h | 67 +++++++++++++---------- data/inputs/ptcombust.cti | 3 +- 14 files changed, 167 insertions(+), 72 deletions(-) create mode 100644 Cantera/python/examples/dustygas.py diff --git a/Cantera/clib/src/ct.cpp b/Cantera/clib/src/ct.cpp index f181d3a45..1fe47a655 100755 --- a/Cantera/clib/src/ct.cpp +++ b/Cantera/clib/src/ct.cpp @@ -124,7 +124,7 @@ extern "C" { int DLL_EXPORT phase_getMassFractions(int n, int leny, double* y) { ThermoPhase* p = ph(n); if (leny >= p->nSpecies()) { - p->getMassFractions(leny, y); + p->getMassFractions(y); return 0; } else @@ -678,11 +678,13 @@ extern "C" { catch (CanteraError) {return -1;} } - int DLL_EXPORT kin_getRevRateConstants(int n, int len, double* krev) { + int DLL_EXPORT kin_getRevRateConstants(int n, int doIrreversible, int len, double* krev) { try { Kinetics* k = kin(n); + bool doirrev = false; + if (doIrreversible != 0) doirrev = true; if (len >= k->nReactions()) { - k->getRevRateConstants(krev); + k->getRevRateConstants(krev, doirrev); return 0; } else diff --git a/Cantera/clib/src/ct.h b/Cantera/clib/src/ct.h index 3c49faeb7..979d630e8 100755 --- a/Cantera/clib/src/ct.h +++ b/Cantera/clib/src/ct.h @@ -103,7 +103,7 @@ extern "C" { int DLL_IMPORT kin_getEquilibriumConstants(int n, int len, double* kc); int DLL_IMPORT kin_getFwdRateConstants(int n, int len, double* kfwd); - int DLL_IMPORT kin_getRevRateConstants(int n, int len, double* krev); + int DLL_IMPORT kin_getRevRateConstants(int n, int doIrreversible, int len, double* krev); int DLL_IMPORT kin_getActivationEnergies(int n, int len, double* E); int DLL_IMPORT kin_getCreationRates(int n, int len, double* cdot); diff --git a/Cantera/python/Cantera/Kinetics.py b/Cantera/python/Cantera/Kinetics.py index 3f67eddda..a8afe475a 100755 --- a/Cantera/python/Cantera/Kinetics.py +++ b/Cantera/python/Cantera/Kinetics.py @@ -190,17 +190,22 @@ class Kinetics: return _cantera.kin_getarray(self.ckin,30) def equilibriumConstants(self): + """Equilibrium constants in concentration units for all reactions.""" return _cantera.kin_getarray(self.ckin,40) def activationEnergies(self): + """Activation energies in Kelvin for all reactions.""" return _cantera.kin_getarray(self.ckin,32) def fwdRateConstants(self): return _cantera.kin_getarray(self.ckin,34) - def revRateConstants(self): - return _cantera.kin_getarray(self.ckin,36) - + def revRateConstants(self, doIrreversible = 0): + if doIrreversible: + return _cantera.kin_getarray(self.ckin,35) + else: + return _cantera.kin_getarray(self.ckin,36) + def creationRates(self, phase = None): c = _cantera.kin_getarray(self.ckin,50) if phase: diff --git a/Cantera/python/examples/catcomb.py b/Cantera/python/examples/catcomb.py index adf9861b6..e4ccd2663 100644 --- a/Cantera/python/examples/catcomb.py +++ b/Cantera/python/examples/catcomb.py @@ -1,4 +1,4 @@ -# CATCOMB -- Catalytic combustion on platinum. +# CATCOMB -- Catalytic combustion of methane on platinum. # # This script solves a catalytic combustion problem. A stagnation flow # is set up, with a gas inlet 10 cm from a platinum surface at 900 @@ -82,27 +82,42 @@ surf_phase.setTemperature(tsurf) # coverages. surf_phase.advanceCoverages(1.0) +# create the object that simulates the stagnation flow, and specify an +# initial grid sim = StagnationFlow(gas = gas, surfchem = surf_phase, grid = initial_grid) +# Objects of class StagnationFlow have members that represent the gas inlet ('inlet') and the surface ('surface'). Set some parameters of these objects. sim.inlet.set(mdot = mdot, T = tinlet, X = comp1) sim.surface.set(T = tsurf) +# Set error tolerances sim.set(tol = tol_ss, tol_time = tol_ts) +# Method 'init' must be called before beginning a simulation sim.init() + +# Show the initial solution estimate sim.showSolution() -# start with the energy equation on + + +# Solving problems with stiff chemistry coulpled to flow can require +# a sequential approach where solutions are first obtained for +# simpler problems and used as the initial guess for more difficult +# problems. + +# start with the energy equation on (default is 'off') sim.set(energy = 'on') # disable the surface coverage equations, and turn off all gas and -# surface chemistry +# surface chemistry. sim.surface.setCoverageEqs('off') surf_phase.setMultiplier(0.0); gas.setMultiplier(0.0); -# solve the problem, refining the grid if needed +# solve the problem, refining the grid if needed, to determine the +# non-reacting velocity and temperature distributions sim.solve(loglevel, refine_grid) # now turn on the surface coverage equations, and turn the @@ -119,7 +134,7 @@ for iter in range(6): # problem. sim.showSolution() -#Now switch the inlet to the methane/air composition. +# Now switch the inlet to the methane/air composition. sim.inlet.set(X = comp2) # set more stringent grid refinement criteria @@ -138,16 +153,19 @@ sim.save("catcomb.xml", "soln1") # save selected solution components in a CSV file for plotting in # Excel or MATLAB. +# These methods return arrays containing the values at all grid points z = sim.flow.grid() T = sim.T() u = sim.u() V = sim.V() + f = open('catcomb.csv','w') -writeCSV(f, ['z (m)', 'u (m/s)', 'V (1/s)', 'T (K)'] +writeCSV(f, ['z (m)', 'u (m/s)', 'V (1/s)', 'T (K)', 'rho (kg/m3'] + list(gas.speciesNames())) for n in range(sim.flow.nPoints()): sim.setGasState(n) - writeCSV(f, [z[n], u[n], V[n], T[n]]+list(gas.moleFractions())) + writeCSV(f, [z[n], u[n], V[n], T[n], gas.density()] + +list(gas.moleFractions())) # write the surface coverages to the CSV file cov = sim.coverages() @@ -159,4 +177,5 @@ f.close() print 'solution saved to catcomb.csv' +# show some statistics sim.showStats() diff --git a/Cantera/python/examples/dustygas.py b/Cantera/python/examples/dustygas.py new file mode 100644 index 000000000..70949a840 --- /dev/null +++ b/Cantera/python/examples/dustygas.py @@ -0,0 +1,30 @@ +""" + Dusty Gas transport model. + + The Dusty Gas model is a mulicomponent transport model for gas + transport through the pores of a stationary porous medium. This + example shows how to create a transport manager that implements the + Dusty Gas model and use it to compute the multicomponent diffusion + coefficients. + + """ + +from Cantera import * +from Cantera.DustyGasTransport import * + +# create a gas-phase object to represent the gas in the pores +g = importPhase('h2o2.cti') + +# set the gas state +g.setState_TPX(500.0, OneAtm, "OH:1, H:2, O2:3") + +# create a Dusty Gas transport manager for this phase +d = DustyGasTransport(g) + +# set its parameters +d.set(porosity = 0.2, tortuosity = 4.0, + pore_radius = 1.5e-7, diameter = 1.5e-6) # lengths in meters + +# print the multicomponent diffusion coefficients +print d.multiDiffCoeffs() + diff --git a/Cantera/python/examples/flame1.py b/Cantera/python/examples/flame1.py index c7fa58b64..5af347037 100755 --- a/Cantera/python/examples/flame1.py +++ b/Cantera/python/examples/flame1.py @@ -73,11 +73,12 @@ T = f.T() u = f.u() V = f.V() fcsv = open('flame1.csv','w') -writeCSV(fcsv, ['z (m)', 'u (m/s)', 'V (1/s)', 'T (K)'] +writeCSV(fcsv, ['z (m)', 'u (m/s)', 'V (1/s)', 'T (K)', 'rho (kg/m3)'] + list(gas.speciesNames())) for n in range(f.flame.nPoints()): f.setGasState(n) - writeCSV(fcsv, [z[n], u[n], V[n], T[n]]+list(gas.moleFractions())) + writeCSV(fcsv, [z[n], u[n], V[n], T[n], gas.density()] + +list(gas.moleFractions())) fcsv.close() print 'solution saved to flame1.csv' diff --git a/Cantera/python/examples/flame2.py b/Cantera/python/examples/flame2.py index 7b45893b3..da17b04f6 100755 --- a/Cantera/python/examples/flame2.py +++ b/Cantera/python/examples/flame2.py @@ -1,8 +1,6 @@ # -# FLAME1 - A burner-stabilized flat flame -# -# This script simulates a burner-stablized lean hydrogen-oxygen flame -# at low pressure. +# FLAME2 - A burner-stabilized, premixed methane/air flat flame +# with multicomponent transport properties # from Cantera import * from Cantera.OneD import * @@ -35,9 +33,10 @@ refine_grid = 1 # 1 to enable refinement, 0 to ################ create the gas object ######################## # -# This object will be used to evaluate all thermodynamic, kinetic, -# and transport properties -# +# This object will be used to evaluate all thermodynamic, kinetic, and +# transport properties. It is created with two transport managers, to +# enable switching from mixture-averaged to multicomponent transport +# on the last solution. gas = GRI30('Mix') gas.addTransportModel('Multi') @@ -71,19 +70,20 @@ gas.switchTransportModel('Multi') f.flame.setTransportModel(gas) f.solve(loglevel, refine_grid) f.save('ch4_flame1.xml','energy_multi', - 'solution with the energy equation enabled') + 'solution with the energy equation enabled and multicomponent transport') -# write the velocity, temperature, and mole fractions to a CSV file +# write the velocity, temperature, density, and mole fractions to a CSV file z = f.flame.grid() T = f.T() u = f.u() V = f.V() fcsv = open('flame2.csv','w') -writeCSV(fcsv, ['z (m)', 'u (m/s)', 'V (1/s)', 'T (K)'] +writeCSV(fcsv, ['z (m)', 'u (m/s)', 'V (1/s)', 'T (K)', 'rho (kg/m3)'] + list(gas.speciesNames())) for n in range(f.flame.nPoints()): f.setGasState(n) - writeCSV(fcsv, [z[n], u[n], V[n], T[n]]+list(gas.moleFractions())) + writeCSV(fcsv, [z[n], u[n], V[n], T[n], gas.density()] + +list(gas.moleFractions())) fcsv.close() print 'solution saved to flame2.csv' diff --git a/Cantera/python/src/ctkinetics_methods.cpp b/Cantera/python/src/ctkinetics_methods.cpp index 77186cc61..a5209b51d 100644 --- a/Cantera/python/src/ctkinetics_methods.cpp +++ b/Cantera/python/src/ctkinetics_methods.cpp @@ -157,8 +157,10 @@ kin_getarray(PyObject *self, PyObject *args) case 34: iok = kin_getFwdRateConstants(kin, nrxns, xd); break; + case 35: + iok = kin_getRevRateConstants(kin, 1, nrxns, xd); case 36: - iok = kin_getRevRateConstants(kin, nrxns, xd); + iok = kin_getRevRateConstants(kin, 0, nrxns, xd); break; case 40: iok = kin_getEquilibriumConstants(kin, nrxns, xd); diff --git a/Cantera/python/src/methods.h b/Cantera/python/src/methods.h index 420e2100d..8843915fa 100644 --- a/Cantera/python/src/methods.h +++ b/Cantera/python/src/methods.h @@ -238,6 +238,10 @@ static PyMethodDef ct_methods[] = { {"func_del", py_func_del, METH_VARARGS}, {"func_value", py_func_value, METH_VARARGS}, +#ifdef INCL_USER_PYTHON +#include "usermethods.h" +#endif + {NULL, NULL} /* sentinel */ }; diff --git a/Cantera/python/src/pycantera.cpp b/Cantera/python/src/pycantera.cpp index de97be699..3dd882cd9 100644 --- a/Cantera/python/src/pycantera.cpp +++ b/Cantera/python/src/pycantera.cpp @@ -17,7 +17,6 @@ #include "ct.h" #include "ctxml.h" -//#include "ctstagn.h" #include "ctsurf.h" #include "ctbdry.h" #include "ctrpath.h" @@ -39,7 +38,6 @@ static PyObject *ErrorObject; #include "ctkinetics_methods.cpp" #include "cttransport_methods.cpp" #include "ctxml_methods.cpp" -//#include "ctflow_methods.cpp" #include "ctfuncs.cpp" #include "ctsurf_methods.cpp" #include "ctbndry_methods.cpp" @@ -48,6 +46,12 @@ static PyObject *ErrorObject; #include "ctfunc_methods.cpp" #include "ctonedim_methods.cpp" +#ifdef INCL_USER_PYTHON +#include "ctuser.h" +#include "ctuser_methods.cpp" +#endif + + #include "methods.h" extern "C" { diff --git a/Cantera/src/InterfaceKinetics.cpp b/Cantera/src/InterfaceKinetics.cpp index 97bb8e40e..d92c7edeb 100644 --- a/Cantera/src/InterfaceKinetics.cpp +++ b/Cantera/src/InterfaceKinetics.cpp @@ -313,12 +313,22 @@ namespace Cantera { * Update the rates of progress of the reactions in the reaciton * mechanism. This routine operates on internal data. */ - void InterfaceKinetics::getRevRateConstants(doublereal* krev) { + void InterfaceKinetics::getRevRateConstants(doublereal* krev, bool doIrreversible) { getFwdRateConstants(krev); - const vector_fp& rkc = m_kdata->m_rkcn; - multiply_each(krev, krev + nReactions(), rkc.begin()); + if (doIrreversible) { + doublereal *tmpKc = m_kdata->m_ropnet.begin(); + getEquilibriumConstants(tmpKc); + for (int i = 0; i < m_ii; i++) { + krev[i] /= tmpKc[i]; + } + } + else { + const vector_fp& rkc = m_kdata->m_rkcn; + multiply_each(krev, krev + nReactions(), rkc.begin()); + } } + void InterfaceKinetics::getActivationEnergies(doublereal *E) { copy(m_E.begin(), m_E.end(), E); } diff --git a/Cantera/src/InterfaceKinetics.h b/Cantera/src/InterfaceKinetics.h index c8df81f45..19d743a56 100644 --- a/Cantera/src/InterfaceKinetics.h +++ b/Cantera/src/InterfaceKinetics.h @@ -23,7 +23,6 @@ #include "utilities.h" #include "RateCoeffMgr.h" #include "ReactionStoichMgr.h" -//#include "StoichManager.h" namespace Cantera { @@ -35,8 +34,9 @@ namespace Cantera { class SurfPhase; class ImplicitSurfChem; + /** - * Holds mechanism-specific data. + * Holds mechanism-specific data. */ class InterfaceKineticsData { public: @@ -57,6 +57,10 @@ namespace Cantera { }; + /// + /// A kinetics manager for heterogeneous reaction mechanisms. The reactions are + /// assumed to occur at a 2D interface between two 3D phases. + /// class InterfaceKinetics : public Kinetics { public: @@ -73,6 +77,7 @@ namespace Cantera { */ InterfaceKinetics(thermo_t* thermo = 0); + /// Destructor. virtual ~InterfaceKinetics(); @@ -83,13 +88,14 @@ namespace Cantera { virtual int ID() { return cInterfaceKinetics; } /** - * Identifies the subclass of the Kinetics manager type. + * Identifies the subclass of the kinetics manager type. * These are listed in mix_defs.h. */ virtual int type() { return cInterfaceKinetics; } /** - * Set the electric potential in the nth phase + * Set the electric potential in the nth phase. + * @deprecated * * @param n phase Index in this kinetics object. * @param V Electric potential (volts) @@ -100,11 +106,14 @@ namespace Cantera { } //@} - /** - * @name Reaction Rates Of Progress - */ + + + /// + /// @name Reaction Rates Of Progress + /// //@{ + /** * Forward rates of progress. * Return the forward rates of progress in array fwdROP, which @@ -117,6 +126,7 @@ namespace Cantera { copy(m_kdata->m_ropf.begin(), m_kdata->m_ropf.end(), fwdROP); } + /** * Reverse rates of progress. * Return the reverse rates of progress in array revROP, which @@ -266,7 +276,7 @@ namespace Cantera { virtual void getFwdRateConstants(doublereal* kfwd); - virtual void getRevRateConstants(doublereal* krev); + virtual void getRevRateConstants(doublereal* krev, bool doIrreversible = false); virtual void getActivationEnergies(doublereal *E); //@} diff --git a/Cantera/src/Kinetics.h b/Cantera/src/Kinetics.h index 638ec4f0a..d7c229cdf 100755 --- a/Cantera/src/Kinetics.h +++ b/Cantera/src/Kinetics.h @@ -6,7 +6,7 @@ * $Date$ */ -// Copyright 2001 California Institute of Technology +// Copyright 2001-2004 California Institute of Technology /** * @defgroup kineticsGroup Kinetics @@ -21,14 +21,16 @@ namespace Cantera { + // forward references class ReactionData; - /** - * Public interface for kinetics managers. This class serves as a - * base class to derive 'kinetics managers', which are classes - * that manage homogeneous chemistry within one phase, or - * heterogeneous chemistry at one interface. - */ + + /// + /// Public interface for kinetics managers. This class serves as a + /// base class to derive 'kinetics managers', which are classes + /// that manage homogeneous chemistry within one phase, or + /// heterogeneous chemistry at one interface. + /// class Kinetics { public: @@ -46,26 +48,22 @@ namespace Cantera { /** * This Constructor initializes with a starting phase. - * All of the appropriate entries that addPhase() (below) - * sets up are also done here. */ Kinetics(thermo_t* thermo) : m_ii(0), m_index(-1), m_surfphase(-1) { if (thermo) { addPhase(*thermo); -// m_start.push_back(0); -// if (thermo->eosType() == cSurf) m_surfphase = nPhases(); -// m_thermo.push_back(thermo); -// m_phaseindex[m_thermo.back()->id()] = nPhases(); } } /// Destructor. Does nothing. virtual ~Kinetics() {} + /// For internal use. int index(){ return m_index; } void setIndex(int index) { m_index = index; } + /** * Identifies the subclass of the Kinetics manager type. * These are listed in mix_defs.h. @@ -78,13 +76,14 @@ namespace Cantera { //@} + /** * @name Information/Lookup Functions about Phases and Species */ //@{ /** - * Return the number of phases defined within the kinetics + * The number of phases defined within the kinetics * object. */ int nPhases() const { return m_thermo.size(); } @@ -101,6 +100,7 @@ namespace Cantera { * the Kinetics object. * (HKM -> unfound object will create another entry in the * map, suggest rewriting this function) + * @todo rewrite. */ int phaseIndex(string ph) { return m_phaseindex[ph] - 1; } @@ -133,10 +133,10 @@ namespace Cantera { const thermo_t& phase(int n=0) const { return *m_thermo[n]; } /** - * Returns the total number of species in all phases - * participating in the kinetics mechanism. This is useful to - * dimension arrays for use in calls to methods that return - * the species production rates, for example. + * The total number of species in all phases participating in + * the kinetics mechanism. This is useful to dimension arrays + * for use in calls to methods that return the species + * production rates, for example. */ int nTotalSpecies() const { int n=0, np; @@ -155,9 +155,25 @@ namespace Cantera { */ int start(int n) { return m_start[n]; } + /** - * This method returns the index of a species in the source - * term vector for this kinetics object. + * The location of species k of phase n in species arrays. + * Kinetics manager classes return species production rates in + * flat arrays, with the species of each phases following one + * another, in the order the phases were added. This method + * is useful to find the value for a particular species of a + * particular phase in arrrays returned from methods like + * getCreationRates that return an array of species-specific + * quantities. + * + * Example: suppose a heterogeneous mechanism involves three + * phases. The first contains 12 species, the second 26, and + * the third 3. Then species arrays must have size at least + * 41, and positions 0 - 11 are the values for the species in + * the first phase, positions 12 - 37 are the values for the + * species in the second phase, etc. Then + * kineticsSpeciesIndex(7, 0) = 7, kineticsSpeciesIndex(4, 1) + * = 16, and kineticsSpeciesIndex(2, 2) = 40. * * @param k species index * @param n phase index for the species @@ -244,16 +260,7 @@ namespace Cantera { * manager) and returns the species' owning ThermoPhase object. */ thermo_t& speciesPhase(int k) { - return thermo(speciesPhaseIndex(k)); -// int np = m_start.size(); -// for (int n = np-1; n >= 0; n--) { -// if (k >= m_start[n]) { -// return thermo(n); -// } -// } -// throw CanteraError("speciesPhase", -// "illegal species index: "+int2str(k)); - + return thermo(speciesPhaseIndex(k)); } /** diff --git a/data/inputs/ptcombust.cti b/data/inputs/ptcombust.cti index 7de703be4..4cfd1a643 100644 --- a/data/inputs/ptcombust.cti +++ b/data/inputs/ptcombust.cti @@ -249,7 +249,8 @@ surface_reaction( "H(S) + OH(S) <=> H2O(S) + PT(S)", [3.70000E+21, 0, 17400]) surface_reaction( "OH(S) + OH(S) <=> H2O(S) + O(S)", [3.70000E+21, 0, 48200]) # Reaction 15 -surface_reaction( "CO + PT(S) => CO(S)", [1.61800E+20, 0.5, 0], order = "PT(S):2") +surface_reaction( "CO + PT(S) => CO(S)", [1.61800E+20, 0.5, 0], + order = "PT(S):2") # Reaction 16 surface_reaction( "CO(S) => CO + PT(S)", [1.00000E+13, 0, 125500])