*** empty log message ***

This commit is contained in:
Dave Goodwin 2004-03-12 05:59:53 +00:00
parent 2351fc7fea
commit d298e3b4b0
22 changed files with 237 additions and 109 deletions

View file

@ -539,21 +539,23 @@ extern "C" {
catch (CanteraError) { return DERR; }
}
int DLL_EXPORT th_setState_satLiquid(int n) {
int DLL_EXPORT th_setState_Psat(int n, double p, double x) {
try {
purefluid(n)->setState_satLiquid();
purefluid(n)->setState_Psat(p, x);
return 0;
}
catch (CanteraError) { return -1; }
}
int DLL_EXPORT th_setState_satVapor(int n) {
int DLL_EXPORT th_setState_Tsat(int n, double t, double x) {
try {
purefluid(n)->setState_satVapor();
purefluid(n)->setState_Tsat(t, x);
return 0;
}
catch (CanteraError) { return -1; }
}
//-------------- Kinetics ------------------//

View file

@ -83,8 +83,8 @@ extern "C" {
double DLL_IMPORT th_vaporFraction(int n);
double DLL_IMPORT th_satTemperature(int n, double p);
double DLL_IMPORT th_satPressure(int n, double t);
int DLL_IMPORT th_setState_satLiquid(int n);
int DLL_IMPORT th_setState_satVapor(int n);
int DLL_IMPORT th_setState_Psat(int n, double p, double x);
int DLL_IMPORT th_setState_Tsat(int n, double t, double x);
#endif
int DLL_IMPORT newKineticsFromXML(int mxml, int iphase,

View file

@ -4,12 +4,12 @@ function s = GRI30(tr)
%
% GRI-Mech 3.0 is a widely-used reaction mechanism for natural gas
% combustion. It contains 53 species composed of the elements H,
% C, O, N, and/or Ar, and contains 325 reactions, most of which
% are reversible. GRI-Mech 3.0, like most combustion mechanisms,
% is designed for use at pressures where the ideal gas law holds.
% C, O, N, and/or Ar, and 325 reactions, most of which are
% reversible. GRI-Mech 3.0, like most combustion mechanisms, is
% designed for use at pressures where the ideal gas law holds.
%
% Function GRI30 creates the solution according to the
% specifications in file gri30.xml. The ideal gas equation of
% specifications in file gri30.cti. The ideal gas equation of
% state is used. Transport property evaluation is disabled by
% default. To enable transport properties, supply the name of
% the transport model to use.

View file

@ -7,5 +7,8 @@ function n = Hydrogen()
% equation of state is taken from W. C. Reynolds, "Thermodynamic
% Properties in SI."
%
% For more details, see classes Cantera::PureFluid and tpx::hydrogen in the
% Cantera C++ source code documentation.
%
n = importPhase('purefluids.cti','hydrogen');

View file

@ -9,7 +9,7 @@ function m = MassFlowController(upstream, downstream)
% mass flow controller that maintains a constant mass flow rate
% independent of upstream or downstream conditions. If two reactor
% objects are supplied as arguments, the controller is installed
% between the two reactors.
% between the two reactors.
%
% see also: FlowDevice
%

View file

@ -7,5 +7,8 @@ function w = Water()
% equation of state is taken from W. C. Reynolds, "Thermodynamic
% Properties in SI."
%
% For more details, see classes Cantera::PureFluid and tpx::water in the
% Cantera C++ source code documentation.
%
w = importPhase('purefluids.cti','water');

View file

@ -29,9 +29,9 @@ static void thermoset( int nlhs, mxArray *plhs[],
case 1:
ierr = th_setPressure(th,*ptr); break;
case 2:
ierr = th_setState_satLiquid(th); break;
ierr = th_setState_Psat(th,ptr[0],ptr[1]); break;
case 3:
ierr = th_setState_satVapor(th); break;
ierr = th_setState_Tsat(th,ptr[0],ptr[1]); break;
default:
mexErrMsgTxt("unknown attribute.");
}
@ -62,7 +62,8 @@ static void thermoset( int nlhs, mxArray *plhs[],
else if (job == 50) {
int xy = int(*ptr);
ierr = th_equil(th, xy);
}
}
if (ierr < 0) reportError();
plhs[0] = mxCreateNumericMatrix(1,1,mxDOUBLE_CLASS,mxREAL);
double *h = mxGetPr(plhs[0]);

View file

@ -262,11 +262,11 @@ class ThermoPhase(Phase):
"""Vapor fraction."""
return _cantera.thermo_getfp(self._phase_id,53)
def setState_satLiquid(self):
_cantera.thermo_setfp(self._phase_id,7,0.0,0.0)
def setState_Psat(self, p, vaporFraction):
_cantera.thermo_setfp(self._phase_id,8, p, vaporFraction)
def setState_satVapor(self):
_cantera.thermo_setfp(self._phase_id,8,0.0,0.0)
def setState_Tsat(self, t, vaporFraction):
_cantera.thermo_setfp(self._phase_id,7, t, vaporFraction)
def thermophase(self):

View file

@ -1012,9 +1012,15 @@ class ideal_gas(phase):
def is_ideal_gas(self):
return 1
class pure_solid(phase):
"""A pure solid."""
class stoichiometric_solid(phase):
"""A solid compound or pure element.Stoichiometric solid phases
contain exactly one species, which always has unit activity. The
solid is assumed to have constant density. Therefore the rates of
reactions involving these phases do not contain any concentration
terms for the (one) species in the phase, since the concentration
is always the same. """
def __init__(self,
name = '',
elements = '',
@ -1038,7 +1044,7 @@ class pure_solid(phase):
def build(self, p):
ph = phase.build(self, p)
e = ph.addChild("thermo")
e['model'] = 'SolidCompound'
e['model'] = 'StoichCompound'
addFloat(e, 'density', self._dens, defunits = _umass+'/'+_ulen+'3')
if self._tr:
t = ph.addChild('transport')
@ -1047,6 +1053,50 @@ class pure_solid(phase):
k['model'] = 'none'
class stoichiometric_liquid(stoichiometric_solid):
"""A stoichiometric liquid. Currently, there is no distinction
between stoichiometric liquids and solids."""
def __init__(self,
name = '',
elements = '',
species = '',
density = -1.0,
transport = 'None',
initial_state = None,
options = []):
stoichiometric_solid.__init__(self, name, 3, elements,
species, 'none',
initial_state, options)
self._dens = density
self._pure = 1
if self._dens < 0.0:
raise 'density must be specified.'
self._tr = transport
class pure_solid(stoichiometric_solid):
"""Deprecated. Use stoichiometric_solid"""
def __init__(self,
name = '',
elements = '',
species = '',
density = -1.0,
transport = 'None',
initial_state = None,
options = []):
stoichiometric_solid.__init__(self, name, 3, elements,
species, 'none',
initial_state, options)
self._dens = density
self._pure = 1
if self._dens < 0.0:
raise 'density must be specified.'
self._tr = transport
print 'WARNING: entry type pure_solid is deprecated.'
print 'Use stoichiometric_solid instead.'
class metal(phase):
"""A metal."""
def __init__(self,
@ -1062,8 +1112,6 @@ class metal(phase):
initial_state, options)
self._dens = density
self._pure = 0
#if self._dens < 0.0:
# raise 'density must be specified.'
self._tr = transport
def conc_dim(self):
@ -1115,8 +1163,13 @@ class incompressible_solid(phase):
k['model'] = 'none'
class pure_fluid(phase):
"""A pure fluid."""
class liquid_vapor(phase):
"""A fluid with a complete liquid/vapor equation of state.
This entry type selects one of a set of predefined fluids with
built-in liquid/vapor equations of state. The substance_flag
parameter selects the fluid. See purefluids.py for the usage
of this entry type."""
def __init__(self,
name = '',
elements = '',
@ -1310,7 +1363,10 @@ validate()
# $Revision$
# $Date$
# $Log$
# Revision 1.30 2004-02-08 13:25:21 dggoodwin
# Revision 1.31 2004-03-12 05:59:59 dggoodwin
# *** empty log message ***
#
# Revision 1.30 2004/02/08 13:25:21 dggoodwin
# *** empty log message ***
#
# Revision 1.29 2004/02/08 13:22:31 dggoodwin

View file

@ -1,12 +1,11 @@
#
# A CVD example. This example computes the growth rate of a diamond film according to
# a simplified version of a particular published growth mechanism (see file diamond.cti
# for details). Only the surface coverage equations are solved here; the gas composition
# is fixed. (For an example of coupled gas-phase and surface, see catcomb.py.)
#
# Atomic hydrogen plays an important role in diamond CVD, and this
# example computes the growth rate and surface coverages as a function
# of [H] at the surface for fixed temperature and [CH3].
# A CVD example. This example computes the growth rate of a diamond
#film according to a simplified version of a particular published
#growth mechanism (see file diamond.cti for details). Only the surface
#coverage equations are solved here; the gas composition is
#fixed. (For an example of coupled gas-phase and surface, see
#catcomb.py.) Atomic hydrogen plays an important role in diamond CVD,
#and this example computes the growth rate and surface coverages as a
#function of [H] at the surface for fixed temperature and [CH3].
from Cantera import *
import math

View file

@ -1,37 +1,83 @@
#
# an ideal Rankine cycle
# an Rankine cycle
#
from Cantera import *
from Cantera.pureFluids import Water
from Cantera.liquidvapor import Water
# parameters
eta_pump = 0.6 # pump isentropic efficiency
et_turbine = 0.8 # turbine isentropic efficiency
pmax = 8.0e5 # maximum pressure
########################################################
#
# some useful functions
#
def pump(fluid, pfinal, eta):
"""Adiabatically pump a fluid to pressure pfinal, using
a pump with isentropic efficiency eta."""
h0 = fluid.enthalpy_mass()
s0 = fluid.entropy_mass()
fluid.setState_SP(s0, pfinal)
h1s = fluid.enthalpy_mass()
isentropic_work = h1s - h0
actual_work = isentropic_work / eta
h1 = h0 + actual_work
fluid.setState_HP(h1, pfinal)
return actual_work
def expand(fluid, pfinal, eta):
"""Adiabatically expand a fluid to pressure pfinal, using
a turbine with isentropic efficiency eta."""
h0 = fluid.enthalpy_mass()
s0 = fluid.entropy_mass()
fluid.setState_SP(s0, pfinal)
h1s = fluid.enthalpy_mass()
isentropic_work = h0 - h1s
actual_work = isentropic_work * eta
h1 = h0 - actual_work
fluid.setState_HP(h1, pfinal)
return actual_work
###############################################################
# create an object representing water
w = Water()
# start with saturated liquid water at 300 K
w.setTemperature(300.0)
w.setState_satLiquid()
h1 = w.enthalpy_mass()
p1 = w.pressure()
w.setState_Tsat(0.0)
hf = w.enthalpy_mass()
print w
w.setState_Tsat(1.0)
hv = w.enthalpy_mass()
print hv - hf
# pump it isentropically to 10 MPa
w.setState_SP(w.entropy_mass(), 1.0e7)
h2 = w.enthalpy_mass()
print w
pump_work = h2 - h1
# pump it adiabatically to pmax
pump_work = pump(w, pmax, eta_pump)
print pump_work
# heat at constant pressure to 1500 K
w.setState_TP(1500.0, w.pressure())
h3 = w.enthalpy_mass()
# heat it at constant pressure until it reaches the
# saturated vapor state at this pressure
#w.setState_Psat(1.0)
#print w
heat_in = h3 - h2
# expand isentropically back to 300 K
w.setState_SP(w.entropy_mass(), p1)
h4 = w.enthalpy_mass()
work_out = h3 - h4
heat_out = h4 - h1
efficiency = (work_out - pump_work)/heat_in
print 'efficiency = ',efficiency
w.setTemperature(273.16)
w.setState_Tsat(0.0)
h0 = w.enthalpy_mass()
for t in [300.0, 350.0, 400.0, 450.0, 500.0]:
w.setTemperature(t)
w.setState_Tsat(0.0)
hf = w.enthalpy_mass()
w.setState_Tsat(1.0)
hv = w.enthalpy_mass()
print t, 0.001*(hf - h0), 0.001*(hv - h0), 0.001*(hv - hf)
for t in [750.0, 800.0, 850.0, 1150.0]:
w.setState_TP(t, 2.0e4)
print t, w.enthalpy_mass() - h0

View file

@ -159,9 +159,9 @@ thermo_setfp(PyObject *self, PyObject *args)
case 6:
iok = th_setElectricPotential(th, v[0]); break;
case 7:
iok = th_setState_satLiquid(th); break;
iok = th_setState_Tsat(th, v1, v2); break;
case 8:
iok = th_setState_satVapor(th); break;
iok = th_setState_Psat(th, v1, v2); break;
default:
iok = -10;
}

View file

@ -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 SolidCompound.o SpeciesThermoFactory.o ThermoFactory.o
THERMO = $(BASE) ThermoPhase.o IdealGasPhase.o ConstDensityThermo.o StoichSubstance.o SpeciesThermoFactory.o ThermoFactory.o
# homogeneous kinetics
KINETICS = GRI_30_Kinetics.o KineticsFactory.o GasKinetics.o FalloffFactory.o \

View file

@ -38,6 +38,10 @@ namespace Cantera {
}
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()) {
@ -197,7 +201,11 @@ namespace Cantera {
virtual doublereal satPressure(doublereal t) const {
doublereal tsv = m_sub->Temp();
doublereal vsv = m_sub->v();
m_sub->Set(tpx::TP, t, 0.5*m_sub->Pcrit());
if (t < 0.0)
m_sub->Set(tpx::TP, temperature(), 0.5*m_sub->Pcrit());
else
m_sub->Set(tpx::TP, t, 0.5*m_sub->Pcrit());
doublereal ps = m_sub->Ps();
m_sub->Set(tpx::TV,tsv,vsv);
check();
@ -211,19 +219,20 @@ namespace Cantera {
return x;
}
virtual void setState_satLiquid() {
virtual void setState_Tsat(doublereal t, doublereal x) {
setTemperature(t);
setTPXState();
m_sub->Set(tpx::TX, temperature(), 0.0);
m_sub->Set(tpx::TX, t, x);
setDensity(1.0/m_sub->v());
check();
}
virtual void setState_satVapor() {
virtual void setState_Psat(doublereal p, doublereal x) {
setTPXState();
m_sub->Set(tpx::TX, temperature(), 1.0);
setDensity(1.0/m_sub->v());
check();
}
m_sub->Set(tpx::PX, p, x);
setTemperature(m_sub->Temp());
setDensity(1.0/m_sub->v());
}
protected:

View file

@ -1,6 +1,6 @@
/**
*
* @file IdealGasPhase.cpp
* @file StoichSubstance.cpp
*
*/
@ -11,16 +11,16 @@
#include "ct_defs.h"
#include "mix_defs.h"
#include "SolidCompound.h"
#include "StoichSubstance.h"
#include "SpeciesThermo.h"
namespace Cantera {
void SolidCompound::initThermo() {
void StoichSubstance::initThermo() {
m_kk = nSpecies();
if (m_kk > 1) {
throw CanteraError("initThermo",
"solid compounds may only contain one species.");
"stoichiometric substances may only contain one species.");
}
doublereal tmin = m_spthermo->minTemp();
doublereal tmax = m_spthermo->maxTemp();
@ -35,7 +35,7 @@ namespace Cantera {
}
void SolidCompound::_updateThermo() const {
void StoichSubstance::_updateThermo() const {
doublereal tnow = temperature();
if (m_tlast != tnow) {
m_spthermo->update(tnow, m_cp0_R.begin(), m_h0_RT.begin(),

View file

@ -223,6 +223,7 @@ namespace Cantera {
bool ready() const { return (m_kk > 0); }
protected:
/**
@ -240,6 +241,11 @@ namespace Cantera {
*/
int m_kk;
void setMolecularWeight(int k, double mw) {
m_molwts[k] = mw;
m_rmolwts[k] = 1.0/mw;
}
private:
/**

View file

@ -24,7 +24,8 @@
#include "SurfPhase.h"
#include "EdgePhase.h"
#include "MetalPhase.h"
#include "SolidCompound.h"
//#include "SolidCompound.h"
#include "StoichSubstance.h"
#include "importCTML.h"
namespace Cantera {
@ -33,11 +34,11 @@ namespace Cantera {
static int ntypes = 7;
static string _types[] = {"IdealGas", "Incompressible",
"Surface", "Edge", "Metal", "SolidCompound",
"Surface", "Edge", "Metal", "StoichSubstance",
"PureFluid"};
static int _itypes[] = {cIdealGas, cIncompressible,
cSurf, cEdge, cMetal, cSolidCompound,
cSurf, cEdge, cMetal, cStoichSubstance,
cPureFluid};
ThermoPhase* ThermoFactory::newThermoPhase(string model) {
@ -72,8 +73,8 @@ namespace Cantera {
th = new MetalPhase;
break;
case cSolidCompound:
th = new SolidCompound;
case cStoichSubstance:
th = new StoichSubstance;
break;
#ifdef INCL_PURE_FLUIDS

View file

@ -96,7 +96,7 @@ namespace Cantera {
doublereal tol) {
doublereal dt;
setPressure(p);
for (int n = 0; n < 20; n++) {
for (int n = 0; n < 50; n++) {
dt = (h - enthalpy_mass())/cp_mass();
if (dt > 100.0) dt = 100.0;
else if (dt < -100.0) dt = -100.0;
@ -112,7 +112,7 @@ namespace Cantera {
doublereal tol) {
doublereal dt;
setDensity(1.0/v);
for (int n = 0; n < 20; n++) {
for (int n = 0; n < 50; n++) {
dt = (u - intEnergy_mass())/cv_mass();
if (dt > 100.0) dt = 100.0;
else if (dt < -100.0) dt = -100.0;
@ -128,7 +128,7 @@ namespace Cantera {
doublereal tol) {
doublereal dt;
setPressure(p);
for (int n = 0; n < 20; n++) {
for (int n = 0; n < 50; n++) {
dt = (s - entropy_mass())*temperature()/cp_mass();
if (dt > 100.0) dt = 100.0;
else if (dt < -100.0) dt = -100.0;
@ -144,10 +144,8 @@ namespace Cantera {
doublereal tol) {
doublereal dt;
setDensity(1.0/v);
for (int n = 0; n < 20; n++) {
cout << "n = " << n << endl;
for (int n = 0; n < 50; n++) {
dt = (s - entropy_mass())*temperature()/cv_mass();
cout << "dt = " << dt << endl;
if (dt > 100.0) dt = 100.0;
else if (dt < -100.0) dt = -100.0;
setTemperature(temperature() + dt);

View file

@ -453,6 +453,9 @@ namespace Cantera {
* Specific entropy. Units: J/kg/K.
*/
doublereal entropy_mass() const {
//cout << "entropy_mass. " << endl;
//cout << "entropy_mole = " << entropy_mole() << endl;
//cout << "meanMolecularWeight = " << meanMolecularWeight() << endl;
return entropy_mole()/meanMolecularWeight();
}
@ -601,14 +604,15 @@ namespace Cantera {
err("vaprFraction"); return -1.0;
}
virtual void setState_satLiquid() {
err("setState_satLiquid");
virtual void setState_Tsat(doublereal t, doublereal x) {
err("setState_sat");
}
virtual void setState_satVapor() {
err("setState_satVapor");
virtual void setState_Psat(doublereal p, doublereal x) {
err("setState_sat");
}
/**
* @internal Install a species thermodynamic property
* manager. The species thermodynamic property manager

View file

@ -42,7 +42,8 @@
using namespace ctml;
//#include <stdio.h>
//#include <stdio.h>
// these are all used to check for duplicate reactions
vector< map<int, doublereal> > _reactiondata;
@ -818,8 +819,8 @@ namespace Cantera {
eoserror = true;
}
}
else if (eos["model"] == "SolidCompound") {
if (th->eosType() == cSolidCompound) {
else if (eos["model"] == "StoichSubstance") {
if (th->eosType() == cStoichSubstance) {
doublereal rho = getFloat(eos, "density", "-");
th->setDensity(rho);
}
@ -855,16 +856,9 @@ namespace Cantera {
else if (eos["model"] == "PureFluid") {
if (th->eosType() == cPureFluid) {
subflag = atoi(eos["fluid_type"].c_str());
//doublereal h0 = getFloat(eos, "h0", "-");
//doublereal s0 = getFloat(eos, "s0", "-");
if (subflag < 0)
throw CanteraError("importCTML",
"missing fluid type flag");
//doublereal c[3];
//c[0] = doublereal(subflag);
//c[1] = h0;
//c[2] = s0;
//th->setParameters(3, c);
}
else {
eoserror = true;

View file

@ -34,7 +34,8 @@ namespace Cantera {
const int cIncompressible = 2; // ConstDensityThermo in ConstDensityThermo.h
const int cSurf = 3; // SurfPhase in SurfPhase.h
const int cMetal = 4; // MetalPhase in MetalPhase.h
const int cSolidCompound = 5; // SolidCompound in SolidCompound.h
// const int cSolidCompound = 5; // SolidCompound in SolidCompound.h
const int cStoichSubstance = 5; // StoichSubstance.h
// pure fluids with liquid/vapor eqs of state
const int cPureFluid = 10;

View file

@ -108,13 +108,13 @@ namespace tpx {
// absolute tolerances
double TolAbsH = 0.01; // J/kg
double TolAbsU = 0.01;
double TolAbsS = 1.e-5;
double TolAbsH = 0.0001; // J/kg
double TolAbsU = 0.0001;
double TolAbsS = 1.e-6;
double TolAbsP = 0.000; // Pa
double TolAbsV = 1.e-7;
double TolAbsT = 1.e-3;
double TolRel = 3.e-6;
double TolRel = 3.e-7;
void Substance::Set(int XY, double x0, double y0) {
double temp;
@ -324,6 +324,8 @@ namespace tpx {
}
double Substance::vprop(int ijob) {
//cout << "vprop: T, Rho = " << T << " " << Rho << endl;
//cout << "entropy = " << sp() << endl;
switch (ijob) {
case EvalH: return hp();
case EvalS: return sp();
@ -413,6 +415,9 @@ namespace tpx {
y_here = prop(ify);
err_x = fabs(X - x_here);
err_y = fabs(Y - y_here);
//cout << x_here << " " << y_here << endl;
//cout << err_x << " " << err_y << endl;
//cout << X << " " << Y << endl;
if ((err_x < atx + rtx*Xa) && (err_y < aty + rty*Ya)) break;