One-dimensional simulations. More...
#include <Sim1D.h>
One-dimensional simulations.
Class Sim1D extends class OneDim by storing the solution vector, and by adding a hybrid Newton/time-stepping solver.
Public Member Functions | |
Sim1D () | |
Default constructor. | |
Sim1D (vector< shared_ptr< Domain1D > > &domains) | |
Standard constructor. | |
void | setTimeStep (double stepsize, size_t n, const int *tsteps) |
Set the number of time steps to try when the steady Newton solver is unsuccessful. | |
void | solve (int loglevel=0, bool refine_grid=true) |
Performs the hybrid Newton steady/time-stepping solution. | |
void | eval (double rdt=-1.0, int count=1) |
void | getResidual (double rdt, double *resid) |
Evaluate the governing equations and return the vector of residuals. | |
int | refine (int loglevel=0) |
Refine the grid in all domains. | |
int | setFixedTemperature (double t) |
Add node for fixed temperature point of freely propagating flame. | |
double | fixedTemperature () |
Return temperature at the point used to fix the flame location. | |
double | fixedTemperatureLocation () |
Return location of the point where temperature is fixed. | |
void | setLeftControlPoint (double temperature) |
Set the left control point location using the specified temperature. | |
void | setRightControlPoint (double temperature) |
Set the right control point location using the specified temperature. | |
void | setRefineCriteria (int dom=-1, double ratio=10.0, double slope=0.8, double curve=0.8, double prune=-0.1) |
Set grid refinement criteria. | |
vector< double > | getRefineCriteria (int dom) |
Get the grid refinement criteria. | |
void | setMaxGridPoints (int dom, int npoints) |
Set the maximum number of grid points in the domain. | |
size_t | maxGridPoints (size_t dom) |
Get the maximum number of grid points in this domain. | |
void | setGridMin (int dom, double gridmin) |
Set the minimum grid spacing in the specified domain(s). | |
void | restoreTimeSteppingSolution () |
Set the current solution vector to the last successful time-stepping solution. | |
void | restoreSteadySolution () |
Set the current solution vector and grid to the last successful steady- state solution. | |
void | getInitialSoln () |
Get the initial value of the system state from each domain in the simulation. | |
double | jacobian (int i, int j) |
Get the Jacobian element \( J_{ij} = \partial f_i / \partial x_j \). | |
void | evalSSJacobian () |
Evaluate the Jacobian in steady-state mode. | |
void | solveAdjoint (const double *b, double *lambda) |
Solve the equation \( J^T \lambda = b \). | |
void | resize () override |
Call after one or more grids has changed size, for example after being refined. | |
void | setSteadyCallback (Func1 *callback) |
Set a function that will be called after each successful steady-state solve, before regridding. | |
Setting initial values | |
These methods are used to set the initial values of solution components. | |
void | setInitialGuess (const string &component, vector< double > &locs, vector< double > &vals) |
Set initial guess for one component for all domains. | |
void | setValue (size_t dom, size_t comp, size_t localPoint, double value) |
Set a single value in the solution vector. | |
double | value (size_t dom, size_t comp, size_t localPoint) const |
Get one entry in the solution vector. | |
double | workValue (size_t dom, size_t comp, size_t localPoint) const |
Get an entry in the work vector, which may contain either a new system state or the current residual of the system. | |
void | setProfile (size_t dom, size_t comp, const vector< double > &pos, const vector< double > &values) |
Specify a profile for one component of one domain. | |
void | setFlatProfile (size_t dom, size_t comp, double v) |
Set component 'comp' of domain 'dom' to value 'v' at all points. | |
Logging, saving and restoring of solutions | |
void | show (std::ostream &s) |
Output information on current solution for all domains to stream. | |
void | show () |
Show logging information on current solution for all domains. | |
void | save (const string &fname, const string &name, const string &desc, bool overwrite=false, int compression=0, const string &basis="") |
Save current simulation data to a container file or CSV format. | |
void | saveResidual (const string &fname, const string &name, const string &desc, bool overwrite=false, int compression=0) |
Save the residual of the current solution to a container file. | |
AnyMap | restore (const string &fname, const string &name) |
Retrieve data and settings from a previously saved simulation. | |
void | clearDebugFile () |
Deletes a debug_sim1d.yaml file if it exists. | |
void | writeDebugInfo (const string &header_suffix, const string &message, int loglevel, int attempt_counter) |
Write solver debugging information to a YAML file based on the specified log level. | |
Public Member Functions inherited from OneDim | |
OneDim () | |
Default constructor. | |
OneDim (vector< shared_ptr< Domain1D > > &domains) | |
Construct a OneDim container for the domains in the list domains. | |
OneDim (const OneDim &)=delete | |
OneDim & | operator= (const OneDim &)=delete |
void | addDomain (shared_ptr< Domain1D > d) |
Add a domain. Domains are added left-to-right. | |
MultiJac & | jacobian () |
Return a reference to the Jacobian evaluator of an OneDim object. | |
MultiNewton & | newton () |
Return a reference to the Newton iterator. | |
int | solve (double *x0, double *x1, int loglevel) |
Solve F(x) = 0, where F(x) is the multi-domain residual function. | |
size_t | nDomains () const |
Number of domains. | |
Domain1D & | domain (size_t i) const |
Return a reference to domain i. | |
size_t | domainIndex (const string &name) |
Get the index of the domain named name . | |
void | checkDomainIndex (size_t n) const |
Check that the specified domain index is in range. | |
void | checkDomainArraySize (size_t nn) const |
Check that an array size is at least nDomains(). | |
size_t | start (size_t i) const |
The index of the start of domain i in the solution vector. | |
size_t | size () const |
Total solution vector length;. | |
Domain1D * | left () |
Pointer to left-most domain (first added). | |
Domain1D * | right () |
Pointer to right-most domain (last added). | |
size_t | nVars (size_t jg) |
Number of solution components at global point jg . | |
size_t | loc (size_t jg) |
Location in the solution vector of the first component of global point jg . | |
std::tuple< string, size_t, string > | component (size_t i) |
Return the domain, local point index, and component name for the i-th component of the global solution vector. | |
size_t | bandwidth () const |
Jacobian bandwidth. | |
void | init () |
Initialize all domains. | |
size_t | points () |
Total number of points. | |
double | ssnorm (double *x, double *r) |
Steady-state max norm (infinity norm) of the residual evaluated using solution x. | |
double | rdt () const |
Reciprocal of the time step. | |
void | initTimeInteg (double dt, double *x) |
Prepare for time stepping beginning with solution x and timestep dt. | |
bool | transient () const |
True if transient mode. | |
bool | steady () const |
True if steady mode. | |
void | setSteadyMode () |
Prepare to solve the steady-state problem. | |
void | eval (size_t j, double *x, double *r, double rdt=-1.0, int count=1) |
Evaluate the multi-domain residual function. | |
Domain1D * | pointDomain (size_t i) |
Return a pointer to the domain global point i belongs to. | |
virtual void | resize () |
Call after one or more grids has changed size, for example after being refined. | |
vector< int > & | transientMask () |
Access the vector indicating which equations contain a transient term. | |
double | timeStep (int nsteps, double dt, double *x, double *r, int loglevel) |
Take time steps using Backward Euler. | |
void | resetBadValues (double *x) |
Reset values such as negative species concentrations in each domain. | |
void | writeStats (int printTime=1) |
Write statistics about the number of iterations and Jacobians at each grid level. | |
void | setJacAge (int ss_age, int ts_age=-1) |
Set the maximum number of steps that can be taken using the same Jacobian before it must be re-evaluated. | |
void | saveStats () |
Save statistics on function and Jacobian evaluation, and reset the counters. | |
void | clearStats () |
Clear saved statistics. | |
const vector< size_t > & | gridSizeStats () |
Return total grid size in each call to solve() | |
const vector< double > & | jacobianTimeStats () |
Return CPU time spent evaluating Jacobians in each call to solve() | |
const vector< double > & | evalTimeStats () |
Return CPU time spent on non-Jacobian function evaluations in each call to solve() | |
const vector< int > & | jacobianCountStats () |
Return number of Jacobian evaluations made in each call to solve() | |
const vector< int > & | evalCountStats () |
Return number of non-Jacobian function evaluations made in each call to solve() | |
const vector< int > & | timeStepStats () |
Return number of time steps taken in each call to solve() | |
void | setInterrupt (Func1 *interrupt) |
Set a function that will be called every time eval is called. | |
void | setTimeStepCallback (Func1 *callback) |
Set a function that will be called after each successful timestep. | |
void | setMinTimeStep (double tmin) |
Set the minimum time step allowed during time stepping. | |
void | setMaxTimeStep (double tmax) |
Set the maximum time step allowed during time stepping. | |
void | setTimeStepFactor (double tfactor) |
Sets a factor by which the time step is reduced if the time stepping fails. | |
void | setMaxTimeStepCount (int nmax) |
Set the maximum number of timeteps allowed before successful steady-state solve. | |
int | maxTimeStepCount () const |
Return the maximum number of timeteps allowed before successful steady-state solve. | |
Protected Attributes | |
vector< double > | m_xlast_ts |
the solution vector after the last successful timestepping | |
vector< double > | m_xlast_ss |
the solution vector after the last successful steady-state solve (stored before grid refinement) | |
vector< vector< double > > | m_grid_last_ss |
the grids for each domain after the last successful steady-state solve (stored before grid refinement) | |
vector< double > | m_xnew |
a work array used to hold the residual or the new solution | |
double | m_tstep |
timestep | |
vector< int > | m_steps |
array of number of steps to take before re-attempting the steady-state solution | |
Func1 * | m_steady_callback |
User-supplied function called after a successful steady-state solve. | |
Protected Attributes inherited from OneDim | |
double | m_tmin = 1e-16 |
minimum timestep size | |
double | m_tmax = 1e+08 |
maximum timestep size | |
double | m_tfactor = 0.5 |
factor time step is multiplied by if time stepping fails ( < 1 ) | |
shared_ptr< vector< double > > | m_state |
Solution vector. | |
unique_ptr< MultiJac > | m_jac |
Jacobian evaluator. | |
unique_ptr< MultiNewton > | m_newt |
Newton iterator. | |
double | m_rdt = 0.0 |
reciprocal of time step | |
bool | m_jac_ok = false |
if true, Jacobian is current | |
size_t | m_bw = 0 |
Jacobian bandwidth. | |
size_t | m_size = 0 |
solution vector size | |
vector< shared_ptr< Domain1D > > | m_dom |
All domains comprising the system. | |
vector< shared_ptr< Domain1D > > | m_connect |
All connector and boundary domains. | |
vector< shared_ptr< Domain1D > > | m_bulk |
All bulk/flow domains. | |
bool | m_init = false |
Indicates whether one-time initialization for each domain has been completed. | |
vector< size_t > | m_nvars |
Number of variables at each point, across all domains. | |
vector< size_t > | m_loc |
Location in the state vector of the first component of each point, across all domains. | |
vector< int > | m_mask |
Transient mask. See transientMask(). | |
size_t | m_pts = 0 |
Total number of points. | |
int | m_ss_jac_age = 20 |
Maximum age of the Jacobian in steady-state mode. | |
int | m_ts_jac_age = 20 |
Maximum age of the Jacobian in time-stepping mode. | |
Func1 * | m_interrupt = nullptr |
Function called at the start of every call to eval. | |
Func1 * | m_time_step_callback = nullptr |
User-supplied function called after each successful timestep. | |
int | m_nsteps = 0 |
Number of time steps taken in the current call to solve() | |
int | m_nsteps_max = 500 |
Maximum number of timesteps allowed per call to solve() | |
Private Member Functions | |
void | finalize () |
Calls method _finalize in each domain. | |
int | newtonSolve (int loglevel) |
Wrapper around the Newton solver. | |
Additional Inherited Members | |
Protected Member Functions inherited from OneDim | |
void | evalSSJacobian (double *x, double *rsd) |
Evaluate the steady-state Jacobian, accessible via jacobian() | |
|
inline |
Standard constructor.
domains | A vector of shared pointers to the domains to be linked together. The domain pointers must be entered in left-to-right order — that is, the pointer to the leftmost domain is domain[0], the pointer to the domain to its right is domain[1], etc. |
void setInitialGuess | ( | const string & | component, |
vector< double > & | locs, | ||
vector< double > & | vals | ||
) |
Set initial guess for one component for all domains.
component | component name |
locs | A vector of relative positions, beginning with 0.0 at the left of the domain, and ending with 1.0 at the right of the domain. |
vals | A vector of values corresponding to the relative position locations. |
void setValue | ( | size_t | dom, |
size_t | comp, | ||
size_t | localPoint, | ||
double | value | ||
) |
double value | ( | size_t | dom, |
size_t | comp, | ||
size_t | localPoint | ||
) | const |
double workValue | ( | size_t | dom, |
size_t | comp, | ||
size_t | localPoint | ||
) | const |
void setProfile | ( | size_t | dom, |
size_t | comp, | ||
const vector< double > & | pos, | ||
const vector< double > & | values | ||
) |
Specify a profile for one component of one domain.
dom | domain number, beginning with 0 for the leftmost domain. |
comp | component number |
pos | A vector of relative positions, beginning with 0.0 at the left of the domain, and ending with 1.0 at the right of the domain. |
values | A vector of values corresponding to the relative position locations. |
Note that the vector pos and values can have lengths different than the number of grid points, but their lengths must be equal. The values at the grid points will be linearly interpolated based on the (pos, values) specification.
void setFlatProfile | ( | size_t | dom, |
size_t | comp, | ||
double | v | ||
) |
void show | ( | std::ostream & | s | ) |
Output information on current solution for all domains to stream.
s | Output stream |
void show | ( | ) |
void save | ( | const string & | fname, |
const string & | name, | ||
const string & | desc, | ||
bool | overwrite = false , |
||
int | compression = 0 , |
||
const string & | basis = "" |
||
) |
Save current simulation data to a container file or CSV format.
In order to save the content of a Sim1D object, individual domains are converted to SolutionArray objects and saved using the SolutionArray::save() method. For HDF and YAML output, all domains are written to a single container file with shared header information. Simulation settings of individual domains are preserved as meta data of the corresponding SolutionArray objects. For CSV files, only state and auxiliary data of the main 1D domain are saved.
The complete state of the current object can be restored from HDF and YAML container files using the restore() method, while individual domains can be loaded using SolutionArray::restore() for further analysis. While CSV do not contain complete information, they can still be used for setting initial states of individual simulation objects for some Cantera API's.
fname | Name of output file (CSV, YAML or HDF) |
name | Identifier of storage location within the container file; this node/group contains header information and multiple subgroups holding domain-specific SolutionArray data (YAML/HDF only) |
desc | Custom comment describing the dataset to be stored (YAML/HDF only) |
overwrite | Force overwrite if file/name exists; optional (default=false) |
compression | Compression level (0-9); optional (default=0; HDF only) |
basis | Output mass ("Y"/"mass") or mole ("X"/"mole") fractions; if not specified (default=""), the native basis of the underlying ThermoPhase manager is used - |
void saveResidual | ( | const string & | fname, |
const string & | name, | ||
const string & | desc, | ||
bool | overwrite = false , |
||
int | compression = 0 |
||
) |
Save the residual of the current solution to a container file.
fname | Name of output container file |
name | Identifier of solution within the container file |
desc | Description of the solution |
overwrite | Force overwrite if name exists; optional (default=false) |
compression | Compression level (optional; HDF only) |
AnyMap restore | ( | const string & | fname, |
const string & | name | ||
) |
Retrieve data and settings from a previously saved simulation.
This method restores a simulation object from YAML or HDF data previously saved using the save() method.
fname | Name of container file (YAML or HDF) |
name | Identifier of location within the container file; this node/group contains header information and subgroups with domain-specific SolutionArray data |
void clearDebugFile | ( | ) |
void writeDebugInfo | ( | const string & | header_suffix, |
const string & | message, | ||
int | loglevel, | ||
int | attempt_counter | ||
) |
Write solver debugging information to a YAML file based on the specified log level.
This method writes solver debug information to a specified YAML file (debug_sim1d.yaml
). The section headers are formatted according to the provided header_suffix
and attempt_counter
arguments. Depending on the log level, the method will save either the solution information or the residual information for each attempted solution.
header_suffix | Header used to construct a unique section in the YAML file where the information will be written to. |
message | A string that is written to the description tag in the YAML file. |
loglevel | Controls the type of output that will be written. A loglevel greater than 6 saves the solution, and a loglevel greater than 7 saves the residual additionally. |
attempt_counter | An integer counter used to uniquely identify the attempt which is included in the file header to differentiate between multiple solution attempts. |
void setTimeStep | ( | double | stepsize, |
size_t | n, | ||
const int * | tsteps | ||
) |
Set the number of time steps to try when the steady Newton solver is unsuccessful.
stepsize | Initial time step size [s] |
n | Length of tsteps array |
tsteps | A sequence of time step counts to take after subsequent failures of the steady-state solver. The last value in tsteps will be used again after further unsuccessful solution attempts. |
void solve | ( | int | loglevel = 0 , |
bool | refine_grid = true |
||
) |
Performs the hybrid Newton steady/time-stepping solution.
The solver attempts to solve the steady-state problem first. If the steady-state solver fails, the time-stepping solver is used to take multiple time steps to move the solution closer to the steady-state solution. The steady-state solver is called again after the timesteps to make further progress towards the steady-state solution. This process is repeated until the steady-state solver converges or the maximum number of timesteps is reached.
At the end of a successful solve, if the refine_grid
flag is set, the grid will be analyzed and refined if necessary. If the grid is refined, the solution process described above is repeated with the new grid. This process is repeated until the grid no longer needs refinement based on the refine criteria.
loglevel | Controls the amount of diagnostic output. |
refine_grid | If true , the grid will be refined |
|
inline |
int refine | ( | int | loglevel = 0 | ) |
Refine the grid in all domains.
int setFixedTemperature | ( | double | t | ) |
double fixedTemperature | ( | ) |
double fixedTemperatureLocation | ( | ) |
void setLeftControlPoint | ( | double | temperature | ) |
Set the left control point location using the specified temperature.
This is used when two-point flame control is active.
The provided temperature will be used to locate the closest grid point to that temperature, which will serve to locate the left control point's coordinate. Starting from the left boundary, the first grid point that is equal to or exceeds the specified temperature will be used to locate the left control point's coordinate.
void setRightControlPoint | ( | double | temperature | ) |
Set the right control point location using the specified temperature.
This is used when two-point flame control is active.
The provided temperature will be used to locate the closest grid point to that temperature, which will serve to locate the right control point's coordinate. Starting from the right boundary, the first grid point that is equal to or exceeds the specified temperature will be used to locate the right control point's coordinate.
void setRefineCriteria | ( | int | dom = -1 , |
double | ratio = 10.0 , |
||
double | slope = 0.8 , |
||
double | curve = 0.8 , |
||
double | prune = -0.1 |
||
) |
Set grid refinement criteria.
If dom >= 0, then the settings apply only to the specified domain. If dom < 0, the settings are applied to each domain.
vector< double > getRefineCriteria | ( | int | dom | ) |
Get the grid refinement criteria.
dom must be greater than or equal to zero (that is, the domain must be specified).
void setMaxGridPoints | ( | int | dom, |
int | npoints | ||
) |
Set the maximum number of grid points in the domain.
If dom >= 0, then the settings apply only to the specified domain. If dom < 0, the settings are applied to each domain.
size_t maxGridPoints | ( | size_t | dom | ) |
Get the maximum number of grid points in this domain.
dom | domain number, beginning with 0 for the leftmost domain. |
void setGridMin | ( | int | dom, |
double | gridmin | ||
) |
void restoreTimeSteppingSolution | ( | ) |
void restoreSteadySolution | ( | ) |
void getInitialSoln | ( | ) |
double jacobian | ( | int | i, |
int | j | ||
) |
void evalSSJacobian | ( | ) |
void solveAdjoint | ( | const double * | b, |
double * | lambda | ||
) |
Solve the equation \( J^T \lambda = b \).
Here, \( J = \partial f/\partial x \) is the Jacobian matrix of the system of equations \( f(x,p)=0 \). This can be used to efficiently solve for the sensitivities of a scalar objective function \( g(x,p) \) to a vector of parameters \( p \) by solving:
\[ J^T \lambda = \left( \frac{\partial g}{\partial x} \right)^T \]
for \( \lambda \) and then computing:
\[ \left.\frac{dg}{dp}\right|_{f=0} = \frac{\partial g}{\partial p} - \lambda^T \frac{\partial f}{\partial p} \]
|
overridevirtual |
|
inline |
|
private |
|
private |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |