Go to the documentation of this file.
20 ImplicitSurfChem::ImplicitSurfChem(
21 vector<InterfaceKinetics*> k,
double rtol,
double atol,
22 double maxStepSize,
size_t maxSteps,
23 size_t maxErrTestFails) :
25 m_numTotalBulkSpecies(0),
29 m_maxstep(maxStepSize),
31 m_maxErrTestFails(maxErrTestFails),
32 m_mediumSpeciesStart(-1),
33 m_bulkSpeciesStart(-1),
34 m_surfSpeciesStart(-1),
35 m_commonTempPressForPhases(true),
39 size_t kinSpIndex = 0;
41 for (
size_t n = 0; n < k.size(); n++) {
44 size_t ns = k[n]->surfacePhaseIndex();
47 "kinetics manager contains no surface phase");
51 size_t nsp =
m_surf.back()->nSpecies();
54 size_t nt = k[n]->nTotalSpecies();
55 ntmax = std::max(nt, ntmax);
56 m_specStartIndex.push_back(kinSpIndex);
58 size_t nPhases = kinPtr->
nPhases();
61 for (
size_t ip = 0; ip < nPhases; ip++) {
67 m_numTotalBulkSpecies += nsp;
70 pLocTmp[ip] = int(imatch);
72 pLocTmp[ip] = -int(n);
75 pLocVec.push_back(pLocTmp);
77 m_numTotalSpecies =
m_nv + m_numTotalBulkSpecies;
79 m_concSpeciesSave.resize(m_numTotalSpecies, 0.0);
81 m_integ.reset(newIntegrator(
"CVODE"));
86 m_integ->setProblemType(DENSE + NOJAC);
90 int ImplicitSurfChem::checkMatch(std::vector<ThermoPhase*> m_vec,
ThermoPhase* thPtr)
93 for (
int i = 0; i < (int) m_vec.size(); i++) {
105 for (
size_t n = 0; n <
m_surf.size(); n++) {
106 m_surf[n]->getCoverages(c + loc);
123 m_integ->setTolerances(m_rtol, m_atol);
144 m_integ->initialize(t0, *
this);
152 m_integ->setMaxStepSize(t1 - t0);
167 for (
size_t n = 0; n <
m_surf.size(); n++) {
168 m_surf[n]->setCoverages(c + loc);
174 doublereal* ydot, doublereal* p)
178 for (
size_t n = 0; n <
m_surf.size(); n++) {
179 double rs0 = 1.0/
m_surf[n]->siteDensity();
183 for (
size_t k = 1; k <
m_nsp[n]; k++) {
184 ydot[k + loc] = m_work[kstart + k] * rs0 *
m_surf[n]->size(k);
193 doublereal timeScaleOverride)
199 doublereal time_scale = timeScaleOverride;
209 if (ifuncOverride >= 0) {
210 ifunc = ifuncOverride;
231 doublereal reltol = 1.0E-6;
232 doublereal atol = 1.0E-20;
237 for (
size_t k = 0; k <
m_nv; k++) {
252 int retn =
m_surfSolver->solveSurfProb(ifunc, time_scale, TKelvin, PGas,
259 retn =
m_surfSolver->solveSurfProb(ifunc, time_scale, TKelvin, PGas,
262 throw CanteraError(
"ImplicitSurfChem::solvePseudoSteadyStateProblem",
263 "solveSP return an error condition!");
271 for (
size_t ip = 0; ip <
m_surf.size(); ip++) {
273 kstart = m_specStartIndex[ip];
287 for (
size_t ip = 0; ip <
m_surf.size(); ip++) {
289 kstart = m_specStartIndex[ip];
302 for (
size_t ip = 0; ip <
m_surf.size(); ip++) {
void solvePseudoSteadyStateProblem(int ifuncOverride=-1, doublereal timeScaleOverride=1.0)
Solve for the pseudo steady-state of the surface problem.
void setCommonState_TP(doublereal TKelvin, doublereal PresPa)
Sets the state variable in all thermodynamic phases (surface and surrounding bulk phases) to the inpu...
size_t m_nv
Total number of surface species in all surface phases.
std::vector< InterfaceKinetics * > m_vecKinPtrs
vector of pointers to InterfaceKinetics objects
size_t nPhases() const
The number of phases participating in the reaction mechanism.
std::vector< size_t > m_surfindex
index of the surface phase in each InterfaceKinetics object
std::vector< SurfPhase * > m_surf
vector of pointers to surface phases.
size_t reactionPhaseIndex() const
Phase where the reactions occur.
virtual void setState_TP(doublereal t, doublereal p)
Set the temperature (K) and pressure (Pa)
void getConcentrations(double *const c) const
Get the species concentrations (kmol/m^3).
virtual void setMaxErrTestFails(size_t maxErrTestFails=7)
int m_ioFlag
Controls the amount of printing from this routine and underlying routines.
const int BULK_ETCH
Etching of a bulk phase is to be expected.
virtual void eval(doublereal t, doublereal *y, doublereal *ydot, doublereal *p)
Evaluate the value of ydot[k] at the current conditions.
@ BDF_Method
Backward Differentiation.
std::unique_ptr< solveSP > m_surfSolver
Pointer to the helper method, Placid, which solves the surface problem.
void integrate0(doublereal t0, doublereal t1)
Integrate from t0 to t1 without reinitializing the integrator.
void setConcSpecies(const doublereal *const vecConcSpecies)
Sets the concentrations within phases that are unknowns in the surface problem.
virtual void setTolerances(double rtol=1.e-7, double atol=1.e-14)
const int SFLUX_RESIDUAL
Need to solve the surface problem in order to calculate the surface fluxes of gas-phase species.
const int SFLUX_INITIALIZE
This assumes that the initial guess supplied to the routine is far from the correct one.
Base class for a phase with thermodynamic properties.
std::unique_ptr< Integrator > m_integ
Pointer to the CVODE integrator.
A kinetics manager for heterogeneous reaction mechanisms.
std::vector< ThermoPhase * > m_bulkPhases
Vector of pointers to bulk phases.
virtual double pressure() const
Return the thermodynamic pressure (Pa).
virtual void initialize(doublereal t0=0.0)
doublereal temperature() const
Temperature (K).
doublereal m_maxstep
max step size
std::vector< size_t > m_nsp
Vector of number of species in each Surface Phase.
size_t nSpecies() const
Returns the number of species in the phase.
Method to solve a pseudo steady state surface problem.
A simple thermodynamic model for a surface phase, assuming an ideal solution model.
void getConcSpecies(doublereal *const vecConcSpecies) const
vector_fp m_concSpecies
Temporary vector - length num species in the Kinetics object.
virtual void getState(doublereal *y)
Get the current state of the solution vector.
thermo_t & thermo(size_t n=0)
This method returns a reference to the nth ThermoPhase object defined in this kinetics mechanism.
void integrate(doublereal t0, doublereal t1)
Integrate from t0 to t1. The integrator is reinitialized first.
size_t m_maxErrTestFails
maximum number of error test failures allowed
virtual void setMaxSteps(size_t maxsteps=20000)
std::vector< int > vector_int
Vector of ints.
Base class for exceptions thrown by Cantera classes.
void updateState(doublereal *y)
Set the mixture to a state consistent with solution vector y.
size_t m_nmax
maximum number of steps allowed
const size_t npos
index returned by functions to indicate "no position"
virtual void setConcentrations(const double *const conc)
Set the concentrations to the specified values within the phase.
virtual void setMaxStepSize(double maxstep=0.0)
Namespace for the Cantera kernel.
bool m_commonTempPressForPhases
If true, a common temperature and pressure for all surface and bulk phases associated with the surfac...