Cantera  2.2.1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
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 (doublereal stepsize, size_t n, integer *tsteps)
 
void solve (int loglevel=0, bool refine_grid=true)
 
void eval (doublereal rdt=-1.0, int count=1)
 
int refine (int loglevel=0)
 Refine the grid in all domains. More...
 
int setFixedTemperature (doublereal t)
 Add node for fixed temperature point of freely propagating flame. More...
 
void setRefineCriteria (int dom=-1, doublereal ratio=10.0, doublereal slope=0.8, doublereal curve=0.8, doublereal prune=-0.1)
 Set grid refinement criteria. More...
 
void setMaxGridPoints (int dom=-1, int npoints=300)
 
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 getInitialSoln ()
 
void setSolution (const doublereal *soln)
 
const doublereal * solution () const
 
doublereal jacobian (int i, int j)
 
void evalSSJacobian ()
 
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 based on equilibrium. 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...
 
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 Throws an exception if n is greater than nDomains()-1. More...
 
void checkDomainArraySize (size_t nn) const
 Check that an array size is at least nDomains() Throws an exception if nn is less than 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...
 
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...
 
void resize ()
 Call after one or more grids has been refined. More...
 
vector_inttransientMask ()
 
double timeStep (int nsteps, double dt, double *x, double *r, int loglevel)
 
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)
 
void setMinTimeStep (doublereal tmin)
 
void setMaxTimeStep (doublereal tmax)
 
void setTimeStepFactor (doublereal tfactor)
 
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...
 
void setInterrupt (Func1 *interrupt)
 Set a function that will be called every time eval is called. More...
 

Protected Attributes

vector_fp m_x
 the solution vector 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...
 
- Protected Attributes inherited from OneDim
doublereal m_tmin
 
doublereal m_tmax
 
doublereal m_tfactor
 
MultiJacm_jac
 
MultiNewtonm_newt
 
doublereal m_rdt
 
bool m_jac_ok
 
size_t m_nd
 number of domains More...
 
size_t m_bw
 
size_t m_size
 
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
 
doublereal m_solve_time
 
int m_ss_jac_age
 
int m_ts_jac_age
 
Func1m_interrupt
 Function called at the start of every call to eval. More...
 

Private Member Functions

void finalize ()
 Calls method _finalize in each domain. More...
 
int newtonSolve (int loglevel)
 

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 18 of file Sim1D.h.

Constructor & Destructor Documentation

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 28 of file Sim1D.h.

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 — i.e., the pointer to the leftmost domain is domain[0], the pointer to the domain to its right is domain[1], etc.

Definition at line 18 of file Sim1D.cpp.

References Domain1D::_getInitialSoln(), DATA_PTR, OneDim::domain(), OneDim::m_nd, Sim1D::m_steps, Sim1D::m_tstep, Sim1D::m_x, Sim1D::m_xnew, OneDim::size(), and OneDim::start().

Member Function Documentation

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

Set initial guess based on equilibrium.

Definition at line 38 of file Sim1D.cpp.

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

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 51 of file Sim1D.cpp.

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

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

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 60 of file Sim1D.cpp.

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

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

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().

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 149 of file Sim1D.cpp.

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

void showSolution ( std::ostream &  s)

Print to stream s the current solution for all domains.

Definition at line 158 of file Sim1D.cpp.

References DATA_PTR, OneDim::domain(), OneDim::m_nd, Sim1D::m_x, and OneDim::start().

int refine ( int  loglevel = 0)
int setFixedTemperature ( doublereal  t)
void setRefineCriteria ( int  dom = -1,
doublereal  ratio = 10.0,
doublereal  slope = 0.8,
doublereal  curve = 0.8,
doublereal  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 514 of file Sim1D.cpp.

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

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 528 of file Sim1D.cpp.

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

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

Calls method _finalize in each domain.

Definition at line 185 of file Sim1D.cpp.

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

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

int newtonSolve ( int  loglevel)
private

Wrapper around the Newton solver.

Returns
0 if successful, -1 on failure

Definition at line 201 of file Sim1D.cpp.

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

Member Data Documentation

vector_fp m_x
protected
vector_fp m_xnew
protected

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

Definition at line 163 of file Sim1D.h.

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

doublereal m_tstep
protected

timestep

Definition at line 166 of file Sim1D.h.

Referenced by Sim1D::Sim1D().

vector_int m_steps
protected

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

Definition at line 170 of file Sim1D.h.

Referenced by Sim1D::Sim1D().


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