Cantera  2.2.1
Domain1D.cpp
Go to the documentation of this file.
1 /**
2  * @file Domain1D.cpp
3  */
4
7 #include "cantera/base/ctml.h"
8
9 #include <cstdio>
10
11 using namespace std;
12
13 namespace Cantera
14 {
15
16 void Domain1D::setTransientTolerances(doublereal rtol, doublereal atol, size_t n)
17 {
18  if (n == npos) {
19  for (n = 0; n < m_nv; n++) {
20  m_rtol_ts[n] = rtol;
21  m_atol_ts[n] = atol;
22  }
23  } else {
24  m_rtol_ts[n] = rtol;
25  m_atol_ts[n] = atol;
26  }
27 }
28
29 void Domain1D::setSteadyTolerances(doublereal rtol, doublereal atol, size_t n)
30 {
31  if (n == npos) {
32  for (n = 0; n < m_nv; n++) {
33  m_rtol_ss[n] = rtol;
34  m_atol_ss[n] = atol;
35  }
36  } else {
37  m_rtol_ss[n] = rtol;
38  m_atol_ss[n] = atol;
39  }
40 }
41
42 void Domain1D::needJacUpdate()
43 {
44  if (m_container) {
45  m_container->jacobian().setAge(10000);
46  m_container->saveStats();
47  }
48 }
49
50 void Domain1D::eval(size_t jg, doublereal* xg, doublereal* rg,
52 {
53
54  if (jg != npos && (jg + 1 < firstPoint() || jg > lastPoint() + 1)) {
55  return;
56  }
57
58  // if evaluating a Jacobian, compute the steady-state residual
59  if (jg != npos) {
60  rdt = 0.0;
61  }
62
63  // start of local part of global arrays
64  doublereal* x = xg + loc();
65  doublereal* rsd = rg + loc();
66  integer* diag = mask + loc();
67
68  size_t jmin, jmax, jpt, j, i;
69  jpt = jg - firstPoint();
70
71  if (jg == npos) { // evaluate all points
72  jmin = 0;
73  jmax = m_points - 1;
74  } else { // evaluate points for Jacobian
75  jmin = std::max<size_t>(jpt, 1) - 1;
76  jmax = std::min(jpt+1,m_points-1);
77  }
78
79  for (j = jmin; j <= jmax; j++) {
80  if (j == 0 || j == m_points - 1) {
81  for (i = 0; i < m_nv; i++) {
82  rsd[index(i,j)] = residual(x,i,j);
83  diag[index(i,j)] = 0;
84  }
85  } else {
86  for (i = 0; i < m_nv; i++) {
87  rsd[index(i,j)] = residual(x,i,j)
88  - timeDerivativeFlag(i)*rdt*(value(x,i,j) - prevSoln(i,j));
89  diag[index(i,j)] = timeDerivativeFlag(i);
90  }
91  }
92  }
93 }
94
95 XML_Node& Domain1D::save(XML_Node& o, const doublereal* const sol)
96 {
105  return d;
106 }
107
108 void Domain1D::restore(const XML_Node& dom, doublereal* soln, int loglevel)
109 {
110  vector_fp values;
111  vector<XML_Node*> nodes = dom.getChildren("floatArray");
112  for (size_t i = 0; i < nodes.size(); i++) {
113  string title = nodes[i]->attrib("title");
114  getFloatArray(*nodes[i], values, false);
115  if (values.size() != nComponents()) {
116  if (loglevel > 0) {
117  writelog("Warning: Domain1D::restore: Got an array of length " +
118  int2str(values.size()) + " when one of length " +
119  int2str(nComponents()) + " was expected. " +
120  "Tolerances for individual species may not be preserved.\n");
121  }
122  // The number of components will differ when restoring from a
123  // mechanism with a different number of species. Assuming that
124  // tolerances are the same for all species, we can just copy the
125  // tolerance from the last species.
126  if (!values.empty()) {
127  values.resize(nComponents(), values[values.size()-1]);
128  } else {
129  // If the tolerance vector is empty, just leave the defaults
130  // in place.
131  continue;
132  }
133  }
134  if (title == "abstol_transient") {
135  m_atol_ts = values;
136  } else if (title == "reltol_transient") {
137  m_rtol_ts = values;
138  } else if (title == "abstol_steady") {
139  m_atol_ss = values;
140  } else if (title == "reltol_steady") {
141  m_rtol_ss = values;
142  } else {
143  throw CanteraError("Domain1D::restore",
144  "Got an unexpected array, '" + title + "'");
145  }
146  }
147 }
148
149 void Domain1D::setupGrid(size_t n, const doublereal* z)
150 {
151  if (n > 1) {
152  resize(m_nv, n);
153  for (size_t j = 0; j < m_points; j++) {
154  m_z[j] = z[j];
155  }
156  }
157 }
158
159 void Domain1D::showSolution(const doublereal* x)
160 {
161  size_t nn = m_nv/5;
162  size_t i, j, n;
163  char buf[100];
164  doublereal v;
165  for (i = 0; i < nn; i++) {
166  writeline('-', 79, false, true);
167  sprintf(buf, "\n z ");
168  writelog(buf);
169  for (n = 0; n < 5; n++) {
170  sprintf(buf, " %10s ",componentName(i*5 + n).c_str());
171  writelog(buf);
172  }
173  writeline('-', 79, false, true);
174  for (j = 0; j < m_points; j++) {
175  sprintf(buf, "\n %10.4g ",m_z[j]);
176  writelog(buf);
177  for (n = 0; n < 5; n++) {
178  v = value(x, i*5+n, j);
179  sprintf(buf, " %10.4g ",v);
180  writelog(buf);
181  }
182  }
183  writelog("\n");
184  }
185  size_t nrem = m_nv - 5*nn;
186  writeline('-', 79, false, true);
187  sprintf(buf, "\n z ");
188  writelog(buf);
189  for (n = 0; n < nrem; n++) {
190  sprintf(buf, " %10s ", componentName(nn*5 + n).c_str());
191  writelog(buf);
192  }
193  writeline('-', 79, false, true);
194  for (j = 0; j < m_points; j++) {
195  sprintf(buf, "\n %10.4g ",m_z[j]);
196  writelog(buf);
197  for (n = 0; n < nrem; n++) {
198  v = value(x, nn*5+n, j);
199  sprintf(buf, " %10.4g ", v);
200  writelog(buf);
201  }
202  }
203  writelog("\n");
204 }
205
206 void Domain1D::_getInitialSoln(doublereal* x)
207 {
208  for (size_t j = 0; j < m_points; j++) {
209  for (size_t n = 0; n < m_nv; n++) {
210  x[index(n,j)] = initialValue(n,j);
211  }
212  }
213 }
214
215 doublereal Domain1D::initialValue(size_t n, size_t j)
216 {
217  throw CanteraError("Domain1D::initialValue",
218  "base class method called!");
219  return 0.0;
220 }
221
222 } // namespace
void addFloatArray(XML_Node &node, const std::string &title, const size_t n, const doublereal *const vals, const std::string &units, const std::string &type, const doublereal minval, const doublereal maxval)
This function adds a child node with the name, "floatArray", with a value consisting of a comma separ...
Definition: ctml.cpp:51
std::string int2str(const int n, const std::string &fmt)
Convert an int to a string using a format converter.
Definition: stringUtils.cpp:39
CTML ("Cantera Markup Language") is the variant of XML that Cantera uses to store data...
const size_t npos
index returned by functions to indicate "no position"
Definition: ct_defs.h:165
Class XML_Node is a tree-based representation of the contents of an XML file.
Definition: xml.h:100
Base class for exceptions thrown by Cantera classes.
Definition: ctexceptions.h:99
void addAttribute(const std::string &attrib, const std::string &value)
Add or modify an attribute of the current node.
Definition: xml.cpp:501
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
size_t getFloatArray(const XML_Node &node, std::vector< doublereal > &v, const bool convert, const std::string &unitsString, const std::string &nodeName)
This function reads the current node or a child node of the current node with the default name...
Definition: ctml.cpp:323
void writelog(const std::string &msg)
Write a message to the screen.
Definition: global.cpp:33
void getChildren(const std::string &name, std::vector< XML_Node * > &children) const
Get a vector of pointers to XML_Node containing all of the children of the current node which matches...
Definition: xml.cpp:915