diff --git a/Cantera/src/ConstDensityThermo.cpp b/Cantera/src/ConstDensityThermo.cpp index b817734b5..d7b476721 100755 --- a/Cantera/src/ConstDensityThermo.cpp +++ b/Cantera/src/ConstDensityThermo.cpp @@ -111,6 +111,13 @@ namespace Cantera { m_tlast = tnow; } } + + void ConstDensityThermo::setParametersFromXML(const XML_Node& eosdata) { + eosdata.require("model","Incompressible"); + doublereal rho = getFloat(eosdata, "density", "-"); + setDensity(rho); + } + } diff --git a/Cantera/src/ConstDensityThermo.h b/Cantera/src/ConstDensityThermo.h index b02aabdf7..1a53493e7 100755 --- a/Cantera/src/ConstDensityThermo.h +++ b/Cantera/src/ConstDensityThermo.h @@ -170,6 +170,7 @@ namespace Cantera { setDensity(c[0]); } + virtual void setParametersFromXML(const XML_Node& eosdata); protected: diff --git a/Cantera/src/Constituents.cpp b/Cantera/src/Constituents.cpp index 492f24fd0..27b9dd178 100755 --- a/Cantera/src/Constituents.cpp +++ b/Cantera/src/Constituents.cpp @@ -216,6 +216,10 @@ namespace Cantera { m_Elements->addUniqueElement(e); } + void Constituents::addElementsFromXML(const XML_Node& phase) { + m_Elements->addElementsFromXML(phase); + } + /******************************************************************* * * freezeElements() diff --git a/Cantera/src/Constituents.h b/Cantera/src/Constituents.h index 9a53f5737..82e4e7f55 100755 --- a/Cantera/src/Constituents.h +++ b/Cantera/src/Constituents.h @@ -8,7 +8,10 @@ * $Revision$ * * $Log$ - * Revision 1.5 2003-11-12 18:58:17 dggoodwin + * Revision 1.6 2004-06-09 00:59:24 dggoodwin + * moved functions to read XML input files out of importCTML.cpp and into specific classes + * + * Revision 1.5 2003/11/12 18:58:17 dggoodwin * *** empty log message *** * * Revision 1.4 2003/09/03 18:15:50 hkmoffa @@ -134,6 +137,7 @@ namespace Cantera { void addUniqueElement(const XML_Node& e); + void addElementsFromXML(const XML_Node& phase); /// Prohibit addition of more elements, and prepare to add /// species. diff --git a/Cantera/src/EdgePhase.h b/Cantera/src/EdgePhase.h index 3cb665a0f..13ddd2eda 100644 --- a/Cantera/src/EdgePhase.h +++ b/Cantera/src/EdgePhase.h @@ -28,7 +28,7 @@ namespace Cantera { EdgePhase(doublereal n0 = 0.0); virtual ~EdgePhase() {} virtual int eosType() const { return cEdge; } - + virtual void setParametersFromXML(const XML_Node& eosdata); }; } diff --git a/Cantera/src/Elements.h b/Cantera/src/Elements.h index f5a0daf46..de3af4f2e 100644 --- a/Cantera/src/Elements.h +++ b/Cantera/src/Elements.h @@ -102,6 +102,8 @@ namespace Cantera { doublereal weight = -12345.0); void addUniqueElement(const XML_Node& e); + void addElementsFromXML(const XML_Node& phase); + /** * Prohibit addition of more elements, and prepare to add * species. diff --git a/Cantera/src/Makefile.in b/Cantera/src/Makefile.in index a550d31f6..50fe913e8 100755 --- a/Cantera/src/Makefile.in +++ b/Cantera/src/Makefile.in @@ -26,7 +26,7 @@ BASE = State.o Elements.o Constituents.o stringUtils.o misc.o importCTML.o xml.o Phase.o DenseMatrix.o ctml.o funcs.o ctvector.o phasereport.o ct2ctml.o # thermodynamic properties -THERMO = $(BASE) ThermoPhase.o IdealGasPhase.o ConstDensityThermo.o StoichSubstance.o SpeciesThermoFactory.o ThermoFactory.o +THERMO = $(BASE) ThermoPhase.o IdealGasPhase.o ConstDensityThermo.o StoichSubstance.o PureFluidPhase.o SpeciesThermoFactory.o ThermoFactory.o # homogeneous kinetics KINETICS = GRI_30_Kinetics.o KineticsFactory.o GasKinetics.o FalloffFactory.o \ @@ -55,7 +55,7 @@ EVERYTHING = $(KINETICS) $(HETEROKIN) $(ELECTROCHEM) $(EQUIL) $(CK) \ $(TRANSPORT) $(REACTOR) $(RPATH) $(SOLVERS) $(FLOW1D) -PCH = ct_defs.h.gch utilities.h.gch ThermoPhase.h.gch +PCH = ct_defs.h.gch utilities.h.gch ThermoPhase.h.gch Kinetic.h.gch all: config.h $(PCH) @KERNEL@ lib @@ -71,6 +71,11 @@ else @echo 'skipping precompiling header file $*.h' endif +ct_defs.h.gch: ct_defs.h +utilities.h.gch: utilities.h +ThermoPhase.h.gch: ThermoPhase.h +Kinetics.h.gch: Kinetics.h + base: $(BASE) thermo: $(THERMO) diff --git a/Cantera/src/PureFluidPhase.h b/Cantera/src/PureFluidPhase.h index de07c5e8b..71e022826 100644 --- a/Cantera/src/PureFluidPhase.h +++ b/Cantera/src/PureFluidPhase.h @@ -18,6 +18,7 @@ #ifdef INCL_PURE_FLUIDS +#include "mix_defs.h" #include "../../ext/tpx/Sub.h" #include "../../ext/tpx/utils.h" @@ -28,42 +29,11 @@ namespace Cantera { public: - PureFluid() : ThermoPhase(), m_sub(0) {} + PureFluid() : ThermoPhase(), m_sub(0), m_subflag(0), + m_mw(-1.0), m_verbose(true) {} virtual ~PureFluid() { delete m_sub; } - virtual void setParameters(int n, doublereal* c) { - if (n == 1) { - int subflag = int(c[0]); - if (m_sub) delete m_sub; - m_sub = tpx::GetSub(subflag); - if (m_sub == 0) { - throw CanteraError("PureFluid::setParameters", - "could not create new substance object."); - } - m_subflag = subflag; - m_mw = m_sub->MolWt(); - m_weight[0] = m_mw; - setMolecularWeight(0,m_mw); - double one = 1.0; - setMoleFractions(&one); - double cp0_R, h0_RT, s0_R, T0, p; - T0 = 298.15; - if (T0 < m_sub->Tcrit()) { - m_sub->Set(tpx::TX, T0, 1.0); - p = 0.01*m_sub->P(); - } - else { - p = 0.001*m_sub->Pcrit(); - } - m_sub->Set(tpx::TP, T0, p); - - m_spthermo->update_one(0, T0, &cp0_R, &h0_RT, &s0_R); - double s_R = s0_R - log(p/refPressure()); - m_sub->setStdState(h0_RT*GasConstant*298.15/m_mw, - s_R*GasConstant/m_mw, T0, p); - } - } virtual int eosType() const { return cPureFluid; } @@ -247,6 +217,8 @@ namespace Cantera { check(); } + virtual void initThermo(); + virtual void setParametersFromXML(const XML_Node& eosdata); protected: @@ -265,6 +237,7 @@ private: mutable tpx::Substance* m_sub; int m_subflag; doublereal m_mw; + bool m_verbose; }; } diff --git a/Cantera/src/SpeciesThermoFactory.cpp b/Cantera/src/SpeciesThermoFactory.cpp index 373ccbd66..4456a3da6 100755 --- a/Cantera/src/SpeciesThermoFactory.cpp +++ b/Cantera/src/SpeciesThermoFactory.cpp @@ -27,6 +27,8 @@ #include "speciesThermoTypes.h" #include "xml.h" +#include "ctml.h" +using namespace ctml; namespace Cantera { @@ -161,4 +163,182 @@ namespace Cantera { } } + + /** + * Install a NASA polynomial thermodynamic property + * parameterization for species k into a SpeciesThermo instance. + */ + static void installNasaThermoFromXML(SpeciesThermo& sp, int k, + const XML_Node* f0ptr, const XML_Node* f1ptr) { + doublereal tmin0, tmax0, tmin1, tmax1, tmin, tmid, tmax; + + const XML_Node& f0 = *f0ptr; + bool dualRange = false; + if (f1ptr) {dualRange = true;} + tmin0 = fpValue(f0["Tmin"]); + tmax0 = fpValue(f0["Tmax"]); + tmin1 = tmax0; + tmax1 = tmin1 + 0.0001; + if (dualRange) { + tmin1 = fpValue((*f1ptr)["Tmin"]); + tmax1 = fpValue((*f1ptr)["Tmax"]); + } + + vector_fp c0, c1; + if (fabs(tmax0 - tmin1) < 0.01) { + tmin = tmin0; + tmid = tmax0; + tmax = tmax1; + getFloatArray(f0.child("floatArray"), c0, false); + if (dualRange) + getFloatArray(f1ptr->child("floatArray"), c1, false); + else + c1.resize(7,0.0); + } + else if (fabs(tmax1 - tmin0) < 0.01) { + tmin = tmin1; + tmid = tmax1; + tmax = tmax0; + getFloatArray(f1ptr->child("floatArray"), c0, false); + getFloatArray(f0.child("floatArray"), c1, false); + } + else { + throw CanteraError("installNasaThermo", + "non-continuous temperature ranges."); + } + array_fp c(15); + c[0] = tmid; + doublereal p0 = OneAtm; + c[1] = c0[5]; + c[2] = c0[6]; + copy(c0.begin(), c0.begin()+5, c.begin() + 3); + c[8] = c1[5]; + c[9] = c1[6]; + copy(c1.begin(), c1.begin()+5, c.begin() + 10); + sp.install(k, NASA, c.begin(), tmin, tmax, p0); + } + + + /** + * Install a NASA polynomial thermodynamic property + * parameterization for species k. + */ + static void installShomateThermoFromXML(SpeciesThermo& sp, int k, + const XML_Node* f0ptr, const XML_Node* f1ptr) { + doublereal tmin0, tmax0, tmin1, tmax1, tmin, tmid, tmax; + + const XML_Node& f0 = *f0ptr; + bool dualRange = false; + if (f1ptr) {dualRange = true;} + tmin0 = fpValue(f0["Tmin"]); + tmax0 = fpValue(f0["Tmax"]); + tmin1 = tmax0; + tmax1 = tmin1 + 0.0001; + if (dualRange) { + tmin1 = fpValue((*f1ptr)["Tmin"]); + tmax1 = fpValue((*f1ptr)["Tmax"]); + } + + vector_fp c0, c1; + if (fabs(tmax0 - tmin1) < 0.01) { + tmin = tmin0; + tmid = tmax0; + tmax = tmax1; + getFloatArray(f0.child("floatArray"), c0, false); + if (dualRange) + getFloatArray(f1ptr->child("floatArray"), c1, false); + else + c1.resize(7,0.0); + } + else if (fabs(tmax1 - tmin0) < 0.01) { + tmin = tmin1; + tmid = tmax1; + tmax = tmax0; + getFloatArray(f1ptr->child("floatArray"), c0, false); + getFloatArray(f0.child("floatArray"), c1, false); + } + else { + throw CanteraError("installShomateThermo", + "non-continuous temperature ranges."); + } + array_fp c(15); + c[0] = tmid; + doublereal p0 = OneAtm; + copy(c0.begin(), c0.begin()+7, c.begin() + 1); + copy(c1.begin(), c1.begin()+7, c.begin() + 8); + sp.install(k, SHOMATE, c.begin(), tmin, tmax, p0); + } + + + + /** + * Install a constant-cp thermodynamic property + * parameterization for species k. + */ + static void installSimpleThermoFromXML(SpeciesThermo& sp, int k, + const XML_Node& f) { + doublereal tmin, tmax; + tmin = fpValue(f["Tmin"]); + tmax = fpValue(f["Tmax"]); + if (tmax == 0.0) tmax = 1.0e30; + + vector_fp c(4); + c[0] = getFloat(f, "t0", "-"); + c[1] = getFloat(f, "h0", "-"); + c[2] = getFloat(f, "s0", "-"); + c[3] = getFloat(f, "cp0", "-"); + doublereal p0 = OneAtm; + sp.install(k, SIMPLE, c.begin(), tmin, tmax, p0); + } + + + /** + * Install a species thermodynamic property parameterization + * for one species into a species thermo manager. + * @param k species number + * @param s XML node specifying species + * @param spthermo species thermo manager + */ + void SpeciesThermoFactory:: + installThermoForSpecies(int k, const XML_Node& s, + SpeciesThermo& spthermo) { + + const XML_Node& thermo = s.child("thermo"); + const vector& tp = thermo.children(); + int nc = tp.size(); + if (nc == 1) { + const XML_Node* f = tp[0]; + if (f->name() == "Shomate") { + installShomateThermoFromXML(spthermo, k, f, 0); + } + else if (f->name() == "const_cp") { + installSimpleThermoFromXML(spthermo, k, *f); + } + else if (f->name() == "NASA") { + installNasaThermoFromXML(spthermo, k, f, 0); + } + else { + UnknownSpeciesThermoModel("installSpecies", s["name"], f->name()); + } + } + else if (nc == 2) { + const XML_Node* f0 = tp[0]; + const XML_Node* f1 = tp[1]; + if (f0->name() == "NASA" && f1->name() == "NASA") { + installNasaThermoFromXML(spthermo, k, f0, f1); + } + else if (f0->name() == "Shomate" && f1->name() == "Shomate") { + installShomateThermoFromXML(spthermo, k, f0, f1); + } + else { + UnknownSpeciesThermoModel("installSpecies", s["name"], + f0->name() + " and " + f1->name()); + } + } + else { + UnknownSpeciesThermoModel("installSpecies", s["name"], + "multiple"); + } + } + } diff --git a/Cantera/src/SpeciesThermoFactory.h b/Cantera/src/SpeciesThermoFactory.h index 0d9eeb842..5c2ca32dc 100755 --- a/Cantera/src/SpeciesThermoFactory.h +++ b/Cantera/src/SpeciesThermoFactory.h @@ -73,6 +73,10 @@ namespace Cantera { virtual SpeciesThermo* newSpeciesThermo(vector nodes); virtual SpeciesThermo* newSpeciesThermoOpt(vector nodes); + + virtual void installThermoForSpecies(int k, const XML_Node& s, + SpeciesThermo& spthermo); + private: static SpeciesThermoFactory* __factory; SpeciesThermoFactory(){} diff --git a/Cantera/src/StoichSubstance.cpp b/Cantera/src/StoichSubstance.cpp index c7f7ca9dd..93f93c316 100644 --- a/Cantera/src/StoichSubstance.cpp +++ b/Cantera/src/StoichSubstance.cpp @@ -43,6 +43,13 @@ namespace Cantera { m_tlast = tnow; } } + + void StoichSubstance::setParametersFromXML(const XML_Node& eosdata) { + eosdata.require("model","Incompressible"); + doublereal rho = getFloat(eosdata, "density", "-"); + setDensity(rho); + } + } diff --git a/Cantera/src/StoichSubstance.h b/Cantera/src/StoichSubstance.h index 553678b51..72f4743ee 100644 --- a/Cantera/src/StoichSubstance.h +++ b/Cantera/src/StoichSubstance.h @@ -189,6 +189,7 @@ namespace Cantera { virtual void initThermo(); + virtual void setParametersFromXML(const XML_Node& eosdata); protected: diff --git a/Cantera/src/SurfPhase.cpp b/Cantera/src/SurfPhase.cpp index 5d50faaa6..3f94bf40b 100644 --- a/Cantera/src/SurfPhase.cpp +++ b/Cantera/src/SurfPhase.cpp @@ -92,6 +92,8 @@ namespace Cantera { return m_logn0 - m_logsize[k]; } + + /// The only parameter that can be set is the site density. void SurfPhase:: setParameters(int n, doublereal* c) { m_n0 = c[0]; @@ -153,6 +155,7 @@ namespace Cantera { // _updateThermo(true); //} + /** * Set the coverage fractions to a specified * state. This routine converts to concentrations @@ -236,10 +239,46 @@ namespace Cantera { } } + void SurfPhase:: + setParametersFromXML(const XML_Node& eosdata) { + eosdata.require("model","Surface"); + doublereal n = getFloat(eosdata, "site_density", "-"); + if (n <= 0.0) + throw CanteraError("SurfPhase::setParametersFromXML", + "missing or negative site density"); + m_n0 = n; + m_logn0 = log(m_n0); + } + + + void SurfPhase::setStateFromXML(const XML_Node& state) { + + if (state.hasChild("temperature")) { + double t = getFloat(state, "temperature", "temperature"); + setTemperature(t); + } + + if (state.hasChild("coverages")) { + string comp = getString(state,"coverages"); + setCoveragesByName(comp); + } + } EdgePhase::EdgePhase(doublereal n0) : SurfPhase(n0) { setNDim(1); } + void EdgePhase:: + setParametersFromXML(const XML_Node& eosdata) { + eosdata.require("model","Edge"); + doublereal n = getFloat(eosdata, "site_density", "-"); + if (n <= 0.0) + throw CanteraError("EdgePhase::setParametersFromXML", + "missing or negative site density"); + m_n0 = n; + m_logn0 = log(m_n0); + } + + } diff --git a/Cantera/src/SurfPhase.h b/Cantera/src/SurfPhase.h index 54c6ac6e8..3d345731c 100644 --- a/Cantera/src/SurfPhase.h +++ b/Cantera/src/SurfPhase.h @@ -19,6 +19,7 @@ #include "mix_defs.h" #include "ThermoPhase.h" + namespace Cantera { @@ -45,7 +46,9 @@ namespace Cantera { virtual doublereal standardConcentration(int k = 0) const; virtual doublereal logStandardConc(int k=0) const; virtual void setParameters(int n, doublereal* c); + virtual void setParametersFromXML(const XML_Node& eosdata); virtual void initThermo(); + virtual void setStateFromXML(const XML_Node& state); doublereal siteDensity(){ return m_n0; } void setPotentialEnergy(int k, doublereal pe); doublereal potentialEnergy(int k) {return m_pe[k];} diff --git a/Cantera/src/ThermoFactory.cpp b/Cantera/src/ThermoFactory.cpp index 551ad6010..d7315c097 100644 --- a/Cantera/src/ThermoFactory.cpp +++ b/Cantera/src/ThermoFactory.cpp @@ -89,4 +89,77 @@ namespace Cantera { } return th; } + + +// void setEOSParameters(const XML_Node& xmlphase, ThermoPhase* th) { + +// // if no thermo model is specified for the phase, simply +// // return +// if (!phase.hasChild("thermo")) return; + +// const XML_Node& eos = phase.child("thermo"); + +// // set the parameters for the particular equation of state type, +// // and +// if (eos["model"] == "Incompressible") { +// if (th->eosType() == cIncompressible) { +// doublereal rho = getFloat(eos, "density", "-"); +// th->setParameters(1, &rho); +// } +// else { +// eoserror = true; +// } +// } +// else if (eos["model"] == "StoichSubstance") { +// if (th->eosType() == cStoichSubstance) { +// doublereal rho = getFloat(eos, "density", "-"); +// th->setDensity(rho); +// } +// else { +// eoserror = true; +// } +// } +// else if (eos["model"] == "Surface") { +// if (th->eosType() == cSurf) { +// doublereal n = getFloat(eos, "site_density", "-"); +// if (n <= 0.0) +// throw CanteraError("importCTML", +// "missing or negative site density"); +// th->setParameters(1, &n); +// } +// else { +// eoserror = true; +// } +// } +// else if (eos["model"] == "Edge") { +// if (th->eosType() == cEdge) { +// doublereal n = getFloat(eos, "site_density", "-"); +// if (n <= 0.0) +// throw CanteraError("importCTML", +// "missing or negative site density"); +// th->setParameters(1, &n); +// } +// else { +// eoserror = true; +// } +// } +// #ifdef INCL_PURE_FLUIDS +// else if (eos["model"] == "PureFluid") { +// if (th->eosType() == cPureFluid) { +// subflag = atoi(eos["fluid_type"].c_str()); +// if (subflag < 0) +// throw CanteraError("importCTML", +// "missing fluid type flag"); +// } +// else { +// eoserror = true; +// } +// } +// #endif +// if (eoserror) { +// string msg = "Wrong equation of state type for phase "+phase["id"]+"\n"; +// msg += eos["model"]+" is not consistent with eos type "+int2str(th->eosType()); +// throw CanteraError("importCTML",msg); +// } + } diff --git a/Cantera/src/ThermoPhase.cpp b/Cantera/src/ThermoPhase.cpp index 0292ff50d..a7618d0df 100644 --- a/Cantera/src/ThermoPhase.cpp +++ b/Cantera/src/ThermoPhase.cpp @@ -197,6 +197,35 @@ namespace Cantera { } } + + + /** + * Set the thermodynamic state. + */ + void ThermoPhase::setStateFromXML(const XML_Node& state) { + + string comp = getString(state,"moleFractions"); + if (comp != "") + setMoleFractionsByName(comp); + else { + comp = getString(state,"massFractions"); + if (comp != "") + setMassFractionsByName(comp); + } + if (state.hasChild("temperature")) { + double t = getFloat(state, "temperature", "temperature"); + setTemperature(t); + } + if (state.hasChild("pressure")) { + double p = getFloat(state, "pressure", "pressure"); + setPressure(p); + } + if (state.hasChild("density")) { + double rho = getFloat(state, "density", "density"); + setDensity(rho); + } + } + } diff --git a/Cantera/src/ThermoPhase.h b/Cantera/src/ThermoPhase.h index 2848b4690..9da24ea51 100755 --- a/Cantera/src/ThermoPhase.h +++ b/Cantera/src/ThermoPhase.h @@ -31,15 +31,15 @@ namespace Cantera { */ /** - * A phase with thermodynamic properties. - * Extends class Phase by adding methods that compute - * thermodynamic properties. + * A phase with thermodynamic properties. Extends class Phase by + * adding methods that compute thermodynamic properties that + * require knowledge of the equation of state. * * Class ThermoPhase is the base class for the family of classes * that represent phases of matter with particular equations of * state. Instances of subclasses of ThermoPhase should be created * using the factory class ThermoFactory, not by calling the - * constructor directly. + * constructor directly. * * To implement a new equation of state, derive a class from * ThermoPhase and overload the virtual methods in @@ -69,60 +69,6 @@ namespace Cantera { * @{ */ - - /** - * @internal - * Index number. This method can be used to identify the - * location of a phase object in a list, and is used by the - * interface library (clib) routines for this purpose. - */ - int index() { return m_index; } - - - /** - * @internal Set the index number. The Cantera interface - * library uses this method to set the index number to the - * location of the pointer to this object in the pointer array - * it maintains. Using this method for any other purpose will - * lead to unpredictable results if used in conjunction with - * the interface library. - */ - void setIndex(int m) { m_index = m; } - - - /// used to access data needed to construct transport manager - /// later. - void saveSpeciesData(const XML_Node* data) { - m_speciesData = data; - } - - const XML_Node* speciesData() { - if (m_speciesData) - return m_speciesData; - else { - throw CanteraError("ThermoPhase::speciesData", - "m_speciesData is NULL"); - return 0; - } - } - - - /** - * @internal Initialize. This method is provided to allow - * subclasses to perform any initialization required after all - * species have been added. For example, it might be used to - * resize internal work arrays that must have an entry for - * each species. The base class implementation does nothing, - * and subclasses that do not require initialization do not - * need to overload this method. When importing a CTML phase - * description, this method is called just prior to returning - * from function importPhase. - * - * @see importCTML.cpp - */ - virtual void initThermo() {} - - /** * Equation of state type flag. The base class returns * zero. Subclasses should define this to return a unique @@ -132,57 +78,43 @@ namespace Cantera { virtual int eosType() const { return 0; } - /** * @} * @name Molar Thermodynamic Properties * @{ */ - - /** - * Molar enthalpy. Units: J/kmol. - */ + /// Molar enthalpy. Units: J/kmol. virtual doublereal enthalpy_mole() const { return err("enthalpy_mole"); } - /** - * Molar internal energy. Units: J/kmol. - */ + /// Molar internal energy. Units: J/kmol. virtual doublereal intEnergy_mole() const { return err("intEnergy_mole"); } - /** - * Molar entropy. Units: J/kmol/K. - */ + /// Molar entropy. Units: J/kmol/K. virtual doublereal entropy_mole() const { return err("entropy_mole"); } - /** - * Molar Gibbs function. Units: J/kmol. - */ + /// Molar Gibbs function. Units: J/kmol. virtual doublereal gibbs_mole() const { return err("gibbs_mole"); } - /** - * Molar heat capacity at constant pressure. Units: J/kmol/K. - */ + /// Molar heat capacity at constant pressure. Units: J/kmol/K. virtual doublereal cp_mole() const { return err("cp_mole"); } - /** - * Molar heat capacity at constant volume. Units: J/kmol/K. - */ + /// Molar heat capacity at constant volume. Units: J/kmol/K. virtual doublereal cv_mole() const { return err("cv_mole"); } @@ -246,11 +178,18 @@ namespace Cantera { return err("potentialEnergy"); } + /** + * Set the electric potential of this phase (V). + * This is used by classes InterfaceKinetics and EdgeKinetics to + * compute the rates of charge-transfer reactions, and in computing + * the electrochemical potentials of the species. + */ void setElectricPotential(doublereal v) { m_phi = v; } - doublereal electricPotential() { return m_phi; } + /// The electric potential of this phase (V). + doublereal electricPotential() const { return m_phi; } /** @@ -271,8 +210,8 @@ namespace Cantera { * C^0_k, \f$ where \f$ C^0_k \f$ is a standard concentration * defined below. These generalized concentrations are used * by kinetics manager classes to compute the forward and - * reverse rates of elementary reactions. - * + * reverse rates of elementary reactions. + * * @param c Array of generalized concentrations. The * units depend upon the implementation of the * reaction rate expressions within the phase. @@ -298,7 +237,6 @@ namespace Cantera { } /** - * * Returns the natural logarithm of the standard * concentration of the kth species */ @@ -315,10 +253,10 @@ namespace Cantera { } /** - * Returns the units of the standard and general concentrations - * Note they have the same units, as their divisor is - * defined to be equal to the activity of the kth species - * in the solution, which is unitless. + * Returns the units of the standard and generalized + * concentrations Note they have the same units, as their + * ratio is defined to be equal to the activity of the kth + * species in the solution, which is unitless. * * This routine is used in print out applications where the * units are needed. Usually, MKS units are assumed throughout @@ -350,6 +288,19 @@ namespace Cantera { err("getChemPotentials_RT"); } + /** + * Get the species electrochemical potentials. Units: J/kmol. + * This method adds a term \f$ Fz_k \phi_k$ to the + * to each chemical potential. + */ + void getElectrochemPotentials(doublereal* mu) const { + getChemPotentials(mu); + double ve = Faraday * electricPotential(); + for (int k = 0; k < m_kk; k++) { + mu[k] += ve*charge(k); + } + } + //@} /// @name Partial Molar Properties //@{ @@ -568,7 +519,9 @@ namespace Cantera { * @param c array of \i n coefficients */ virtual void setParameters(int n, doublereal* c) {} - + virtual void setParametersFromXML(const XML_Node& eosdata) {} + virtual void setStateFromXML(const XML_Node& state); + virtual doublereal isothermalCompressibility() { err("isothermalCompressibility"); return -1.0; } @@ -579,6 +532,7 @@ namespace Cantera { //--------------------------------------------------------- /// @name Critical state properties. + /// These methods are only implemented by some subclasses. //@{ @@ -599,6 +553,10 @@ namespace Cantera { //@} + /// @name Saturation properties. + /// These methods are only implemented by subclasses that + /// implement full liquid-vapor equations of state. + /// virtual doublereal satTemperature(doublereal p) const { err("satTemperature"); return -1.0; } @@ -619,6 +577,49 @@ namespace Cantera { err("setState_sat"); } + //@} + + + /** + * Returns the reference pressure in Pa. This function is a wrapper + * that calls the species thermo refPressure function. + */ + doublereal refPressure() const { + return m_spthermo->refPressure(); + } + + doublereal minTemp(int k = -1) { + return m_spthermo->minTemp(k); + } + + doublereal maxTemp(int k = -1) { + return m_spthermo->maxTemp(k); + } + + + /// The following methods are used in the process of constructing + /// the phase and setting its parameters from a specification in an + /// input file. They are not normally used in application programs. + /// To see how they are used, see files importCTML.cpp and + /// ThermoFactory.cpp. + + + /// used to access data needed to construct transport manager + /// later. + void saveSpeciesData(const XML_Node* data) { + m_speciesData = data; + } + + const XML_Node* speciesData() { + if (m_speciesData) + return m_speciesData; + else { + throw CanteraError("ThermoPhase::speciesData", + "m_speciesData is NULL"); + return 0; + } + } + /** * @internal Install a species thermodynamic property @@ -638,22 +639,48 @@ namespace Cantera { */ SpeciesThermo& speciesThermo() { return *m_spthermo; } - /** - * Returns the reference pressure in Pa. This function is a wrapper - * that calls the species thermo refPressure function. - */ - doublereal refPressure() const { - return m_spthermo->refPressure(); - } - doublereal minTemp(int k = -1) { - return m_spthermo->minTemp(k); - } + /** + * @internal Initialize. This method is provided to allow + * subclasses to perform any initialization required after all + * species have been added. For example, it might be used to + * resize internal work arrays that must have an entry for + * each species. The base class implementation does nothing, + * and subclasses that do not require initialization do not + * need to overload this method. When importing a CTML phase + * description, this method is called just prior to returning + * from function importPhase. + * + * @see importCTML.cpp + */ + virtual void initThermo() {} + + + + // The following methods are used by the clib interface + // library, and should not be used by application programs. + + /** + * @internal + * Index number. This method can be used to identify the + * location of a phase object in a list, and is used by the + * interface library (clib) routines for this purpose. + */ + int index() { return m_index; } + + + /** + * @internal Set the index number. The Cantera interface + * library uses this method to set the index number to the + * location of the pointer to this object in the pointer array + * it maintains. Using this method for any other purpose will + * lead to unpredictable results if used in conjunction with + * the interface library. + */ + void setIndex(int m) { m_index = m; } + + - doublereal maxTemp(int k = -1) { - return m_spthermo->maxTemp(k); - } - protected: diff --git a/Cantera/src/importCTML.cpp b/Cantera/src/importCTML.cpp index 42e243cee..1e62c234e 100755 --- a/Cantera/src/importCTML.cpp +++ b/Cantera/src/importCTML.cpp @@ -42,8 +42,6 @@ using namespace ctml; -//#include - // these are all used to check for duplicate reactions vector< map > _reactiondata; @@ -55,7 +53,7 @@ vector _rev; namespace Cantera { /* - * First we define a coule of typedef's which will + * First we define a couple of typedefs that will * be used throught this file */ typedef const vector nodeset_t; @@ -64,6 +62,8 @@ namespace Cantera { const doublereal DefaultPref = 1.01325e5; // one atm + /// split a string at a '#' sign. Used to separate a file name + /// from an id string. static void split(const string& src, string& file, string& id) { string::size_type ipound = src.find('#'); if (ipound != string::npos) { @@ -76,6 +76,7 @@ namespace Cantera { } } + /** * This routine will locate an XML node in either the input * XML tree or in another input file specified by the file @@ -100,13 +101,16 @@ namespace Cantera { XML_Node *db, *doc; split(file_ID, fname, idstr); if (fname == "") { - if (!root) throw CanteraError("get_XML_Node","no file name given. file_ID = "+file_ID); + if (!root) throw CanteraError("get_XML_Node", + "no file name given. file_ID = "+file_ID); db = root->findID(idstr, 3); } else { doc = get_XML_File(fname); - if (!doc) throw CanteraError("get_XML_Node", "get_XML_File failed trying to open "+fname); + if (!doc) throw CanteraError("get_XML_Node", + "get_XML_File failed trying to open "+fname); db = doc->findID(idstr, 3); - if (!db) throw CanteraError("get_XML_Node", "id tag "+idstr+" not found."); + if (!db) throw CanteraError("get_XML_Node", + "id tag "+idstr+" not found."); } return db; } @@ -152,161 +156,6 @@ namespace Cantera { return db; } - /** - * Install a NASA polynomial thermodynamic property - * parameterization for species k. - */ - static void installNasaThermo(SpeciesThermo& sp, int k, const XML_Node* f0ptr, - const XML_Node* f1ptr) { - doublereal tmin0, tmax0, tmin1, tmax1, tmin, tmid, tmax; - - const XML_Node& f0 = *f0ptr; - bool dualRange = false; - if (f1ptr) {dualRange = true;} - tmin0 = fpValue(f0["Tmin"]); - tmax0 = fpValue(f0["Tmax"]); - tmin1 = tmax0; - tmax1 = tmin1 + 0.0001; - if (dualRange) { - tmin1 = fpValue((*f1ptr)["Tmin"]); - tmax1 = fpValue((*f1ptr)["Tmax"]); - } - - vector_fp c0, c1; - if (fabs(tmax0 - tmin1) < 0.01) { - tmin = tmin0; - tmid = tmax0; - tmax = tmax1; - getFloatArray(f0.child("floatArray"), c0, false); - if (dualRange) - getFloatArray(f1ptr->child("floatArray"), c1, false); - else - c1.resize(7,0.0); - } - else if (fabs(tmax1 - tmin0) < 0.01) { - tmin = tmin1; - tmid = tmax1; - tmax = tmax0; - getFloatArray(f1ptr->child("floatArray"), c0, false); - getFloatArray(f0.child("floatArray"), c1, false); - } - else { - throw CanteraError("installNasaThermo", - "non-continuous temperature ranges."); - } - array_fp c(15); - c[0] = tmid; - doublereal p0 = OneAtm; - c[1] = c0[5]; - c[2] = c0[6]; - copy(c0.begin(), c0.begin()+5, c.begin() + 3); - c[8] = c1[5]; - c[9] = c1[6]; - copy(c1.begin(), c1.begin()+5, c.begin() + 10); - sp.install(k, NASA, c.begin(), tmin, tmax, p0); - } - - - - /** - * Install a NASA polynomial thermodynamic property - * parameterization for species k. - */ - static void installShomateThermo(SpeciesThermo& sp, int k, const XML_Node* f0ptr, - const XML_Node* f1ptr) { - doublereal tmin0, tmax0, tmin1, tmax1, tmin, tmid, tmax; - - const XML_Node& f0 = *f0ptr; - bool dualRange = false; - if (f1ptr) {dualRange = true;} - tmin0 = fpValue(f0["Tmin"]); - tmax0 = fpValue(f0["Tmax"]); - tmin1 = tmax0; - tmax1 = tmin1 + 0.0001; - if (dualRange) { - tmin1 = fpValue((*f1ptr)["Tmin"]); - tmax1 = fpValue((*f1ptr)["Tmax"]); - } - - vector_fp c0, c1; - if (fabs(tmax0 - tmin1) < 0.01) { - tmin = tmin0; - tmid = tmax0; - tmax = tmax1; - getFloatArray(f0.child("floatArray"), c0, false); - if (dualRange) - getFloatArray(f1ptr->child("floatArray"), c1, false); - else - c1.resize(7,0.0); - } - else if (fabs(tmax1 - tmin0) < 0.01) { - tmin = tmin1; - tmid = tmax1; - tmax = tmax0; - getFloatArray(f1ptr->child("floatArray"), c0, false); - getFloatArray(f0.child("floatArray"), c1, false); - } - else { - throw CanteraError("installShomateThermo", - "non-continuous temperature ranges."); - } - array_fp c(15); - c[0] = tmid; - doublereal p0 = OneAtm; - copy(c0.begin(), c0.begin()+7, c.begin() + 1); - copy(c1.begin(), c1.begin()+7, c.begin() + 8); - sp.install(k, SHOMATE, c.begin(), tmin, tmax, p0); - } - - -// /** -// * Install a Shomate polynomial thermodynamic property -// * parameterization for species k. -// */ -// static void installShomateThermo(SpeciesThermo& sp, int k, const XML_Node& f) { -// doublereal tmin, tmid, tmax; -// tmin = fpValue(f["Tmin"]); -// tmid = fpValue(f["Tmid"]); -// tmax = fpValue(f["Tmax"]); - -// vector fa; -// f.getChildren("floatArray",fa); -// vector_fp c0, c1; -// getFloatArray(*fa[0], c0, false); -// getFloatArray(*fa[1], c1, false); -// array_fp c(15); -// c[0] = tmid; -// doublereal p0 = OneAtm; -// if ((*fa[0])["title"] == "low") { -// copy(c0.begin(), c0.end(), c.begin() + 1); -// copy(c1.begin(), c1.end(), c.begin() + 8); -// } -// else { -// copy(c1.begin(), c1.end(), c.begin() + 1); -// copy(c0.begin(), c0.end(), c.begin() + 8); -// } -// sp.install(k, SHOMATE, c.begin(), tmin, tmax, p0); -// } - - - /** - * Install a constant-cp thermodynamic property - * parameterization for species k. - */ - static void installSimpleThermo(SpeciesThermo& sp, int k, const XML_Node& f) { - doublereal tmin, tmax; - tmin = fpValue(f["Tmin"]); - tmax = fpValue(f["Tmax"]); - if (tmax == 0.0) tmax = 1.0e30; - - vector_fp c(4); - c[0] = getFloat(f, "t0", "-"); - c[1] = getFloat(f, "h0", "-"); - c[2] = getFloat(f, "s0", "-"); - c[3] = getFloat(f, "cp0", "-"); - doublereal p0 = OneAtm; - sp.install(k, SIMPLE, c.begin(), tmin, tmax, p0); - } /** * Install a species into a ThermoPhase object, which defines @@ -329,7 +178,7 @@ namespace Cantera { * an "UnknownSpeciesThermoModel" exception being thrown. */ bool installSpecies(int k, const XML_Node& s, thermo_t& p, - SpeciesThermo& spthermo, int rule) { + SpeciesThermo& spthermo, int rule, SpeciesThermoFactory* factory) { // get the composition of the species const XML_Node& a = s.child("atomArray"); @@ -342,77 +191,46 @@ namespace Cantera { // otherwise, throw an exception map::const_iterator _b = comp.begin(); for (; _b != comp.end(); ++_b) { - if (p.elementIndex(_b->first) < 0) { - if (rule == 0) { - throw CanteraError("installSpecies", - "Species " + s["name"] + - " contains undeclared element " + _b->first); - } - else - return false; - } + if (p.elementIndex(_b->first) < 0) { + if (rule == 0) { + throw CanteraError("installSpecies", + "Species " + s["name"] + + " contains undeclared element " + _b->first); + } + else + return false; + } } + // construct a vector of atom numbers for each + // element in phase p. Elements not declared in the + // species (i.e., not in map comp) will have zero + // entries in the vector. int m, nel = p.nElements(); vector_fp ecomp(nel, 0.0); for (m = 0; m < nel; m++) { - ecomp[m] = atoi(comp[p.elementName(m)].c_str()); + ecomp[m] = atoi(comp[p.elementName(m)].c_str()); } - /* - * Define a map and get all of the floats in the - * current XML species block - */ + + // get the species charge, if any. Note that the charge need + // not be explicitly specified if special element 'E' + // (electron) is one of the elements. doublereal chrg = 0.0; if (s.hasChild("charge")) chrg = getFloat(s, "charge"); + + // get the species size, if any. (This is used by surface + // phases to represent how many sites a species occupies.) doublereal sz = 1.0; if (s.hasChild("size")) sz = getFloat(s, "size"); + // add the species to phase p. p.addUniqueSpecies(s["name"], ecomp.begin(), chrg, sz); - // get thermo. We currently only support single-range Shomate - // and const_cp, and dual-range NASA - if (!s.hasChild("thermo")) { - throw - UnknownSpeciesThermoModel("installSpecies", s["name"], "missing"); - - } - const XML_Node& thermo = s.child("thermo"); - const vector& tp = thermo.children(); - int nc = tp.size(); - if (nc == 1) { - const XML_Node* f = tp[0]; - if (f->name() == "Shomate") { - installShomateThermo(spthermo, k, f, 0); - } - else if (f->name() == "const_cp") { - installSimpleThermo(spthermo, k, *f); - } - else if (f->name() == "NASA") { - installNasaThermo(spthermo, k, f, 0); - } - else { - UnknownSpeciesThermoModel("installSpecies", s["name"], f->name()); - } - } - else if (nc == 2) { - const XML_Node* f0 = tp[0]; - const XML_Node* f1 = tp[1]; - if (f0->name() == "NASA" && f1->name() == "NASA") { - installNasaThermo(spthermo, k, f0, f1); - } - else if (f0->name() == "Shomate" && f1->name() == "Shomate") { - installShomateThermo(spthermo, k, f0, f1); - } - else { - UnknownSpeciesThermoModel("installSpecies", s["name"], - f0->name() + " and " + f1->name()); - } - } - else { - UnknownSpeciesThermoModel("installSpecies", s["name"], - "multiple"); - } + // install the thermo parameterization for this species into + // the species thermo manager for phase p. + factory->installThermoForSpecies(k, s, spthermo); + return true; } @@ -815,15 +633,15 @@ namespace Cantera { /** - * Create a new ThermoPhase object and initializes it according - * to the XML tree database. - * This routine first looks up the identity of the model for the - * solution thermodynamics in the model attribute of the thermo - * child of the xml phase node. Then, it does a string lookup on - * the model to figure out what ThermoPhase derived class is - * assigned. It mallocs a new instance of that class, and then - * calls importPhase() to populate that class with the correct - * parameters from the XML tree. + * Create a new ThermoPhase object and initializes it according to + * the XML tree database. This routine first looks up the + * identity of the model for the solution thermodynamics in the + * model attribute of the thermo child of the xml phase + * node. Then, it does a string lookup on the model to figure out + * what ThermoPhase derived class is assigned. It creates a new + * instance of that class, and then calls importPhase() to + * populate that class with the correct parameters from the XML + * tree. */ ThermoPhase* newPhase(XML_Node& xmlphase) { const XML_Node& th = xmlphase.child("thermo"); @@ -843,44 +661,44 @@ namespace Cantera { return 0; } - /** - * Set the thermodynamic state. - */ - static void setState(const XML_Node& phase, ThermoPhase* th) { - if (!phase.hasChild("state")) return; - const XML_Node state = phase.child("state"); - doublereal t, p, rho; - string comp = getString(state,"moleFractions"); - if (comp != "") - th->setMoleFractionsByName(comp); - else { - comp = getString(state,"massFractions"); - if (comp != "") - th->setMassFractionsByName(comp); - } - if (state.hasChild("temperature")) { - t = getFloat(state, "temperature", "temperature"); - th->setTemperature(t); - } - if (state.hasChild("pressure")) { - p = getFloat(state, "pressure", "pressure"); - th->setPressure(p); - } - if (state.hasChild("density")) { - rho = getFloat(state, "density", "density"); - th->setDensity(rho); - } - if (th->eosType() == cSurf && state.hasChild("coverages")) { - comp = getString(state,"coverages"); - SurfPhase* s = (SurfPhase*)th; - s->setCoveragesByName(comp); - } - if (th->eosType() == cEdge && state.hasChild("coverages")) { - comp = getString(state,"coverages"); - EdgePhase* s = (EdgePhase*)th; - s->setCoveragesByName(comp); - } - } +// /** +// * Set the thermodynamic state. +// */ +// static void setState(const XML_Node& phase, ThermoPhase* th) { +// if (!phase.hasChild("state")) return; +// const XML_Node state = phase.child("state"); +// doublereal t, p, rho; +// string comp = getString(state,"moleFractions"); +// if (comp != "") +// th->setMoleFractionsByName(comp); +// else { +// comp = getString(state,"massFractions"); +// if (comp != "") +// th->setMassFractionsByName(comp); +// } +// if (state.hasChild("temperature")) { +// t = getFloat(state, "temperature", "temperature"); +// th->setTemperature(t); +// } +// if (state.hasChild("pressure")) { +// p = getFloat(state, "pressure", "pressure"); +// th->setPressure(p); +// } +// if (state.hasChild("density")) { +// rho = getFloat(state, "density", "density"); +// th->setDensity(rho); +// } +// if (th->eosType() == cSurf && state.hasChild("coverages")) { +// comp = getString(state,"coverages"); +// SurfPhase* s = (SurfPhase*)th; +// s->setCoveragesByName(comp); +// } +// if (th->eosType() == cEdge && state.hasChild("coverages")) { +// comp = getString(state,"coverages"); +// EdgePhase* s = (EdgePhase*)th; +// s->setCoveragesByName(comp); +// } +// } /** * Import a phase specification. @@ -907,15 +725,23 @@ namespace Cantera { * here, especially for those objects which are * part of the Cantera Kernel. */ - bool importPhase(XML_Node& phase, ThermoPhase* th) { - - int subflag = -1; + bool importPhase(XML_Node& phase, ThermoPhase* th, + SpeciesThermoFactory* spfactory) { + // Check the the supplied XML node in fact represents a + // phase. if (phase.name() != "phase") throw CanteraError("importPhase", "Current const XML_Node is not a phase element."); - th->setID(phase.id()); // set the phase id + // if no species thermo factory was supplied, + // use the default one. + if (!spfactory) + spfactory = SpeciesThermoFactory::factory(); + + // set the id attribute of the phase to the 'id' attribute + // in the XML tree. + th->setID(phase.id()); // Number of spatial dimensions. Defaults to 3 (bulk phase) if (phase.hasAttrib("dim")) { @@ -928,116 +754,30 @@ namespace Cantera { else th->setNDim(3); // default - /** - * Equation of State: We initialize the ThermoPhase objects that - * we know about here, with additional parameters obtained from - * the xml tree. EOS's that we don't know about don't create an - * error condition. - */ - bool eoserror = false; + + + // set equation of state parameters. The parameters are + // specific to each subclass of ThermoPhase, so this is done + // by method setParametersFromXML in each subclass. if (phase.hasChild("thermo")) { const XML_Node& eos = phase.child("thermo"); - if (eos["model"] == "Incompressible") { - if (th->eosType() == cIncompressible) { - doublereal rho = getFloat(eos, "density", "-"); - //doublereal rho = d["density"]; - th->setParameters(1, &rho); - } - else { - eoserror = true; - } - } - else if (eos["model"] == "StoichSubstance") { - if (th->eosType() == cStoichSubstance) { - doublereal rho = getFloat(eos, "density", "-"); - th->setDensity(rho); - } - else { - eoserror = true; - } - } - else if (eos["model"] == "Surface") { - if (th->eosType() == cSurf) { - doublereal n = getFloat(eos, "site_density", "-"); - if (n <= 0.0) - throw CanteraError("importCTML", - "missing or negative site density"); - th->setParameters(1, &n); - } - else { - eoserror = true; - } - } - else if (eos["model"] == "Edge") { - if (th->eosType() == cEdge) { - doublereal n = getFloat(eos, "site_density", "-"); - if (n <= 0.0) - throw CanteraError("importCTML", - "missing or negative site density"); - th->setParameters(1, &n); - } - else { - eoserror = true; - } - } -#ifdef INCL_PURE_FLUIDS - else if (eos["model"] == "PureFluid") { - if (th->eosType() == cPureFluid) { - subflag = atoi(eos["fluid_type"].c_str()); - if (subflag < 0) - throw CanteraError("importCTML", - "missing fluid type flag"); - } - else { - eoserror = true; - } - } -#endif - if (eoserror) { - string msg = "Wrong equation of state type for phase "+phase["id"]+"\n"; - msg += eos["model"]+" is not consistent with eos type "+int2str(th->eosType()); - throw CanteraError("importCTML",msg); - } - } - - - /************************************************* - * Add elements. - ************************************************/ - - - // get the declared element names - XML_Node& elements = phase.child("elementArray"); - vector enames; - getStringArray(elements, enames); - - // // element database defaults to elements.xml - string element_database = "elements.xml"; - if (elements.hasAttrib("datasrc")) - element_database = elements["datasrc"]; - XML_Node* doc = get_XML_File(element_database); - XML_Node* dbe = &doc->child("ctml/elementData"); - - int nel = enames.size(); - int i; - string enm; - for (i = 0; i < nel; i++) { - XML_Node* e = dbe->findByAttr("name",enames[i]); - if (e) { - th->addUniqueElement(*e); - } - else { - throw CanteraError("importPhase","no data for element " - +enames[i]); - } + th->setParametersFromXML(eos); } /*************************************************************** - * Add the species. First get the speciesArray element, then - * the species database. + * Add the elements. ***************************************************************/ + th->addElementsFromXML(phase); + + /*************************************************************** + * Add the species. + * + * Species definitions may be imported from multiple + * sources. For each one, a speciesArray element must be + * present. + ***************************************************************/ XML_Node* db = 0; vector sparrays; phase.getChildren("speciesArray", sparrays); @@ -1045,10 +785,16 @@ namespace Cantera { vector dbases; vector_int sprule(nspa,0); + // loop over the speciesArray elements for (jsp = 0; jsp < nspa; jsp++) { const XML_Node& species = *sparrays[jsp]; + // If the speciesArray element has a child element + // + // then set sprule[jsp] to 1, so + // that any species with an undeclared element will be + // quietly skipped when importing species. if (species.hasChild("skip")) { const XML_Node& sk = species.child("skip"); string eskip = sk["element"]; @@ -1056,42 +802,54 @@ namespace Cantera { sprule[jsp] = 1; } } + string fname, idstr; - + + // get a pointer to the node containing the species + // definitions for the species declared in this + // speciesArray element. This may be in the local file + // containing the phase element, or may be in another + // file. db = get_XML_Node(species["datasrc"], &phase.root()); - //db = find_XML(species["datasrc"], &phase.root(), species["idRef"], - // "","speciesData"); + + // add this node to the list of species database nodes. dbases.push_back(db); } - /******************************************************* - * Set the species thermo manager. - * Function 'newSpeciesThermoMgr' looks at the species - * in the database to see what thermodynamic property - * parameterizations are used, and selects a class - * that can handle the parameterizations found. - ******************************************************/ - + // if the phase has a species thermo manager already installed, + // delete it since we are adding new species. delete &th->speciesThermo(); + + // create a new species thermo manager. Function + // 'newSpeciesThermoMgr' looks at the species in the database + // to see what thermodynamic property parameterizations are + // used, and selects a class that can handle the + // parameterizations found. SpeciesThermo* spth = newSpeciesThermoMgr(dbases); + + // install it in the phase object th->setSpeciesThermo(spth); SpeciesThermo& spthermo = th->speciesThermo(); + // used to check that each species is declared only once map declared; - int k = 0; + + int i, k = 0; + + // loop over the species arrays for (jsp = 0; jsp < nspa; jsp++) { const XML_Node& species = *sparrays[jsp]; db = dbases[jsp]; - /* - * Get the array of species name strings. - */ + // Get the array of species name strings. vector spnames; getStringArray(species, spnames); int nsp = spnames.size(); + // if 'all' is specified, then add all species + // defined in this database to the phase if (nsp == 1 && spnames[0] == "all") { vector allsp; db->getChildren("species",allsp); @@ -1113,12 +871,11 @@ namespace Cantera { } declared[name] = true; - /* - * Find the species in the database by name. - */ + // Find the species in the database by name. XML_Node* s = db->findByAttr("name",spnames[i]); if (s) { - if (installSpecies(k, *s, *th, spthermo, sprule[jsp])) + if (installSpecies(k, *s, *th, spthermo, sprule[jsp], + spfactory)) ++k; } else { @@ -1127,20 +884,26 @@ namespace Cantera { } } } + + // done adding species. th->freezeSpecies(); + + // perform any required subclass-specific initialization. th->initThermo(); th->saveSpeciesData(db); - if (th->eosType() == cPureFluid) { - doublereal dsub = doublereal(subflag); - th->setParameters(1, &dsub); + // set the state of the phase from the XML specification + if (phase.hasChild("state")) { + XML_Node& state = phase.child("state"); + th->setStateFromXML(state); } - setState(phase, th); return true; } + + /** * This function returns true if two reactions are duplicates of * one another, and false otherwise. The input arguments are two diff --git a/Cantera/src/importCTML.h b/Cantera/src/importCTML.h index 802378160..f8e9011b5 100755 --- a/Cantera/src/importCTML.h +++ b/Cantera/src/importCTML.h @@ -22,6 +22,7 @@ using namespace std; namespace Cantera { class Kinetics; + class SpeciesThermoFactory; //class ThermoPhase; class XML_Node; @@ -72,10 +73,11 @@ namespace Cantera { XML_Node* get_XML_NameID(const string& nameTarget, const string& file_ID, XML_Node* root); - bool installSpecies(int k, const XML_Node& s, thermo_t& p, - SpeciesThermo& spthermo, int rule); + //bool installSpecies(int k, const XML_Node& s, thermo_t& p, + // SpeciesThermo& spthermo, int rule); - bool importPhase(XML_Node& phase, ThermoPhase* th); + bool importPhase(XML_Node& phase, ThermoPhase* th, + SpeciesThermoFactory* spfactory = 0); /** * This function returns true if two reactions are duplicates of diff --git a/Cantera/src/xml.cpp b/Cantera/src/xml.cpp index af98c7f53..05a57d430 100755 --- a/Cantera/src/xml.cpp +++ b/Cantera/src/xml.cpp @@ -696,82 +696,87 @@ namespace Cantera { } } - //const XML_Node* XML_Node::getRef() const { - // if (!hasAttrib("idRef")) return this; - // XML_Node& node = *this; - // return find_XML(node["src"], &root(), node["idRef"]); - //} -#ifdef FIND_XML - /* - * Find a particular XML element by a fairly complicated hierarchal - * search objective. - * - * HKM -Note: Right now this routine contains a memory leak. - * A "new" operation is conditionally carried out and - * the pointer may or may not be returned to the calling - * program. Therefore, it can't be deleted in the - * calling program. This - * eventually needs to be fixed by extracting the xml - * malloc and build operation from the search operation. - */ - XML_Node* find_XML(string src, XML_Node* root, string id, string loc, - string name) { - string file, id2; - split(src, file, id2); - src = file; - if (id2 != "") id = id2; - - XML_Node *doc = 0, *r = 0; - if (src != "") { - doc = new XML_Node("doc"); - string spath = findInputFile(src); - ifstream fin(spath.c_str()); - if (!fin) - throw CanteraError("find_XML","could not open file "+src+ - " for input."); - doc->build(fin); - root = 0; - } - else if (root) { - doc = root; - } - else { - throw CanteraError("find_XML", - "either root or src must be specified."); - } - - try { - if (id != "") - r = doc->findID(id); - else if (loc != "") - r = &doc->child(loc); - else if (name != "") - r = doc->findByName(name); - if (!r) { - string opt = " src="+src+", loc="+loc+", id=" - +id+", name="+name; - throw CanteraError("find_XML", "XML element with "+opt+ - " not found."); - } - return r; - } - catch (CanteraError) { - - // root was used, but element was not found. Try src. - if (root && src != "") { - return find_XML(src, 0, id, loc, name); - } - else { - string opt = " src="+src+", loc="+loc+", id=" - +id+", name="+name; - throw CanteraError("find_XML", "XML element with "+opt+ - " not found."); - return 0; - } + void XML_Node::require(string a, string v) const { + if (hasAttrib(a)) { + if (attrib(a) == v) return; } + string msg="XML_Node "+name()+" is required to have the value " + "\""+v+"\", but instead is \""+attrib(a); + throw CanteraError("XML_Node::require",msg); } -#endif + + +// #ifdef FIND_XML +// /* +// * Find a particular XML element by a fairly complicated hierarchal +// * search objective. +// * +// * HKM -Note: Right now this routine contains a memory leak. +// * A "new" operation is conditionally carried out and +// * the pointer may or may not be returned to the calling +// * program. Therefore, it can't be deleted in the +// * calling program. This +// * eventually needs to be fixed by extracting the xml +// * malloc and build operation from the search operation. +// */ +// XML_Node* find_XML(string src, XML_Node* root, string id, string loc, +// string name) { +// string file, id2; +// split(src, file, id2); +// src = file; +// if (id2 != "") id = id2; + +// XML_Node *doc = 0, *r = 0; +// if (src != "") { +// doc = new XML_Node("doc"); +// string spath = findInputFile(src); +// ifstream fin(spath.c_str()); +// if (!fin) +// throw CanteraError("find_XML","could not open file "+src+ +// " for input."); +// doc->build(fin); +// root = 0; +// } +// else if (root) { +// doc = root; +// } +// else { +// throw CanteraError("find_XML", +// "either root or src must be specified."); +// } + +// try { +// if (id != "") +// r = doc->findID(id); +// else if (loc != "") +// r = &doc->child(loc); +// else if (name != "") +// r = doc->findByName(name); +// if (!r) { +// string opt = " src="+src+", loc="+loc+", id=" +// +id+", name="+name; +// throw CanteraError("find_XML", "XML element with "+opt+ +// " not found."); +// } +// return r; +// } +// catch (CanteraError) { + +// // root was used, but element was not found. Try src. +// if (root && src != "") { +// return find_XML(src, 0, id, loc, name); +// } +// else { +// string opt = " src="+src+", loc="+loc+", id=" +// +id+", name="+name; +// throw CanteraError("find_XML", "XML element with "+opt+ +// " not found."); +// return 0; +// } +// } +// } +// #endif XML_Node * findXMLPhase(XML_Node *root, diff --git a/Cantera/src/xml.h b/Cantera/src/xml.h index d248703a8..5f0654141 100755 --- a/Cantera/src/xml.h +++ b/Cantera/src/xml.h @@ -129,7 +129,7 @@ namespace Cantera { int nChildren() const { return m_nchildren; } void build(istream& f); - + void require(string a, string v) const; /** * This routine carries out a search for an XML node based * on both the xml element name and the attribute ID. diff --git a/Cantera/src/zeroD/Wall.h b/Cantera/src/zeroD/Wall.h index 46fecfba9..a1f822288 100644 --- a/Cantera/src/zeroD/Wall.h +++ b/Cantera/src/zeroD/Wall.h @@ -1,12 +1,15 @@ /** - * @file FlowDevice.h - * - * $Author$ + * @file Wall.h + * Header file for class Wall. + */ + +/* $Author$ * $Date$ * $Revision$ */ -// Copyright 2001 California Institute of Technology +// Copyright 2001-2004 California Institute of Technology + #ifndef CT_WALL_H #define CT_WALL_H @@ -17,17 +20,21 @@ #endif #include "../ct_defs.h" +#include "../ctexceptions.h" #include "../Func1.h" namespace Cantera { - class ReactorBase; // forward reference + // forward references + class ReactorBase; class Kinetics; class Func1; class SurfPhase; - const int Rigid_Type = 1; - const int Flexible_Type = 2; + // const int Rigid_Type = 1; + // const int Flexible_Type = 2; + + class Wall { @@ -36,14 +43,17 @@ namespace Cantera { /// Constructor Wall(); - /// Destructor + /// Destructor. Since Wall instances do not allocate memory, + /// the destructor does nothing. virtual ~Wall() {} - /** - * Rate of volume change (kg/s). Positive value increases - * volume of reactor on left, and decreases volume on right. - */ + + /// Rate of volume change (kg/s). Positive value increases + /// volume of reactor on left, and decreases volume on right. virtual doublereal vdot(doublereal t); + + /// Heat flow rate through the wall (W). Positive values + /// denote a flux from left to right. virtual doublereal Q(doublereal t); /// Area in m^2. @@ -51,16 +61,18 @@ namespace Cantera { /// Set the area [m^2]. void setArea(doublereal a) { m_area = a; } - + void setThermalResistance(doublereal Rth) { m_rrth = 1.0/Rth; } /// Set the overall heat transfer coefficient [W/m^2/K]. void setHeatTransferCoeff(doublereal U) { m_rrth = U; } - void setEmissivity(doublereal epsilon) { m_emiss = epsilon; } - - // /** Set the rate of volume change to a specified function.*/ - // void setExpansionRate(Func1* f=0) {if (f) m_vf = f;} + /// Set the emissivity. + void setEmissivity(doublereal epsilon) { + if (epsilon > 1.0 || epsilon < 0.0) + throw CanteraError("Wall::setEmissivity", + "emissivity must be between 0.0 and 1.0"); + m_emiss = epsilon; } /** Set the piston velocity to a specified function. */ void setVelocity(Func1* f=0) {if (f) m_vf = f;} @@ -70,31 +82,39 @@ namespace Cantera { */ void setExpansionRateCoeff(doublereal k) {m_k = k;} - /** - * Specify the heat flux q(t). - */ + + /// Specify the heat flux function \f$ q_0(t) \f$. void setHeatFlux(Func1* q) { m_qf = q;} - bool install(ReactorBase& in, ReactorBase& out); + /// Install the wall between two reactors or reservoirs + bool install(ReactorBase& leftReactor, ReactorBase& rightReactor); + + /// True if the wall is correctly configured and ready to use. virtual bool ready() { return (m_left != 0 && m_right != 0); } - int type() { return 0; } + // int type() { return 0; } - /// Return a reference to the left reactor. + + /// Return a reference to the reactor or reservoir to the left + /// of the wall. ReactorBase& left() const { return *m_left; } - /// Return a reference to the right-hand reactor. + /// Return a reference to the reactor or reservoir to the + /// right of the wall. const ReactorBase& right() { return *m_right; } - /// set parameters - virtual void setParameters(int n, doublereal* coeffs) { - m_coeffs.resize(n); - copy(coeffs, coeffs + n, m_coeffs.begin()); - } + // /// Set wall parameters. + //virtual void setParameters(int n, doublereal* coeffs) { + // m_coeffs.resize(n); + // copy(coeffs, coeffs + n, m_coeffs.begin()); + //} - void setKinetics(Kinetics* left = 0, - Kinetics* right = 0); + // Specify the heterogeneous reaction mechanisms for each side + // of the wall. + void setKinetics(Kinetics* leftMechanism, Kinetics* rightMechanism); + /// Return a pointer to the surface phase object for the left + /// or right wall surface. SurfPhase* surface(int leftright) { return m_surf[leftright]; } @@ -103,16 +123,22 @@ namespace Cantera { return m_chem[leftright]; } + /// Set the surface coverages on the left or right surface to + /// the values in array 'cov'. void setCoverages(int leftright, const doublereal* cov); + /// Write the coverages of the left or right surface into + /// array cov. void getCoverages(int leftright, doublereal* cov); + /// Set the coverages in the surface phase object to the + /// values for this wall surface. void syncCoverages(int leftright); protected: - vector_fp m_coeffs; + //vector_fp m_coeffs; ReactorBase* m_left; ReactorBase* m_right; @@ -128,53 +154,6 @@ namespace Cantera { private: }; - - -// class Piston : public Wall { -// public: -// Piston() -// : m_omega(2.0*3.1415926*freq), Wall() { -// //m_vdisp = stroke * Pi * bore * bore / 4.0; -// //m_vclear = tdc * bore; -// //m_ra = 1.0/crankradius; -// } -// ~Piston() {} -// virtual doublereal vdot(double t) { -// doublereal theta = m_omega * t; -// doublereal sinth = sin(theta); -// return 0.0; -// //return m_vclear + 0.5*m_vdist*(1.0 + m_ra - m_omega*sin(theta) -// // - sqrt(m_ra * m_ra - sinth*sinth)); -// } -// protected: -// doublereal m_tdc, m_bdc, m_stroke, m_bore, m_rodlen, -// m_radius, m_omega; -// }; - - - class Piston : public Wall { - public: - Piston(doublereal freq, - doublereal tdc, doublereal bdc, - doublereal stroke, doublereal bore, - doublereal rodlen, doublereal crankradius) - : Wall(), m_omega(2.0*3.1415926*freq) { - //m_vdisp = stroke * Pi * bore * bore / 4.0; - //m_vclear = tdc * bore; - //m_ra = 1.0/crankradius; - } - virtual ~Piston() {} - virtual doublereal vdot(double t) { - // doublereal theta = m_omega * t; - //doublereal sinth = sin(theta); - return 0.0; - //return m_vclear + 0.5*m_vdist*(1.0 + m_ra - m_omega*sin(theta) - // - sqrt(m_ra * m_ra - sinth*sinth)); - } - protected: - doublereal m_tdc, m_bdc, m_stroke, m_bore, m_rodlen, - m_radius, m_omega; - }; }