Reaktoro 
A unified framework for modeling chemically reactive systems
ChemicalVector.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/Common/ChemicalScalar.hpp>
22 #include <Reaktoro/Common/ThermoScalar.hpp>
23 #include <Reaktoro/Common/ThermoVector.hpp>
24 #include <Reaktoro/Math/Matrix.hpp>
25 
26 namespace Reaktoro {
27 
28 // Forward declaration
29 template<typename V, typename T, typename P, typename N>
31 
40 
43 template<typename V, typename T, typename P, typename N>
45 {
46 public:
48  V val;
49 
51  T ddT;
52 
54  P ddP;
55 
57  N ddn;
58 
61 
64  explicit ChemicalVectorBase(Index nspecies)
65  : ChemicalVectorBase(nspecies, nspecies) {}
66 
70  ChemicalVectorBase(Index nrows, Index nspecies)
71  : ChemicalVectorBase(zeros(nrows), zeros(nrows), zeros(nrows), zeros(nrows, nspecies)) {}
72 
77  ChemicalVectorBase(Index nrows, Index nspecies, double val)
78  : ChemicalVectorBase(constants(nrows, val), zeros(nrows), zeros(nrows), zeros(nrows, nspecies)) {}
79 
85  ChemicalVectorBase(const V& val, const T& ddT, const P& ddP, const N& ddn)
86  : val(val), ddT(ddT), ddP(ddP), ddn(ddn) {}
87 
89  template<typename VR, typename TR, typename PR, typename NR>
91  : val(other.val), ddT(other.ddT), ddP(other.ddP), ddn(other.ddn) {}
92 
94  auto size() const -> Index
95  {
96  return val.size();
97  }
98 
102  auto resize(Index nrows, Index nspecies) -> void
103  {
104  val.resize(nrows);
105  ddT.resize(nrows);
106  ddP.resize(nrows);
107  ddn.resize(nrows, nspecies);
108  }
109 
112  auto resize(Index nspecies) -> void
113  {
114  resize(nspecies, nspecies);
115  }
116 
118  template<typename VR, typename NR>
119  auto fill(const ChemicalScalarBase<VR, NR>& other) -> void
120  {
121  val.fill(other.val);
122  ddT.fill(other.ddT);
123  ddP.fill(other.ddP);
124  for(auto i = 0; i < ddn.rows(); ++i) ddn.row(i) = other.ddn;
125  }
126 
128  template<typename VR>
129  auto fill(const ThermoScalarBase<VR>& other) -> void
130  {
131  val.fill(other.val);
132  ddT.fill(other.ddT);
133  ddP.fill(other.ddP);
134  ddn.fill(0.0);
135  }
136 
138  auto fill(double value) -> void
139  {
140  val.fill(value);
141  ddT.fill(0.0);
142  ddP.fill(0.0);
143  ddn.fill(0.0);
144  }
145 
147  template<typename VR, typename TR, typename PR, typename NR>
149  {
150  val = other.val;
151  ddT = other.ddT;
152  ddP = other.ddP;
153  ddn = other.ddn;
154  return *this;
155  }
156 
158  template<typename VR, typename NR>
160  {
161  val.fill(other.val);
162  ddT.fill(other.ddT);
163  ddP.fill(other.ddP);
164  for(auto i = 0; i < ddn.rows(); ++i) ddn.row(i) = other.ddn;
165  return *this;
166  }
167 
169  template<typename VR>
171  {
172  val.fill(other.val);
173  ddT.fill(other.ddT);
174  ddP.fill(other.ddP);
175  ddn.fill(0.0);
176  return *this;
177  }
178 
180  auto operator=(double other) -> ChemicalVectorBase&
181  {
182  val.fill(other);
183  ddT.fill(0.0);
184  ddP.fill(0.0);
185  ddn.fill(0.0);
186  return *this;
187  }
188 
190  template<typename VR, typename TR, typename PR, typename NR>
192  {
193  val += other.val;
194  ddT += other.ddT;
195  ddP += other.ddP;
196  ddn += other.ddn;
197  return *this;
198  }
199 
201  template<typename VR, typename TR, typename PR>
203  {
204  val += other.val;
205  ddT += other.ddT;
206  ddP += other.ddP;
207  return *this;
208  }
209 
211  template<typename VR>
213  {
214  val.array() += other.val;
215  ddT.array() += other.ddT;
216  ddP.array() += other.ddP;
217  return *this;
218  }
219 
221  auto operator+=(double other) -> ChemicalVectorBase&
222  {
223  val.array() += other;
224  return *this;
225  }
226 
228  template<typename VR, typename TR, typename PR, typename NR>
230  {
231  val -= other.val;
232  ddT -= other.ddT;
233  ddP -= other.ddP;
234  ddn -= other.ddn;
235  return *this;
236  }
237 
239  template<typename VR, typename TR, typename PR>
241  {
242  val -= other.val;
243  ddT -= other.ddT;
244  ddP -= other.ddP;
245  return *this;
246  }
247 
249  template<typename VR>
251  {
252  val.array() -= other.val;
253  ddT.array() -= other.ddT;
254  ddP.array() -= other.ddP;
255  return *this;
256  }
257 
259  auto operator-=(double other) -> ChemicalVectorBase&
260  {
261  val.array() -= other;
262  return *this;
263  }
264 
266  template<typename VR, typename TR, typename PR, typename NR>
268  {
269  ddT = diag(val) * other.ddT + diag(other.val) * ddT;
270  ddP = diag(val) * other.ddP + diag(other.val) * ddP;
271  ddn = diag(val) * other.ddn + diag(other.val) * ddn;
272  val = diag(val) * other.val;
273  return *this;
274  }
275 
277  auto operator*=(double other) -> ChemicalVectorBase&
278  {
279  val *= other;
280  ddT *= other;
281  ddP *= other;
282  ddn *= other;
283  return *this;
284  }
285 
287  auto operator/=(double other) -> ChemicalVectorBase&
288  {
289  *this *= 1.0/other;
290  return *this;
291  }
292 
294  auto operator[](Index irow) -> ChemicalScalarBase<double&, decltype(ddn.row(irow).transpose())>
295  {
296  return {val[irow], ddT[irow], ddP[irow], ddn.row(irow).transpose()};
297  }
298 
300  auto operator[](Index irow) const -> ChemicalScalarBase<const double&, decltype(ddn.row(irow).transpose())>
301  {
302  return {val[irow], ddT[irow], ddP[irow], ddn.row(irow).transpose()};
303  }
304 
306  explicit operator Vector() const
307  {
308  return val;
309  }
310 };
311 
313 template<typename Derived>
314 class Composition : public ChemicalVectorBase<const Eigen::MatrixBase<Derived>&, decltype(zeros(0)), decltype(zeros(0)), decltype(identity(0,0))>
315 {
316 public:
318  using Base = ChemicalVectorBase<const Eigen::MatrixBase<Derived>&, decltype(zeros(0)), decltype(zeros(0)), decltype(identity(0,0))>;
319 
321  Composition(const Eigen::MatrixBase<Derived>& n) : Base(n, zeros(n.rows()), zeros(n.rows()), identity(n.rows(), n.rows())) {}
322 };
323 
325 template<typename Derived>
326 auto composition(const Eigen::MatrixBase<Derived>& n) -> Composition<Derived>
327 {
328  return Composition<Derived>(n);
329 }
330 
331 template<typename V, typename T, typename P, typename N>
332 auto operator<<(std::ostream& out, const ChemicalVectorBase<V,T,P,N>& r) -> std::ostream&
333 {
334  out << r.val;
335  return out;
336 }
337 
338 template<typename V, typename T, typename P, typename N>
340 {
341  return r;
342 }
343 
344 template<typename V, typename T, typename P, typename N>
345 auto operator-(const ChemicalVectorBase<V,T,P,N>& r) -> ChemicalVectorBase<decltype(-r.val), decltype(-r.ddT), decltype(-r.ddP), decltype(-r.ddn)>
346 {
347  return {-r.val, -r.ddT, -r.ddP, -r.ddn};
348 }
349 
350 template<typename VL, typename TL, typename PL, typename NL, typename VR, typename TR, typename PR, typename NR>
351 auto operator+(const ChemicalVectorBase<VL,TL,PL,NL>& l, const ChemicalVectorBase<VR,TR,PR,NR>& r) -> ChemicalVectorBase<decltype(l.val + r.val), decltype(l.ddT + r.ddT), decltype(l.ddP + r.ddP), decltype(l.ddn + r.ddn)>
352 {
353  return {l.val + r.val, l.ddT + r.ddT, l.ddP + r.ddP, l.ddn + r.ddn};
354 }
355 
356 template<typename VL, typename TL, typename PL, typename NL, typename VR, typename TR, typename PR>
357 auto operator+(const ChemicalVectorBase<VL,TL,PL,NL>& l, const ThermoVectorBase<VR,TR,PR>& r) -> ChemicalVectorBase<decltype(l.val + r.val), decltype(l.ddT + r.ddT), decltype(l.ddP + r.ddP), decltype(l.ddn)>
358 {
359  return {l.val + r.val, l.ddT + r.ddT, l.ddP + r.ddP, l.ddn};
360 }
361 
362 template<typename VL, typename TL, typename PL, typename VR, typename TR, typename PR, typename NR>
363 auto operator+(const ThermoVectorBase<VL,TL,PL>& l, const ChemicalVectorBase<VR,TR,PR,NR>& r) -> decltype(r + l)
364 {
365  return r + l;
366 }
367 
368 template<typename VL, typename TL, typename PL, typename NL, typename VR>
369 auto operator+(const ChemicalVectorBase<VL,TL,PL,NL>& l, const ThermoScalarBase<VR>& r) -> ChemicalVectorBase<decltype(l.val + r.val*ones(l.size())), decltype(l.ddT + r.ddT*ones(l.size())), decltype(l.ddP + r.ddP*ones(l.size())), decltype(l.ddn)>
370 {
371  return {l.val + r.val*ones(l.size()), l.ddT + r.ddT*ones(l.size()), l.ddP + r.ddP*ones(l.size()), l.ddn};
372 }
373 
374 template<typename VL, typename VR, typename TR, typename PR, typename NR>
375 auto operator+(const ThermoScalarBase<VL>& l, const ChemicalVectorBase<VR,TR,PR,NR>& r) -> decltype(r + l)
376 {
377  return r + l;
378 }
379 
380 template<typename V, typename T, typename P, typename N>
381 auto operator+(const ChemicalVectorBase<V,T,P,N>& l, const Vector& r) -> ChemicalVectorBase<decltype(l.val + r),T,P,N>
382 {
383  return {l.val + r, l.ddT, l.ddP, l.ddn};
384 }
385 
386 template<typename V, typename T, typename P, typename N>
387 auto operator+(const Vector& l, const ChemicalVectorBase<V,T,P,N>& r) -> decltype(r + l)
388 {
389  return r + l;
390 }
391 
392 template<typename VL, typename TL, typename PL, typename NL, typename VR, typename TR, typename PR, typename NR>
393 auto operator-(const ChemicalVectorBase<VL,TL,PL,NL>& l, const ChemicalVectorBase<VR,TR,PR,NR>& r) -> ChemicalVectorBase<decltype(l.val - r.val), decltype(l.ddT - r.ddT), decltype(l.ddP - r.ddP), decltype(l.ddn - r.ddn)>
394 {
395  return {l.val - r.val, l.ddT - r.ddT, l.ddP - r.ddP, l.ddn - r.ddn};
396 }
397 
398 template<typename VL, typename TL, typename PL, typename NL, typename VR, typename TR, typename PR>
399 auto operator-(const ChemicalVectorBase<VL,TL,PL,NL>& l, const ThermoVectorBase<VR,TR,PR>& r) -> ChemicalVectorBase<decltype(l.val - r.val), decltype(l.ddT - r.ddT), decltype(l.ddP - r.ddP), decltype(l.ddn)>
400 {
401  return {l.val - r.val, l.ddT - r.ddT, l.ddP - r.ddP, l.ddn};
402 }
403 
404 template<typename VL, typename TL, typename PL, typename VR, typename TR, typename PR, typename NR>
405 auto operator-(const ThermoVectorBase<VL,TL,PL>& l, const ChemicalVectorBase<VR,TR,PR,NR>& r) -> decltype(-(r - l))
406 {
407  return -(r - l);
408 }
409 
410 template<typename VL, typename TL, typename PL, typename NL, typename VR>
411 auto operator-(const ChemicalVectorBase<VL,TL,PL,NL>& l, const ThermoScalarBase<VR>& r) -> ChemicalVectorBase<decltype(l.val - r.val*ones(l.size())), decltype(l.ddT - r.ddT*ones(l.size())), decltype(l.ddP - r.ddP*ones(l.size())), decltype(l.ddn)>
412 {
413  return {l.val - r.val*ones(l.size()), l.ddT - r.ddT*ones(l.size()), l.ddP - r.ddP*ones(l.size()), l.ddn};
414 }
415 
416 template<typename VL, typename VR, typename TR, typename PR, typename NR>
417 auto operator-(const ThermoScalarBase<VL>& l, const ChemicalVectorBase<VR,TR,PR,NR>& r) -> decltype(r + l)
418 {
419  return -(r - l);
420 }
421 
422 template<typename V, typename T, typename P, typename N>
423 auto operator-(const ChemicalVectorBase<V,T,P,N>& l, const Vector& r) -> ChemicalVectorBase<decltype(l.val - r),T,P,N>
424 {
425  return {l.val - r, l.ddT, l.ddP, l.ddn};
426 }
427 
428 template<typename V, typename T, typename P, typename N>
429 auto operator-(const Vector& l, const ChemicalVectorBase<V,T,P,N>& r) -> decltype(-(r - l))
430 {
431  return -(r - l);
432 }
433 
434 template<typename V, typename T, typename P, typename N>
435 auto operator*(double l, const ChemicalVectorBase<V,T,P,N>& r) -> ChemicalVectorBase<decltype(l * r.val), decltype(l * r.ddT), decltype(l * r.ddP), decltype(l * r.ddn)>
436 {
437  return {l * r.val, l * r.ddT, l * r.ddP, l * r.ddn};
438 }
439 
440 template<typename V, typename T, typename P, typename N>
441 auto operator*(const ChemicalVectorBase<V,T,P,N>& l, double r) -> decltype(r * l)
442 {
443  return r * l;
444 }
445 
446 template<typename VL, typename VR, typename T, typename P, typename N>
447 auto operator*(const ThermoScalarBase<VL>& l, const ChemicalVectorBase<VR,T,P,N>& r) -> ChemicalVectorBase<decltype(l.val * r.val), decltype(l.val * r.ddT + l.ddT * r.val), decltype(l.val * r.ddP + l.ddP * r.val), decltype(l.val * r.ddn)>
448 {
449  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};
450 }
451 
452 template<typename VL, typename VR, typename T, typename P, typename N>
453 auto operator*(const ChemicalVectorBase<VL,T,P,N>& l, const ThermoScalarBase<VR>& r) -> decltype(r * l)
454 {
455  return r * l;
456 }
457 
458 template<typename VL, typename TL, typename PL, typename NL, typename VR, typename NR>
459 auto operator*(const ChemicalVectorBase<VL,TL,PL,NL>& l, const ChemicalScalarBase<VR,NR>& r) -> ChemicalVectorBase<decltype(l.val * r.val), decltype(l.val * r.ddT + l.ddT * r.val), decltype(l.val * r.ddP + l.ddP * r.val), decltype(l.val * tr(r.ddn) + l.ddn * r.val)>
460 {
461  return {l.val * r.val, l.val * r.ddT + l.ddT * r.val, l.val * r.ddP + l.ddP * r.val, l.val * tr(r.ddn) + l.ddn * r.val};
462 }
463 
464 template<typename VL, typename NL, typename VR, typename TR, typename PR, typename NR>
465 auto operator*(const ChemicalScalarBase<VL,NL>& l, const ChemicalVectorBase<VR,TR,PR,NR>& r) -> decltype(r * l)
466 {
467  return r * l;
468 }
469 
470 template<typename VL, typename TL, typename PL, typename NL, typename VR, typename TR, typename PR, typename NR>
471 auto operator%(const ChemicalVectorBase<VL,TL,PL,NL>& l, const ChemicalVectorBase<VR,TR,PR,NR>& r) -> ChemicalVectorBase<decltype(diag(l.val) * r.val), decltype(diag(l.val) * r.ddT + diag(r.val) * l.ddT), decltype(diag(l.val) * r.ddP + diag(r.val) * l.ddP), decltype(diag(l.val) * r.ddn + diag(r.val) * l.ddn)>
472 {
473  return {diag(l.val) * r.val,
474  diag(l.val) * r.ddT + diag(r.val) * l.ddT,
475  diag(l.val) * r.ddP + diag(r.val) * l.ddP,
476  diag(l.val) * r.ddn + diag(r.val) * l.ddn};
477 }
478 
479 template<typename VL, typename TL, typename PL, typename VR, typename TR, typename PR, typename NR>
480 auto operator%(const ThermoVectorBase<VL,TL,PL>& l, const ChemicalVectorBase<VR,TR,PR,NR>& r) -> ChemicalVectorBase<decltype(diag(l.val) * r.val), decltype(diag(l.val) * r.ddT + diag(r.val) * l.ddT), decltype(diag(l.val) * r.ddP + diag(r.val) * l.ddP), decltype(diag(l.val) * r.ddn)>
481 {
482  return {diag(l.val) * r.val,
483  diag(l.val) * r.ddT + diag(r.val) * l.ddT,
484  diag(l.val) * r.ddP + diag(r.val) * l.ddP,
485  diag(l.val) * r.ddn};
486 }
487 
488 template<typename VL, typename TL, typename PL, typename NL, typename VR, typename TR, typename PR>
489 auto operator%(const ChemicalVectorBase<VL,TL,PL,NL>& l, const ThermoVectorBase<VR,TR,PR>& r) -> decltype(r % l)
490 {
491  return r % l;
492 }
493 
494 template<typename V, typename T, typename P, typename N>
495 auto operator%(const Vector& l, const ChemicalVectorBase<V,T,P,N>& r) -> ChemicalVectorBase<decltype(diag(l) * r.val), decltype(diag(l) * r.ddT), decltype(diag(l) * r.ddP), decltype(diag(l) * r.ddn)>
496 {
497  return {diag(l) * r.val, diag(l) * r.ddT, diag(l) * r.ddP, diag(l) * r.ddn};
498 }
499 
500 template<typename VL, typename TL, typename PL, typename NL>
501 auto operator%(const ChemicalVectorBase<VL,TL,PL,NL>& l, const Vector& r) -> decltype(r % l)
502 {
503  return r % l;
504 }
505 
506 template<typename VL, typename TL, typename PL, typename NL, typename VR, typename TR, typename PR, typename NR>
508 {
509  const Vector tmp = 1.0/(r.val % r.val);
510  return {l.val/r.val,
511  diag(tmp) * (diag(r.val) * l.ddT - diag(l.val) * r.ddT),
512  diag(tmp) * (diag(r.val) * l.ddP - diag(l.val) * r.ddP),
513  diag(tmp) * (diag(r.val) * l.ddn - diag(l.val) * r.ddn)};
514 }
515 
516 template<typename VL, typename TL, typename PL, typename NL, typename VR, typename TR, typename PR, typename NR>
518 {
519  const Vector tmp = 1.0/(r.val % r.val);
520  return {l.val/r.val,
521  diag(tmp) * (diag(r.val) * l.ddT - diag(l.val) * r.ddT),
522  diag(tmp) * (diag(r.val) * l.ddP - diag(l.val) * r.ddP),
523  diag(tmp) * (diag(r.val) * l.ddn)};
524 }
525 
526 template<typename VL, typename TL, typename PL, typename NL, typename VR, typename NR>
528 {
529  const double tmp = 1.0/(r.val * r.val);
530  return {l.val/r.val,
531  tmp * (r.val * l.ddT - l.val * r.ddT),
532  tmp * (r.val * l.ddP - l.val * r.ddP),
533  tmp * (r.val * l.ddn - l.val * tr(r.ddn))};
534 }
535 
536 template<typename VL, typename VR, typename T, typename P, typename N>
537 auto operator/(const ChemicalVectorBase<VL,T,P,N>& l, const ThermoScalarBase<VR>& r) -> ChemicalVectorBase<decltype(l.val/r.val), decltype(double() * (l.ddT * r.val - l.val * r.ddT)), decltype(double() * (l.ddP * r.val - l.val * r.ddP)), decltype(double() * (l.ddn * r.val))>
538 {
539  const double tmp = 1.0/(r.val * r.val);
540  return {l.val/r.val,
541  tmp * (l.ddT * r.val - l.val * r.ddT),
542  tmp * (l.ddP * r.val - l.val * r.ddP),
543  tmp * (l.ddn * r.val)};
544 }
545 
546 template<typename V, typename T, typename P, typename N>
547 auto operator/(const ChemicalVectorBase<V,T,P,N>& l, double r) -> decltype((1.0/r) * l)
548 {
549  return (1.0/r) * l;
550 }
551 
552 template<typename VL, typename TL, typename PL, typename NL, typename VR, typename TR, typename PR, typename NR>
553 auto operator<(const ChemicalVectorBase<VL,TL,PL,NL>& l, const ChemicalVectorBase<VR,TR,PR,NR>& r) -> bool
554 {
555  return l.val < r.val;
556 }
557 
558 template<typename VL, typename TL, typename PL, typename NL, typename VR, typename TR, typename PR, typename NR>
559 auto operator<=(const ChemicalVectorBase<VL,TL,PL,NL>& l, const ChemicalVectorBase<VR,TR,PR,NR>& r) -> bool
560 {
561  return l.val <= r.val;
562 }
563 
564 template<typename VL, typename TL, typename PL, typename NL, typename VR, typename TR, typename PR, typename NR>
565 auto operator>(const ChemicalVectorBase<VL,TL,PL,NL>& l, const ChemicalVectorBase<VR,TR,PR,NR>& r) -> bool
566 {
567  return l.val > r.val;
568 }
569 
570 template<typename VL, typename TL, typename PL, typename NL, typename VR, typename TR, typename PR, typename NR>
571 auto operator>=(const ChemicalVectorBase<VL,TL,PL,NL>& l, const ChemicalVectorBase<VR,TR,PR,NR>& r) -> bool
572 {
573  return l.val >= r.val;
574 }
575 
576 template<typename VL, typename TL, typename PL, typename NL, typename VR, typename TR, typename PR, typename NR>
577 auto operator==(const ChemicalVectorBase<VL,TL,PL,NL>& l, const ChemicalVectorBase<VR,TR,PR,NR>& r) -> bool
578 {
579  return l.val == r.val;
580 }
581 
582 template<typename VL, typename TL, typename PL, typename NL, typename VR, typename TR, typename PR, typename NR>
583 auto operator!=(const ChemicalVectorBase<VL,TL,PL,NL>& l, const ChemicalVectorBase<VR,TR,PR,NR>& r) -> bool
584 {
585  return l.val == r.val;
586 }
587 
589 template<typename V, typename T, typename P, typename N>
590 auto row(ChemicalVectorBase<V,T,P,N>& vec, Index irow) -> ChemicalScalarBase<double&, decltype(tr(vec.ddn.row(irow)))>
591 {
592  return {vec.val[irow], vec.ddT[irow], vec.ddP[irow], tr(vec.ddn.row(irow))};
593 }
594 
596 template<typename V, typename T, typename P, typename N>
597 auto row(const ChemicalVectorBase<V,T,P,N>& vec, Index irow) -> ChemicalScalarBase<const double&, decltype(tr(vec.ddn.row(irow)))>
598 {
599  return {vec.val[irow], vec.ddT[irow], vec.ddP[irow], tr(vec.ddn.row(irow))};
600 }
601 
603 template<typename V, typename T, typename P, typename N>
604 auto row(ChemicalVectorBase<V,T,P,N>& vec, Index irow, Index icol, Index ncols) -> ChemicalScalarBase<double&, decltype(vec.ddn.row(irow).segment(icol, ncols))>
605 {
606  return {vec.val[irow], vec.ddT[irow], vec.ddP[irow], vec.ddn.row(irow).segment(icol, ncols)};
607 }
608 
610 template<typename V, typename T, typename P, typename N>
611 auto row(const ChemicalVectorBase<V,T,P,N>& vec, Index irow, Index icol, Index ncols) -> ChemicalScalarBase<const double&, decltype(vec.ddn.row(irow).segment(icol, ncols))>
612 {
613  return {vec.val[irow], vec.ddT[irow], vec.ddP[irow], vec.ddn.row(irow).segment(icol, ncols)};
614 }
615 
617 template<typename V, typename T, typename P, typename N>
618 auto rows(ChemicalVectorBase<V,T,P,N>& vec, Index irow, Index nrows) -> ChemicalVectorBase<decltype(rows(vec.val, irow, nrows)), decltype(rows(vec.ddT, irow, nrows)), decltype(rows(vec.ddP, irow, nrows)), decltype(rows(vec.ddn, irow, nrows))>
619 {
620  return {rows(vec.val, irow, nrows), rows(vec.ddT, irow, nrows), rows(vec.ddP, irow, nrows), rows(vec.ddn, irow, nrows)};
621 }
622 
624 template<typename V, typename T, typename P, typename N>
625 auto rows(const ChemicalVectorBase<V,T,P,N>& vec, Index irow, Index nrows) -> ChemicalVectorBase<decltype(rows(vec.val, irow, nrows)), decltype(rows(vec.ddT, irow, nrows)), decltype(rows(vec.ddP, irow, nrows)), decltype(rows(vec.ddn, irow, nrows))>
626 {
627  return {rows(vec.val, irow, nrows), rows(vec.ddT, irow, nrows), rows(vec.ddP, irow, nrows), rows(vec.ddn, irow, nrows)};
628 }
629 
631 template<typename V, typename T, typename P, typename N>
632 auto rows(ChemicalVectorBase<V,T,P,N>& vec, Index irow, Index icol, Index nrows, Index ncols) -> ChemicalVectorBase<decltype(rows(vec.val, irow, nrows)), decltype(rows(vec.ddT, irow, nrows)), decltype(rows(vec.ddP, irow, nrows)), decltype(block(vec.ddn, irow, icol, nrows, ncols))>
633 {
634  return {rows(vec.val, irow, nrows), rows(vec.ddT, irow, nrows), rows(vec.ddP, irow, nrows), block(vec.ddn, irow, icol, nrows, ncols)};
635 }
636 
638 template<typename V, typename T, typename P, typename N>
639 auto rows(const ChemicalVectorBase<V,T,P,N>& vec, Index irow, Index icol, Index nrows, Index ncols) -> ChemicalVectorBase<decltype(rows(vec.val, irow, nrows)), decltype(rows(vec.ddT, irow, nrows)), decltype(rows(vec.ddP, irow, nrows)), decltype(block(vec.ddn, irow, icol, nrows, ncols))>
640 {
641  return {rows(vec.val, irow, nrows), rows(vec.ddT, irow, nrows), rows(vec.ddP, irow, nrows), block(vec.ddn, irow, icol, nrows, ncols)};
642 }
643 
645 template<typename V, typename T, typename P, typename N>
647 {
648  return {rows(vec.val, irows), rows(vec.ddT, irows), rows(vec.ddP, irows), rows(vec.ddn, irows)};
649 }
650 
652 template<typename V, typename T, typename P, typename N>
653 auto rows(const ChemicalVectorBase<V,T,P,N>& vec, const Indices& irows) -> ChemicalVector
654 {
655  return {rows(vec.val, irows), rows(vec.ddT, irows), rows(vec.ddP, irows), rows(vec.ddn, irows)};
656 }
657 
659 template<typename V, typename T, typename P, typename N>
660 auto rows(ChemicalVectorBase<V,T,P,N>& vec, const Indices& irows, const Indices& icols) -> ChemicalVector
661 {
662  return {rows(vec.val, irows), rows(vec.ddT, irows), rows(vec.ddP, irows), submatrix(vec.ddn, irows, icols)};
663 }
664 
666 template<typename V, typename T, typename P, typename N>
667 auto rows(const ChemicalVectorBase<V,T,P,N>& vec, const Indices& irows, const Indices& icols) -> ChemicalVector
668 {
669  return {rows(vec.val, irows), rows(vec.ddT, irows), rows(vec.ddP, irows), submatrix(vec.ddn, irows, icols)};
670 }
671 
672 template<typename V, typename T, typename P, typename N>
673 auto abs(const ChemicalVectorBase<V,T,P,N>& l) -> ChemicalVector
674 {
675  const Vector tmp1 = abs(l.val);
676  const Vector tmp2 = l.val/tmp1;
677  return {tmp1, diag(tmp2) * l.ddT, diag(tmp2) * l.ddP, diag(tmp2) * l.ddn};
678 }
679 
680 template<typename V, typename T, typename P, typename N>
681 auto sqrt(const ChemicalVectorBase<V,T,P,N>& l) -> ChemicalVector
682 {
683  const Vector tmp1 = sqrt(l.val);
684  const Vector tmp2 = 0.5 * tmp1/l.val;
685  return {tmp1, diag(tmp2) * l.ddT, diag(tmp2) * l.ddP, diag(tmp2) * l.ddn};
686 }
687 
688 template<typename V, typename T, typename P, typename N>
689 auto pow(const ChemicalVectorBase<V,T,P,N>& l, double power) -> ChemicalVector
690 {
691  const Vector tmp1 = pow(l.val, power);
692  const Vector tmp2 = power * tmp1/l.val;
693  return {tmp1, diag(tmp2) * l.ddT, diag(tmp2) * l.ddP, diag(tmp2) * l.ddn};
694 }
695 
696 template<typename V, typename T, typename P, typename N>
697 auto exp(const ChemicalVectorBase<V,T,P,N>& l) -> ChemicalVector
698 {
699  const Vector tmp = exp(l.val);
700  return {tmp, diag(tmp) * l.ddT, diag(tmp) * l.ddP, diag(tmp) * l.ddn};
701 }
702 
703 template<typename V, typename T, typename P, typename N>
704 auto log(const ChemicalVectorBase<V,T,P,N>& l) -> ChemicalVector
705 {
706  const Vector tmp1 = log(l.val);
707  const Vector tmp2 = 1.0/l.val;
708  return {tmp1, diag(tmp2) * l.ddT, diag(tmp2) * l.ddP, diag(tmp2) * l.ddn};
709 }
710 
711 template<typename V, typename T, typename P, typename N>
712 auto log10(const ChemicalVectorBase<V,T,P,N>& l) -> decltype(log(l)/double())
713 {
714  const double ln10 = 2.302585092994046;
715  return log(l)/ln10;
716 }
717 
718 template<typename V, typename T, typename P, typename N>
719 auto min(const ChemicalVectorBase<V,T,P,N>& r) -> double
720 {
721  return min(r.val);
722 }
723 
724 template<typename V, typename T, typename P, typename N>
725 auto max(const ChemicalVectorBase<V,T,P,N>& r) -> double
726 {
727  return max(r.val);
728 }
729 
730 template<typename V, typename T, typename P, typename N>
731 auto sum(const ChemicalVectorBase<V,T,P,N>& r) -> ChemicalScalarBase<double, decltype(tr(r.ddn).rowwise().sum())>
732 {
733  return {r.val.sum(), r.ddT.sum(), r.ddP.sum(), tr(r.ddn).rowwise().sum()};
734 }
735 
736 template<typename VL, typename TL, typename PL, typename NL, typename VR, typename TR, typename PR, typename NR>
737 auto dot(const ChemicalVectorBase<VL,TL,PL,NL>& l, const ChemicalVectorBase<VR,TR,PR,NR>& r) -> ChemicalScalarBase<double, decltype(tr(l.val) * tr(r.ddn) + tr(r.val) * tr(l.ddn))>
738 {
739  return {tr(l.val) * r.val, tr(l.val) * r.ddT + tr(r.val) * l.ddT, tr(l.val) * r.ddP + tr(r.val) * l.ddP, tr(l.val) * tr(r.ddn) + tr(r.val) * tr(l.ddn)};
740 }
741 
742 template<typename VL, typename TL, typename PL, typename NL, typename VR, typename TR, typename PR, typename NR>
743 auto dot(const ThermoVectorBase<VL,TL,PL>& l, const ChemicalVectorBase<VR,TR,PR,NR>& r) -> ChemicalScalarBase<double, decltype(tr(l.val) * tr(r.ddn))>
744 {
745  return {tr(l.val) * r.val, tr(l.val) * r.ddT + tr(r.val) * l.ddT, tr(l.val) * r.ddP + tr(r.val) * l.ddP, tr(l.val) * tr(r.ddn)};
746 }
747 
748 template<typename VL, typename TL, typename PL, typename NL, typename VR, typename TR, typename PR, typename NR>
749 auto dot(const ChemicalVectorBase<VR,TR,PR,NR>& l, const ThermoVectorBase<VL,TL,PL>& r) -> decltype(dot(r, l))
750 {
751  return dot(r, l);
752 }
753 
754 template<typename V, typename T, typename P, typename N>
755 auto dot(const Vector& l, const ChemicalVectorBase<V,T,P,N>& r) -> ChemicalScalarBase<double, decltype(tr(l) * tr(r.ddn))>
756 {
757  return {tr(l) * r.val, tr(l) * r.ddT, tr(l) * r.ddP, tr(l) * tr(r.ddn)};
758 }
759 
760 template<typename V, typename T, typename P, typename N>
761 auto dot(const ChemicalVectorBase<V,T,P,N>& l, const Vector& r) -> decltype(dot(r, l))
762 {
763  return dot(r, l);
764 }
765 
766 } // namespace Reaktoro
auto row(ChemicalVectorBase< V, T, P, N > &vec, Index irow) -> ChemicalScalarBase< double &, decltype(tr(vec.ddn.row(irow)))>
Return a ChemicalScalarBase with reference to the chemical scalar in a given row. ...
Definition: ChemicalVector.hpp:590
V val
The vector of values of the thermodynamic properties.
Definition: ThermoVector.hpp:45
auto fill(const ChemicalScalarBase< VR, NR > &other) -> void
Assign a ChemicalScalarBase instance to this.
Definition: ChemicalVector.hpp:119
V val
The value of the thermodynamic property.
Definition: ThermoScalar.hpp:49
auto operator+=(const ThermoVectorBase< VR, TR, PR > &other) -> ChemicalVectorBase &
Assign-addition of a ThermoVectorBase instance to this.
Definition: ChemicalVector.hpp:202
auto operator-=(double other) -> ChemicalVectorBase &
Assign-subtraction of a scalar to this.
Definition: ChemicalVector.hpp:259
A template base class to represent a vector of chemical scalars and their partial derivatives...
Definition: ChemicalVector.hpp:30
auto resize(Index nspecies) -> void
Resize this ChemicalVectorBase instance with number of rows equal the number of species.
Definition: ChemicalVector.hpp:112
ChemicalVectorBase(Index nspecies)
Construct a ChemicalVectorBase instance with number of rows equal to given number of species...
Definition: ChemicalVector.hpp:64
std::vector< Index > Indices
Define a type that represents a collection of indices.
Definition: Index.hpp:29
V val
The value of the chemical scalar.
Definition: ChemicalScalar.hpp:46
ChemicalVectorBase()
Construct a default ChemicalVectorBase instance.
Definition: ChemicalVector.hpp:60
ChemicalVectorBase(Index nrows, Index nspecies, double val)
Construct a ChemicalVectorBase instance with given number of rows and species, and a constant value...
Definition: ChemicalVector.hpp:77
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
auto diag(const Eigen::MatrixBase< Derived > &vec) -> decltype(vec.asDiagonal())
Return a diagonal matrix representation of a vector.
Definition: Matrix.hxx:175
auto operator-=(const ThermoVectorBase< VR, TR, PR > &other) -> ChemicalVectorBase &
Assign-subtraction of a ThermoVectorBase instance to this.
Definition: ChemicalVector.hpp:240
Composition(const Eigen::MatrixBase< Derived > &n)
Construct a Composition instance with given composition vector.
Definition: ChemicalVector.hpp:321
ChemicalVectorBase(const ChemicalVectorBase< VR, TR, PR, NR > &other)
Construct a ChemicalVectorBase instance from another.
Definition: ChemicalVector.hpp:90
T ddT
The vector of partial temperature derivatives of the chemical scalars.
Definition: ChemicalVector.hpp:51
T ddT
The vector of partial temperature derivatives of the thermodynamic properties.
Definition: ThermoVector.hpp:48
auto operator=(double other) -> ChemicalVectorBase &
Assign a scalar to this.
Definition: ChemicalVector.hpp:180
P ddP
The vector of partial pressure derivatives of the thermodynamic properties.
Definition: ThermoVector.hpp:51
auto operator+=(double other) -> ChemicalVectorBase &
Assign-addition of a scalar to this.
Definition: ChemicalVector.hpp:221
V ddT
The partial temperature derivative of the thermodynamic property.
Definition: ThermoScalar.hpp:52
auto operator*=(double other) -> ChemicalVectorBase &
Assign-multiplication of a scalar to this.
Definition: ChemicalVector.hpp:277
auto rows(ChemicalVectorBase< V, T, P, N > &vec, Index irow, Index nrows) -> ChemicalVectorBase< decltype(rows(vec.val, irow, nrows)), decltype(rows(vec.ddT, irow, nrows)), decltype(rows(vec.ddP, irow, nrows)), decltype(rows(vec.ddn, irow, nrows))>
Return a reference of a sequence of rows of this ChemicalVectorBase instance.
Definition: ChemicalVector.hpp:618
auto operator[](Index irow) const -> ChemicalScalarBase< const double &, decltype(ddn.row(irow).transpose())>
Return a ChemicalScalarBase with const reference to the chemical scalar in a given row...
Definition: ChemicalVector.hpp:300
A type that describes temperature in units of K.
Definition: ChemicalVector.hpp:314
V ddT
The partial temperature derivative of the chemical scalar.
Definition: ChemicalScalar.hpp:49
auto operator+=(const ChemicalVectorBase< VR, TR, PR, NR > &other) -> ChemicalVectorBase &
Assign-addition of a ChemicalVectorBase instance to this.
Definition: ChemicalVector.hpp:191
auto operator[](Index irow) -> ChemicalScalarBase< double &, decltype(ddn.row(irow).transpose())>
Return a ChemicalScalarBase with reference to the chemical scalar in a given row. ...
Definition: ChemicalVector.hpp:294
V val
The vector of chemical scalars.
Definition: ChemicalVector.hpp:48
Eigen::VectorXd Vector
Define an alias to the vector type of the Eigen library.
Definition: Matrix.hpp:384
A template base class to represent a thermodynamic scalar and its partial derivatives.
Definition: ScalarTypes.hpp:29
auto operator+=(const ThermoScalarBase< VR > &other) -> ChemicalVectorBase &
Assign-addition of a ThermoScalarBase instance to this.
Definition: ChemicalVector.hpp:212
A template base class to represent a vector of thermodynamic scalars and their partial derivatives...
Definition: ScalarTypes.hpp:32
auto block(Eigen::MatrixBase< Derived > &mat, Index irow, Index icol, Index nrows, Index ncols) -> decltype(mat.block(irow, icol, nrows, ncols))
Return a view of some rows and columns of a matrix.
Definition: Matrix.hxx:133
auto fill(const ThermoScalarBase< VR > &other) -> void
Assign a ThermoScalarBase instance to this.
Definition: ChemicalVector.hpp:129
N ddn
The partial mole derivatives of the chemical scalar.
Definition: ChemicalScalar.hpp:55
auto identity(Index rows, Index cols) -> decltype(Matrix::Identity(rows, cols))
Return an expression of an identity matrix.
Definition: Matrix.hxx:67
auto operator=(const ChemicalVectorBase< VR, TR, PR, NR > &other) -> ChemicalVectorBase &
Assign another ChemicalVectorBase instance to this.
Definition: ChemicalVector.hpp:148
auto ones(Index rows) -> decltype(Vector::Ones(rows))
Return an expression of a vector with entries equal to one.
Definition: Matrix.hxx:27
N ddn
The matrix of partial mole derivatives of the chemical scalars.
Definition: ChemicalVector.hpp:57
auto operator/=(double other) -> ChemicalVectorBase &
Assign-division of a scalar to this.
Definition: ChemicalVector.hpp:287
auto operator=(const ThermoScalarBase< VR > &other) -> ChemicalVectorBase &
Assign a ThermoScalarBase instance to this.
Definition: ChemicalVector.hpp:170
std::size_t Index
Define a type that represents an index.
Definition: Index.hpp:26
auto operator*=(const ChemicalVectorBase< VR, TR, PR, NR > &other) -> ChemicalVectorBase &
Assign-multiplication of a ChemicalVectorBase instance to this.
Definition: ChemicalVector.hpp:267
ChemicalVectorBase(Index nrows, Index nspecies)
Construct a ChemicalVectorBase instance with given number of rows and species.
Definition: ChemicalVector.hpp:70
P ddP
The vector of partial pressure derivatives of the chemical scalars.
Definition: ChemicalVector.hpp:54
ChemicalVectorBase(const V &val, const T &ddT, const P &ddP, const N &ddn)
Construct a ChemicalVectorBase instance with given values and derivatives.
Definition: ChemicalVector.hpp:85
auto tr(Eigen::MatrixBase< Derived > &mat) -> decltype(mat.transpose())
Return the transpose of the matrix.
Definition: Matrix.hxx:157
auto operator-=(const ChemicalVectorBase< VR, TR, PR, NR > &other) -> ChemicalVectorBase &
Assign-subtraction of a ChemicalVectorBase instance to this.
Definition: ChemicalVector.hpp:229
auto operator=(const ChemicalScalarBase< VR, NR > &other) -> ChemicalVectorBase &
Assign a ChemicalScalarBase instance to this.
Definition: ChemicalVector.hpp:159
auto composition(const Eigen::MatrixBase< Derived > &n) -> Composition< Derived >
Return a Composition instance for a given MatrixBase instance.
Definition: ChemicalVector.hpp:326
The namespace containing all components of the Reaktoro library.
Definition: ChemicalScalar.hpp:24
auto fill(double value) -> void
Assign a scalarsto this.
Definition: ChemicalVector.hpp:138
auto submatrix(Eigen::MatrixBase< Derived > &mat, const Indices &irows, const Indices &icols) -> Eigen::MatrixSubView< Derived, Indices >
Return a view of some rows and columns of a matrix.
Definition: Matrix.hxx:145
auto operator-=(const ThermoScalarBase< VR > &other) -> ChemicalVectorBase &
Assign-subtraction of a ThermoScalarBase instance to this.
Definition: ChemicalVector.hpp:250
auto size() const -> Index
Return the number of rows in this ChemicalVectorBase instance.
Definition: ChemicalVector.hpp:94
A template base class to represent a chemical scalar and its partial derivatives. ...
Definition: ChemicalScalar.hpp:28
auto zeros(Index rows) -> decltype(Vector::Zero(rows))
Return an expression of a zero vector.
Definition: Matrix.hxx:22
auto resize(Index nrows, Index nspecies) -> void
Resize this ChemicalVectorBase instance with new number of rows and number of species.
Definition: ChemicalVector.hpp:102