Cantera  2.3.0
Sim1D.cpp
Go to the documentation of this file.
1 /**
2  * @file Sim1D.cpp
3  */
4 
5 // This file is part of Cantera. See License.txt in the top-level directory or
6 // at http://www.cantera.org/license.txt for license and copyright information.
7 
8 #include "cantera/oneD/Sim1D.h"
10 #include "cantera/oneD/StFlow.h"
12 #include "cantera/numerics/funcs.h"
13 #include "cantera/base/xml.h"
14 #include "cantera/numerics/Func1.h"
15 
16 using namespace std;
17 
18 namespace Cantera
19 {
20 
21 Sim1D::Sim1D(vector<Domain1D*>& domains) :
22  OneDim(domains),
23  m_steady_callback(0)
24 {
25  // resize the internal solution vector and the work array, and perform
26  // domain-specific initialization of the solution vector.
27  resize();
28  for (size_t n = 0; n < nDomains(); n++) {
30  }
31 
32  // set some defaults
33  m_tstep = 1.0e-5;
34  m_steps = { 10 };
35 }
36 
37 void Sim1D::setInitialGuess(const std::string& component, vector_fp& locs, vector_fp& vals)
38 {
39  for (size_t dom=0; dom<nDomains(); dom++) {
40  Domain1D& d = domain(dom);
41  size_t ncomp = d.nComponents();
42  for (size_t comp=0; comp<ncomp; comp++) {
43  if (d.componentName(comp)==component) {
44  setProfile(dom,comp,locs,vals);
45  }
46  }
47  }
48 }
49 
50 void Sim1D::setValue(size_t dom, size_t comp, size_t localPoint, doublereal value)
51 {
52  size_t iloc = domain(dom).loc() + domain(dom).index(comp, localPoint);
53  AssertThrowMsg(iloc < m_x.size(), "Sim1D::setValue",
54  "Index out of bounds: {} > {}", iloc, m_x.size());
55  m_x[iloc] = value;
56 }
57 
58 doublereal Sim1D::value(size_t dom, size_t comp, size_t localPoint) const
59 {
60  size_t iloc = domain(dom).loc() + domain(dom).index(comp, localPoint);
61  AssertThrowMsg(iloc < m_x.size(), "Sim1D::value",
62  "Index out of bounds: {} > {}", iloc, m_x.size());
63  return m_x[iloc];
64 }
65 
66 doublereal Sim1D::workValue(size_t dom, size_t comp, size_t localPoint) const
67 {
68  size_t iloc = domain(dom).loc() + domain(dom).index(comp, localPoint);
69  AssertThrowMsg(iloc < m_x.size(), "Sim1D::workValue",
70  "Index out of bounds: {} > {}", iloc, m_x.size());
71  return m_xnew[iloc];
72 }
73 
74 void Sim1D::setProfile(size_t dom, size_t comp,
75  const vector_fp& pos, const vector_fp& values)
76 {
77  if (pos.front() != 0.0 || pos.back() != 1.0) {
78  throw CanteraError("Sim1D::setProfile",
79  "`pos` vector must span the range [0, 1]. Got a vector spanning "
80  "[{}, {}] instead.", pos.front(), pos.back());
81  }
82  Domain1D& d = domain(dom);
83  doublereal z0 = d.zmin();
84  doublereal z1 = d.zmax();
85  for (size_t n = 0; n < d.nPoints(); n++) {
86  double zpt = d.z(n);
87  double frac = (zpt - z0)/(z1 - z0);
88  double v = linearInterp(frac, pos, values);
89  setValue(dom, comp, n, v);
90  }
91 }
92 
93 void Sim1D::save(const std::string& fname, const std::string& id,
94  const std::string& desc, int loglevel)
95 {
96  OneDim::save(fname, id, desc, m_x.data(), loglevel);
97 }
98 
99 void Sim1D::saveResidual(const std::string& fname, const std::string& id,
100  const std::string& desc, int loglevel)
101 {
102  vector_fp res(m_x.size(), -999);
103  OneDim::eval(npos, &m_x[0], &res[0], 0.0);
104  OneDim::save(fname, id, desc, &res[0], loglevel);
105 }
106 
107 void Sim1D::restore(const std::string& fname, const std::string& id,
108  int loglevel)
109 {
110  XML_Node root;
111  root.build(fname);
112 
113  XML_Node* f = root.findID(id);
114  if (!f) {
115  throw CanteraError("Sim1D::restore","No solution with id = "+id);
116  }
117 
118  vector<XML_Node*> xd = f->getChildren("domain");
119  if (xd.size() != nDomains()) {
120  throw CanteraError("Sim1D::restore", "Solution does not contain the "
121  " correct number of domains. Found {} expected {}.\n",
122  xd.size(), nDomains());
123  }
124  for (size_t m = 0; m < nDomains(); m++) {
125  Domain1D& dom = domain(m);
126  if (loglevel > 0 && xd[m]->attrib("id") != dom.id()) {
127  writelog("Warning: domain names do not match: '" +
128  (*xd[m])["id"] + + "' and '" + dom.id() + "'\n");
129  }
130  dom.resize(domain(m).nComponents(), intValue((*xd[m])["points"]));
131  }
132  resize();
133  m_xlast_ts.clear();
134  for (size_t m = 0; m < nDomains(); m++) {
135  domain(m).restore(*xd[m], &m_x[domain(m).loc()], loglevel);
136  }
137  finalize();
138 }
139 
140 void Sim1D::setFlatProfile(size_t dom, size_t comp, doublereal v)
141 {
142  size_t np = domain(dom).nPoints();
143  for (size_t n = 0; n < np; n++) {
144  setValue(dom, comp, n, v);
145  }
146 }
147 
148 void Sim1D::showSolution(ostream& s)
149 {
150  for (size_t n = 0; n < nDomains(); n++) {
151  if (domain(n).domainType() != cEmptyType) {
152  domain(n).showSolution_s(s, &m_x[start(n)]);
153  }
154  }
155 }
156 
157 void Sim1D::showSolution()
158 {
159  for (size_t n = 0; n < nDomains(); n++) {
160  if (domain(n).domainType() != cEmptyType) {
161  writelog("\n\n>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> "+domain(n).id()
162  +" <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<\n\n");
163  domain(n).showSolution(&m_x[start(n)]);
164  }
165  }
166 }
167 
169 {
170  if (m_xlast_ts.empty()) {
171  throw CanteraError("Sim1D::restoreTimeSteppingSolution",
172  "No successful time steps taken on this grid.");
173  }
174  m_x = m_xlast_ts;
175 }
176 
178 {
179  if (m_xlast_ss.empty()) {
180  throw CanteraError("Sim1D::restoreSteadySolution",
181  "No successful steady state solution");
182  }
183  m_x = m_xlast_ss;
184  for (size_t n = 0; n < nDomains(); n++) {
185  vector_fp& z = m_grid_last_ss[n];
186  domain(n).setupGrid(z.size(), z.data());
187  }
188 }
189 
190 void Sim1D::getInitialSoln()
191 {
192  for (size_t n = 0; n < nDomains(); n++) {
193  domain(n)._getInitialSoln(&m_x[start(n)]);
194  }
195 }
196 
198 {
199  for (size_t n = 0; n < nDomains(); n++) {
200  domain(n)._finalize(&m_x[start(n)]);
201  }
202 }
203 
204 void Sim1D::setTimeStep(double stepsize, size_t n, const int* tsteps)
205 {
206  m_tstep = stepsize;
207  m_steps.resize(n);
208  for (size_t i = 0; i < n; i++) {
209  m_steps[i] = tsteps[i];
210  }
211 }
212 
213 int Sim1D::newtonSolve(int loglevel)
214 {
215  int m = OneDim::solve(m_x.data(), m_xnew.data(), loglevel);
216  if (m >= 0) {
217  m_x = m_xnew;
218  return 0;
219  } else if (m > -10) {
220  return -1;
221  } else {
222  throw CanteraError("Sim1D::newtonSolve",
223  "ERROR: OneDim::solve returned m = {}", m);
224  }
225 }
226 
227 void Sim1D::solve(int loglevel, bool refine_grid)
228 {
229  int new_points = 1;
230  doublereal dt = m_tstep;
231  m_nsteps = 0;
232  int soln_number = -1;
233  finalize();
234 
235  while (new_points > 0) {
236  size_t istep = 0;
237  int nsteps = m_steps[istep];
238 
239  bool ok = false;
240  if (loglevel > 0) {
241  writeline('.', 78, true, true);
242  }
243  while (!ok) {
244  // Attempt to solve the steady problem
245  setSteadyMode();
246  newton().setOptions(m_ss_jac_age);
247  debuglog("Attempt Newton solution of steady-state problem...", loglevel);
248  int status = newtonSolve(loglevel-1);
249 
250  if (status == 0) {
251  if (loglevel > 0) {
252  writelog(" success.\n\n");
253  writelog("Problem solved on [");
254  for (size_t mm = 1; mm < nDomains(); mm+=2) {
255  writelog("{}", domain(mm).nPoints());
256  if (mm + 2 < nDomains()) {
257  writelog(", ");
258  }
259  }
260  writelog("] point grid(s).\n");
261  }
262  if (m_steady_callback) {
264  }
265 
266  if (loglevel > 6) {
267  save("debug_sim1d.xml", "debug",
268  "After successful Newton solve");
269  }
270  if (loglevel > 7) {
271  saveResidual("debug_sim1d.xml", "residual",
272  "After successful Newton solve");
273  }
274  ok = true;
275  soln_number++;
276  } else {
277  debuglog(" failure. \n", loglevel);
278  if (loglevel > 6) {
279  save("debug_sim1d.xml", "debug",
280  "After unsuccessful Newton solve");
281  }
282  if (loglevel > 7) {
283  saveResidual("debug_sim1d.xml", "residual",
284  "After unsuccessful Newton solve");
285  }
286  if (loglevel > 0) {
287  writelog("Take {} timesteps ", nsteps);
288  }
289  dt = timeStep(nsteps, dt, m_x.data(), m_xnew.data(),
290  loglevel-1);
291  m_xlast_ts = m_x;
292  if (loglevel > 6) {
293  save("debug_sim1d.xml", "debug", "After timestepping");
294  }
295  if (loglevel > 7) {
296  saveResidual("debug_sim1d.xml", "residual",
297  "After timestepping");
298  }
299 
300  if (loglevel == 1) {
301  writelog(" {:10.4g} {:10.4g}\n", dt,
302  log10(ssnorm(m_x.data(), m_xnew.data())));
303  }
304  istep++;
305  if (istep >= m_steps.size()) {
306  nsteps = m_steps.back();
307  } else {
308  nsteps = m_steps[istep];
309  }
310  dt = std::min(dt, m_tmax);
311  }
312  }
313  if (loglevel > 0) {
314  writeline('.', 78, true, true);
315  }
316  if (loglevel > 2) {
317  showSolution();
318  }
319 
320  if (refine_grid) {
321  new_points = refine(loglevel);
322  if (new_points) {
323  // If the grid has changed, preemptively reduce the timestep
324  // to avoid multiple successive failed time steps.
325  dt = m_tstep;
326  }
327  if (new_points && loglevel > 6) {
328  save("debug_sim1d.xml", "debug", "After regridding");
329  }
330  if (new_points && loglevel > 7) {
331  saveResidual("debug_sim1d.xml", "residual",
332  "After regridding");
333  }
334  } else {
335  debuglog("grid refinement disabled.\n", loglevel);
336  new_points = 0;
337  }
338  }
339 }
340 
341 int Sim1D::refine(int loglevel)
342 {
343  int ianalyze, np = 0;
344  vector_fp znew, xnew;
345  std::vector<size_t> dsize;
346 
347  m_xlast_ss = m_x;
348  m_grid_last_ss.clear();
349 
350  for (size_t n = 0; n < nDomains(); n++) {
351  Domain1D& d = domain(n);
352  Refiner& r = d.refiner();
353 
354  // Save the old grid corresponding to the converged solution
355  m_grid_last_ss.push_back(d.grid());
356 
357  // determine where new points are needed
358  ianalyze = r.analyze(d.grid().size(), d.grid().data(), &m_x[start(n)]);
359  if (ianalyze < 0) {
360  return ianalyze;
361  }
362 
363  if (loglevel > 0) {
364  r.show();
365  }
366 
367  np += r.nNewPoints();
368  size_t comp = d.nComponents();
369 
370  // loop over points in the current grid
371  size_t npnow = d.nPoints();
372  size_t nstart = znew.size();
373  for (size_t m = 0; m < npnow; m++) {
374  if (r.keepPoint(m)) {
375  // add the current grid point to the new grid
376  znew.push_back(d.grid(m));
377 
378  // do the same for the solution at this point
379  for (size_t i = 0; i < comp; i++) {
380  xnew.push_back(value(n, i, m));
381  }
382 
383  // now check whether a new point is needed in the interval to
384  // the right of point m, and if so, add entries to znew and xnew
385  // for this new point
386  if (r.newPointNeeded(m) && m + 1 < npnow) {
387  // add new point at midpoint
388  double zmid = 0.5*(d.grid(m) + d.grid(m+1));
389  znew.push_back(zmid);
390  np++;
391 
392  // for each component, linearly interpolate
393  // the solution to this point
394  for (size_t i = 0; i < comp; i++) {
395  double xmid = 0.5*(value(n, i, m) + value(n, i, m+1));
396  xnew.push_back(xmid);
397  }
398  }
399  } else {
400  if (loglevel > 0) {
401  writelog("refine: discarding point at {}\n", d.grid(m));
402  }
403  }
404  }
405  dsize.push_back(znew.size() - nstart);
406  }
407 
408  // At this point, the new grid znew and the new solution vector xnew have
409  // been constructed, but the domains themselves have not yet been modified.
410  // Now update each domain with the new grid.
411 
412  size_t gridstart = 0, gridsize;
413  for (size_t n = 0; n < nDomains(); n++) {
414  Domain1D& d = domain(n);
415  gridsize = dsize[n];
416  d.setupGrid(gridsize, &znew[gridstart]);
417  gridstart += gridsize;
418  }
419 
420  // Replace the current solution vector with the new one
421  m_x = xnew;
422  resize();
423  finalize();
424  return np;
425 }
426 
427 int Sim1D::setFixedTemperature(doublereal t)
428 {
429  int np = 0;
430  vector_fp znew, xnew;
431  doublereal zfixed;
432  doublereal z1 = 0.0, z2 = 0.0, t1,t2;
433  size_t m1 = 0;
434  std::vector<size_t> dsize;
435 
436  for (size_t n = 0; n < nDomains(); n++) {
437  bool addnewpt=false;
438  Domain1D& d = domain(n);
439  size_t comp = d.nComponents();
440 
441  // loop over points in the current grid to determine where new point is
442  // needed.
443  FreeFlame* d_free = dynamic_cast<FreeFlame*>(&domain(n));
444  size_t npnow = d.nPoints();
445  size_t nstart = znew.size();
446  if (d_free) {
447  for (size_t m = 0; m < npnow-1; m++) {
448  if (value(n,2,m) == t) {
449  zfixed = d.grid(m);
450  d_free->m_zfixed = zfixed;
451  d_free->m_tfixed = t;
452  addnewpt = false;
453  break;
454  } else if ((value(n,2,m)<t) && (value(n,2,m+1)>t)) {
455  z1 = d.grid(m);
456  m1 = m;
457  z2 = d.grid(m+1);
458  t1 = value(n,2,m);
459  t2 = value(n,2,m+1);
460 
461  zfixed = (z1-z2)/(t1-t2)*(t-t2)+z2;
462  d_free->m_zfixed = zfixed;
463  d_free->m_tfixed = t;
464  addnewpt = true;
465  break;
466  //copy solution domain and push back values
467  }
468  }
469  }
470 
471  for (size_t m = 0; m < npnow; m++) {
472  // add the current grid point to the new grid
473  znew.push_back(d.grid(m));
474 
475  // do the same for the solution at this point
476  for (size_t i = 0; i < comp; i++) {
477  xnew.push_back(value(n, i, m));
478  }
479  if (m==m1 && addnewpt) {
480  //add new point at zfixed
481  znew.push_back(zfixed);
482  np++;
483  double interp_factor = (zfixed-z2) / (z1-z2);
484  // for each component, linearly interpolate
485  // the solution to this point
486  for (size_t i = 0; i < comp; i++) {
487  double xmid = interp_factor*(value(n, i, m) - value(n, i, m+1)) + value(n,i,m+1);
488  xnew.push_back(xmid);
489  }
490  }
491  }
492  dsize.push_back(znew.size() - nstart);
493  }
494 
495  // At this point, the new grid znew and the new solution vector xnew have
496  // been constructed, but the domains themselves have not yet been modified.
497  // Now update each domain with the new grid.
498  size_t gridstart = 0;
499  for (size_t n = 0; n < nDomains(); n++) {
500  Domain1D& d = domain(n);
501  size_t gridsize = dsize[n];
502  d.setupGrid(gridsize, &znew[gridstart]);
503  gridstart += gridsize;
504  }
505 
506  // Replace the current solution vector with the new one
507  m_x = xnew;
508 
509  resize();
510  finalize();
511  return np;
512 }
513 
514 void Sim1D::setRefineCriteria(int dom, doublereal ratio,
515  doublereal slope, doublereal curve, doublereal prune)
516 {
517  if (dom >= 0) {
518  Refiner& r = domain(dom).refiner();
519  r.setCriteria(ratio, slope, curve, prune);
520  } else {
521  for (size_t n = 0; n < nDomains(); n++) {
522  Refiner& r = domain(n).refiner();
523  r.setCriteria(ratio, slope, curve, prune);
524  }
525  }
526 }
527 
528 void Sim1D::setGridMin(int dom, double gridmin)
529 {
530  if (dom >= 0) {
531  Refiner& r = domain(dom).refiner();
532  r.setGridMin(gridmin);
533  } else {
534  for (size_t n = 0; n < nDomains(); n++) {
535  Refiner& r = domain(n).refiner();
536  r.setGridMin(gridmin);
537  }
538  }
539 }
540 
541 void Sim1D::setMaxGridPoints(int dom, int npoints)
542 {
543  if (dom >= 0) {
544  Refiner& r = domain(dom).refiner();
545  r.setMaxPoints(npoints);
546  } else {
547  for (size_t n = 0; n < nDomains(); n++) {
548  Refiner& r = domain(n).refiner();
549  r.setMaxPoints(npoints);
550  }
551  }
552 }
553 
554 size_t Sim1D::maxGridPoints(size_t dom)
555 {
556  Refiner& r = domain(dom).refiner();
557  return r.maxPoints();
558 }
559 
560 doublereal Sim1D::jacobian(int i, int j)
561 {
562  return OneDim::jacobian().value(i,j);
563 }
564 
565 void Sim1D::evalSSJacobian()
566 {
567  OneDim::evalSSJacobian(m_x.data(), m_xnew.data());
568 }
569 
570 void Sim1D::solveAdjoint(const double* b, double* lambda)
571 {
572  for (auto& D : m_dom) {
573  D->forceFullUpdate(true);
574  }
575  evalSSJacobian();
576  for (auto& D : m_dom) {
577  D->forceFullUpdate(false);
578  }
579 
580  // Form J^T
581  size_t bw = bandwidth();
582  BandMatrix Jt(size(), bw, bw);
583  for (size_t i = 0; i < size(); i++) {
584  size_t j1 = (i > bw) ? i - bw : 0;
585  size_t j2 = (i + bw >= size()) ? size() - 1: i + bw;
586  for (size_t j = j1; j <= j2; j++) {
587  Jt(j,i) = m_jac->value(i,j);
588  }
589  }
590 
591  Jt.solve(b, lambda);
592 }
593 
595 {
596  OneDim::resize();
597  m_x.resize(size(), 0.0);
598  m_xnew.resize(size(), 0.0);
599 }
600 
601 }
Container class for multiple-domain 1D problems.
Definition: OneDim.h:25
std::vector< XML_Node * > getChildren(const std::string &name) const
Get a vector of pointers to XML_Node containing all of the children of the current node which match t...
Definition: xml.cpp:864
std::vector< vector_fp > m_grid_last_ss
the grids for each domain after the last successful steady-state solve (stored before grid refinement...
Definition: Sim1D.h:228
void showSolution(std::ostream &s)
Print to stream s the current solution for all domains.
Definition: Sim1D.cpp:148
virtual void _getInitialSoln(doublereal *x)
Writes some or all initial solution values into the global solution array, beginning at the location ...
Definition: Domain1D.cpp:293
size_t maxGridPoints(size_t dom)
Get the maximum number of grid points in this domain.
Definition: Sim1D.cpp:554
void restore(const std::string &fname, const std::string &id, int loglevel=2)
Initialize the solution with a previously-saved solution.
Definition: Sim1D.cpp:107
double timeStep(int nsteps, double dt, double *x, double *r, int loglevel)
Definition: OneDim.cpp:348
vector_fp m_x
the solution vector
Definition: Sim1D.h:217
vector_fp m_xlast_ts
the solution vector after the last successful timestepping
Definition: Sim1D.h:220
virtual void resize()
Call after one or more grids has changed size, e.g. after being refined.
Definition: OneDim.cpp:168
doublereal ssnorm(doublereal *x, doublereal *r)
Steady-state max norm (infinity norm) of the residual evaluated using solution x. ...
Definition: OneDim.cpp:290
void eval(size_t j, double *x, double *r, doublereal rdt=-1.0, int count=1)
Evaluate the multi-domain residual function.
Definition: OneDim.cpp:258
vector_fp m_xlast_ss
the solution vector after the last successful steady-state solve (stored before grid refinement) ...
Definition: Sim1D.h:224
size_t size() const
Total solution vector length;.
Definition: OneDim.h:86
MultiJac & jacobian()
Return a reference to the Jacobian evaluator.
Definition: OneDim.cpp:101
std::unique_ptr< MultiJac > m_jac
Jacobian evaluator.
Definition: OneDim.h:324
const size_t npos
index returned by functions to indicate "no position"
Definition: ct_defs.h:165
void setValue(size_t dom, size_t comp, size_t localPoint, doublereal value)
Set a single value in the solution vector.
Definition: Sim1D.cpp:50
int solve(const doublereal *const b, doublereal *const x)
Solve the matrix problem Ax = b.
Definition: BandMatrix.cpp:261
void writelog(const std::string &fmt, const Args &... args)
Write a formatted message to the screen.
Definition: global.h:179
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.
Definition: Sim1D.cpp:514
int solve(doublereal *x0, doublereal *x1, int loglevel)
Solve F(x) = 0, where F(x) is the multi-domain residual function.
Definition: OneDim.cpp:224
Class XML_Node is a tree-based representation of the contents of an XML file.
Definition: xml.h:97
void setInitialGuess(const std::string &component, vector_fp &locs, vector_fp &vals)
Set initial guess for one component for all domains.
Definition: Sim1D.cpp:37
STL namespace.
virtual doublereal eval(doublereal t) const
Evaluate the function.
Definition: Func1.cpp:60
virtual void restore(const XML_Node &dom, doublereal *soln, int loglevel)
Restore the solution for this domain from an XML_Node.
Definition: Domain1D.cpp:173
void solveAdjoint(const double *b, double *lambda)
Solve the equation .
Definition: Sim1D.cpp:570
size_t maxPoints() const
Returns the maximum number of points allowed in the domain.
Definition: refine.h:49
size_t loc(size_t jg)
Location in the solution vector of the first component of global point jg.
Definition: OneDim.h:107
MultiNewton & newton()
Return a reference to the Newton iterator.
Definition: OneDim.cpp:105
vector_int m_steps
array of number of steps to take before re-attempting the steady-state solution
Definition: Sim1D.h:238
int newtonSolve(int loglevel)
Wrapper around the Newton solver.
Definition: Sim1D.cpp:213
A class for freely-propagating premixed flames.
Definition: StFlow.h:457
size_t nDomains() const
Number of domains.
Definition: OneDim.h:52
doublereal m_tfixed
Temperature at the point used to fix the flame location.
Definition: StFlow.h:481
Domain1D & domain(size_t i) const
Return a reference to domain i.
Definition: OneDim.h:57
Base class for one-dimensional domains.
Definition: Domain1D.h:36
virtual void resize(size_t nv, size_t np)
Definition: Domain1D.cpp:34
int setFixedTemperature(doublereal t)
Add node for fixed temperature point of freely propagating flame.
Definition: Sim1D.cpp:427
void setCriteria(doublereal ratio=10.0, doublereal slope=0.8, doublereal curve=0.8, doublereal prune=-0.1)
Set grid refinement criteria.
Definition: refine.cpp:23
void setMaxPoints(int npmax)
Set the maximum number of points allowed in the domain.
Definition: refine.h:44
void finalize()
Calls method _finalize in each domain.
Definition: Sim1D.cpp:197
virtual void _finalize(const doublereal *x)
In some cases, a domain may need to set parameters that depend on the initial solution estimate...
Definition: Domain1D.h:540
void restoreTimeSteppingSolution()
Set the current solution vector to the last successful time-stepping solution.
Definition: Sim1D.cpp:168
size_t bandwidth() const
Jacobian bandwidth.
Definition: OneDim.h:116
Classes providing support for XML data files.
virtual void resize()
Call after one or more grids has changed size, e.g. after being refined.
Definition: Sim1D.cpp:594
void setSteadyMode()
Definition: OneDim.cpp:319
virtual void showSolution(const doublereal *x)
Print the solution.
Definition: Domain1D.cpp:244
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 solutio...
Definition: OneDim.cpp:65
Func1 * m_steady_callback
User-supplied function called after a successful steady-state solve.
Definition: Sim1D.h:241
Base class for exceptions thrown by Cantera classes.
Definition: ctexceptions.h:65
doublereal m_tstep
timestep
Definition: Sim1D.h:234
void build(const std::string &filename)
Populate the XML tree from an input file.
Definition: xml.cpp:716
doublereal & value(size_t i, size_t j)
Return a changeable reference to element (i,j).
Definition: BandMatrix.cpp:146
size_t start(size_t i) const
The index of the start of domain i in the solution vector.
Definition: OneDim.h:81
doublereal m_zfixed
Location of the point where temperature is fixed.
Definition: StFlow.h:478
void setGridMin(double gridmin)
Set the minimum allowable spacing between adjacent grid points [m].
Definition: refine.h:54
void setFlatProfile(size_t dom, size_t comp, doublereal v)
Set component &#39;comp&#39; of domain &#39;dom&#39; to value &#39;v&#39; at all points.
Definition: Sim1D.cpp:140
virtual std::string componentName(size_t n) const
Name of the nth component. May be overloaded.
Definition: Domain1D.cpp:58
Refiner & refiner()
Return a reference to the grid refiner.
Definition: Domain1D.h:125
int domainType()
Domain type flag.
Definition: Domain1D.h:50
void debuglog(const std::string &msg, int loglevel)
Write a message to the log only if loglevel > 0.
Definition: global.h:161
int intValue(const std::string &val)
Translate a string into one integer value.
vector_fp m_xnew
a work array used to hold the residual or the new solution
Definition: Sim1D.h:231
size_t nComponents() const
Number of components at each grid point.
Definition: Domain1D.h:130
#define AssertThrowMsg(expr, procedure,...)
Assertion must be true or an error is thrown.
Definition: ctexceptions.h:270
doublereal value(size_t dom, size_t comp, size_t localPoint) const
Get one entry in the solution vector.
Definition: Sim1D.cpp:58
std::vector< double > vector_fp
Turn on the use of stl vectors for the basic array type within cantera Vector of doubles.
Definition: ct_defs.h:157
int refine(int loglevel=0)
Refine the grid in all domains.
Definition: Sim1D.cpp:341
void restoreSteadySolution()
Set the current solution vector and grid to the last successful steady- state solution.
Definition: Sim1D.cpp:177
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.
Definition: Sim1D.cpp:74
Refine Domain1D grids so that profiles satisfy adaptation tolerances.
Definition: refine.h:16
void setMaxGridPoints(int dom, int npoints)
Set the maximum number of grid points in the domain.
Definition: Sim1D.cpp:541
size_t nPoints() const
Number of grid points in this domain.
Definition: Domain1D.h:152
Namespace for the Cantera kernel.
Definition: application.cpp:29
Header for a file containing miscellaneous numerical functions.
int m_nsteps
Number of time steps taken in the current call to solve()
Definition: OneDim.h:351
void setGridMin(int dom, double gridmin)
Set the minimum grid spacing in the specified domain(s).
Definition: Sim1D.cpp:528
doublereal m_tmax
maximum timestep size
Definition: OneDim.h:319
void setOptions(int maxJacAge=5)
Set options.
Definition: MultiNewton.h:66
virtual void setupGrid(size_t n, const doublereal *z)
called to set up initial grid, and after grid refinement
Definition: Domain1D.cpp:234
virtual size_t loc(size_t j=0) const
Location of the start of the local solution vector in the global solution vector,.
Definition: Domain1D.h:403
XML_Node * findID(const std::string &id, const int depth=100) const
This routine carries out a recursive search for an XML node based on the XML element attribute...
Definition: xml.cpp:645
doublereal linearInterp(doublereal x, const vector_fp &xpts, const vector_fp &fpts)
Linearly interpolate a function defined on a discrete grid.
Definition: funcs.cpp:13
A class for banded matrices, involving matrix inversion processes.
Definition: BandMatrix.h:41