Reaktoro 
A unified framework for modeling chemically reactive systems
ChemicalScalar.hpp
1 // Reaktoro is a unified framework for modeling chemically reactive systems.
2 //
3 // Copyright (C) 2014-2015 Allan Leal
4 //
5 // This program is free software: you can redistribute it and/or modify
6 // it under the terms of the GNU General Public License as published by
7 // the Free Software Foundation, either version 3 of the License, or
8 // (at your option) any later version.
9 //
10 // This program is distributed in the hope that it will be useful,
11 // but WITHOUT ANY WARRANTY; without even the implied warranty of
12 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13 // GNU General Public License for more details.
14 //
15 // You should have received a copy of the GNU General Public License
16 // along with this program. If not, see <http://www.gnu.org/licenses/>.
17 
18 #pragma once
19 
20 // Reaktoro includes
21 #include <Reaktoro/Math/Matrix.hpp>
22 #include <Reaktoro/Common/ThermoScalar.hpp>
23 
24 namespace Reaktoro {
25 
26 // Forward declaration
27 template<typename V, typename N>
29 
36 
41 template<typename V, typename N>
43 {
44 public:
46  V val;
47 
49  V ddT;
50 
52  V ddP;
53 
55  N ddn;
56 
59  : val(), ddT(), ddP(), ddn() {}
60 
63  explicit ChemicalScalarBase(Index nspecies)
64  : ChemicalScalarBase(0.0, 0.0, 0.0, zeros(nspecies)) {}
65 
69  ChemicalScalarBase(Index nspecies, double val)
70  : ChemicalScalarBase(val, 0.0, 0.0, zeros(nspecies)) {}
71 
77  ChemicalScalarBase(const V& val, const V& ddT, const V& ddP, const N& ddn)
78  : val(val), ddT(ddT), ddP(ddP), ddn(ddn) {}
79 
81  template<typename VR, typename NR>
83  : val(other.val), ddT(other.ddT), ddP(other.ddP), ddn(other.ddn) {}
84 
86  template<typename VR, typename NR>
88  {
89  val = other.val;
90  ddT = other.ddT;
91  ddP = other.ddP;
92  ddn = other.ddn;
93  return *this;
94  }
95 
97  template<typename VR>
99  {
100  val = other.val;
101  ddT = other.ddT;
102  ddP = other.ddP;
103  return *this;
104  }
105 
107  auto operator=(double other) -> ChemicalScalarBase&
108  {
109  val = other;
110  return *this;
111  }
112 
114  template<typename VR, typename NR>
116  {
117  val += other.val;
118  ddT += other.ddT;
119  ddP += other.ddP;
120  ddn += other.ddn;
121  return *this;
122  }
123 
125  template<typename VR>
127  {
128  val += other.val;
129  ddT += other.ddT;
130  ddP += other.ddP;
131  return *this;
132  }
133 
135  auto operator+=(double other) -> ChemicalScalarBase&
136  {
137  val += other;
138  return *this;
139  }
140 
142  template<typename VR, typename NR>
144  {
145  val -= other.val;
146  ddT -= other.ddT;
147  ddP -= other.ddP;
148  ddn -= other.ddn;
149  return *this;
150  }
151 
153  template<typename VR>
155  {
156  val -= other.val;
157  ddT -= other.ddT;
158  ddP -= other.ddP;
159  return *this;
160  }
161 
163  auto operator-=(double other) -> ChemicalScalarBase&
164  {
165  val -= other;
166  return *this;
167  }
168 
170  template<typename VR, typename NR>
172  {
173  ddT = ddT * other.val + val * other.ddT;
174  ddP = ddP * other.val + val * other.ddP;
175  ddn = ddn * other.val + val * other.ddn;
176  val *= other.val;
177  return *this;
178  }
179 
181  auto operator*=(double other) -> ChemicalScalarBase&
182  {
183  val *= other;
184  ddT *= other;
185  ddP *= other;
186  ddn *= other;
187  return *this;
188  }
189 
191  auto operator/=(double other) -> ChemicalScalarBase&
192  {
193  *this *= 1.0/other;
194  return *this;
195  }
196 
198  explicit operator double() const
199  {
200  return val;
201  }
202 };
203 
209 {
210  return {value, 0.0, 0.0, unit(size, index)};
211 }
212 
213 template<typename V, typename N>
214 auto operator+(const ChemicalScalarBase<V,N>& l) -> ChemicalScalarBase<V,N>
215 {
216  return l;
217 }
218 
219 template<typename V, typename N>
220 auto operator-(const ChemicalScalarBase<V,N>& l) -> ChemicalScalarBase<double, decltype(-l.ddn)>
221 {
222  return {-l.val, -l.ddT, -l.ddP, -l.ddn};
223 }
224 
225 template<typename VL, typename NL, typename VR, typename NR>
226 auto operator+(const ChemicalScalarBase<VL,NL>& l, const ChemicalScalarBase<VR,NR>& r) -> ChemicalScalarBase<double, decltype(l.ddn + r.ddn)>
227 {
228  return {l.val + r.val, l.ddT + r.ddT, l.ddP + r.ddP, l.ddn + r.ddn};
229 }
230 
231 template<typename VL, typename NL, typename VR>
232 auto operator+(const ChemicalScalarBase<VL,NL>& l, const ThermoScalarBase<VR>& r) -> ChemicalScalarBase<double, decltype(l.ddn)>
233 {
234  return {l.val + r.val, l.ddT + r.ddT, l.ddP + r.ddP, l.ddn};
235 }
236 
237 template<typename VL, typename VR, typename NR>
238 auto operator+(const ThermoScalarBase<VL>& l, const ChemicalScalarBase<VR,NR>& r) -> decltype(r + l)
239 {
240  return r + l;
241 }
242 
243 template<typename V, typename N>
244 auto operator+(const ChemicalScalarBase<V,N>& l, double r) -> ChemicalScalarBase<double, decltype(l.ddn)>
245 {
246  return {l.val + r, l.ddT, l.ddP, l.ddn};
247 }
248 
249 template<typename V, typename N>
250 auto operator+(double l, const ChemicalScalarBase<V,N>& r) -> decltype(r + l)
251 {
252  return r + l;
253 }
254 
255 template<typename VL, typename NL, typename VR, typename NR>
256 auto operator-(const ChemicalScalarBase<VL,NL>& l, const ChemicalScalarBase<VR,NR>& r) -> ChemicalScalarBase<double, decltype(l.ddn - r.ddn)>
257 {
258  return {l.val - r.val, l.ddT - r.ddT, l.ddP - r.ddP, l.ddn - r.ddn};
259 }
260 
261 template<typename VL, typename NL, typename VR>
262 auto operator-(const ChemicalScalarBase<VL,NL>& l, const ThermoScalarBase<VR>& r) -> ChemicalScalarBase<double, decltype(l.ddn)>
263 {
264  return {l.val - r.val, l.ddT - r.ddT, l.ddP - r.ddP, l.ddn};
265 }
266 
267 template<typename VL, typename VR, typename NR>
268 auto operator-(const ThermoScalarBase<VL>& l, const ChemicalScalarBase<VR,NR>& r) -> ChemicalScalarBase<double, decltype(-r.ddn)>
269 {
270  return {l.val - r.val, l.ddT - r.ddT, l.ddP - r.ddP, -r.ddn};
271 }
272 
273 template<typename V, typename N>
274 auto operator-(const ChemicalScalarBase<V,N>& l, double r) -> ChemicalScalarBase<double, decltype(l.ddn)>
275 {
276  return {l.val - r, l.ddT, l.ddP, l.ddn};
277 }
278 
279 template<typename V, typename N>
280 auto operator-(double l, const ChemicalScalarBase<V,N>& r) -> ChemicalScalarBase<double, decltype(-r.ddn)>
281 {
282  return {l - r.val, -r.ddT, -r.ddP, -r.ddn};
283 }
284 
285 template<typename V, typename N>
286 auto operator*(double l, const ChemicalScalarBase<V,N>& r) -> ChemicalScalarBase<double, decltype(l * r.ddn)>
287 {
288  return {l * r.val, l * r.ddT, l * r.ddP, l * r.ddn};
289 }
290 
291 template<typename V, typename N>
292 auto operator*(const ChemicalScalarBase<V,N>& l, double r) -> decltype(r * l)
293 {
294  return r * l;
295 }
296 
297 template<typename VL, typename VR, typename NR>
298 auto operator*(const ThermoScalarBase<VL>& l, const ChemicalScalarBase<VR,NR>& r) -> ChemicalScalarBase<double, decltype(l.val * r.ddn)>
299 {
300  return {l.val * r.val, l.val * r.ddT + l.ddT * r.val, l.val * r.ddP + l.ddP * r.val, l.val * r.ddn};
301 }
302 
303 template<typename VL, typename NL, typename VR>
304 auto operator*(const ChemicalScalarBase<VL,NL>& l, const ThermoScalarBase<VR>& r) -> decltype(r * l)
305 {
306  return r * l;
307 }
308 
309 template<typename VL, typename NL, typename VR, typename NR>
310 auto operator*(const ChemicalScalarBase<VL,NL>& l, const ChemicalScalarBase<VR,NR>& r) -> ChemicalScalarBase<double, decltype(l.val * r.ddn + l.ddn * r.val)>
311 {
312  return {l.val * r.val, l.val * r.ddT + l.ddT * r.val, l.val * r.ddP + l.ddP * r.val, l.val * r.ddn + l.ddn * r.val};
313 }
314 
315 template<typename V, typename N>
316 auto operator/(double l, const ChemicalScalarBase<V,N>& r) -> ChemicalScalarBase<double, decltype(double() * r.ddn)>
317 {
318  const double tmp1 = 1.0/r.val;
319  const double tmp2 = -l * tmp1 * tmp1;
320  return {tmp1 * l, tmp2 * r.ddT, tmp2 * r.ddP, tmp2 * r.ddn};
321 }
322 
323 template<typename V, typename N>
324 auto operator/(const ChemicalScalarBase<V,N>& l, double r) -> decltype((1.0/r) * l)
325 {
326  return (1.0/r) * l;
327 }
328 
329 template<typename VL, typename VR, typename NR>
330 auto operator/(const ThermoScalarBase<VL>& l, const ChemicalScalarBase<VR,NR>& r) -> ChemicalScalarBase<double, decltype(-(l.val * r.ddn) * double())>
331 {
332  const double tmp = 1.0/(r.val * r.val);
333  return {l.val/r.val, (l.ddT * r.val - l.val * r.ddT) * tmp, (l.ddP * r.val - l.val * r.ddP) * tmp, -(l.val * r.ddn) * tmp};
334 }
335 
336 template<typename VL, typename NL, typename VR>
337 auto operator/(const ChemicalScalarBase<VL,NL>& l, const ThermoScalarBase<VR>& r) -> ChemicalScalarBase<double, decltype(l.ddn/r.val)>
338 {
339  const double tmp = 1.0/(r.val * r.val);
340  return {l.val/r.val, (l.ddT * r.val - l.val * r.ddT) * tmp, (l.ddP * r.val - l.val * r.ddP) * tmp, l.ddn/r.val};
341 }
342 
343 template<typename VL, typename NL, typename VR, typename NR>
344 auto operator/(const ChemicalScalarBase<VL,NL>& l, const ChemicalScalarBase<VR,NR>& r) -> ChemicalScalarBase<double, decltype((l.ddn * r.val - l.val * r.ddn) * double())>
345 {
346  const double tmp = 1.0/(r.val * r.val);
347  return {l.val/r.val, (l.ddT * r.val - l.val * r.ddT) * tmp, (l.ddP * r.val - l.val * r.ddP) * tmp, (l.ddn * r.val - l.val * r.ddn) * tmp};
348 }
349 
350 template<typename VL, typename NL, typename VR, typename NR>
351 auto operator<(const ChemicalScalarBase<VL,NL>& l, const ChemicalScalarBase<VR,NR>& r) -> bool
352 {
353  return l.val < r.val;
354 }
355 
356 template<typename V, typename N>
357 auto operator<(const ChemicalScalarBase<V,N>& l, double r) -> bool
358 {
359  return l.val < r;
360 }
361 
362 template<typename V, typename N>
363 auto operator<(double l, const ChemicalScalarBase<V,N>& r) -> bool
364 {
365  return l < r.val;
366 }
367 
368 template<typename VL, typename NL, typename VR, typename NR>
369 auto operator<=(const ChemicalScalarBase<VL,NL>& l, const ChemicalScalarBase<VR,NR>& r) -> bool
370 {
371  return l.val <= r.val;
372 }
373 
374 template<typename V, typename N>
375 auto operator<=(const ChemicalScalarBase<V,N>& l, double r) -> bool
376 {
377  return l.val <= r;
378 }
379 
380 template<typename V, typename N>
381 auto operator<=(double l, const ChemicalScalarBase<V,N>& r) -> bool
382 {
383  return l <= r.val;
384 }
385 
386 template<typename VL, typename NL, typename VR, typename NR>
387 auto operator>(const ChemicalScalarBase<VL,NL>& l, const ChemicalScalarBase<VR,NR>& r) -> bool
388 {
389  return l.val > r.val;
390 }
391 
392 template<typename V, typename N>
393 auto operator>(const ChemicalScalarBase<V,N>& l, double r) -> bool
394 {
395  return l.val > r;
396 }
397 
398 template<typename V, typename N>
399 auto operator>(double l, const ChemicalScalarBase<V,N>& r) -> bool
400 {
401  return l > r.val;
402 }
403 
404 template<typename VL, typename NL, typename VR, typename NR>
405 auto operator>=(const ChemicalScalarBase<VL,NL>& l, const ChemicalScalarBase<VR,NR>& r) -> bool
406 {
407  return l.val >= r.val;
408 }
409 
410 template<typename V, typename N>
411 auto operator>=(const ChemicalScalarBase<V,N>& l, double r) -> bool
412 {
413  return l.val >= r;
414 }
415 
416 template<typename V, typename N>
417 auto operator>=(double l, const ChemicalScalarBase<V,N>& r) -> bool
418 {
419  return l >= r.val;
420 }
421 
422 template<typename VL, typename NL, typename VR, typename NR>
423 auto operator==(const ChemicalScalarBase<VL,NL>& l, const ChemicalScalarBase<VR,NR>& r) -> bool
424 {
425  return l.val == r.val;
426 }
427 
428 template<typename V, typename N>
429 auto operator==(const ChemicalScalarBase<V,N>& l, double r) -> bool
430 {
431  return l.val == r;
432 }
433 
434 template<typename V, typename N>
435 auto operator==(double l, const ChemicalScalarBase<V,N>& r) -> bool
436 {
437  return l == r.val;
438 }
439 
440 template<typename VL, typename NL, typename VR, typename NR>
441 auto operator!=(const ChemicalScalarBase<VL,NL>& l, const ChemicalScalarBase<VR,NR>& r) -> bool
442 {
443  return l.val == r.val;
444 }
445 
446 template<typename V, typename N>
447 auto operator!=(const ChemicalScalarBase<V,N>& l, double r) -> bool
448 {
449  return l.val != r;
450 }
451 
452 template<typename V, typename N>
453 auto operator!=(double l, const ChemicalScalarBase<V,N>& r) -> bool
454 {
455  return l != r.val;
456 }
457 
458 inline auto operator<<(std::ostream& out, const ChemicalScalar& scalar) -> std::ostream&
459 {
460  out << scalar.val;
461  return out;
462 }
463 
464 template<typename V, typename N>
465 auto abs(const ChemicalScalarBase<V,N>& l) -> ChemicalScalarBase<double, decltype(double() * l.ddn)>
466 {
467  const double tmp1 = std::abs(l.val);
468  const double tmp2 = l.val/tmp1;
469  return {tmp1, tmp2 * l.ddT, tmp2 * l.ddP, tmp2 * l.ddn};
470 }
471 
472 template<typename V, typename N>
473 auto sqrt(const ChemicalScalarBase<V,N>& l) -> ChemicalScalarBase<double, decltype(double() * l.ddn)>
474 {
475  const double tmp1 = std::sqrt(l.val);
476  const double tmp2 = 0.5 * tmp1/l.val;
477  return {tmp1, tmp2 * l.ddT, tmp2 * l.ddP, tmp2 * l.ddn};
478 }
479 
480 template<typename V, typename N>
481 auto pow(const ChemicalScalarBase<V,N>& l, double power) -> ChemicalScalarBase<double, decltype(double() * l.ddn)>
482 {
483  const double tmp1 = std::pow(l.val, power);
484  const double tmp2 = power * tmp1/l.val;
485  return {tmp1, tmp2 * l.ddT, tmp2 * l.ddP, tmp2 * l.ddn};
486 }
487 
488 template<typename V, typename N>
489 auto exp(const ChemicalScalarBase<V,N>& l) -> ChemicalScalarBase<double, decltype(double() * l.ddn)>
490 {
491  const double tmp1 = std::exp(l.val);
492  return {tmp1, tmp1 * l.ddT, tmp1 * l.ddP, tmp1 * l.ddn};
493 }
494 
495 template<typename V, typename N>
496 auto log(const ChemicalScalarBase<V,N>& l) -> ChemicalScalarBase<double, decltype(double() * l.ddn)>
497 {
498  const double tmp1 = std::log(l.val);
499  const double tmp2 = 1.0/l.val;
500  return {tmp1, tmp2 * l.ddT, tmp2 * l.ddP, tmp2 * l.ddn};
501 }
502 
503 template<typename V, typename N>
504 auto log10(const ChemicalScalarBase<V,N>& l) -> decltype(log(l)/1.0)
505 {
506  const double ln10 = 2.302585092994046;
507  return log(l)/ln10;
508 }
509 
510 } // namespace Reaktoro
auto operator=(const ChemicalScalarBase< VR, NR > &other) -> ChemicalScalarBase &
Assign another ChemicalScalarBase instance to this.
Definition: ChemicalScalar.hpp:87
auto amount(double value, Index size, Index index) -> ChemicalScalarBase< double, decltype(unit(size, index))>
Return a ChemicalScalar representation of a mole amount of a species.
Definition: ChemicalScalar.hpp:208
V val
The value of the thermodynamic property.
Definition: ThermoScalar.hpp:49
auto operator+=(double other) -> ChemicalScalarBase &
Assign-addition of a scalar to this.
Definition: ChemicalScalar.hpp:135
auto index(const T &value, const std::vector< T > &values) -> Index
Find the index of a value in a container of values.
Definition: SetUtils.hxx:21
V val
The value of the chemical scalar.
Definition: ChemicalScalar.hpp:46
V ddP
The partial pressure derivative of the thermodynamic property.
Definition: ThermoScalar.hpp:55
V ddP
The partial pressure derivative of the chemical scalar.
Definition: ChemicalScalar.hpp:52
ChemicalScalarBase()
Construct a default ChemicalScalarBase instance.
Definition: ChemicalScalar.hpp:58
auto operator=(double other) -> ChemicalScalarBase &
Assign a scalar to this.
Definition: ChemicalScalar.hpp:107
V ddT
The partial temperature derivative of the thermodynamic property.
Definition: ThermoScalar.hpp:52
auto operator+=(const ThermoScalarBase< VR > &other) -> ChemicalScalarBase &
Assign-addition of a ThermoScalarBase instance to this.
Definition: ChemicalScalar.hpp:126
auto operator*=(double other) -> ChemicalScalarBase &
Assign-multiplication of a scalar to this.
Definition: ChemicalScalar.hpp:181
V ddT
The partial temperature derivative of the chemical scalar.
Definition: ChemicalScalar.hpp:49
auto operator-=(double other) -> ChemicalScalarBase &
Assign-subtraction of a scalar to this.
Definition: ChemicalScalar.hpp:163
auto operator+=(const ChemicalScalarBase< VR, NR > &other) -> ChemicalScalarBase &
Assign-addition of a ChemicalScalarBase instance to this.
Definition: ChemicalScalar.hpp:115
A template base class to represent a thermodynamic scalar and its partial derivatives.
Definition: ScalarTypes.hpp:29
N ddn
The partial mole derivatives of the chemical scalar.
Definition: ChemicalScalar.hpp:55
auto operator*=(const ChemicalScalarBase< VR, NR > &other) -> ChemicalScalarBase &
Assign-multiplication of a ChemicalScalarBase instance to this.
Definition: ChemicalScalar.hpp:171
auto unit(Index rows, Index i) -> decltype(Vector::Unit(rows, i))
Return an expression of a unit vector.
Definition: Matrix.hxx:47
ChemicalScalarBase(Index nspecies, double val)
Construct a ChemicalScalarBase instance with given number of species and a constant value...
Definition: ChemicalScalar.hpp:69
std::size_t Index
Define a type that represents an index.
Definition: Index.hpp:26
auto operator-=(const ThermoScalarBase< VR > &other) -> ChemicalScalarBase &
Assign-subtraction of a ThermoScalarBase instance to this.
Definition: ChemicalScalar.hpp:154
ChemicalScalarBase(Index nspecies)
Construct a ChemicalScalarBase instance with given number of species.
Definition: ChemicalScalar.hpp:63
auto operator=(const ThermoScalarBase< VR > &other) -> ChemicalScalarBase &
Assign a ThermoScalarBase instance to this.
Definition: ChemicalScalar.hpp:98
The namespace containing all components of the Reaktoro library.
Definition: ChemicalScalar.hpp:24
auto operator/=(double other) -> ChemicalScalarBase &
Assign-division of a scalar to this.
Definition: ChemicalScalar.hpp:191
A template base class to represent a chemical scalar and its partial derivatives. ...
Definition: ChemicalScalar.hpp:28
ChemicalScalarBase(const V &val, const V &ddT, const V &ddP, const N &ddn)
Construct a ChemicalScalarBase instance with given values and derivatives.
Definition: ChemicalScalar.hpp:77
auto operator-=(const ChemicalScalarBase< VR, NR > &other) -> ChemicalScalarBase &
Assign-subtraction of a ChemicalScalarBase instance to this.
Definition: ChemicalScalar.hpp:143
ChemicalScalarBase(const ChemicalScalarBase< VR, NR > &other)
Construct a ChemicalScalarBase instance from another.
Definition: ChemicalScalar.hpp:82
auto zeros(Index rows) -> decltype(Vector::Zero(rows))
Return an expression of a zero vector.
Definition: Matrix.hxx:22