 |
Cantera
3.1.0a1
|
Method to solve a pseudo steady state surface problem. More...
#include <solveSP.h>
Method to solve a pseudo steady state surface problem.
The following class handles solving the surface problem. The calculation uses Newton's method to obtain the surface fractions of the surface and bulk species by requiring that the surface species production rate = 0 and that the either the bulk fractions are proportional to their production rates or they are constants.
Currently, the bulk mole fractions are treated as constants. Implementation of their being added to the unknown solution vector is delayed.
Lets introduce the unknown vector for the "surface problem". The surface problem is defined as the evaluation of the surface site fractions for multiple surface phases. The unknown vector will consist of the vector of surface concentrations for each species in each surface vector. Species are grouped first by their surface phases
The unknown solution vector is defined as follows:
| C_i_j | kindexSP |
|---|---|
| C_0_0 | 0 |
| C_1_0 | 1 |
| C_2_0 | 2 |
| . . . | ... |
| C_N0-1_0 | N0-1 |
| C_0_1 | N0 |
| C_1_1 | N0+1 |
| C_2_1 | N0+2 |
| . . . | ... |
| C_N1-1_1 | NO+N1-1 |
Note there are a couple of different types of species indices floating around in the formulation of this object.
kindexSP: This is the species index in the contiguous vector of unknowns for the surface problem.
Note, in the future, BULK_DEPOSITION systems will be added, and the solveSP unknown vector will get more complicated. It will include the mole fraction and growth rates of specified bulk phases
Indices which relate to individual kinetics objects use the suffix KSI (kinetics species index).
This routine is typically used within a residual calculation in a large code. It's typically invoked millions of times for large calculations, and it must work every time. Therefore, requirements demand that it be robust but also efficient.
The solution methodology is largely determined by the ifunc parameter, that is input to the solution object. This parameter may have one of the values defined in Surface Problem Solver Methods.
The time step is determined from sdot[], so that the time step doesn't ever change the value of a variable by more than 100%.
This algorithm does use a damped Newton's method to relax the equations. Damping is based on a "delta damping" technique. The solution unknowns are not allowed to vary too much between iterations.
EXTRA_ACCURACY: A constant that is the ratio of the required update norm in this Newton iteration compared to that in the nonlinear solver. A value of 0.1 is used so surface species are safely overconverged.
Functions called:
ct_dgetrf – First half of LAPACK direct solve of a full Matrixct_dgetrs – Second half of LAPACK direct solve of a full matrix. Returns solution vector in the right-hand-side vector, resid. Public Member Functions | |
| solveSP (ImplicitSurfChem *surfChemPtr, int bulkFunc=BULK_ETCH) | |
| Constructor for the object. | |
| ~solveSP ()=default | |
| Destructor. | |
| solveSP (const solveSP &)=delete | |
| solveSP & | operator= (const solveSP &)=delete |
| int | solveSurfProb (int ifunc, double time_scale, double TKelvin, double PGas, double reltol, double abstol) |
| Main routine that actually calculates the pseudo steady state of the surface problem. | |
Public Attributes | |
| int | m_ioflag = 0 |
Private Member Functions | |
| void | print_header (int ioflag, int ifunc, double time_scale, int damping, double reltol, double abstol) |
| Printing routine that optionally gets called at the start of every invocation. | |
| void | printIteration (int ioflag, double damp, int label_d, int label_t, double inv_t, double t_real, size_t iter, double update_norm, double resid_norm, bool do_time, bool final=false) |
| Printing routine that gets called after every iteration. | |
| double | calc_t (double netProdRateSolnSP[], double XMolSolnSP[], int *label, int *label_old, double *label_factor, int ioflag) |
| Calculate a conservative delta T to use in a pseudo-steady state algorithm. | |
| void | calcWeights (double wtSpecies[], double wtResid[], const Array2D &Jac, const double CSolnSP[], const double abstol, const double reltol) |
| Calculate the solution and residual weights. | |
| void | updateState (const double *cSurfSpec) |
| Update the surface states of the surface phases. | |
| void | updateMFSolnSP (double *XMolSolnSP) |
| Update mole fraction vector consisting of unknowns in surface problem. | |
| void | updateMFKinSpecies (double *XMolKinSp, int isp) |
| Update the mole fraction vector for a specific kinetic species vector corresponding to one InterfaceKinetics object. | |
| void | evalSurfLarge (const double *CSolnSP) |
| Update the vector that keeps track of the largest species in each surface phase. | |
| void | fun_eval (double *resid, const double *CSolnSP, const double *CSolnOldSP, const bool do_time, const double deltaT) |
| Main Function evaluation. | |
| void | resjac_eval (DenseMatrix &jac, double *resid, double *CSolnSP, const double *CSolnSPOld, const bool do_time, const double deltaT) |
| Main routine that calculates the current residual and Jacobian. | |
Private Attributes | |
| vector< InterfaceKinetics * > & | m_objects |
| Vector of interface kinetics objects. | |
| size_t | m_neq = 0 |
| Total number of equations to solve in the implicit problem. | |
| int | m_bulkFunc |
| This variable determines how the bulk phases are to be handled. | |
| size_t | m_numSurfPhases = 0 |
| Number of surface phases in the surface problem. | |
| size_t | m_numTotSurfSpecies = 0 |
| Total number of surface species in all surface phases. | |
| vector< size_t > | m_indexKinObjSurfPhase |
| Mapping between the surface phases and the InterfaceKinetics objects. | |
| vector< size_t > | m_nSpeciesSurfPhase |
| Vector of length number of surface phases containing the number of surface species in each phase. | |
| vector< SurfPhase * > | m_ptrsSurfPhase |
| Vector of surface phase pointers. | |
| vector< size_t > | m_eqnIndexStartSolnPhase |
| Index of the start of the unknowns for each solution phase. | |
| size_t | m_numBulkPhasesSS = 0 |
| Total number of volumetric condensed phases included in the steady state problem handled by this routine. | |
| vector< size_t > | m_numBulkSpecies |
| Vector of number of species in the m_numBulkPhases phases. | |
| size_t | m_numTotBulkSpeciesSS = 0 |
| Total number of species in all bulk phases. | |
| vector< ThermoPhase * > | m_bulkPhasePtrs |
| Vector of bulk phase pointers, length is equal to m_numBulkPhases. | |
| vector< size_t > | m_kinSpecIndex |
| Index between the equation index and the position in the kinetic species array for the appropriate kinetics operator. | |
| vector< size_t > | m_kinObjIndex |
| Index between the equation index and the index of the InterfaceKinetics object. | |
| vector< size_t > | m_spSurfLarge |
| Vector containing the indices of the largest species in each surface phase. | |
| size_t | m_maxTotSpecies = 0 |
| Maximum number of species in any single kinetics operator -> also maxed wrt the total # of solution species. | |
| vector< double > | m_netProductionRatesSave |
| Temporary vector with length equal to max m_maxTotSpecies. | |
| vector< double > | m_numEqn1 |
| Temporary vector with length equal to max m_maxTotSpecies. | |
| vector< double > | m_numEqn2 |
| Temporary vector with length equal to max m_maxTotSpecies. | |
| vector< double > | m_CSolnSave |
| Temporary vector with length equal to max m_maxTotSpecies. | |
| vector< double > | m_CSolnSP |
| Solution vector. length MAX(1, m_neq) | |
| vector< double > | m_CSolnSPInit |
| Saved initial solution vector. length MAX(1, m_neq) | |
| vector< double > | m_CSolnSPOld |
| Saved solution vector at the old time step. length MAX(1, m_neq) | |
| vector< double > | m_wtResid |
| Weights for the residual norm calculation. length MAX(1, m_neq) | |
| vector< double > | m_wtSpecies |
| Weights for the species concentrations norm calculation. | |
| vector< double > | m_resid |
| Residual for the surface problem. | |
| vector< double > | m_XMolKinSpecies |
| Vector of mole fractions. length m_maxTotSpecies. | |
| DenseMatrix | m_Jac |
| Jacobian. | |
| solveSP | ( | ImplicitSurfChem * | surfChemPtr, |
| int | bulkFunc = BULK_ETCH |
||
| ) |
Constructor for the object.
| surfChemPtr | Pointer to the ImplicitSurfChem object that defines the surface problem to be solved. |
| bulkFunc | Integer representing how the bulk phases should be handled. See Surface Problem Bulk Phase Mode. Currently, only the default value of BULK_ETCH is supported. |
Definition at line 22 of file solveSP.cpp.
|
default |
Destructor.
| int solveSurfProb | ( | int | ifunc, |
| double | time_scale, | ||
| double | TKelvin, | ||
| double | PGas, | ||
| double | reltol, | ||
| double | abstol | ||
| ) |
Main routine that actually calculates the pseudo steady state of the surface problem.
The actual converged solution is returned as part of the internal state of the InterfaceKinetics objects.
Uses Newton's method to get the surface fractions of the surface and bulk species by requiring that the surface species production rate = 0 and that the bulk fractions are proportional to their production rates.
| ifunc | Determines the type of solution algorithm to be used. See Surface Problem Solver Methods for possible values. |
| time_scale | Time over which to integrate the surface equations, where applicable |
| TKelvin | Temperature (kelvin) |
| PGas | Pressure (pascals) |
| reltol | Relative tolerance to use |
| abstol | absolute tolerance. |
Definition at line 89 of file solveSP.cpp.
|
private |
Printing routine that optionally gets called at the start of every invocation.
Definition at line 602 of file solveSP.cpp.
|
private |
Printing routine that gets called after every iteration.
Definition at line 651 of file solveSP.cpp.
|
private |
Calculate a conservative delta T to use in a pseudo-steady state algorithm.
This routine calculates a pretty conservative 1/del_t based on MAX_i(sdot_i/(X_i*SDen0)). This probably guarantees diagonal dominance.
Small surface fractions are allowed to intervene in the del_t determination, no matter how small. This may be changed. Now minimum changed to 1.0e-12,
Maximum time step set to time_scale.
| netProdRateSolnSP | Output variable. Net production rate of all of the species in the solution vector. |
| XMolSolnSP | output variable. Mole fraction of all of the species in the solution vector |
| label | Output variable. Pointer to the value of the species index (kindexSP) that is controlling the time step |
| label_old | Output variable. Pointer to the value of the species index (kindexSP) that controlled the time step at the previous iteration |
| label_factor | Output variable. Pointer to the current factor that is used to indicate the same species is controlling the time step. |
| ioflag | Level of the output requested. |
Definition at line 563 of file solveSP.cpp.
|
private |
Calculate the solution and residual weights.
| wtSpecies | Weights to use for the soln unknowns. These are in concentration units |
| wtResid | Weights to sue for the residual unknowns. |
| Jac | Jacobian. Row sum scaling is used for the Jacobian |
| CSolnSP | Solution vector for the surface problem |
| abstol | Absolute error tolerance |
| reltol | Relative error tolerance |
Definition at line 530 of file solveSP.cpp.
|
private |
Update the surface states of the surface phases.
Definition at line 272 of file solveSP.cpp.
|
private |
Update mole fraction vector consisting of unknowns in surface problem.
| XMolSolnSP | Vector of mole fractions for the unknowns in the surface problem. |
Definition at line 286 of file solveSP.cpp.
|
private |
Update the mole fraction vector for a specific kinetic species vector corresponding to one InterfaceKinetics object.
| XMolKinSp | Mole fraction vector corresponding to a particular kinetic species for a single InterfaceKinetics Object This is a vector over all the species in all of the phases in the InterfaceKinetics object |
| isp | ID of the InterfaceKinetics Object. |
Definition at line 294 of file solveSP.cpp.
|
private |
Update the vector that keeps track of the largest species in each surface phase.
| CSolnSP | Vector of the current values of the surface concentrations in all of the surface species. |
Definition at line 303 of file solveSP.cpp.
|
private |
Main Function evaluation.
| resid | output Vector of residuals, length = m_neq |
| CSolnSP | Vector of species concentrations, unknowns in the problem, length = m_neq |
| CSolnOldSP | Old Vector of species concentrations, unknowns in the problem, length = m_neq |
| do_time | Calculate a time dependent residual |
| deltaT | Delta time for time dependent problem. |
Definition at line 319 of file solveSP.cpp.
|
private |
Main routine that calculates the current residual and Jacobian.
| jac | Jacobian to be evaluated. |
| resid | output Vector of residuals, length = m_neq |
| CSolnSP | Vector of species concentrations, unknowns in the problem, length = m_neq. These are tweaked in order to derive the columns of the Jacobian. |
| CSolnSPOld | Old Vector of species concentrations, unknowns in the problem, length = m_neq |
| do_time | Calculate a time dependent residual |
| deltaT | Delta time for time dependent problem. |
Definition at line 418 of file solveSP.cpp.
|
private |
|
private |
|
private |
This variable determines how the bulk phases are to be handled.
Possible values are given in Surface Problem Bulk Phase Mode.
|
private |
Number of surface phases in the surface problem.
This number is equal to the number of InterfaceKinetics objects in the problem. (until further noted)
|
private |
|
private |
Mapping between the surface phases and the InterfaceKinetics objects.
Currently this is defined to be a 1-1 mapping (and probably assumed in some places) m_surfKinObjID[i] = i
|
private |
|
private |
|
private |
Index of the start of the unknowns for each solution phase.
i_eqn = m_eqnIndexStartPhase[isp]
isp is the phase id in the list of phases solved by the surface problem.
i_eqn is the equation number of the first unknown in the solution vector corresponding to isp'th phase.
|
private |
Total number of volumetric condensed phases included in the steady state problem handled by this routine.
This is equal to or less than the total number of volumetric phases in all of the InterfaceKinetics objects. We usually do not include bulk phases. Bulk phases are only included in the calculation when their domain isn't included in the underlying continuum model conservation equation system.
This is equal to 0, for the time being
|
private |
|
private |
|
private |
|
private |
|
private |
Index between the equation index and the index of the InterfaceKinetics object.
Length m_neq
|
private |
Vector containing the indices of the largest species in each surface phase.
k = m_spSurfLarge[i] where k is the local species index, that is, it varies from 0 to (num species in phase - 1) and i is the surface phase index in the problem. Length is equal to m_numSurfPhases.
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
Residual for the surface problem.
The residual vector of length "dim" that, that has the value of "sdot" for surface species. The residuals for the bulk species are a function of the sdots for all species in the bulk phase. The last residual of each phase enforces {Sum(fractions) = 1}. After linear solve (dgetrf_ & dgetrs_), resid holds the update vector.
length MAX(1, m_neq)
|
private |
|
private |
1.9.7