My Project
Evaluation7.hpp
Go to the documentation of this file.
1 // -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2 // vi: set et ts=4 sw=4 sts=4:
3 /*
4  This file is part of the Open Porous Media project (OPM).
5 
6  OPM is free software: you can redistribute it and/or modify
7  it under the terms of the GNU General Public License as published by
8  the Free Software Foundation, either version 2 of the License, or
9  (at your option) any later version.
10 
11  OPM is distributed in the hope that it will be useful,
12  but WITHOUT ANY WARRANTY; without even the implied warranty of
13  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  GNU General Public License for more details.
15 
16  You should have received a copy of the GNU General Public License
17  along with OPM. If not, see <http://www.gnu.org/licenses/>.
18 
19  Consult the COPYING file in the top-level source directory of this
20  module for the precise wording of the license and the list of
21  copyright holders.
22 */
31 #ifndef OPM_DENSEAD_EVALUATION7_HPP
32 #define OPM_DENSEAD_EVALUATION7_HPP
33 
34 #include "Evaluation.hpp"
35 #include "Math.hpp"
36 
38 
39 #include <array>
40 #include <cmath>
41 #include <cassert>
42 #include <cstring>
43 #include <iostream>
44 #include <algorithm>
45 
46 namespace Opm {
47 namespace DenseAd {
48 
49 template <class ValueT>
50 class Evaluation<ValueT, 7>
51 {
52 public:
55  static const int numVars = 7;
56 
58  typedef ValueT ValueType;
59 
61  constexpr int size() const
62  { return 7; };
63 
64 protected:
66  constexpr int length_() const
67  { return size() + 1; }
68 
69 
71  constexpr int valuepos_() const
72  { return 0; }
74  constexpr int dstart_() const
75  { return 1; }
77  constexpr int dend_() const
78  { return length_(); }
79 
82  void checkDefined_() const
83  {
84 #ifndef NDEBUG
85  for (const auto& v: data_)
86  Valgrind::CheckDefined(v);
87 #endif
88  }
89 
90 public:
92  Evaluation() : data_()
93  {}
94 
96  Evaluation(const Evaluation& other) = default;
97 
98 
99  // create an evaluation which represents a constant function
100  //
101  // i.e., f(x) = c. this implies an evaluation with the given value and all
102  // derivatives being zero.
103  template <class RhsValueType>
104  Evaluation(const RhsValueType& c)
105  {
106  setValue(c);
107  clearDerivatives();
108 
109  checkDefined_();
110  }
111 
112  // create an evaluation which represents a constant function
113  //
114  // i.e., f(x) = c. this implies an evaluation with the given value and all
115  // derivatives being zero.
116  template <class RhsValueType>
117  Evaluation(const RhsValueType& c, int varPos)
118  {
119  // The variable position must be in represented by the given variable descriptor
120  assert(0 <= varPos && varPos < size());
121 
122  setValue( c );
123  clearDerivatives();
124 
125  data_[varPos + dstart_()] = 1.0;
126 
127  checkDefined_();
128  }
129 
130  // set all derivatives to zero
131  void clearDerivatives()
132  {
133  data_[1] = 0.0;
134  data_[2] = 0.0;
135  data_[3] = 0.0;
136  data_[4] = 0.0;
137  data_[5] = 0.0;
138  data_[6] = 0.0;
139  data_[7] = 0.0;
140  }
141 
142  // create an uninitialized Evaluation object that is compatible with the
143  // argument, but not initialized
144  //
145  // This basically boils down to the copy constructor without copying
146  // anything. If the number of derivatives is known at compile time, this
147  // is equivalent to creating an uninitialized object using the default
148  // constructor, while for dynamic evaluations, it creates an Evaluation
149  // object which exhibits the same number of derivatives as the argument.
150  static Evaluation createBlank(const Evaluation&)
151  { return Evaluation(); }
152 
153  // create an Evaluation with value and all the derivatives to be zero
154  static Evaluation createConstantZero(const Evaluation&)
155  { return Evaluation(0.); }
156 
157  // create an Evaluation with value to be one and all the derivatives to be zero
158  static Evaluation createConstantOne(const Evaluation&)
159  { return Evaluation(1.); }
160 
161  // create a function evaluation for a "naked" depending variable (i.e., f(x) = x)
162  template <class RhsValueType>
163  static Evaluation createVariable(const RhsValueType& value, int varPos)
164  {
165  // copy function value and set all derivatives to 0, except for the variable
166  // which is represented by the value (which is set to 1.0)
167  return Evaluation(value, varPos);
168  }
169 
170  template <class RhsValueType>
171  static Evaluation createVariable(int nVars, const RhsValueType& value, int varPos)
172  {
173  if (nVars != 7)
174  throw std::logic_error("This statically-sized evaluation can only represent objects"
175  " with 7 derivatives");
176 
177  // copy function value and set all derivatives to 0, except for the variable
178  // which is represented by the value (which is set to 1.0)
179  return Evaluation(nVars, value, varPos);
180  }
181 
182  template <class RhsValueType>
183  static Evaluation createVariable(const Evaluation&, const RhsValueType& value, int varPos)
184  {
185  // copy function value and set all derivatives to 0, except for the variable
186  // which is represented by the value (which is set to 1.0)
187  return Evaluation(value, varPos);
188  }
189 
190 
191  // "evaluate" a constant function (i.e. a function that does not depend on the set of
192  // relevant variables, f(x) = c).
193  template <class RhsValueType>
194  static Evaluation createConstant(int nVars, const RhsValueType& value)
195  {
196  if (nVars != 7)
197  throw std::logic_error("This statically-sized evaluation can only represent objects"
198  " with 7 derivatives");
199  return Evaluation(value);
200  }
201 
202  // "evaluate" a constant function (i.e. a function that does not depend on the set of
203  // relevant variables, f(x) = c).
204  template <class RhsValueType>
205  static Evaluation createConstant(const RhsValueType& value)
206  {
207  return Evaluation(value);
208  }
209 
210  // "evaluate" a constant function (i.e. a function that does not depend on the set of
211  // relevant variables, f(x) = c).
212  template <class RhsValueType>
213  static Evaluation createConstant(const Evaluation&, const RhsValueType& value)
214  {
215  return Evaluation(value);
216  }
217 
218  // print the value and the derivatives of the function evaluation
219  void print(std::ostream& os = std::cout) const
220  {
221  // print value
222  os << "v: " << value() << " / d:";
223 
224  // print derivatives
225  for (int varIdx = 0; varIdx < size(); ++varIdx) {
226  os << " " << derivative(varIdx);
227  }
228  }
229 
230  // copy all derivatives from other
231  void copyDerivatives(const Evaluation& other)
232  {
233  assert(size() == other.size());
234 
235  data_[1] = other.data_[1];
236  data_[2] = other.data_[2];
237  data_[3] = other.data_[3];
238  data_[4] = other.data_[4];
239  data_[5] = other.data_[5];
240  data_[6] = other.data_[6];
241  data_[7] = other.data_[7];
242  }
243 
244 
245  // add value and derivatives from other to this values and derivatives
246  Evaluation& operator+=(const Evaluation& other)
247  {
248  assert(size() == other.size());
249 
250  data_[0] += other.data_[0];
251  data_[1] += other.data_[1];
252  data_[2] += other.data_[2];
253  data_[3] += other.data_[3];
254  data_[4] += other.data_[4];
255  data_[5] += other.data_[5];
256  data_[6] += other.data_[6];
257  data_[7] += other.data_[7];
258 
259  return *this;
260  }
261 
262  // add value from other to this values
263  template <class RhsValueType>
264  Evaluation& operator+=(const RhsValueType& other)
265  {
266  // value is added, derivatives stay the same
267  data_[valuepos_()] += other;
268 
269  return *this;
270  }
271 
272  // subtract other's value and derivatives from this values
273  Evaluation& operator-=(const Evaluation& other)
274  {
275  assert(size() == other.size());
276 
277  data_[0] -= other.data_[0];
278  data_[1] -= other.data_[1];
279  data_[2] -= other.data_[2];
280  data_[3] -= other.data_[3];
281  data_[4] -= other.data_[4];
282  data_[5] -= other.data_[5];
283  data_[6] -= other.data_[6];
284  data_[7] -= other.data_[7];
285 
286  return *this;
287  }
288 
289  // subtract other's value from this values
290  template <class RhsValueType>
291  Evaluation& operator-=(const RhsValueType& other)
292  {
293  // for constants, values are subtracted, derivatives stay the same
294  data_[valuepos_()] -= other;
295 
296  return *this;
297  }
298 
299  // multiply values and apply chain rule to derivatives: (u*v)' = (v'u + u'v)
300  Evaluation& operator*=(const Evaluation& other)
301  {
302  assert(size() == other.size());
303 
304  // while the values are multiplied, the derivatives follow the product rule,
305  // i.e., (u*v)' = (v'u + u'v).
306  const ValueType u = this->value();
307  const ValueType v = other.value();
308 
309  // value
310  data_[valuepos_()] *= v ;
311 
312  // derivatives
313  data_[1] = data_[1] * v + other.data_[1] * u;
314  data_[2] = data_[2] * v + other.data_[2] * u;
315  data_[3] = data_[3] * v + other.data_[3] * u;
316  data_[4] = data_[4] * v + other.data_[4] * u;
317  data_[5] = data_[5] * v + other.data_[5] * u;
318  data_[6] = data_[6] * v + other.data_[6] * u;
319  data_[7] = data_[7] * v + other.data_[7] * u;
320 
321  return *this;
322  }
323 
324  // m(c*u)' = c*u'
325  template <class RhsValueType>
326  Evaluation& operator*=(const RhsValueType& other)
327  {
328  data_[0] *= other;
329  data_[1] *= other;
330  data_[2] *= other;
331  data_[3] *= other;
332  data_[4] *= other;
333  data_[5] *= other;
334  data_[6] *= other;
335  data_[7] *= other;
336 
337  return *this;
338  }
339 
340  // m(u*v)' = (vu' - uv')/v^2
341  Evaluation& operator/=(const Evaluation& other)
342  {
343  assert(size() == other.size());
344 
345  // values are divided, derivatives follow the rule for division, i.e., (u/v)' = (v'u -
346  // u'v)/v^2.
347  ValueType& u = data_[valuepos_()];
348  const ValueType& v = other.value();
349  data_[1] = (v*data_[1] - u*other.data_[1])/(v*v);
350  data_[2] = (v*data_[2] - u*other.data_[2])/(v*v);
351  data_[3] = (v*data_[3] - u*other.data_[3])/(v*v);
352  data_[4] = (v*data_[4] - u*other.data_[4])/(v*v);
353  data_[5] = (v*data_[5] - u*other.data_[5])/(v*v);
354  data_[6] = (v*data_[6] - u*other.data_[6])/(v*v);
355  data_[7] = (v*data_[7] - u*other.data_[7])/(v*v);
356  u /= v;
357 
358  return *this;
359  }
360 
361  // divide value and derivatives by value of other
362  template <class RhsValueType>
363  Evaluation& operator/=(const RhsValueType& other)
364  {
365  const ValueType tmp = 1.0/other;
366 
367  data_[0] *= tmp;
368  data_[1] *= tmp;
369  data_[2] *= tmp;
370  data_[3] *= tmp;
371  data_[4] *= tmp;
372  data_[5] *= tmp;
373  data_[6] *= tmp;
374  data_[7] *= tmp;
375 
376  return *this;
377  }
378 
379  // add two evaluation objects
380  Evaluation operator+(const Evaluation& other) const
381  {
382  assert(size() == other.size());
383 
384  Evaluation result(*this);
385 
386  result += other;
387 
388  return result;
389  }
390 
391  // add constant to this object
392  template <class RhsValueType>
393  Evaluation operator+(const RhsValueType& other) const
394  {
395  Evaluation result(*this);
396 
397  result += other;
398 
399  return result;
400  }
401 
402  // subtract two evaluation objects
403  Evaluation operator-(const Evaluation& other) const
404  {
405  assert(size() == other.size());
406 
407  Evaluation result(*this);
408 
409  result -= other;
410 
411  return result;
412  }
413 
414  // subtract constant from evaluation object
415  template <class RhsValueType>
416  Evaluation operator-(const RhsValueType& other) const
417  {
418  Evaluation result(*this);
419 
420  result -= other;
421 
422  return result;
423  }
424 
425  // negation (unary minus) operator
426  Evaluation operator-() const
427  {
428  Evaluation result;
429 
430  // set value and derivatives to negative
431  result.data_[0] = - data_[0];
432  result.data_[1] = - data_[1];
433  result.data_[2] = - data_[2];
434  result.data_[3] = - data_[3];
435  result.data_[4] = - data_[4];
436  result.data_[5] = - data_[5];
437  result.data_[6] = - data_[6];
438  result.data_[7] = - data_[7];
439 
440  return result;
441  }
442 
443  Evaluation operator*(const Evaluation& other) const
444  {
445  assert(size() == other.size());
446 
447  Evaluation result(*this);
448 
449  result *= other;
450 
451  return result;
452  }
453 
454  template <class RhsValueType>
455  Evaluation operator*(const RhsValueType& other) const
456  {
457  Evaluation result(*this);
458 
459  result *= other;
460 
461  return result;
462  }
463 
464  Evaluation operator/(const Evaluation& other) const
465  {
466  assert(size() == other.size());
467 
468  Evaluation result(*this);
469 
470  result /= other;
471 
472  return result;
473  }
474 
475  template <class RhsValueType>
476  Evaluation operator/(const RhsValueType& other) const
477  {
478  Evaluation result(*this);
479 
480  result /= other;
481 
482  return result;
483  }
484 
485  template <class RhsValueType>
486  Evaluation& operator=(const RhsValueType& other)
487  {
488  setValue( other );
489  clearDerivatives();
490 
491  return *this;
492  }
493 
494  // copy assignment from evaluation
495  Evaluation& operator=(const Evaluation& other) = default;
496 
497  template <class RhsValueType>
498  bool operator==(const RhsValueType& other) const
499  { return value() == other; }
500 
501  bool operator==(const Evaluation& other) const
502  {
503  assert(size() == other.size());
504 
505  for (int idx = 0; idx < length_(); ++idx) {
506  if (data_[idx] != other.data_[idx]) {
507  return false;
508  }
509  }
510  return true;
511  }
512 
513  bool operator!=(const Evaluation& other) const
514  { return !operator==(other); }
515 
516  template <class RhsValueType>
517  bool operator!=(const RhsValueType& other) const
518  { return !operator==(other); }
519 
520  template <class RhsValueType>
521  bool operator>(RhsValueType other) const
522  { return value() > other; }
523 
524  bool operator>(const Evaluation& other) const
525  {
526  assert(size() == other.size());
527 
528  return value() > other.value();
529  }
530 
531  template <class RhsValueType>
532  bool operator<(RhsValueType other) const
533  { return value() < other; }
534 
535  bool operator<(const Evaluation& other) const
536  {
537  assert(size() == other.size());
538 
539  return value() < other.value();
540  }
541 
542  template <class RhsValueType>
543  bool operator>=(RhsValueType other) const
544  { return value() >= other; }
545 
546  bool operator>=(const Evaluation& other) const
547  {
548  assert(size() == other.size());
549 
550  return value() >= other.value();
551  }
552 
553  template <class RhsValueType>
554  bool operator<=(RhsValueType other) const
555  { return value() <= other; }
556 
557  bool operator<=(const Evaluation& other) const
558  {
559  assert(size() == other.size());
560 
561  return value() <= other.value();
562  }
563 
564  // return value of variable
565  const ValueType& value() const
566  { return data_[valuepos_()]; }
567 
568  // set value of variable
569  template <class RhsValueType>
570  void setValue(const RhsValueType& val)
571  { data_[valuepos_()] = val; }
572 
573  // return varIdx'th derivative
574  const ValueType& derivative(int varIdx) const
575  {
576  assert(0 <= varIdx && varIdx < size());
577 
578  return data_[dstart_() + varIdx];
579  }
580 
581  // set derivative at position varIdx
582  void setDerivative(int varIdx, const ValueType& derVal)
583  {
584  assert(0 <= varIdx && varIdx < size());
585 
586  data_[dstart_() + varIdx] = derVal;
587  }
588 
589 private:
590 
591  std::array<ValueT, 8> data_;
592 };
593 
594 } // namespace DenseAd
595 } // namespace Opm
596 
597 #endif // OPM_DENSEAD_EVALUATION7_HPP
Representation of an evaluation of a function and its derivatives w.r.t.
A number of commonly used algebraic functions for the localized OPM automatic differentiation (AD) fr...
Some templates to wrap the valgrind client request macros.
constexpr int size() const
number of derivatives
Definition: Evaluation7.hpp:61
Evaluation(const Evaluation &other)=default
copy other function evaluation
Evaluation()
default constructor
Definition: Evaluation7.hpp:92
constexpr int length_() const
length of internal data vector
Definition: Evaluation7.hpp:66
ValueT ValueType
field type
Definition: Evaluation7.hpp:58
constexpr int valuepos_() const
position index for value
Definition: Evaluation7.hpp:71
void checkDefined_() const
instruct valgrind to check that the value and all derivatives of the Evaluation object are well-defin...
Definition: Evaluation7.hpp:82
constexpr int dend_() const
end+1 index for derivatives
Definition: Evaluation7.hpp:77
constexpr int dstart_() const
start index for derivatives
Definition: Evaluation7.hpp:74
Represents a function evaluation and its derivatives w.r.t.
Definition: Evaluation.hpp:59
Evaluation()
default constructor
Definition: Evaluation.hpp:100
ValueT ValueType
field type
Definition: Evaluation.hpp:66
void checkDefined_() const
instruct valgrind to check that the value and all derivatives of the Evaluation object are well-defin...
Definition: Evaluation.hpp:90
static const int numVars
the template argument which specifies the number of derivatives (-1 == "DynamicSize" means runtime de...
Definition: Evaluation.hpp:63
constexpr int size() const
number of derivatives
Definition: Evaluation.hpp:69
constexpr int valuepos_() const
position index for value
Definition: Evaluation.hpp:79
constexpr int length_() const
length of internal data vector
Definition: Evaluation.hpp:74
constexpr int dstart_() const
start index for derivatives
Definition: Evaluation.hpp:82