Cantera 2.6.0
Public Member Functions | Protected Attributes | Private Member Functions | List of all members
Sim1D Class Reference

One-dimensional simulations. More...

#include <Sim1D.h>

Inheritance diagram for Sim1D:
[legend]
Collaboration diagram for Sim1D:
[legend]

Public Member Functions

 Sim1D ()
 Default constructor. More...
 
 Sim1D (std::vector< Domain1D * > &domains)
 Standard constructor. More...
 
void save (const std::string &fname, const std::string &id, const std::string &desc, int loglevel=1)
 
void saveResidual (const std::string &fname, const std::string &id, const std::string &desc, int loglevel=1)
 
void showSolution (std::ostream &s)
 Print to stream s the current solution for all domains. More...
 
void showSolution ()
 
const doublereal * solution ()
 
void setTimeStep (double stepsize, size_t n, const int *tsteps)
 
void solve (int loglevel=0, bool refine_grid=true)
 
void eval (doublereal rdt=-1.0, int count=1)
 
void getResidual (double rdt, double *resid)
 
int refine (int loglevel=0)
 Refine the grid in all domains. More...
 
int setFixedTemperature (double t)
 Add node for fixed temperature point of freely propagating flame. More...
 
double fixedTemperature ()
 Return temperature at the point used to fix the flame location. More...
 
double fixedTemperatureLocation ()
 Return location of the point where temperature is fixed. More...
 
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. More...
 
vector_fp getRefineCriteria (int dom)
 Get the grid refinement criteria. More...
 
void setMaxGridPoints (int dom, int npoints)
 Set the maximum number of grid points in the domain. More...
 
size_t maxGridPoints (size_t dom)
 Get the maximum number of grid points in this domain. More...
 
void setGridMin (int dom, double gridmin)
 Set the minimum grid spacing in the specified domain(s). More...
 
void restore (const std::string &fname, const std::string &id, int loglevel=2)
 Initialize the solution with a previously-saved solution. More...
 
void restoreTimeSteppingSolution ()
 Set the current solution vector to the last successful time-stepping solution. More...
 
void restoreSteadySolution ()
 Set the current solution vector and grid to the last successful steady- state solution. More...
 
void getInitialSoln ()
 
void setSolution (const doublereal *soln)
 
const doublereal * solution () const
 
doublereal jacobian (int i, int j)
 
void evalSSJacobian ()
 
void solveAdjoint (const double *b, double *lambda)
 Solve the equation \( J^T \lambda = b \). More...
 
virtual void resize ()
 Call after one or more grids has changed size, for example after being refined. More...
 
void setSteadyCallback (Func1 *callback)
 Set a function that will be called after each successful steady-state solve, before regridding. More...
 
Setting initial values

These methods are used to set the initial values of solution components.

void setInitialGuess (const std::string &component, vector_fp &locs, vector_fp &vals)
 Set initial guess for one component for all domains. More...
 
void setValue (size_t dom, size_t comp, size_t localPoint, doublereal value)
 Set a single value in the solution vector. More...
 
doublereal value (size_t dom, size_t comp, size_t localPoint) const
 Get one entry in the solution vector. More...
 
doublereal workValue (size_t dom, size_t comp, size_t localPoint) const
 
void setProfile (size_t dom, size_t comp, const vector_fp &pos, const vector_fp &values)
 Specify a profile for one component of one domain. More...
 
void setFlatProfile (size_t dom, size_t comp, doublereal v)
 Set component 'comp' of domain 'dom' to value 'v' at all points. More...
 
- Public Member Functions inherited from OneDim
 OneDim (std::vector< Domain1D * > domains)
 Construct a OneDim container for the domains in the list domains. More...
 
 OneDim (const OneDim &)=delete
 
OneDimoperator= (const OneDim &)=delete
 
void addDomain (Domain1D *d)
 Add a domain. Domains are added left-to-right. More...
 
MultiJacjacobian ()
 Return a reference to the Jacobian evaluator. More...
 
MultiNewtonnewton ()
 Return a reference to the Newton iterator. More...
 
int solve (doublereal *x0, doublereal *x1, int loglevel)
 Solve F(x) = 0, where F(x) is the multi-domain residual function. More...
 
size_t nDomains () const
 Number of domains. More...
 
Domain1Ddomain (size_t i) const
 Return a reference to domain i. More...
 
size_t domainIndex (const std::string &name)
 
void checkDomainIndex (size_t n) const
 Check that the specified domain index is in range. More...
 
void checkDomainArraySize (size_t nn) const
 Check that an array size is at least nDomains(). More...
 
size_t start (size_t i) const
 The index of the start of domain i in the solution vector. More...
 
size_t size () const
 Total solution vector length;. More...
 
Domain1Dleft ()
 Pointer to left-most domain (first added). More...
 
Domain1Dright ()
 Pointer to right-most domain (last added). More...
 
size_t nVars (size_t jg)
 Number of solution components at global point jg. More...
 
size_t loc (size_t jg)
 Location in the solution vector of the first component of global point jg. More...
 
std::tuple< std::string, size_t, std::string > component (size_t i)
 Return the domain, local point index, and component name for the i-th component of the global solution vector. More...
 
size_t bandwidth () const
 Jacobian bandwidth. More...
 
void init ()
 
size_t points ()
 Total number of points. More...
 
doublereal ssnorm (doublereal *x, doublereal *r)
 Steady-state max norm (infinity norm) of the residual evaluated using solution x. More...
 
doublereal rdt () const
 Reciprocal of the time step. More...
 
void initTimeInteg (doublereal dt, doublereal *x)
 Prepare for time stepping beginning with solution x and timestep dt. More...
 
bool transient () const
 True if transient mode. More...
 
bool steady () const
 True if steady mode. More...
 
void setSteadyMode ()
 
void eval (size_t j, double *x, double *r, doublereal rdt=-1.0, int count=1)
 Evaluate the multi-domain residual function. More...
 
Domain1DpointDomain (size_t i)
 Return a pointer to the domain global point i belongs to. More...
 
vector_inttransientMask ()
 
double timeStep (int nsteps, double dt, double *x, double *r, int loglevel)
 
void resetBadValues (double *x)
 
void writeStats (int printTime=1)
 Write statistics about the number of iterations and Jacobians at each grid level. More...
 
void save (const std::string &fname, std::string id, const std::string &desc, doublereal *sol, int loglevel)
 
AnyMap serialize (const double *soln) const
 
void setMinTimeStep (doublereal tmin)
 
void setMaxTimeStep (doublereal tmax)
 
void setTimeStepFactor (doublereal tfactor)
 
void setMaxTimeStepCount (int nmax)
 Set the maximum number of timeteps allowed before successful steady-state solve. More...
 
int maxTimeStepCount () const
 Return the maximum number of timeteps allowed before successful steady-state solve. More...
 
void setJacAge (int ss_age, int ts_age=-1)
 
void saveStats ()
 Save statistics on function and Jacobian evaluation, and reset the counters. More...
 
void clearStats ()
 Clear saved statistics. More...
 
const std::vector< size_t > & gridSizeStats ()
 Return total grid size in each call to solve() More...
 
const vector_fpjacobianTimeStats ()
 Return CPU time spent evaluating Jacobians in each call to solve() More...
 
const vector_fpevalTimeStats ()
 Return CPU time spent on non-Jacobian function evaluations in each call to solve() More...
 
const vector_intjacobianCountStats ()
 Return number of Jacobian evaluations made in each call to solve() More...
 
const vector_intevalCountStats ()
 Return number of non-Jacobian function evaluations made in each call to solve() More...
 
const vector_inttimeStepStats ()
 Return number of time steps taken in each call to solve() More...
 
void setInterrupt (Func1 *interrupt)
 Set a function that will be called every time eval is called. More...
 
void setTimeStepCallback (Func1 *callback)
 Set a function that will be called after each successful timestep. More...
 

Protected Attributes

vector_fp m_x
 the solution vector More...
 
vector_fp m_xlast_ts
 the solution vector after the last successful timestepping More...
 
vector_fp m_xlast_ss
 the solution vector after the last successful steady-state solve (stored before grid refinement) More...
 
std::vector< vector_fpm_grid_last_ss
 the grids for each domain after the last successful steady-state solve (stored before grid refinement) More...
 
vector_fp m_xnew
 a work array used to hold the residual or the new solution More...
 
doublereal m_tstep
 timestep More...
 
vector_int m_steps
 array of number of steps to take before re-attempting the steady-state solution More...
 
Func1m_steady_callback
 User-supplied function called after a successful steady-state solve. More...
 
- Protected Attributes inherited from OneDim
doublereal m_tmin
 minimum timestep size More...
 
doublereal m_tmax
 maximum timestep size More...
 
doublereal m_tfactor
 factor time step is multiplied by if time stepping fails ( < 1 ) More...
 
std::unique_ptr< MultiJacm_jac
 Jacobian evaluator. More...
 
std::unique_ptr< MultiNewtonm_newt
 Newton iterator. More...
 
doublereal m_rdt
 reciprocal of time step More...
 
bool m_jac_ok
 if true, Jacobian is current More...
 
size_t m_bw
 Jacobian bandwidth. More...
 
size_t m_size
 solution vector size More...
 
std::vector< Domain1D * > m_dom
 
std::vector< Domain1D * > m_connect
 
std::vector< Domain1D * > m_bulk
 
bool m_init
 
std::vector< size_t > m_nvars
 
std::vector< size_t > m_loc
 
vector_int m_mask
 
size_t m_pts
 
int m_ss_jac_age
 
int m_ts_jac_age
 
Func1m_interrupt
 Function called at the start of every call to eval. More...
 
Func1m_time_step_callback
 User-supplied function called after each successful timestep. More...
 
int m_nsteps
 Number of time steps taken in the current call to solve() More...
 
int m_nsteps_max
 Maximum number of timesteps allowed per call to solve() More...
 

Private Member Functions

void finalize ()
 Calls method _finalize in each domain. More...
 
int newtonSolve (int loglevel)
 Wrapper around the Newton solver. More...
 

Additional Inherited Members

- Protected Member Functions inherited from OneDim
void evalSSJacobian (doublereal *x, doublereal *xnew)
 

Detailed Description

One-dimensional simulations.

Class Sim1D extends class OneDim by storing the solution vector, and by adding a hybrid Newton/time-stepping solver.

Definition at line 21 of file Sim1D.h.

Constructor & Destructor Documentation

◆ Sim1D() [1/2]

Sim1D ( )
inline

Default constructor.

This constructor is provided to make the class default-constructible, but is not meant to be used in most applications. Use the next constructor

Definition at line 29 of file Sim1D.h.

◆ Sim1D() [2/2]

Sim1D ( std::vector< Domain1D * > &  domains)

Standard constructor.

Parameters
domainsA vector of 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.

Definition at line 25 of file Sim1D.cpp.

References Domain1D::_getInitialSoln(), OneDim::domain(), Sim1D::m_steps, Sim1D::m_tstep, Sim1D::m_x, OneDim::nDomains(), Sim1D::resize(), and OneDim::start().

Member Function Documentation

◆ setInitialGuess()

void setInitialGuess ( const std::string &  component,
vector_fp locs,
vector_fp vals 
)

Set initial guess for one component for all domains.

Parameters
componentcomponent name
locsA 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.
valsA vector of values corresponding to the relative position locations.

Definition at line 41 of file Sim1D.cpp.

References OneDim::component(), Domain1D::componentName(), OneDim::domain(), Domain1D::nComponents(), OneDim::nDomains(), and Sim1D::setProfile().

◆ setValue()

void setValue ( size_t  dom,
size_t  comp,
size_t  localPoint,
doublereal  value 
)

Set a single value in the solution vector.

Parameters
domdomain number, beginning with 0 for the leftmost domain.
compcomponent number
localPointgrid point within the domain, beginning with 0 for the leftmost grid point in the domain.
valuethe value.

Definition at line 54 of file Sim1D.cpp.

References AssertThrowMsg, OneDim::domain(), Domain1D::loc(), Sim1D::m_x, and Sim1D::value().

Referenced by Sim1D::setFlatProfile(), and Sim1D::setProfile().

◆ value()

doublereal value ( size_t  dom,
size_t  comp,
size_t  localPoint 
) const

Get one entry in the solution vector.

Parameters
domdomain number, beginning with 0 for the leftmost domain.
compcomponent number
localPointgrid point within the domain, beginning with 0 for the leftmost grid point in the domain.

Definition at line 62 of file Sim1D.cpp.

References AssertThrowMsg, OneDim::domain(), Domain1D::loc(), and Sim1D::m_x.

Referenced by Sim1D::refine(), and Sim1D::setValue().

◆ workValue()

doublereal workValue ( size_t  dom,
size_t  comp,
size_t  localPoint 
) const

Definition at line 70 of file Sim1D.cpp.

◆ setProfile()

void setProfile ( size_t  dom,
size_t  comp,
const vector_fp pos,
const vector_fp values 
)

Specify a profile for one component of one domain.

Parameters
domdomain number, beginning with 0 for the leftmost domain.
compcomponent number
posA 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.
valuesA 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.

Definition at line 78 of file Sim1D.cpp.

References OneDim::domain(), Cantera::linearInterp(), Domain1D::nPoints(), and Sim1D::setValue().

Referenced by Sim1D::setInitialGuess().

◆ setFlatProfile()

void setFlatProfile ( size_t  dom,
size_t  comp,
doublereal  v 
)

Set component 'comp' of domain 'dom' to value 'v' at all points.

Definition at line 220 of file Sim1D.cpp.

References OneDim::domain(), Domain1D::nPoints(), and Sim1D::setValue().

◆ save()

void save ( const std::string &  fname,
const std::string &  id,
const std::string &  desc,
int  loglevel = 1 
)

Definition at line 97 of file Sim1D.cpp.

◆ saveResidual()

void saveResidual ( const std::string &  fname,
const std::string &  id,
const std::string &  desc,
int  loglevel = 1 
)

Definition at line 150 of file Sim1D.cpp.

◆ showSolution() [1/2]

void showSolution ( std::ostream &  s)

Print to stream s the current solution for all domains.

Definition at line 228 of file Sim1D.cpp.

References OneDim::domain(), Domain1D::domainType(), and OneDim::nDomains().

◆ showSolution() [2/2]

void showSolution ( )

Definition at line 237 of file Sim1D.cpp.

◆ solution() [1/2]

const doublereal * solution ( )
inline

Definition at line 111 of file Sim1D.h.

◆ setTimeStep()

void setTimeStep ( double  stepsize,
size_t  n,
const int *  tsteps 
)

Definition at line 284 of file Sim1D.cpp.

◆ solve()

void solve ( int  loglevel = 0,
bool  refine_grid = true 
)

Definition at line 307 of file Sim1D.cpp.

◆ eval()

void eval ( doublereal  rdt = -1.0,
int  count = 1 
)
inline

Definition at line 119 of file Sim1D.h.

◆ getResidual()

void getResidual ( double  rdt,
double *  resid 
)
inline

Definition at line 124 of file Sim1D.h.

◆ refine()

int refine ( int  loglevel = 0)

◆ setFixedTemperature()

int setFixedTemperature ( double  t)

Add node for fixed temperature point of freely propagating flame.

Definition at line 507 of file Sim1D.cpp.

References OneDim::domain(), Domain1D::domainType(), Domain1D::nComponents(), OneDim::nDomains(), Domain1D::nPoints(), and Cantera::npos.

◆ fixedTemperature()

double fixedTemperature ( )

Return temperature at the point used to fix the flame location.

Definition at line 597 of file Sim1D.cpp.

References OneDim::domain(), Domain1D::domainType(), and OneDim::nDomains().

◆ fixedTemperatureLocation()

double fixedTemperatureLocation ( )

Return location of the point where temperature is fixed.

Definition at line 610 of file Sim1D.cpp.

References OneDim::domain(), Domain1D::domainType(), and OneDim::nDomains().

◆ setRefineCriteria()

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.

See also
Refiner::setCriteria.

Definition at line 623 of file Sim1D.cpp.

References OneDim::domain(), OneDim::nDomains(), Domain1D::refiner(), and Refiner::setCriteria().

◆ getRefineCriteria()

vector_fp getRefineCriteria ( int  dom)

Get the grid refinement criteria.

dom must be greater than or equal to zero (that is, the domain must be specified).

See also
Refiner::getCriteria

Definition at line 637 of file Sim1D.cpp.

References OneDim::domain(), Refiner::getCriteria(), and Domain1D::refiner().

◆ setMaxGridPoints()

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.

See also
Refiner::setMaxPoints.

Definition at line 661 of file Sim1D.cpp.

References OneDim::domain(), OneDim::nDomains(), Domain1D::refiner(), and Refiner::setMaxPoints().

◆ maxGridPoints()

size_t maxGridPoints ( size_t  dom)

Get the maximum number of grid points in this domain.

See also
Refiner::maxPoints
Parameters
domdomain number, beginning with 0 for the leftmost domain.

Definition at line 674 of file Sim1D.cpp.

References OneDim::domain(), Refiner::maxPoints(), and Domain1D::refiner().

◆ setGridMin()

void setGridMin ( int  dom,
double  gridmin 
)

Set the minimum grid spacing in the specified domain(s).

Parameters
domDomain index. If dom == -1, the specified spacing is applied to all domains.
gridminThe minimum allowable grid spacing [m]

Definition at line 648 of file Sim1D.cpp.

References OneDim::domain(), OneDim::nDomains(), Domain1D::refiner(), and Refiner::setGridMin().

◆ restore()

void restore ( const std::string &  fname,
const std::string &  id,
int  loglevel = 2 
)

◆ restoreTimeSteppingSolution()

void restoreTimeSteppingSolution ( )

Set the current solution vector to the last successful time-stepping solution.

This can be used to examine the solver progress after a failed integration.

Definition at line 248 of file Sim1D.cpp.

References Sim1D::m_x, and Sim1D::m_xlast_ts.

◆ restoreSteadySolution()

void restoreSteadySolution ( )

Set the current solution vector and grid to the last successful steady- state solution.

This can be used to examine the solver progress after a failure during grid refinement.

Definition at line 257 of file Sim1D.cpp.

References OneDim::domain(), Sim1D::m_grid_last_ss, Sim1D::m_x, Sim1D::m_xlast_ss, OneDim::nDomains(), and Domain1D::setupGrid().

◆ getInitialSoln()

void getInitialSoln ( )

Definition at line 270 of file Sim1D.cpp.

◆ setSolution()

void setSolution ( const doublereal *  soln)
inline

Definition at line 193 of file Sim1D.h.

◆ solution() [2/2]

const doublereal * solution ( ) const
inline

Definition at line 197 of file Sim1D.h.

◆ jacobian()

doublereal jacobian ( int  i,
int  j 
)

Definition at line 680 of file Sim1D.cpp.

◆ evalSSJacobian()

void evalSSJacobian ( )

Definition at line 685 of file Sim1D.cpp.

◆ solveAdjoint()

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} \]

Definition at line 690 of file Sim1D.cpp.

◆ resize()

void resize ( )
virtual

Call after one or more grids has changed size, for example after being refined.

Reimplemented from OneDim.

Definition at line 714 of file Sim1D.cpp.

References Sim1D::m_x, Sim1D::m_xnew, OneDim::resize(), and OneDim::size().

Referenced by Sim1D::refine(), Sim1D::restore(), and Sim1D::Sim1D().

◆ setSteadyCallback()

void setSteadyCallback ( Func1 callback)
inline

Set a function that will be called after each successful steady-state solve, before regridding.

Intended to be used for observing solver progress for debugging purposes.

Definition at line 225 of file Sim1D.h.

References Sim1D::m_steady_callback.

◆ finalize()

void finalize ( )
private

Calls method _finalize in each domain.

Definition at line 277 of file Sim1D.cpp.

References Domain1D::_finalize(), OneDim::domain(), Sim1D::m_x, OneDim::nDomains(), and OneDim::start().

Referenced by Sim1D::refine(), and Sim1D::restore().

◆ newtonSolve()

int newtonSolve ( int  loglevel)
private

Wrapper around the Newton solver.

Returns
0 if successful, -1 on failure

Definition at line 293 of file Sim1D.cpp.

References Sim1D::m_x, Sim1D::m_xnew, and OneDim::solve().

Member Data Documentation

◆ m_x

vector_fp m_x
protected

◆ m_xlast_ts

vector_fp m_xlast_ts
protected

the solution vector after the last successful timestepping

Definition at line 234 of file Sim1D.h.

Referenced by Sim1D::restore(), and Sim1D::restoreTimeSteppingSolution().

◆ m_xlast_ss

vector_fp m_xlast_ss
protected

the solution vector after the last successful steady-state solve (stored before grid refinement)

Definition at line 238 of file Sim1D.h.

Referenced by Sim1D::refine(), and Sim1D::restoreSteadySolution().

◆ m_grid_last_ss

std::vector<vector_fp> m_grid_last_ss
protected

the grids for each domain after the last successful steady-state solve (stored before grid refinement)

Definition at line 242 of file Sim1D.h.

Referenced by Sim1D::refine(), and Sim1D::restoreSteadySolution().

◆ m_xnew

vector_fp m_xnew
protected

a work array used to hold the residual or the new solution

Definition at line 245 of file Sim1D.h.

Referenced by Sim1D::newtonSolve(), and Sim1D::resize().

◆ m_tstep

doublereal m_tstep
protected

timestep

Definition at line 248 of file Sim1D.h.

Referenced by Sim1D::Sim1D().

◆ m_steps

vector_int m_steps
protected

array of number of steps to take before re-attempting the steady-state solution

Definition at line 252 of file Sim1D.h.

Referenced by Sim1D::Sim1D().

◆ m_steady_callback

Func1* m_steady_callback
protected

User-supplied function called after a successful steady-state solve.

Definition at line 255 of file Sim1D.h.

Referenced by Sim1D::setSteadyCallback().


The documentation for this class was generated from the following files: