My Project
Evaluation6.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_EVALUATION6_HPP
32 #define OPM_DENSEAD_EVALUATION6_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, 6>
51 {
52 public:
55  static const int numVars = 6;
56 
58  typedef ValueT ValueType;
59 
61  constexpr int size() const
62  { return 6; };
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  }
140 
141  // create an uninitialized Evaluation object that is compatible with the
142  // argument, but not initialized
143  //
144  // This basically boils down to the copy constructor without copying
145  // anything. If the number of derivatives is known at compile time, this
146  // is equivalent to creating an uninitialized object using the default
147  // constructor, while for dynamic evaluations, it creates an Evaluation
148  // object which exhibits the same number of derivatives as the argument.
149  static Evaluation createBlank(const Evaluation&)
150  { return Evaluation(); }
151 
152  // create an Evaluation with value and all the derivatives to be zero
153  static Evaluation createConstantZero(const Evaluation&)
154  { return Evaluation(0.); }
155 
156  // create an Evaluation with value to be one and all the derivatives to be zero
157  static Evaluation createConstantOne(const Evaluation&)
158  { return Evaluation(1.); }
159 
160  // create a function evaluation for a "naked" depending variable (i.e., f(x) = x)
161  template <class RhsValueType>
162  static Evaluation createVariable(const RhsValueType& value, int varPos)
163  {
164  // copy function value and set all derivatives to 0, except for the variable
165  // which is represented by the value (which is set to 1.0)
166  return Evaluation(value, varPos);
167  }
168 
169  template <class RhsValueType>
170  static Evaluation createVariable(int nVars, const RhsValueType& value, int varPos)
171  {
172  if (nVars != 6)
173  throw std::logic_error("This statically-sized evaluation can only represent objects"
174  " with 6 derivatives");
175 
176  // copy function value and set all derivatives to 0, except for the variable
177  // which is represented by the value (which is set to 1.0)
178  return Evaluation(nVars, value, varPos);
179  }
180 
181  template <class RhsValueType>
182  static Evaluation createVariable(const Evaluation&, const RhsValueType& value, int varPos)
183  {
184  // copy function value and set all derivatives to 0, except for the variable
185  // which is represented by the value (which is set to 1.0)
186  return Evaluation(value, varPos);
187  }
188 
189 
190  // "evaluate" a constant function (i.e. a function that does not depend on the set of
191  // relevant variables, f(x) = c).
192  template <class RhsValueType>
193  static Evaluation createConstant(int nVars, const RhsValueType& value)
194  {
195  if (nVars != 6)
196  throw std::logic_error("This statically-sized evaluation can only represent objects"
197  " with 6 derivatives");
198  return Evaluation(value);
199  }
200 
201  // "evaluate" a constant function (i.e. a function that does not depend on the set of
202  // relevant variables, f(x) = c).
203  template <class RhsValueType>
204  static Evaluation createConstant(const RhsValueType& value)
205  {
206  return Evaluation(value);
207  }
208 
209  // "evaluate" a constant function (i.e. a function that does not depend on the set of
210  // relevant variables, f(x) = c).
211  template <class RhsValueType>
212  static Evaluation createConstant(const Evaluation&, const RhsValueType& value)
213  {
214  return Evaluation(value);
215  }
216 
217  // print the value and the derivatives of the function evaluation
218  void print(std::ostream& os = std::cout) const
219  {
220  // print value
221  os << "v: " << value() << " / d:";
222 
223  // print derivatives
224  for (int varIdx = 0; varIdx < size(); ++varIdx) {
225  os << " " << derivative(varIdx);
226  }
227  }
228 
229  // copy all derivatives from other
230  void copyDerivatives(const Evaluation& other)
231  {
232  assert(size() == other.size());
233 
234  data_[1] = other.data_[1];
235  data_[2] = other.data_[2];
236  data_[3] = other.data_[3];
237  data_[4] = other.data_[4];
238  data_[5] = other.data_[5];
239  data_[6] = other.data_[6];
240  }
241 
242 
243  // add value and derivatives from other to this values and derivatives
244  Evaluation& operator+=(const Evaluation& other)
245  {
246  assert(size() == other.size());
247 
248  data_[0] += other.data_[0];
249  data_[1] += other.data_[1];
250  data_[2] += other.data_[2];
251  data_[3] += other.data_[3];
252  data_[4] += other.data_[4];
253  data_[5] += other.data_[5];
254  data_[6] += other.data_[6];
255 
256  return *this;
257  }
258 
259  // add value from other to this values
260  template <class RhsValueType>
261  Evaluation& operator+=(const RhsValueType& other)
262  {
263  // value is added, derivatives stay the same
264  data_[valuepos_()] += other;
265 
266  return *this;
267  }
268 
269  // subtract other's value and derivatives from this values
270  Evaluation& operator-=(const Evaluation& other)
271  {
272  assert(size() == other.size());
273 
274  data_[0] -= other.data_[0];
275  data_[1] -= other.data_[1];
276  data_[2] -= other.data_[2];
277  data_[3] -= other.data_[3];
278  data_[4] -= other.data_[4];
279  data_[5] -= other.data_[5];
280  data_[6] -= other.data_[6];
281 
282  return *this;
283  }
284 
285  // subtract other's value from this values
286  template <class RhsValueType>
287  Evaluation& operator-=(const RhsValueType& other)
288  {
289  // for constants, values are subtracted, derivatives stay the same
290  data_[valuepos_()] -= other;
291 
292  return *this;
293  }
294 
295  // multiply values and apply chain rule to derivatives: (u*v)' = (v'u + u'v)
296  Evaluation& operator*=(const Evaluation& other)
297  {
298  assert(size() == other.size());
299 
300  // while the values are multiplied, the derivatives follow the product rule,
301  // i.e., (u*v)' = (v'u + u'v).
302  const ValueType u = this->value();
303  const ValueType v = other.value();
304 
305  // value
306  data_[valuepos_()] *= v ;
307 
308  // derivatives
309  data_[1] = data_[1] * v + other.data_[1] * u;
310  data_[2] = data_[2] * v + other.data_[2] * u;
311  data_[3] = data_[3] * v + other.data_[3] * u;
312  data_[4] = data_[4] * v + other.data_[4] * u;
313  data_[5] = data_[5] * v + other.data_[5] * u;
314  data_[6] = data_[6] * v + other.data_[6] * u;
315 
316  return *this;
317  }
318 
319  // m(c*u)' = c*u'
320  template <class RhsValueType>
321  Evaluation& operator*=(const RhsValueType& other)
322  {
323  data_[0] *= other;
324  data_[1] *= other;
325  data_[2] *= other;
326  data_[3] *= other;
327  data_[4] *= other;
328  data_[5] *= other;
329  data_[6] *= other;
330 
331  return *this;
332  }
333 
334  // m(u*v)' = (vu' - uv')/v^2
335  Evaluation& operator/=(const Evaluation& other)
336  {
337  assert(size() == other.size());
338 
339  // values are divided, derivatives follow the rule for division, i.e., (u/v)' = (v'u -
340  // u'v)/v^2.
341  ValueType& u = data_[valuepos_()];
342  const ValueType& v = other.value();
343  data_[1] = (v*data_[1] - u*other.data_[1])/(v*v);
344  data_[2] = (v*data_[2] - u*other.data_[2])/(v*v);
345  data_[3] = (v*data_[3] - u*other.data_[3])/(v*v);
346  data_[4] = (v*data_[4] - u*other.data_[4])/(v*v);
347  data_[5] = (v*data_[5] - u*other.data_[5])/(v*v);
348  data_[6] = (v*data_[6] - u*other.data_[6])/(v*v);
349  u /= v;
350 
351  return *this;
352  }
353 
354  // divide value and derivatives by value of other
355  template <class RhsValueType>
356  Evaluation& operator/=(const RhsValueType& other)
357  {
358  const ValueType tmp = 1.0/other;
359 
360  data_[0] *= tmp;
361  data_[1] *= tmp;
362  data_[2] *= tmp;
363  data_[3] *= tmp;
364  data_[4] *= tmp;
365  data_[5] *= tmp;
366  data_[6] *= tmp;
367 
368  return *this;
369  }
370 
371  // add two evaluation objects
372  Evaluation operator+(const Evaluation& other) const
373  {
374  assert(size() == other.size());
375 
376  Evaluation result(*this);
377 
378  result += other;
379 
380  return result;
381  }
382 
383  // add constant to this object
384  template <class RhsValueType>
385  Evaluation operator+(const RhsValueType& other) const
386  {
387  Evaluation result(*this);
388 
389  result += other;
390 
391  return result;
392  }
393 
394  // subtract two evaluation objects
395  Evaluation operator-(const Evaluation& other) const
396  {
397  assert(size() == other.size());
398 
399  Evaluation result(*this);
400 
401  result -= other;
402 
403  return result;
404  }
405 
406  // subtract constant from evaluation object
407  template <class RhsValueType>
408  Evaluation operator-(const RhsValueType& other) const
409  {
410  Evaluation result(*this);
411 
412  result -= other;
413 
414  return result;
415  }
416 
417  // negation (unary minus) operator
418  Evaluation operator-() const
419  {
420  Evaluation result;
421 
422  // set value and derivatives to negative
423  result.data_[0] = - data_[0];
424  result.data_[1] = - data_[1];
425  result.data_[2] = - data_[2];
426  result.data_[3] = - data_[3];
427  result.data_[4] = - data_[4];
428  result.data_[5] = - data_[5];
429  result.data_[6] = - data_[6];
430 
431  return result;
432  }
433 
434  Evaluation operator*(const Evaluation& other) const
435  {
436  assert(size() == other.size());
437 
438  Evaluation result(*this);
439 
440  result *= other;
441 
442  return result;
443  }
444 
445  template <class RhsValueType>
446  Evaluation operator*(const RhsValueType& other) const
447  {
448  Evaluation result(*this);
449 
450  result *= other;
451 
452  return result;
453  }
454 
455  Evaluation operator/(const Evaluation& other) const
456  {
457  assert(size() == other.size());
458 
459  Evaluation result(*this);
460 
461  result /= other;
462 
463  return result;
464  }
465 
466  template <class RhsValueType>
467  Evaluation operator/(const RhsValueType& other) const
468  {
469  Evaluation result(*this);
470 
471  result /= other;
472 
473  return result;
474  }
475 
476  template <class RhsValueType>
477  Evaluation& operator=(const RhsValueType& other)
478  {
479  setValue( other );
480  clearDerivatives();
481 
482  return *this;
483  }
484 
485  // copy assignment from evaluation
486  Evaluation& operator=(const Evaluation& other) = default;
487 
488  template <class RhsValueType>
489  bool operator==(const RhsValueType& other) const
490  { return value() == other; }
491 
492  bool operator==(const Evaluation& other) const
493  {
494  assert(size() == other.size());
495 
496  for (int idx = 0; idx < length_(); ++idx) {
497  if (data_[idx] != other.data_[idx]) {
498  return false;
499  }
500  }
501  return true;
502  }
503 
504  bool operator!=(const Evaluation& other) const
505  { return !operator==(other); }
506 
507  template <class RhsValueType>
508  bool operator!=(const RhsValueType& other) const
509  { return !operator==(other); }
510 
511  template <class RhsValueType>
512  bool operator>(RhsValueType other) const
513  { return value() > other; }
514 
515  bool operator>(const Evaluation& other) const
516  {
517  assert(size() == other.size());
518 
519  return value() > other.value();
520  }
521 
522  template <class RhsValueType>
523  bool operator<(RhsValueType other) const
524  { return value() < other; }
525 
526  bool operator<(const Evaluation& other) const
527  {
528  assert(size() == other.size());
529 
530  return value() < other.value();
531  }
532 
533  template <class RhsValueType>
534  bool operator>=(RhsValueType other) const
535  { return value() >= other; }
536 
537  bool operator>=(const Evaluation& other) const
538  {
539  assert(size() == other.size());
540 
541  return value() >= other.value();
542  }
543 
544  template <class RhsValueType>
545  bool operator<=(RhsValueType other) const
546  { return value() <= other; }
547 
548  bool operator<=(const Evaluation& other) const
549  {
550  assert(size() == other.size());
551 
552  return value() <= other.value();
553  }
554 
555  // return value of variable
556  const ValueType& value() const
557  { return data_[valuepos_()]; }
558 
559  // set value of variable
560  template <class RhsValueType>
561  void setValue(const RhsValueType& val)
562  { data_[valuepos_()] = val; }
563 
564  // return varIdx'th derivative
565  const ValueType& derivative(int varIdx) const
566  {
567  assert(0 <= varIdx && varIdx < size());
568 
569  return data_[dstart_() + varIdx];
570  }
571 
572  // set derivative at position varIdx
573  void setDerivative(int varIdx, const ValueType& derVal)
574  {
575  assert(0 <= varIdx && varIdx < size());
576 
577  data_[dstart_() + varIdx] = derVal;
578  }
579 
580 private:
581 
582  std::array<ValueT, 7> data_;
583 };
584 
585 } // namespace DenseAd
586 } // namespace Opm
587 
588 #endif // OPM_DENSEAD_EVALUATION6_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.
Evaluation(const Evaluation &other)=default
copy other function evaluation
constexpr int size() const
number of derivatives
Definition: Evaluation6.hpp:61
constexpr int dend_() const
end+1 index for derivatives
Definition: Evaluation6.hpp:77
Evaluation()
default constructor
Definition: Evaluation6.hpp:92
void checkDefined_() const
instruct valgrind to check that the value and all derivatives of the Evaluation object are well-defin...
Definition: Evaluation6.hpp:82
constexpr int valuepos_() const
position index for value
Definition: Evaluation6.hpp:71
constexpr int length_() const
length of internal data vector
Definition: Evaluation6.hpp:66
constexpr int dstart_() const
start index for derivatives
Definition: Evaluation6.hpp:74
ValueT ValueType
field type
Definition: Evaluation6.hpp:58
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