Added in routines that modify/fix the kinetics operator for heterogeneous

reactions. Now the kinetics operator carries around information
that indicates whether a phase exists or not. If the phase exists
then the species in it can be consumed by the kinetics heterogeneous
operator. If it doesn't exist then species within the phase can
only be created by the kinetics operator.

Worked on Cantera's version control. 

In the config.h file, we now write several forms of version
control information as preprocessor defines

CANTERA_VERSION  = double quoted string id'ing the branch number
CANTERA_VERSION_MAJORNUMBER  major number of Cantera as an int

Defines that are specific to the branch and that are assigned a 1. e.g.:
#define CANTERA_VERSION_18       1
#define CANTERA_VERSION_18_XXX   1
#define CANTERA_VERSION_18_LTD   1

In the Cantera.mak file, we write


CANTERA_VERSION =  double quoted string id'ing the branch number
CANTERA_DEFINES = -DCANTERA_VERSION=@ctversion@
This commit is contained in:
Harry Moffat 2010-05-06 20:25:50 +00:00
parent e8d545fba9
commit 565e57d814
8 changed files with 245 additions and 21 deletions

View file

@ -23,6 +23,8 @@
#
in_CanteraBuildTree = 0
CANTERA_VERSION=@ctversion@
###############################################################################
# CANTERA CORE
###############################################################################
@ -122,6 +124,10 @@ endif
#
CANTERA_TOTAL_INCLUDES= $(CANTERA_CORE_INCLUDES) $(CANTERA_BOOST_INCLUDES) $(CANTERA_CVODE_INCLUDE)
#
# You can add this into the compilation environment to identify the version number
#
CANTERA_DEFINES = -DCANTERA_VERSION=@ctversion@
#
# LIBS and LIBS should be the same ...
#
CANTERA_TOTAL_LIBS2 = -L$(CANTERA_LIBSDIR) @LOCAL_LIBS@

View file

@ -23,6 +23,8 @@
#
in_CanteraBuildTree = 0
CANTERA_VERSION=@ctversion@
###############################################################################
# CANTERA CORE
###############################################################################
@ -117,6 +119,10 @@ endif
####################################################################
CANTERA_TOTAL_INCLUDES= $(CANTERA_CORE_INCLUDES) $(CANTERA_BOOST_INCLUDES) $(CANTERA_CVODE_INCLUDE)
#
# You can add this into the compilation environment to identify the version number
#
CANTERA_DEFINES = -DCANTERA_VERSION=@ctversion@
CANTERA_TOTAL_LIBS2 = @LOCAL_LIB_DIRS@ @LOCAL_LIBS@

View file

@ -52,6 +52,10 @@ namespace Cantera {
m_has_coverage_dependence(false),
m_has_electrochem_rxns(false),
m_has_exchange_current_density_formulation(false),
m_phaseExistsCheck(false),
m_phaseExists(0),
m_rxnPhaseIsReactant(0),
m_rxnPhaseIsProduct(0),
m_ioFlag(0)
{
if (thermo != 0) addPhase(*thermo);
@ -67,6 +71,10 @@ namespace Cantera {
if (m_integrator) {
delete m_integrator;
}
for (int i = 0; i < m_ii; i++) {
delete m_rxnPhaseIsReactant[i];
delete m_rxnPhaseIsProduct[i];
}
}
//====================================================================================================================
// Copy Constructor for the %InterfaceKinetics object.
@ -92,6 +100,10 @@ namespace Cantera {
m_has_coverage_dependence(false),
m_has_electrochem_rxns(false),
m_has_exchange_current_density_formulation(false),
m_phaseExistsCheck(false),
m_phaseExists(0),
m_rxnPhaseIsReactant(0),
m_rxnPhaseIsProduct(0),
m_ioFlag(0)
{
m_kdata = new InterfaceKineticsData;
@ -111,13 +123,21 @@ namespace Cantera {
*/
InterfaceKinetics& InterfaceKinetics::
operator=(const InterfaceKinetics &right) {
int i;
/*
* Check for self assignment.
*/
if (this == &right) return *this;
for (i = 0; i < m_ii; i++) {
delete (m_rxnPhaseIsReactant[i]);
delete (m_rxnPhaseIsProduct[i]);
}
Kinetics::operator=(right);
m_kk = right.m_kk;
m_revindex = right.m_revindex;
m_rates = right.m_rates;
@ -138,7 +158,6 @@ namespace Cantera {
m_E = right.m_E;
m_surf = right.m_surf; //DANGER - shallow copy
m_integrator = right.m_integrator; //DANGER - shallow copy
m_beta = right.m_beta;
m_ctrxn = right.m_ctrxn;
m_ctrxn_ecdf = right.m_ctrxn_ecdf;
@ -149,6 +168,24 @@ namespace Cantera {
m_has_coverage_dependence = right.m_has_coverage_dependence;
m_has_electrochem_rxns = right.m_has_electrochem_rxns;
m_has_exchange_current_density_formulation = right.m_has_exchange_current_density_formulation;
m_phaseExistsCheck = right.m_phaseExistsCheck;
m_phaseExists = right.m_phaseExists;
m_rxnPhaseIsReactant.resize(m_ii, 0);
m_rxnPhaseIsProduct.resize(m_ii, 0);
int np = nPhases();
for (i = 0; i < m_ii; i++) {
m_rxnPhaseIsReactant[i] = new bool[np];
m_rxnPhaseIsProduct[i] = new bool[np];
for (int p = 0; p < np; p++) {
m_rxnPhaseIsReactant[i][p] = right.m_rxnPhaseIsReactant[i][p];
m_rxnPhaseIsProduct[i][p] = right.m_rxnPhaseIsProduct[i][p];
}
}
m_rxnPhaseIsProduct = right.m_rxnPhaseIsProduct;
m_ioFlag = right.m_ioFlag;
return *this;
@ -639,7 +676,57 @@ namespace Cantera {
for (int j = 0; j != m_ii; ++j) {
ropnet[j] = ropf[j] - ropr[j];
}
}
/*
* For reactions involving multiple phases, we must check that the phase
* being consumed actually exists. This is particularly important for
* phases that are stoichiometric phases containing one species with a unity activity
*/
if (m_phaseExistsCheck) {
for (int j = 0; j != m_ii; ++j) {
if ((ropr[j] > ropf[j]) && (ropr[j] > 0.0)) {
for (int p = 0; p < nPhases(); p++) {
if (m_rxnPhaseIsProduct[j][p]) {
if (! m_phaseExists[p]) {
ropnet[j] = 0.0;
ropr[j] = ropf[j];
if (ropf[j] > 0.0) {
for (int rp = 0; rp < nPhases(); rp++) {
if (m_rxnPhaseIsReactant[j][rp]) {
if (! m_phaseExists[rp]) {
ropnet[j] = 0.0;
ropr[j] = ropf[j] = 0.0;;
}
}
}
}
}
}
}
} else if ((ropf[j] > ropr[j]) && (ropf[j] > 0.0)) {
for (int p = 0; p < nPhases(); p++) {
if (m_rxnPhaseIsReactant[j][p]) {
if (! m_phaseExists[p]) {
ropnet[j] = 0.0;
ropf[j] = ropr[j];
if (ropf[j] > 0.0) {
for (int rp = 0; rp < nPhases(); rp++) {
if (m_rxnPhaseIsProduct[j][rp]) {
if (! m_phaseExists[rp]) {
ropnet[j] = 0.0;
ropf[j] = ropr[j] = 0.0;
}
}
}
}
}
}
}
}
}
}
m_kdata->m_ROP_ok = true;
}
@ -851,11 +938,37 @@ namespace Cantera {
addElementaryReaction(r);
// operations common to all reaction types
installReagents( r );
installReagents(r);
//installGroups(reactionNumber(), r.rgroups, r.pgroups);
incrementRxnCount();
m_rxneqn.push_back(r.equation);
}
m_rxnPhaseIsReactant.resize(m_ii, 0);
m_rxnPhaseIsProduct.resize(m_ii, 0);
int np = nPhases();
int i = m_ii -1;
m_rxnPhaseIsReactant[i] = new bool[np];
m_rxnPhaseIsProduct[i] = new bool[np];
for (int p = 0; p < np; p++) {
m_rxnPhaseIsReactant[i][p] = false;
m_rxnPhaseIsProduct[i][p] = false;
}
const vector_int& vr = reactants(i);
for (int ik = 0; ik < (int) vr.size(); ik++) {
int k = vr[ik];
int p = speciesPhaseIndex(k);
m_rxnPhaseIsReactant[i][p] = true;
}
const vector_int& vp = products(i);
for (int ik = 0; ik < (int) vp.size(); ik++) {
int k = vp[ik];
int p = speciesPhaseIndex(k);
m_rxnPhaseIsProduct[i][p] = true;
}
}
//====================================================================================================================
void InterfaceKinetics::addElementaryReaction(const ReactionData& r) {
int iloc;
@ -1008,7 +1121,7 @@ namespace Cantera {
* calculates rates of species production from reaction rates of
* progress.
*/
m_rxnstoich.add( reactionNumber(), r);
m_rxnstoich.add(reactionNumber(), r);
/*
* register reaction in lists of reversible and irreversible rxns.
*/
@ -1016,11 +1129,16 @@ namespace Cantera {
m_revindex.push_back(reactionNumber());
m_nrev++;
} else {
m_irrev.push_back( reactionNumber() );
m_irrev.push_back(reactionNumber());
m_nirrev++;
}
}
//===============================================================================================
void InterfaceKinetics::addPhase(thermo_t &thermo) {
Kinetics::addPhase(thermo);
m_phaseExists.push_back(true);
}
//================================================================================================
/**
* Prepare the class for the addition of reactions. This function
* must be called after instantiation of the class, but before
@ -1044,7 +1162,7 @@ namespace Cantera {
m_pot.resize(m_kk, 0.0);
m_phi.resize(np, 0.0);
}
//================================================================================================
/**
* Finish adding reactions and prepare for use. This function
* must be called after all reactions are entered into the mechanism
@ -1071,10 +1189,13 @@ namespace Cantera {
m_deltaG0.resize(m_ii, 0.0);
m_ProdStanConcReac.resize(m_ii, 0.0);
if (m_thermo.size() != m_phaseExists.size()) {
throw CanteraError("InterfaceKinetics::finalize", "internal error");
}
m_finalized = true;
}
doublereal InterfaceKinetics::electrochem_beta(int irxn) const{
int n = m_ctrxn.size();
for (int i = 0; i < n; i++) {
@ -1089,7 +1210,7 @@ namespace Cantera {
bool InterfaceKinetics::ready() const {
return (m_finalized);
}
//================================================================================================
// Advance the surface coverages in time
/*
* @param tstep Time value to advance the surface coverages
@ -1106,7 +1227,7 @@ namespace Cantera {
delete m_integrator;
m_integrator = 0;
}
//================================================================================================
// Solve for the pseudo steady-state of the surface problem
/*
* Solve for the steady state of the surface problem.
@ -1135,7 +1256,25 @@ namespace Cantera {
*/
m_integrator->solvePseudoSteadyStateProblem(ifuncOverride, timeScaleOverride);
}
//================================================================================================
void InterfaceKinetics::setPhaseExistence(const int iphase, const bool exists) {
if (iphase < 0 || iphase >= (int) m_thermo.size()) {
throw CanteraError("InterfaceKinetics:setPhaseExistence", "out of bounds");
}
if (exists) {
if (!m_phaseExists[iphase]) {
m_phaseExistsCheck--;
m_phaseExists[iphase] = true;
}
} else {
if (m_phaseExists[iphase]) {
m_phaseExistsCheck++;
m_phaseExists[iphase] = false;
}
}
}
//================================================================================================
void EdgeKinetics::finalize() {
m_rwork.resize(nReactions());
int ks = reactionPhaseIndex();
@ -1147,7 +1286,8 @@ namespace Cantera {
"expected interface dimension = 1, but got dimension = "
+int2str(m_surf->nDim()));
m_finalized = true;
}
}
//================================================================================================
}

View file

@ -374,6 +374,19 @@ namespace Cantera {
*/
//@{
//! Add a phase to the kinetics manager object.
/*!
* This must be done before the function init() is called or
* before any reactions are input.
*
* This function calls the Kinetics operator addPhase.
* It also sets the following functions
*
* m_phaseExists[]
*
* @param thermo Reference to the ThermoPhase to be added.
*/
virtual void addPhase(thermo_t& thermo);
//! Prepare the class for the addition of reactions.
/*!
@ -403,7 +416,10 @@ namespace Cantera {
virtual bool ready() const;
//! Internal routine that updates the Rates of Progress of the reactions
/*!
* This is actually the guts of the functionality of the object
*/
void updateROP();
@ -507,6 +523,16 @@ namespace Cantera {
*/
void applyExchangeCurrentDensityFormulation(doublereal* const kfwd);
//! Set the existence of a phase in the reaction object
/*!
* Tell the kinetics object whether a phase in the object exists.
* This is actually an extrinsic specification that must be carried out on top of the
* intrinsic calculation of the reaction rate
*
* @param iphase Index of the phase. This is the order within the internal thermo vector object
* @param exists Boolean indicating whether the phase exists or not
*/
void setPhaseExistence(const int iphase, const bool exists);
protected:
@ -515,7 +541,7 @@ namespace Cantera {
//! m_kk is the number of species in all of the phases
//! that participate in this kinetics mechanism.
int m_kk;
int m_kk;
//! List of reactions numbers which are reversible reactions
/*!
@ -626,7 +652,7 @@ namespace Cantera {
*/
vector_fp m_mu0;
//! Vector of phase potentials
//! Vector of phase electric potentials
/*!
* Temporary vector containing the potential of each phase
* in the kinetics object
@ -720,6 +746,27 @@ namespace Cantera {
*/
bool m_has_exchange_current_density_formulation;
//! Int flag to indicate that some phases in the kinetics mechanism are
//! non-existent.
/*!
* We change the ROP vectors to make sure that non-existent phases are treated
* correctly in the kinetics operator. The value of this is equal to the number
* of phases which don't exist.
*/
int m_phaseExistsCheck;
//! Vector of booleans indicating whether phases exist or not
/*!
* Vector of booleans indicating whether a phase exists or not.
* We use this to set the ROP's so that unphysical things don't happen
*
* length = number of phases in the object
* By default all phases exist.
*/
std::vector<bool> m_phaseExists;
std::vector<bool *> m_rxnPhaseIsReactant;
std::vector<bool *> m_rxnPhaseIsProduct;
int m_ioFlag;
private:

View file

@ -224,7 +224,8 @@ namespace Cantera {
throw CanteraError("speciesPhase", "unknown species "+nm);
}
/**
//==============================================================================================
/*
* This function takes as an argument the kineticsSpecies index
* (i.e., the list index in the list of species in the kinetics
* manager) and returns the index of the phase owning the

View file

@ -330,6 +330,7 @@ namespace Cantera {
*/
IdealGasPhase(const IdealGasPhase &right);
//! Asignment operator
/*!
* Assignment operator for the object. Constructed

View file

@ -4,6 +4,19 @@
#ifndef CT_CONFIG_H
#define CT_CONFIG_H
//---------------------------- Version Flags ------------------//
// Cantera version -> this will be a double-quoted string value
// refering to branch number within svn
#undef CANTERA_VERSION
// Integer for major number of Cantera
#define CANTERA_VERSION_MAJORNUMBER 18
// Flag indicating it's part of major version 18
#define CANTERA_VERSION_18 1
// Flag indicating it's a development version
#define CANTERA_VERSION_18_XXX 1
// Flag indictaing that its part of 1.8_LiquidTransportDevelop branch
#define CANTERA_VERSION_18_LTD 1
//------------------------ Development flags ------------------//
//
@ -64,9 +77,6 @@ typedef int ftnlen; // Fortran hidden string length type
#undef DARWIN
#undef HAS_SSTREAM
// Cantera version
#undef CANTERA_VERSION
// Identify whether the operating system is cygwin's overlay of
// windows, with gcc being used as the compiler.
#undef CYGWIN

View file

@ -488,8 +488,21 @@ CT_SHARED_LIB=${CT_SHARED_LIB:=clib}
# lowercase 'helvetica'.
RPFONT=${RPFONT:="Helvetica"}
# Don't change this.
CANTERA_VERSION=${CANTERA_VERSION:="1.8.0"}
#
# Version Control:
#
# We define a new string for every branch created in svn.
#
# Currently expected versions:
# "1.7" Old version of Cantera
# "1.8.0" Released version of 1.8
# "1.8.x" Development version of 1.8
# "1.8_liquidTransportDevelop" branched version of 1.8 having to deal with liquid transport
#
# This field gets written into config.h and is also an autoconf variable.
#
#
CANTERA_VERSION=${CANTERA_VERSION:="1.8_liquidTransportDevelop"}
#-----------------------------------------------------------------------
#------------------- don't change anything below!! ---------------------