|
| 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.
|
|
|
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.
|
|
|
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.
|
|
| 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.
|
|
shared_ptr< SystemJacobian > | getJacobian () |
|
MultiNewton & | newton () |
| Return a reference to the Newton iterator.
|
|
void | setLinearSolver (shared_ptr< SystemJacobian > solver) |
| Set the linear solver used to hold the Jacobian matrix and solve linear systems as part of each Newton iteration.
|
|
shared_ptr< SystemJacobian > | linearSolver () const |
| Get the type of the linear solver being used.
|
|
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.
|
|
void | evalJacobian (double *x0) |
| Evaluates the Jacobian at x0 using finite differences.
|
|
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 | setJacobianPerturbation (double relative, double absolute, double threshold) |
| Configure perturbations used to evaluate finite difference Jacobian.
|
|
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.
|
|
|
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.
|
|
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.
|
|
shared_ptr< SystemJacobian > | 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< double > | m_work1 |
| Work arrays used during Jacobian evaluation.
|
|
vector< double > | m_work2 |
|
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()
|
|
double | m_jacobianThreshold = 0.0 |
| Threshold for ignoring small elements in Jacobian.
|
|
double | m_jacobianRelPerturb = 1e-5 |
| Relative perturbation of each component in finite difference Jacobian.
|
|
double | m_jacobianAbsPerturb = 1e-10 |
| Absolute perturbation of each component in finite difference Jacobian.
|
|