From 81189197f16610b96eab00ce68d3c016947bc255 Mon Sep 17 00:00:00 2001 From: Dave Goodwin Date: Wed, 13 Aug 2003 13:56:34 +0000 Subject: [PATCH] *** empty log message *** --- Cantera/clib/src/ctbdry.cpp | 18 +- Cantera/clib/src/ctbdry.h | 1 + Cantera/clib/src/ctonedim.cpp | 17 ++ Cantera/clib/src/ctonedim.h | 2 + Cantera/matlab/cantera/examples/flame.m | 4 +- .../matlab/cantera/private/onedimmethods.cpp | 7 + Cantera/python/Cantera/ctml_writer.py | 87 ++++---- Cantera/python/Cantera/solution.py | 3 +- Cantera/src/importCTML.cpp | 167 +++++++++++---- Cantera/src/importCTML.h | 2 +- Cantera/src/oneD/Inlet1D.h | 73 ++++++- Cantera/src/oneD/Sim1D.cpp | 2 +- Cantera/src/oneD/Solid1D.cpp | 116 +++-------- Cantera/src/oneD/boundaries1D.cpp | 195 +++++++++++++++++- Cantera/src/oneD/refine.cpp | 28 ++- Cantera/src/oneD/refine.h | 8 +- Cantera/src/xml.cpp | 1 + 17 files changed, 543 insertions(+), 188 deletions(-) diff --git a/Cantera/clib/src/ctbdry.cpp b/Cantera/clib/src/ctbdry.cpp index 669aa38a4..ba3a12635 100755 --- a/Cantera/clib/src/ctbdry.cpp +++ b/Cantera/clib/src/ctbdry.cpp @@ -12,7 +12,7 @@ // Cantera includes #include "oneD/OneDim.h" #include "oneD/Inlet1D.h" - +#include "InterfaceKinetics.h" #include "Cabinet.h" #include "Storage.h" @@ -34,6 +34,10 @@ inline ThermoPhase* _thermo(int n) { return Storage::__storage->__thtable[n]; } +inline Kinetics* _kin(int n) { + return Storage::__storage->__ktable[n]; +} + extern "C" { int DLL_EXPORT bndry_new(int itype) { @@ -45,6 +49,8 @@ extern "C" { s = new Symm1D(); break; case 3: s = new Surf1D(); break; + case 4: + s = new ReactingSurf1D(); break; default: return -2; } @@ -114,4 +120,14 @@ extern "C" { catch (CanteraError) {return -1;} return 0; } + + int DLL_EXPORT surf_setkinetics(int i, int j) { + try { + ReactingSurf1D* srf = (ReactingSurf1D*)_bndry(i); + InterfaceKinetics* k = (InterfaceKinetics*)_kin(j); + srf->setKineticsMgr(k); + } + catch (CanteraError) {return -1;} + return 0; + } } diff --git a/Cantera/clib/src/ctbdry.h b/Cantera/clib/src/ctbdry.h index d94ae0244..34e5b7255 100755 --- a/Cantera/clib/src/ctbdry.h +++ b/Cantera/clib/src/ctbdry.h @@ -15,5 +15,6 @@ extern "C" { double DLL_IMPORT bndry_mdot(int i); int DLL_IMPORT bndry_setxin(int i, double* xin); int DLL_IMPORT bndry_setxinbyname(int i, char* xin); + int DLL_IMPORT bndry_setkinetics(int i, int j); } #endif diff --git a/Cantera/clib/src/ctonedim.cpp b/Cantera/clib/src/ctonedim.cpp index 387eef95e..05f8e1b20 100644 --- a/Cantera/clib/src/ctonedim.cpp +++ b/Cantera/clib/src/ctonedim.cpp @@ -179,6 +179,14 @@ extern "C" { catch (CanteraError) { return -1; } } + int DLL_EXPORT reactingsurf_new() { + try { + ReactingSurf1D* i = new ReactingSurf1D(); + return Cabinet::cabinet()->add(i); + } + catch (CanteraError) { return -1; } + } + int DLL_EXPORT symm_new() { try { Symm1D* i = new Symm1D(); @@ -240,6 +248,15 @@ extern "C" { catch (CanteraError) { return -1; } } + int DLL_EXPORT reactingsurf_setkineticsmgr(int i, int j) { + try { + ReactingSurf1D* srf = (ReactingSurf1D*)_bdry(i); + InterfaceKinetics* k = (InterfaceKinetics*)_kinetics(j); + srf->setKineticsMgr(k); + return 0; + } + catch (CanteraError) { return -1; } + } //------------------ stagnation flow domains -------------------- diff --git a/Cantera/clib/src/ctonedim.h b/Cantera/clib/src/ctonedim.h index 1eab0aa77..5ddbdd0eb 100644 --- a/Cantera/clib/src/ctonedim.h +++ b/Cantera/clib/src/ctonedim.h @@ -27,11 +27,13 @@ extern "C" { double DLL_IMPORT bdry_temperature(int i); double DLL_IMPORT bdry_massFraction(int i, int k); double DLL_IMPORT bdry_mdot(int i); + int DLL_IMPORT reactingsurf_setkineticsmgr(int i, int j); int DLL_IMPORT inlet_new(); int DLL_IMPORT outlet_new(); int DLL_IMPORT symm_new(); int DLL_IMPORT surf_new(); + int DLL_IMPORT reactingsurf_new(); int DLL_IMPORT stflow_new(int iph, int ikin, int itr); int DLL_IMPORT stflow_setPressure(int i, double p); diff --git a/Cantera/matlab/cantera/examples/flame.m b/Cantera/matlab/cantera/examples/flame.m index e127c6b54..cc9746f17 100644 --- a/Cantera/matlab/cantera/examples/flame.m +++ b/Cantera/matlab/cantera/examples/flame.m @@ -64,8 +64,8 @@ end zz = z(flow); dz = zz(end) - zz(1); setProfile(f, 2, {'u', 'V'}, [0.0 1.0 - mdot0/rho0 -mdot1/rho0 - 0.0 0.0]); + mdot0/rho0 -mdot1/rho0 + 0.0 0.0]); setProfile(f, 2, 'T', [0.0 z1 1.0; t0 2000.0 t1]); for n = 1:nSpecies(gas) diff --git a/Cantera/matlab/cantera/private/onedimmethods.cpp b/Cantera/matlab/cantera/private/onedimmethods.cpp index 9563892da..41818ab35 100644 --- a/Cantera/matlab/cantera/private/onedimmethods.cpp +++ b/Cantera/matlab/cantera/private/onedimmethods.cpp @@ -57,6 +57,13 @@ void onedimmethods( int nlhs, mxArray *plhs[], indx = outlet_new(); break; + // construct a new ReactingSurf1D instance + case 6: + checkNArgs(4, nrhs); + indx = reactingsurf_new(); + reactingsurf_setkineticsmgr(indx, getInt(prhs[3])); + break; + // construct a new Sim1D instance case 8: checkNArgs(5, nrhs); diff --git a/Cantera/python/Cantera/ctml_writer.py b/Cantera/python/Cantera/ctml_writer.py index 8aa5802d8..a0c90d83a 100644 --- a/Cantera/python/Cantera/ctml_writer.py +++ b/Cantera/python/Cantera/ctml_writer.py @@ -49,6 +49,15 @@ _reactions = [] _atw = {} _mw = {} +_valsp = '' +_valrxn = '' + +def validate(species = 'yes', reactions = 'yes'): + global _valsp + global _valrxn + _valsp = species + _valrxn = reactions + def isnum(a): """True if a is an integer or floating-point number.""" if type(a) == types.IntType or type(a) == types.FloatType: @@ -119,6 +128,10 @@ def ufmt(base, n): def write(): """write the CTML file.""" x = XML_Node("ctml") + v = x.addChild("validate") + v["species"] = _valsp + v["reactions"] = _valrxn + for ph in _phases: ph.build(x) s = species_set(name = _name, species = _species) @@ -264,6 +277,8 @@ class species(writer): global _species _species.append(self) global _speciesnames + if name in _speciesnames: + raise CanteraError('species '+name+' multiply defined.') _speciesnames.append(name) @@ -437,11 +452,16 @@ class reaction(writer): def __init__(self, equation = '', kf = None, - id = '', + id = '', + special = [] ): self._id = id self._e = equation + if type(special) == types.StringType: + self._special = [special] + else: + self._special = special global _reactions self._num = len(_reactions)+1 r = '' @@ -493,37 +513,10 @@ class reaction(writer): if ph.is_ideal_gas(): self._igspecies.append(s) break -## if nm < 0: -## if _handle_undeclared_species == _WARN: -## print 'Warning... skipping reaction '+self._e -## print 'because species '+s+' is undeclared' -## if _handle_undeclared_species <= _WARN: -## return 0 -## else: -## raise CanteraError('undeclared species '+s+' in reaction '+self._e) -## else: mdim += nm*ns ldim += nl*ns -## for s in self._p.keys(): -## ns = self._p[s] -## ok = 0 -## for ph in _phases: -## if ph.has_species(s): -## ok = 1 -## break - -## if not ok: -## if _handle_undeclared_species == _WARN: -## print 'Warning... skipping reaction '+self._e -## print 'because species '+s+' is undeclared' -## if _handle_undeclared_species <= _WARN: -## return 0 -## else: -## raise CanteraError('undeclared species '+s+' in reaction '+self._e) - - p.addComment(" reaction "+id+" ") r = p.addChild('reaction') r['id'] = id @@ -531,6 +524,13 @@ class reaction(writer): r['reversible'] = 'yes' else: r['reversible'] = 'no' + nspecial = len(self._special) + for nss in range(nspecial): + s = self._special[nss] + if s == 'duplicate': + r['duplicate'] = 'yes' + elif s == 'negative_A': + r['negative_A'] = 'yes' ee = self._e.replace('<','[') ee = ee.replace('>',']') @@ -596,15 +596,17 @@ class reaction(writer): #------------------- + class three_body_reaction(reaction): def __init__(self, equation = '', kf = None, efficiencies = '', - id = '', + id = '', + special = [] ): - reaction.__init__(self, equation, kf, id) + reaction.__init__(self, equation, kf, id, special) self._type = 'threeBody' self._effm = 1.0 self._eff = efficiencies @@ -615,8 +617,7 @@ class three_body_reaction(reaction): del self._r[r] for p in self._p.keys(): if p == 'M' or p == 'm': - del self._p[p] - + del self._p[p] def build(self, p): @@ -639,7 +640,8 @@ class falloff_reaction(reaction): kf = None, efficiencies = '', falloff = None, - id = '' + id = '', + special = [] ): kf2 = (kf, kf0) @@ -693,7 +695,8 @@ class surface_reaction(reaction): equation = '', kf = None, stick = None, - id = ''): + id = '', + special = []): if stick: reaction.__init__(self, equation, stick, id) else: @@ -775,26 +778,29 @@ class phase(writer): sptoks = spnames.replace(',',' ').split() for s in sptoks: if self._spmap.has_key(s): - raise 'multiply-defined species '+s+' in phase '+self._name + raise CanteraError('Multiply-declared species '+s+' in phase '+self._name) self._spmap[s] = self._dim self._rxns = reactions # check that species have been declared if len(self._spmap) == 0: - raise 'No species declared for phase '+self._name + raise CanteraError('No species declared for phase '+self._name) # and that only one species is declared if it is a pure phase if self.is_pure() and len(self._spmap) > 1: - raise 'Pure phases may only declare one species, but phase '+self._name+' declares '+`len(self._spmap)`+'.' + raise CanteraError('Stoichiometric phases must declare exactly one species, \n'+ + 'but phase '+self._name+' declares '+`len(self._spmap)`+'.') self._initial = initial_state # add this phase to the global phase list global _phases _phases.append(self) + def is_ideal_gas(self): + """True if the entry represents an ideal gas.""" return 0 def is_pure(self): @@ -853,6 +859,12 @@ class phase(writer): ph = p.addChild('phase') ph['id'] = self._name ph['dim'] = `self._dim` + + # ------- error tests ------- + err = ph.addChild('validation') + err.addChild('duplicateReactions','halt') + err.addChild('thermo','warn') + e = ph.addChild('elementArray',self._el) e['datasrc'] = 'elements.xml' for s in self._sp: @@ -1116,3 +1128,4 @@ class Lindemann: get_atomic_wts() +validate() diff --git a/Cantera/python/Cantera/solution.py b/Cantera/python/Cantera/solution.py index 13c6a382b..46da26cd5 100755 --- a/Cantera/python/Cantera/solution.py +++ b/Cantera/python/Cantera/solution.py @@ -39,8 +39,7 @@ class Solution(ThermoPhase, Kinetics, Transport): id = "" if len(fn) > 1: id = fn[1] - fn = fn[0] - + fn = fn[0] fname = os.path.basename(fn) ff = os.path.splitext(fname) diff --git a/Cantera/src/importCTML.cpp b/Cantera/src/importCTML.cpp index d2a6925a8..097945d88 100755 --- a/Cantera/src/importCTML.cpp +++ b/Cantera/src/importCTML.cpp @@ -50,6 +50,12 @@ using namespace ctml; GasKineticsWriter* writer = 0; +vector< map > _reactiondata; +vector _eqn; +vector_int _dup; +vector _rev; + + namespace Cantera { /* * First we define a coule of typedef's which will @@ -69,7 +75,7 @@ namespace Cantera { * Install a NASA polynomial thermodynamic property * parameterization for species k. */ - void installNasaThermo(SpeciesThermo& sp, int k, XML_Node& f0, + static void installNasaThermo(SpeciesThermo& sp, int k, XML_Node& f0, XML_Node& f1) { doublereal tmin0, tmax0, tmin1, tmax1, tmin, tmid, tmax; @@ -113,7 +119,7 @@ namespace Cantera { * Install a Shomate polynomial thermodynamic property * parameterization for species k. */ - void installShomateThermo(SpeciesThermo& sp, int k, XML_Node& f) { + static void installShomateThermo(SpeciesThermo& sp, int k, XML_Node& f) { doublereal tmin, tmid, tmax; tmin = fpValue(f["Tmin"]); tmid = fpValue(f["Tmid"]); @@ -139,17 +145,15 @@ namespace Cantera { } /** - * Install a Shomate polynomial thermodynamic property + * Install a constant-cp thermodynamic property * parameterization for species k. */ - void installSimpleThermo(SpeciesThermo& sp, int k, XML_Node& f) { + static void installSimpleThermo(SpeciesThermo& sp, int k, XML_Node& f) { doublereal tmin, tmax; tmin = fpValue(f["Tmin"]); tmax = fpValue(f["Tmax"]); if (tmax == 0.0) tmax = 1.0e30; - //map fp; - //getFloats(f, fp); vector_fp c(4); c[0] = getFloat(f, "t0", "-"); c[1] = getFloat(f, "h0", "-"); @@ -163,7 +167,7 @@ namespace Cantera { * Install a species into a ThermoPhase object, which defines * the phase thermodynamics and speciation */ - bool installSpecies(int k, XML_Node& s, thermo_t& p, + static bool installSpecies(int k, XML_Node& s, thermo_t& p, SpeciesThermo& spthermo, int rule) { // get the composition of the species @@ -179,9 +183,8 @@ namespace Cantera { if (rule == 0) throw CanteraError("installSpecies", - "species " + s["name"] + - " contains undeclared element " + - _b->first); + "Species " + s["name"] + + " contains undeclared element " + _b->first); else return false; } @@ -202,8 +205,7 @@ namespace Cantera { doublereal sz = 1.0; if (s.hasChild("size")) sz = getFloat(s, "size"); - p.addUniqueSpecies(s["name"], ecomp.begin(), - chrg, sz); + p.addUniqueSpecies(s["name"], ecomp.begin(), chrg, sz); // get thermo XML_Node& thermo = s.child("thermo"); @@ -211,9 +213,6 @@ namespace Cantera { int nc = tp.size(); if (nc == 1) { XML_Node& f = *tp[0]; - //if (f.name() == "NASA") { - // installNasaThermo(spthermo, k, f); - //} if (f.name() == "Shomate") { installShomateThermo(spthermo, k, f); } @@ -268,7 +267,7 @@ namespace Cantera { * rule = If we fail to find a species, we will throw an error * if rule != 1. */ - bool getReagents(XML_Node& rxn, kinetics_t& kin, int rp, + static bool getReagents(XML_Node& rxn, kinetics_t& kin, int rp, string default_phase, vector_int& spnum, vector_int& stoich, vector_fp& order, int rule) { @@ -312,7 +311,7 @@ namespace Cantera { return false; else { throw CanteraError("getReagents", - "undeclared reactant or product species "+sp); + "Undeclared reactant or product species "+sp); return false; } } @@ -337,7 +336,7 @@ namespace Cantera { * Arrhenius expression is * k = A T^(b) exp (-Ea / RT). */ - void getArrhenius(XML_Node& node, int& highlow, doublereal& A, + static void getArrhenius(XML_Node& node, int& highlow, doublereal& A, doublereal& b, doublereal& E) { if (node["name"] == "k0") @@ -366,9 +365,7 @@ namespace Cantera { const ThermoPhase& p = kin.thermo(np); klocal = p.speciesIndex(kin.kineticsSpeciesName(k)); if (p.eosType() == cSurf) { - cout << "f before = " << f << endl; f /= pow(p.standardConcentration(klocal),ns); - cout << "f, k, klocal = " << f << " " << k << " " << klocal << endl; } else not_surf++; @@ -383,14 +380,13 @@ namespace Cantera { b = getFloat(node, "b") + 0.5; E = getFloat(node, "E", "actEnergy"); E /= GasConstant; - cout << "A, cbar, f, b, E: " << A << " " << cbar << " " << f << " " << b << " " << E << endl; } /** * Get falloff parameters for a reaction. */ - void getFalloff(node_t& f, ReactionData& rdata) { + static void getFalloff(node_t& f, ReactionData& rdata) { string type = f["type"]; vector p; getStringArray(f,p); @@ -415,7 +411,7 @@ namespace Cantera { * reaction mechanism is homogeneous, so that all species belong * to phase(0) of 'kin'. */ - void getEfficiencies(node_t& eff, kinetics_t& kin, ReactionData& rdata) { + static void getEfficiencies(node_t& eff, kinetics_t& kin, ReactionData& rdata) { // set the default collision efficiency rdata.default_3b_eff = fpValue(eff["default"]); @@ -440,8 +436,8 @@ namespace Cantera { * This function will fill in more fields in the ReactionData object. * */ - void getRateCoefficient(node_t& kf, kinetics_t& kin, - ReactionData& rdata) { + static void getRateCoefficient(node_t& kf, kinetics_t& kin, + ReactionData& rdata, int negA) { int nc = kf.nChildren(); const nodeset_t& kf_children = kf.children(); @@ -459,6 +455,10 @@ namespace Cantera { || rdata.reactionType == ELEMENTARY_RXN) chigh = coeff; else clow = coeff; + if (coeff[0] <= 0.0 && negA == 0) { + throw CanteraError("getRateCoefficient", + "negative or zero A coefficient for reaction "+int2str(rdata.number)); + } } else if (nm == "Stick") { vector_fp coeff(3); @@ -467,6 +467,10 @@ namespace Cantera { int isp = th.speciesIndex(spname); double mw = th.molecularWeights()[isp]; getStick(c, mw, kin, rdata, coeff[0], coeff[1], coeff[2]); + if (coeff[0] <= 0.0 && negA == 0) { + throw CanteraError("getRateCoefficient", + "negative or zero A coefficient for reaction "+int2str(rdata.number)); + } chigh = coeff; } @@ -761,6 +765,29 @@ namespace Cantera { return true; } + /** + * This function returns true if two reactions are duplicates of + * one another, and false otherwise. The input arguments are two + * maps from species number to stoichiometric coefficient, one for + * each reaction. The reactions are considered duplicates if their + * stoichiometric coefficients have the same ratio for all + * species. + */ + static bool isDuplicateReaction(map& r1, map& r2) { + + map::const_iterator b = r1.begin(), e = r1.end(); + int k1 = b->first; + doublereal ratio = r2[k1]/r1[k1]; + if (r1[k1] == 0.0 || r2[k1] == 0.0) return false; + ++b; + for (; b != e; ++b) { + k1 = b->first; + if (r1[k1] == 0.0 || r2[k1] == 0.0) return false; + if (fabs(r2[k1]/r1[k1] - ratio) > 1.e-8) return false; + } + return true; + } + /** * Install an individual reaction into the kinetics mechanism @@ -773,8 +800,8 @@ namespace Cantera { * rule = Provides a rule for specifying how to handle reactions * which involve missing species. */ - bool installReaction(int i, XML_Node& r, Kinetics* k, - string default_phase, int rule) { + static bool installReaction(int i, XML_Node& r, Kinetics* k, + string default_phase, int rule, bool check_for_duplicates) { Kinetics& kin = *k; /* @@ -790,11 +817,18 @@ namespace Cantera { int nn, eqlen; vector_fp dummy; + int dup = 0; + if (r.hasAttrib("duplicate")) dup = 1; + int negA = 0; + if (r.hasAttrib("negative_A")) negA = 1; + + /* * This seemingly simple expression goes and finds the child element, * "equation". Then it treats all of the contents of the "equation" * as a string, and returns it the variable eqn. We post process - * the string to get rid of [ and ] characters for some reason. + * the string to convert [ and ] characters into < and >, which cannot be + * stored in an XML file. */ if (r.hasChild("equation")) eqn = r("equation"); @@ -821,6 +855,47 @@ namespace Cantera { return false; } + string isrev = r["reversible"]; + if (isrev == "yes" || isrev == "true") + rdata.reversible = true; + + if (check_for_duplicates) { + doublereal c = 0.0; + map rxnstoich; + rxnstoich.clear(); + int nr = rdata.reactants.size(); + for (nn = 0; nn < nr; nn++) { + rxnstoich[rdata.reactants[nn]] -= rdata.rstoich[nn]; + } + int np = rdata.products.size(); + for (nn = 0; nn < np; nn++) { + rxnstoich[rdata.products[nn]] += rdata.pstoich[nn]; + } + int nrxns = _reactiondata.size(); + for (nn = 0; nn < nrxns; nn++) { + c = isDuplicateReaction(rxnstoich, _reactiondata[nn]); + if (c > 0.0 + || (c < 0.0 && rdata.reversible) + || (c < 0.0 && _rev[nn])) { + + if ((!dup || !_dup[nn])) { + string msg = string("Undeclared duplicate reactions detected: \n") + +"Reaction "+int2str(nn+1)+": "+_eqn[nn] + +"\nReaction "+int2str(i+1)+": "+eqn+"\n"; + _reactiondata.clear(); + _eqn.clear(); + _rev.clear(); + _dup.clear(); + throw CanteraError("installReaction",msg); + } + } + } + _dup.push_back(dup); + _rev.push_back(rdata.reversible); + _eqn.push_back(eqn); + _reactiondata.push_back(rxnstoich); + } + rdata.equation = eqn; rdata.reversible = false; rdata.number = i; @@ -849,11 +924,7 @@ namespace Cantera { throw CanteraError("installReaction", "Unknown reaction type: " + typ); - string isrev = r["reversible"]; - if (isrev == "yes" || isrev == "true") - rdata.reversible = true; - - getRateCoefficient(r.child("rateCoeff"), kin, rdata); + getRateCoefficient(r.child("rateCoeff"), kin, rdata, negA); /* * Ok we have read everything in about the reaction. Add it * to the kinetics object by calling the Kinetics member function, @@ -878,7 +949,13 @@ namespace Cantera { * If there is a problem, return false. */ bool installReactionArrays(XML_Node& p, Kinetics& kin, - string default_phase) { + string default_phase, bool check_for_duplicates) { + + _eqn.clear(); + _dup.clear(); + _reactiondata.clear(); + _rev.clear(); + vector rarrays; int itot = 0; /* @@ -944,17 +1021,13 @@ namespace Cantera { XML_Node* r = allrxns[i]; if (r) { if (installReaction(itot, *r, &kin, - default_phase, rxnrule)) ++itot; + default_phase, rxnrule, check_for_duplicates)) ++itot; } } } else { for (int nii = 0; nii < ninc; nii++) { XML_Node& ii = *incl[nii]; - //vector rxn_ids; - //string pref = ii["prefix"]; - //int imin = atoi(ii["min"].c_str()); - //int imax = atoi(ii["max"].c_str()); string imin = ii["min"]; string imax = ii["max"]; for (i = 0; i < nrxns; i++) { @@ -962,8 +1035,6 @@ namespace Cantera { string rxid; if (r) { rxid = (*r)["id"]; - //cout << rxid << " " << imin << " " << imax << endl; - //cout << (rxid >= imin) << " " << (rxid <= imax) << endl; /* * To decide whether the reaction is included or not * we do a lexical min max and operation. This @@ -971,19 +1042,23 @@ namespace Cantera { */ if ((rxid >= imin) && (rxid <= imax)) { if (installReaction(itot, *r, &kin, - default_phase, rxnrule)) ++itot; + default_phase, rxnrule, check_for_duplicates)) ++itot; } } } } } } + /* * Finalize the installation of the kinetics, now that we know * the true number of reactions in the mechanism, itot. */ kin.finalize(); writer = 0; + _eqn.clear(); + _dup.clear(); + _reactiondata.clear(); return true; } @@ -999,6 +1074,12 @@ namespace Cantera { // This phase will be the default one string default_phase = phase["id"]; + bool check_for_duplicates = false; + if (phase.parent()->hasChild("validate")) { + XML_Node& d = phase.parent()->child("validate"); + if (d["reactions"] == "yes") check_for_duplicates = true; + } + // if other phases are involved in the reaction mechanism, // they must be listed in a 'phaseArray' child // element. Homogeneous mechanisms do not need to include a @@ -1049,7 +1130,7 @@ namespace Cantera { kin.init(); // Install the reactions. - return installReactionArrays(phase, kin, default_phase); + return installReactionArrays(phase, kin, default_phase, check_for_duplicates); } /** diff --git a/Cantera/src/importCTML.h b/Cantera/src/importCTML.h index f8208f49c..421fceb08 100755 --- a/Cantera/src/importCTML.h +++ b/Cantera/src/importCTML.h @@ -25,7 +25,7 @@ namespace Cantera { bool importKinetics(XML_Node& phase, vector th, Kinetics* kin); bool installReactionArrays(XML_Node& parent, Kinetics& kin, - string default_phase); + string default_phase, bool check_for_duplicates = false); ThermoPhase* newPhase(XML_Node& phase); bool buildSolutionFromXML(XML_Node& root, string id, string nm, ThermoPhase* th, Kinetics* k); diff --git a/Cantera/src/oneD/Inlet1D.h b/Cantera/src/oneD/Inlet1D.h index 7c9854117..e0f9e1a53 100644 --- a/Cantera/src/oneD/Inlet1D.h +++ b/Cantera/src/oneD/Inlet1D.h @@ -109,7 +109,7 @@ namespace Cantera { * without parameters, a left inlet (facing right) is * constructed). */ - Inlet1D() : m_V0(0.0), m_nsp(0), m_flow(0) { + Inlet1D() : Bdry1D(), m_V0(0.0), m_nsp(0), m_flow(0) { m_type = cInletType; m_xstr = ""; } @@ -200,7 +200,7 @@ namespace Cantera { public: - Symm1D() { + Symm1D() : Bdry1D() { m_type = cSymmType; } virtual ~Symm1D(){} @@ -228,7 +228,7 @@ namespace Cantera { public: - Outlet1D() { + Outlet1D() : Bdry1D() { m_type = cOutletType; } virtual ~Outlet1D(){} @@ -260,7 +260,7 @@ namespace Cantera { public: - Surf1D() { + Surf1D() : Bdry1D() { m_type = cSurfType; } virtual ~Surf1D(){} @@ -299,6 +299,71 @@ namespace Cantera { }; + + /** + * A reacting surface. + * + */ + class ReactingSurf1D : public Bdry1D { + + public: + + ReactingSurf1D() : Bdry1D(), + m_kin(0), m_surfindex(0), m_nsp(0) { + m_type = cSurfType; + } + + void setKineticsMgr(InterfaceKinetics* kin) { + m_surfindex = kin->surfacePhaseIndex(); + m_sphase = (SurfPhase*)&kin->thermo(m_surfindex); + m_nsp = m_sphase->nSpecies(); + m_enabled = true; + } + + virtual ~ReactingSurf1D(){} + + virtual string componentName(int n) const; + + virtual void init(); + + virtual void eval(int jg, doublereal* xg, doublereal* rg, + integer* diagg, doublereal rdt); + + virtual void save(XML_Node& o, doublereal* soln); + virtual void restore(XML_Node& dom, doublereal* soln); + + virtual void _getInitialSoln(doublereal* x) { + x[0] = m_temp; + m_sphase->getCoverages(x+1); + } + + virtual void _finalize(const doublereal* x) { + ; //m_temp = x[0]; + } + + virtual void showSolution(const doublereal* x) { + char buf[80]; + sprintf(buf, " Temperature: %10.4g K \n", x[0]); + writelog(buf); + writelog(" Coverages: \n"); + for (int k = 0; k < m_nsp; k++) { + sprintf(buf, " %20s %10.4g \n", m_sphase->speciesName(k).c_str(), + x[k+1]); + writelog(buf); + } + writelog("\n"); + } + + protected: + + InterfaceKinetics* m_kin; + SurfPhase* m_sphase; + int m_surfindex, m_nsp; + bool m_enabled; + vector_fp m_fixed_cov; + vector_fp m_work; + }; + } // // }; diff --git a/Cantera/src/oneD/Sim1D.cpp b/Cantera/src/oneD/Sim1D.cpp index da35dc12f..7b7ab599a 100644 --- a/Cantera/src/oneD/Sim1D.cpp +++ b/Cantera/src/oneD/Sim1D.cpp @@ -337,7 +337,7 @@ namespace Cantera { } } else { - throw CanteraError("refine","keepPoint is false at m = "+int2str(m)); + ; // throw CanteraError("refine","keepPoint is false at m = "+int2str(m)); } } dsize.push_back(znew.size() - nstart); diff --git a/Cantera/src/oneD/Solid1D.cpp b/Cantera/src/oneD/Solid1D.cpp index 139995ceb..40821b342 100644 --- a/Cantera/src/oneD/Solid1D.cpp +++ b/Cantera/src/oneD/Solid1D.cpp @@ -31,32 +31,23 @@ namespace Cantera { int Solid1D::c_T_loc = 0; - int Solid1D::c_phi_loc = 1; - int Solid1D::c_Y_loc = 2; + int Solid1D::c_C_loc = 1; - - Solid1D::Solid1D(ThermoPhase* ph, int nsp, int points) : - Domain1D(nsp+1, points), - m_nsp(nsp), - m_thermo(0), + Solid1D::Solid1D(ThermoPhase* ph, int points) : + Domain1D(1, points), m_kin(0), m_trans(0), m_jac(0), - m_ok(false), + m_ok(false) { m_type = cSolidType; - m_points = points; m_thermo = ph; - if (ph == 0) return; // used to create a dummy object - - int nsp2 = m_thermo->nSpecies(); - if (nsp2 != m_nsp) { - m_nsp = nsp2; - Domain1D::resize(m_nsp+1, points); - } + if (ph == 0) { m_nsp = 1; return; }// used to create a dummy object + m_nsp = m_thermo->nSpecies(); + Domain1D::resize(m_nsp+1, points); // make a local copy of the species molecular weight vector m_wt = m_thermo->molecularWeights(); @@ -65,13 +56,12 @@ namespace Cantera { // turn off the energy equation at all points m_do_energy.resize(m_points,false); - m_do_gauss.resize(m_points,false); m_do_species.resize(m_nsp,false); m_diff.resize(m_nsp*m_points); m_flux.resize(m_nsp,m_points); m_wdot.resize(m_nsp,m_points, 0.0); - m_ybar.resize(m_nsp); + m_cbar.resize(m_nsp); //-------------- default solution bounds -------------------- @@ -82,11 +72,11 @@ namespace Cantera { vmin[c_T_loc] = 200.0; vmax[c_T_loc]= 1.e9; - // mass fraction bounds + // concentration bounds int k; for (k = 0; k < m_nsp; k++) { - vmin[c_Y_loc + k] = -1.0e-5; - vmax[c_Y_loc + k] = 1.0e5; + vmin[c_C_loc + k] = -1.0e-5; + vmax[c_C_loc + k] = 1.0e5; } setBounds(vmin.size(), vmin.begin(), vmax.size(), vmax.begin()); @@ -117,18 +107,13 @@ namespace Cantera { m_wtm.resize(m_points, 0.0); m_cp.resize(m_points, 0.0); m_tcon.resize(m_points, 0.0); - m_cdens.resize(m_points, 0.0); m_diff.resize(m_nsp*m_points); m_flux.resize(m_nsp,m_points); m_wdot.resize(m_nsp,m_points, 0.0); m_do_energy.resize(m_points,false); - m_do_gauss.resize(m_points,false); - m_do_species.resize(m_nsp,false); - m_fixedtemp.resize(m_points); - m_fixedphi.resize(m_points); m_dz.resize(m_points-1); m_z.resize(m_points); @@ -164,9 +149,8 @@ namespace Cantera { */ void Solid1D::setThermoState(const doublereal* x,int j) { m_thermo->setTemperature(T(x,j)); - m_thermo->setElectricPotential(phi(x,j)); const doublereal* yy = x + m_nv*j + 1; - m_thermo->setMassFractions_NoNorm(yy); + m_thermo->setConcentrations(yy); } @@ -176,12 +160,11 @@ namespace Cantera { */ void Solid1D::setStateAtMidpoint(const doublereal* x,int j) { m_thermo->setTemperature(0.5*(T(x,j)+T(x,j+1))); - m_thermo->setElectricPotential(0.5*(phi(x,j)+phi(x,j+1))); - const doublereal* yyj = x + m_nv*j + 1; - const doublereal* yyjp = x + m_nv*(j+1) + 1; + const doublereal* ccj = x + m_nv*j + 1; + const doublereal* ccjp = x + m_nv*(j+1) + 1; for (int k = 0; k < m_nsp; k++) - m_ybar[k] = 0.5*(yyj[k] + yyjp[k]); - m_thermo->setMassFractions_NoNorm(m_ybar.begin()); + m_ybar[k] = 0.5*(ccj[k] + ccjp[k]); + m_thermo->setConcentrations(m_cbar.begin()); } @@ -221,7 +204,6 @@ namespace Cantera { int j, k; - //----------------------------------------------------- // update properties //----------------------------------------------------- @@ -242,7 +224,6 @@ namespace Cantera { for (j = j0; j <= j1; j++) { setThermoState(j); - m_cdens[j] = m_thermo->chargeDensity(); } //---------------------------------------------------- @@ -259,13 +240,12 @@ namespace Cantera { if (j == 0) { rsd[index(c_T_loc,0)] = T(x,0); - rsd[index(c_phi_loc,0)] = phi(x,0); // The default boundary condition for species is zero // flux. However, the boundary object may modify // this. for (k = 0; k < m_nsp; k++) { - rsd[index(c_Y_loc + k, 0)] = - m_flux(k,0); + rsd[index(c_C_loc + k, 0)] = - m_flux(k,0); } } @@ -278,14 +258,9 @@ namespace Cantera { else if (j == m_points - 1) { rsd[index(c_T_loc,j)] = T(x,j); - rsd[index(c_phi_loc,j)] = phi(x,j); - doublereal sum = 0.0; for (k = 0; k < m_nsp; k++) { - sum += Y(x,k,j); - rsd[index(k+c_Y_loc,j)] = m_flux(k,j-1); + rsd[index(k+c_C_loc,j)] = m_flux(k,j-1); } - rsd[index(c_Y_loc,j)] = 1.0 - sum; - diag[index(c_Y_loc,j)] = 0; } @@ -303,15 +278,14 @@ namespace Cantera { //------------------------------------------------- getWdot(x,j); - doublereal convec, diffus; + doublereal diffus; for (k = 0; k < m_nsp; k++) { diffus = 2.0*(m_flux(k,j) - m_flux(k,j-1)) /(z(j+1) - z(j-1)); - rsd[index(c_Y_loc + k, j)] - = (m_wt[k]*(wdot(k,j) ) - - diffus)/m_rho[j] - - rdt*(Y(x,k,j) - Y_prev(k,j)); - diag[index(c_Y_loc + k, j)] = 1; + rsd[index(c_C_loc + k, j)] + = wdot(k,j) - diffus + - rdt*(C(x,k,j) - C_prev(k,j)); + diag[index(c_C_loc + k, j)] = 1; } //----------------------------------------------- @@ -326,12 +300,6 @@ namespace Cantera { rsd[index(c_T_loc, j)] -= rdt*(T(x,j) - T_prev(j)); diag[index(c_T_loc, j)] = 1; } - - //---------------------------------------------- - // Gauss's equation - //---------------------------------------------- - rsd[index(c_phi_loc, j)] = m_cdens[j] - divDisplCurr(x,j); - } // residual equations if the energy or species equations @@ -341,10 +309,6 @@ namespace Cantera { rsd[index(c_T_loc, j)] = T(x,j) - T_fixed(j); diag[index(c_T_loc, j)] = 0; } - if (!m_do_gauss[j]) { - rsd[index(c_phi_loc, j)] = phi(x,j) - phi_fixed(j); - diag[index(c_phi_loc, j)] = 0; - } } } @@ -424,19 +388,14 @@ namespace Cantera { doublereal dphidz, a1; for (j = j0; j < j1; j++) { sum = 0.0; - wtm = m_wtm[j]; rho = density(j); dz = z(j+1) - z(j); - dphidz = (phi(x,j+1) - phi(x,j))/dz; - a1 = rho*Faraday*dphidz/(GasConstant * T(x,j)); for (k = 0; k < m_nsp; k++) { - m_flux(k,j) = m_wt[k]*(rho*m_diff[k+m_nsp*j]/wtm); - m_flux(k,j) *= (X(x,k,j) - X(x,k,j+1))/dz; - m_flux(k,j) += a1*0.5*(Y(x,k,j) - + Y(x,k,j+1))*m_diff[k+m_nsp*j]*m_charge[k]; + m_flux(k,j) = m_diff[k+m_nsp*j] * + (C(x,k,j) - C(x,k,j+1))/dz; sum -= m_flux(k,j); } - for (k = 0; k < m_nsp; k++) m_flux(k,j) += Y(x,k,j)*sum; + for (k = 0; k < m_nsp; k++) m_flux(k,j) += C(x,k,j)*sum; } break; } @@ -448,14 +407,13 @@ namespace Cantera { s << "TITLE = \"" + title + "\"" << endl; s << "VARIABLES = \"Z (m)\"" << endl; s << "\"T (K)\"" << endl; - s << "\"phi (V)\"" << endl; for (k = 0; k < m_nsp; k++) { s << "\"" << m_thermo->speciesName(k) << "\"" << endl; } s << "ZONE T=\"c" << zone << "\"" << endl; s << " I=" << m_points << ",J=1,K=1,F=POINT" << endl; - s << "DT=(SINGLE SINGLE"; + s << "DT=(SINGLE"; for (k = 0; k < m_nsp; k++) s << " SINGLE"; s << " )" << endl; for (j = 0; j < m_points; j++) { @@ -471,9 +429,8 @@ namespace Cantera { string Solid1D::componentName(int n) const { switch(n) { case c_T_loc: return "T"; - case c_phi_loc: return "phi"; default: - if (n >= (int) 1 && n < (int) (c_Y_loc + m_nsp)) { + if (n >= (int) 1 && n < (int) (c_C_loc + m_nsp)) { return m_thermo->speciesName(n - 1); } else @@ -548,21 +505,13 @@ namespace Cantera { } else goto error; } - else if (nm == "phi") { - writelog("potential "); - if ((int) x.size() == np) { - for (j = 0; j < np; j++) - soln[index(c_phi_loc,j)] = x[j]; - } - else goto error; - } else if (m_thermo->speciesIndex(nm) >= 0) { writelog(nm+" "); if ((int) x.size() == np) { k = m_thermo->speciesIndex(nm); did_species[k] = 1; for (j = 0; j < np; j++) - soln[index(k+c_Y_loc,j)] = x[j]; + soln[index(k+c_C_loc,j)] = x[j]; } } else @@ -615,13 +564,10 @@ namespace Cantera { soln.getRow(c_T_loc,x.begin()); addFloatArray(gv,"T",x.size(),x.begin(),"K","temperature",0.0); - soln.getRow(c_phi_loc,x.begin()); - addFloatArray(gv,"phi",x.size(),x.begin(),"V","potential",0.0); - for (k = 0; k < m_nsp; k++) { - soln.getRow(c_Y_loc+k,x.begin()); + soln.getRow(c_C_loc+k,x.begin()); addFloatArray(gv,m_thermo->speciesName(k), - x.size(),x.begin(),"","massFraction",0.0,1.0); + x.size(),x.begin(),"","concentration",0.0,1.0); } } diff --git a/Cantera/src/oneD/boundaries1D.cpp b/Cantera/src/oneD/boundaries1D.cpp index 6c7977a98..917f61583 100644 --- a/Cantera/src/oneD/boundaries1D.cpp +++ b/Cantera/src/oneD/boundaries1D.cpp @@ -469,10 +469,203 @@ namespace Cantera { soln[0] = x["temperature"]; resize(1,1); } -} + + //----------------------------------------------------------- + // + // ReactingSurf1D + // + //----------------------------------------------------------- + + + + string ReactingSurf1D::componentName(int n) const { + if (n == 0) return "temperature"; + else if (n < m_nsp + 1) + return m_sphase->speciesName(n-1); + else + return ""; + } + + void ReactingSurf1D:: + init() { + _init(m_nsp+1); + m_fixed_cov.resize(m_nsp, 0.0); + m_fixed_cov[0] = 1.0; + m_work.resize(m_kin->nTotalSpecies()); + // set bounds + vector_fp lower(m_nv), upper(m_nv); + lower[0] = 200.0; + upper[0] = 1.e5; + int n; + for (n = 0; n < m_nsp; n++) { + lower[n] = -1.0e-5; + upper[n] = 2.0; + } + setBounds(m_nv, lower.begin(), m_nv, upper.begin()); + + vector_fp rtol(m_nv), atol(m_nv); + for (n = 0; n < m_nv; n++) { + rtol[n] = 1.0e-5; + atol[n] = 1.0e-9; + } + atol[0] = 1.0e-4; + setTolerances(m_nv, rtol.begin(), m_nv, atol.begin()); + } + + + void ReactingSurf1D:: + eval(int jg, doublereal* xg, doublereal* rg, + integer* diagg, doublereal rdt) { + if (jg >= 0 && (jg < firstPoint() - 2 || jg > lastPoint() + 2)) return; + + // start of local part of global arrays + doublereal* x = xg + loc(); + doublereal* r = rg + loc(); + integer* diag = diagg + loc(); + doublereal *xb, *rb; + + // specified surface temp + r[0] = x[0] - m_temp; + + // set the coverages + doublereal sum = 0.0; + int k; + for (k = 0; k < m_nsp; k++) { + m_work[k] = x[k]; + sum += x[k]; + } + m_sphase->setCoverages(m_work.begin()); + + // set the left gas state to the adjacent point + + int leftloc = 0, rightloc = 0; + int pnt = 0; + + if (m_flow_left) { + leftloc = m_flow_left->loc(); + pnt = m_flow_left->nPoints() - 1; + m_flow_left->setGas(xg + leftloc, pnt); + } + + if (m_flow_right) { + rightloc = m_flow_right->loc(); + m_flow_right->setGas(xg + rightloc, 0); + } + + m_kin->getNetProductionRates(m_work.begin()); + doublereal rs0 = 1.0/m_sphase->siteDensity(); + + //scale(m_work.begin(), m_work.end(), m_work.begin(), m_mult[0]); + + bool enabled = true; + int ioffset = m_kin->kineticsSpeciesIndex(0, m_surfindex); + + if (m_enabled) { + doublereal maxx = -1.0; + int imx = -1; + for (k = 0; k < m_nsp; k++) { + r[k] = m_work[k + ioffset] * m_sphase->size(k) * rs0; + r[k] -= rdt*(x[k] - prevSoln(k,0)); + diag[k] = 1; + if (x[k] > maxx) { + maxx = x[k]; + imx = k; + } + } + r[imx] = 1.0 - sum; + diag[imx] = 0; + } + else { + r[k] = x[k] - m_fixed_cov[k]; + diag[k] = 0; + } + + // gas-phase residuals +// doublereal rho; +// if (m_flow_left) { +// rho = m_phase_left->density(); + // doublereal rdz = 2.0/ + // (m_flow_left->z(m_left_points-1) - + // m_flow_left->z(m_left_points - 2)); + + +// for (k = 0; k < m_left_nsp; k++) +// m_work[k + m_start_left] *= m_molwt_left[k]; + +// int ileft = loc() - m_left_nv; + +// // if the energy equation is enabled at this point, +// // set the gas temperature to the surface temperature +// if (m_flow_left->doEnergy(pnt)) { +// rg[ileft + 2] = xg[ileft + 2] - m_sphase->temperature(); +// } + +// for (k = 1; k < m_left_nsp; k++) { +// if (enabled && m_flow_left->doSpecies(k)) { +// rg[ileft + 4 + k] += m_work[k + m_start_left]; +// //+= rdz*m_work[k + m_sp_left]/rho; + +// } +// } +// } + +// if (m_flow_right) { +// for (k = 0; k < m_right_nsp; k++) +// m_work[k + m_start_right] *= m_molwt_right[k]; + +// int iright = loc() + m_nsp; +// rg[iright + 2] -= m_sphase->temperature(); +// //r[iright + 3] = x[iright]; +// for (k = 0; k < m_right_nsp; k++) { +// rg[iright + 4 + k] -= m_work[k + m_start_right]; +// } +// } + + +// diag[0] = 0; +// int nc; + +// if (m_flow_right) { +// rb = r + 1; +// xb = x + 1; +// rb[2] = xb[2] - x[0]; // specified T +// } + +// if (m_flow_left) { +// nc = m_flow_left->nComponents(); +// rb = r - nc; +// xb = x - nc; +// rb[2] = xb[2] - x[0]; // specified T +// } + } + + void ReactingSurf1D:: + save(XML_Node& o, doublereal* soln) { + doublereal* s = soln + loc(); + //XML_Node& inlt = o.addChild("inlet"); + XML_Node& inlt = o.addChild("domain"); + inlt.addAttribute("id",id()); + inlt.addAttribute("points",1); + inlt.addAttribute("type","surface"); + inlt.addAttribute("components",nComponents()); + for (int k = 0; k < nComponents(); k++) { + ctml::addFloat(inlt, componentName(k), s[k], "", "",0.0, 1.0); + } + } + + void ReactingSurf1D:: + restore(XML_Node& dom, doublereal* soln) { + map x; + getFloats(dom, x); + soln[0] = x["temperature"]; + resize(1,1); + } +} + + ///////////////////////////////////////////////////////////// // // surf1D diff --git a/Cantera/src/oneD/refine.cpp b/Cantera/src/oneD/refine.cpp index fa113d5cd..4779bf56f 100644 --- a/Cantera/src/oneD/refine.cpp +++ b/Cantera/src/oneD/refine.cpp @@ -33,8 +33,8 @@ namespace Cantera { Refiner::Refiner(Domain1D& domain) : - m_ratio(10.0), m_slope(0.8), m_curve(0.8), m_min_range(0.001), - m_domain(&domain), m_npmax(200) + m_ratio(10.0), m_slope(0.8), m_curve(0.8), m_prune(0.1), + m_min_range(0.001), m_domain(&domain), m_npmax(300) { m_nv = m_domain->nComponents(); m_active.resize(m_nv, true); @@ -51,6 +51,7 @@ namespace Cantera { m_keep[0] = 1; m_keep[n-1] = 1; + //m_did_analysis = false; if (m_domain->nPoints() <= 1) return 0; @@ -58,10 +59,13 @@ namespace Cantera { m_nv = m_domain->nComponents(); // check consistency - if (n != m_domain->nPoints()) throw CanteraError("analyze","inconsistent"); + if (n != m_domain->nPoints()) + throw CanteraError("analyze","inconsistent"); - if (n >= m_npmax) throw CanteraError("analyze","max points"); + if (n >= m_npmax) { + return -2; // throw CanteraError("analyze","max points"); + } /** * find locations where cell size ratio is too large. @@ -126,12 +130,16 @@ namespace Cantera { if (r > 1.0) { m_loc[j] = 1; m_c[name] = 1; - if (int(m_loc.size()) + n > m_npmax) goto done; + //if (int(m_loc.size()) + n > m_npmax) goto done; } - if (r >= -1.0) { + if (r >= m_prune) { m_keep[j] = 1; m_keep[j+1] = 1; } + else { + if (m_keep[j] == 0) m_keep[j] = -1; + if (m_keep[j+1] == 0) m_keep[j+1] = -1; + } } } @@ -152,17 +160,21 @@ namespace Cantera { m_c[name] = 1; m_loc[j] = 1; m_loc[j+1] = 1; - if (int(m_loc.size()) + n > m_npmax) goto done; + //if (int(m_loc.size()) + n > m_npmax) goto done; } - if (r >= -1.0) { + if (r >= m_prune) { m_keep[j+1] = 1; } + else { + if (m_keep[j+1] == 0) m_keep[j+1] = -1; + } } } } } done: + //m_did_analysis = true; return m_loc.size(); } diff --git a/Cantera/src/oneD/refine.h b/Cantera/src/oneD/refine.h index ec114a25d..a9fb61a7e 100644 --- a/Cantera/src/oneD/refine.h +++ b/Cantera/src/oneD/refine.h @@ -14,8 +14,10 @@ namespace Cantera { void setCriteria(doublereal ratio = 10.0, doublereal slope = 0.8, - doublereal curve = 0.8) { + doublereal curve = 0.8, + doublereal prune = 0.1) { m_ratio = ratio; m_slope = slope; m_curve = curve; + m_prune = prune; } void setActive(int comp, bool state = true) { m_active[comp] = state; } void setMaxPoints(int npmax) { m_npmax = npmax; } @@ -28,7 +30,7 @@ namespace Cantera { return m_loc.find(j) != m_loc.end(); } bool keepPoint(int j) { - return m_keep.find(j) != m_keep.end(); + return m_keep[j] != -1; // m_keep.find(j) != m_keep.end(); } double value(const double* x, int i, int j); @@ -38,7 +40,7 @@ namespace Cantera { map m_keep; map m_c; vector m_active; - doublereal m_ratio, m_slope, m_curve; + doublereal m_ratio, m_slope, m_curve, m_prune; doublereal m_min_range; Domain1D* m_domain; int m_nv, m_npmax; diff --git a/Cantera/src/xml.cpp b/Cantera/src/xml.cpp index 2316bcfb5..b09f8f452 100755 --- a/Cantera/src/xml.cpp +++ b/Cantera/src/xml.cpp @@ -197,6 +197,7 @@ namespace Cantera { if (s[s.size()-1] == '/') { name += "/"; } + // get attributes while (1) { iloc = s.find('=');