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