[Thermo] Make modification of discontinuous NASA polynomials optional

Adding the attribute allow_discontinuities="true" to a phase's
'thermo' node (CTML) or adding 'allow_discontinuous_thermo' to a
phase's 'options' argument (CTI) will disable the automatic adjustment
of the provided polynomials so that enthalpy and entropy are
continuous at the midpoint temperature.
This commit is contained in:
Ray Speth 2013-12-09 01:34:06 +00:00
parent 6b68564245
commit 3e4ccd4027
5 changed files with 36 additions and 22 deletions

View file

@ -269,6 +269,9 @@ Option String Meaning
``'skip_undeclared_third_bodies'`` When importing reactions with third body efficiencies,
ignore any efficiencies for undeclared species, rather
than flagging them as an error.
``'allow_discontinuous_thermo'`` Disable the automatic adjustment of NASA polynomials to
eliminate discontinuities in enthalpy and entropy at the
midpoint temperature.
================================== ================
Using the ``options`` field, it is possible to extract a sub-mechanism from a large

View file

@ -156,7 +156,7 @@ class SpeciesThermo
public:
//! Constructor
SpeciesThermo() {}
SpeciesThermo() : m_allow_discontinuities(false) {}
//! Destructor
virtual ~SpeciesThermo() {}
@ -332,6 +332,8 @@ public:
*/
virtual void modifyOneHf298(const int k, const doublereal Hf298New) = 0;
#endif
bool m_allow_discontinuities;
};
//@}
}

View file

@ -1810,6 +1810,10 @@ class phase(object):
if self._comment:
ph.addChild('note',self._comment)
thermo = ph.addChild('thermo')
if 'allow_discontinuous_thermo' in self._options:
thermo['allow_discontinuities'] = 'true'
return ph
@ -1854,8 +1858,7 @@ class ideal_gas(phase):
def build(self, p):
ph = phase.build(self, p)
e = ph.addChild("thermo")
e['model'] = 'IdealGas'
ph.child('thermo')['model'] = 'IdealGas'
k = ph.addChild("kinetics")
k['model'] = self._kin
t = ph.addChild('transport')
@ -1900,7 +1903,7 @@ class stoichiometric_solid(phase):
def build(self, p):
ph = phase.build(self, p)
e = ph.addChild("thermo")
e = ph.child('thermo')
e['model'] = 'StoichSubstance'
addFloat(e, 'density', self._dens, defunits = _umass+'/'+_ulen+'3')
if self._tr:
@ -1955,7 +1958,7 @@ class metal(phase):
def build(self, p):
ph = phase.build(self, p)
e = ph.addChild("thermo")
e = ph.child("thermo")
e['model'] = 'Metal'
addFloat(e, 'density', self._dens, defunits = _umass+'/'+_ulen+'3')
if self._tr:
@ -1993,7 +1996,7 @@ class semiconductor(phase):
def build(self, p):
ph = phase.build(self, p)
e = ph.addChild("thermo")
e = ph.child("thermo")
e['model'] = 'Semiconductor'
addFloat(e, 'density', self._dens, defunits = _umass+'/'+_ulen+'3')
addFloat(e, 'effectiveMass_e', self._emass, defunits = _umass)
@ -2031,7 +2034,7 @@ class incompressible_solid(phase):
def build(self, p):
ph = phase.build(self, p)
e = ph.addChild("thermo")
e = ph.child("thermo")
e['model'] = 'Incompressible'
addFloat(e, 'density', self._dens, defunits = _umass+'/'+_ulen+'3')
if self._tr:
@ -2071,7 +2074,7 @@ class lattice(phase):
#if visible == 0:
# return
ph = phase.build(self, p)
e = ph.addChild('thermo')
e = ph.child('thermo')
e['model'] = 'Lattice'
addFloat(e, 'site_density', self._n, defunits = _umol+'/'+_ulen+'3')
if self._vac:
@ -2128,7 +2131,7 @@ class lattice_solid(phase):
def build(self, p):
ph = phase.build(self, p)
e = ph.addChild("thermo")
e = ph.child("thermo")
e['model'] = 'LatticeSolid'
if self._lattices:
@ -2171,7 +2174,7 @@ class liquid_vapor(phase):
def build(self, p):
ph = phase.build(self, p)
e = ph.addChild("thermo")
e = ph.child("thermo")
e['model'] = 'PureFluid'
e['fluid_type'] = repr(self._subflag)
k = ph.addChild("kinetics")
@ -2202,7 +2205,7 @@ class RedlichKwongMFTP(phase):
def build(self, p):
ph = phase.build(self,p)
e = ph.addChild("thermo")
e = ph.child("thermo")
e['model'] = 'RedlichKwongMFTP'
if self._activityCoefficients:
a = e.addChild("activityCoefficients")
@ -2250,7 +2253,7 @@ class redlich_kwong(phase):
def build(self, p):
ph = phase.build(self, p)
e = ph.addChild("thermo")
e = ph.child("thermo")
e['model'] = 'PureFluid'
e['fluid_type'] = repr(self._subflag)
addFloat(e, 'Tc', self._tc, defunits = "K")
@ -2300,7 +2303,7 @@ class ideal_interface(phase):
def build(self, p):
ph = phase.build(self, p)
e = ph.addChild("thermo")
e = ph.child("thermo")
e['model'] = 'Surface'
addFloat(e, 'site_density', self._sitedens, defunits = _umol+'/'+_ulen+'2')
k = ph.addChild("kinetics")
@ -2340,7 +2343,7 @@ class edge(phase):
def build(self, p):
ph = phase.build(self, p)
e = ph.addChild("thermo")
e = ph.child("thermo")
e['model'] = 'Edge'
addFloat(e, 'site_density', self._sitedens, defunits = _umol+'/'+_ulen)
k = ph.addChild("kinetics")
@ -2476,7 +2479,7 @@ class edge(phase):
## def build(self, p):
## ph = phase.build(self, p)
## e = ph.addChild("thermo")
## e = ph.child("thermo")
## sc = e.addChild("standardConc")
## sc['model'] = self._stdconc
## e['model'] = 'HMW'

View file

@ -82,12 +82,14 @@ void NasaThermo::install(const std::string& name, size_t index, int type,
vector_fp chigh(c+8, c+15);
vector_fp clow(c+1, c+8);
checkContinuity(name, tmid, &clow[0], &chigh[0]);
//if (maxError > 1e-6) {
// fixDiscontinuities(tlow, tmid, thigh, &clow[0], &chigh[0]);
// AssertThrowMsg(checkContinuity(name, tmid, &clow[0], &chigh[0]) < 1e-12,
// "NasaThermo::install", "Polynomials still not continuous");
//}
if (!m_allow_discontinuities) {
doublereal maxError = checkContinuity(name, tmid, &clow[0], &chigh[0]);
if (maxError > 1e-6) {
fixDiscontinuities(tlow, tmid, thigh, &clow[0], &chigh[0]);
AssertThrowMsg(checkContinuity(name, tmid, &clow[0], &chigh[0]) < 1e-12,
"NasaThermo::install", "Polynomials still not continuous");
}
}
m_high[igrp-1].push_back(NasaPoly1(index, tmid, thigh,
ref_pressure, &chigh[0]));

View file

@ -434,8 +434,8 @@ bool importPhase(XML_Node& phase, ThermoPhase* th,
// Set equation of state parameters. The parameters are
// specific to each subclass of ThermoPhase, so this is done
// by method setParametersFromXML in each subclass.
const XML_Node& eos = phase.child("thermo");
if (phase.hasChild("thermo")) {
const XML_Node& eos = phase.child("thermo");
th->setParametersFromXML(eos);
} else {
throw CanteraError("importPhase",
@ -561,6 +561,10 @@ bool importPhase(XML_Node& phase, ThermoPhase* th,
// used, and selects a class that can handle the
// parameterizations found.
spth = newSpeciesThermoMgr(spDataNodeList);
if (eos["allow_discontinuities"] == "true") {
std::cout << "ALLOWING DISCONTINUOUS THERMO!" << std::endl;
spth->m_allow_discontinuities = true;
}
// install it in the phase object
th->setSpeciesThermo(spth);