All Classes Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
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 
37 #include <opm/common/Valgrind.hpp>
38 
39 #include <dune/common/version.hh>
40 
41 #include <array>
42 #include <cmath>
43 #include <cassert>
44 #include <cstring>
45 #include <iostream>
46 #include <algorithm>
47 
48 namespace Opm {
49 namespace DenseAd {
50 
51 template <class ValueT>
52 class Evaluation<ValueT, 1>
53 {
54 public:
56  typedef ValueT ValueType;
57 
59  static constexpr int size = 1;
60 
61 protected:
63  static constexpr int length_ = size + 1;
64 
66  static constexpr int valuepos_ = 0;
68  static constexpr int dstart_ = 1;
70  static constexpr int dend_ = length_;
71 
72 public:
74  Evaluation() : data_()
75  {}
76 
78  Evaluation(const Evaluation& other) = default;
79 
80  // create an evaluation which represents a constant function
81  //
82  // i.e., f(x) = c. this implies an evaluation with the given value and all
83  // derivatives being zero.
84  template <class RhsValueType>
85  Evaluation(const RhsValueType& c)
86  {
87  setValue( c );
88  clearDerivatives();
89  Valgrind::CheckDefined( data_ );
90  }
91 
92  // create an evaluation which represents a constant function
93  //
94  // i.e., f(x) = c. this implies an evaluation with the given value and all
95  // derivatives being zero.
96  template <class RhsValueType>
97  Evaluation(const RhsValueType& c, int varPos)
98  {
99  // The variable position must be in represented by the given variable descriptor
100  assert(0 <= varPos && varPos < size);
101 
102  setValue( c );
103  clearDerivatives();
104 
105  data_[varPos + dstart_] = 1.0;
106  Valgrind::CheckDefined(data_);
107  }
108 
109  // set all derivatives to zero
110  void clearDerivatives()
111  {
112  data_[1] = 0.0;
113  }
114 
115  // create a function evaluation for a "naked" depending variable (i.e., f(x) = x)
116  template <class RhsValueType>
117  static Evaluation createVariable(const RhsValueType& value, int varPos)
118  {
119  // copy function value and set all derivatives to 0, except for the variable
120  // which is represented by the value (which is set to 1.0)
121  return Evaluation( value, varPos );
122  }
123 
124  // "evaluate" a constant function (i.e. a function that does not depend on the set of
125  // relevant variables, f(x) = c).
126  template <class RhsValueType>
127  static Evaluation createConstant(const RhsValueType& value)
128  {
129  return Evaluation( value );
130  }
131 
132  // print the value and the derivatives of the function evaluation
133  void print(std::ostream& os = std::cout) const
134  {
135  // print value
136  os << "v: " << value() << " / d:";
137 
138  // print derivatives
139  for (int varIdx = 0; varIdx < size; ++varIdx) {
140  os << " " << derivative(varIdx);
141  }
142  }
143 
144  // copy all derivatives from other
145  void copyDerivatives(const Evaluation& other)
146  {
147  data_[1] = other.data_[1];
148  }
149 
150 
151  // add value and derivatives from other to this values and derivatives
152  Evaluation& operator+=(const Evaluation& other)
153  {
154  data_[0] += other.data_[0];
155  data_[1] += other.data_[1];
156 
157  return *this;
158  }
159 
160  // add value from other to this values
161  template <class RhsValueType>
162  Evaluation& operator+=(const RhsValueType& other)
163  {
164  // value is added, derivatives stay the same
165  data_[valuepos_] += other;
166 
167  return *this;
168  }
169 
170  // subtract other's value and derivatives from this values
171  Evaluation& operator-=(const Evaluation& other)
172  {
173  data_[0] -= other.data_[0];
174  data_[1] -= other.data_[1];
175 
176  return *this;
177  }
178 
179  // subtract other's value from this values
180  template <class RhsValueType>
181  Evaluation& operator-=(const RhsValueType& other)
182  {
183  // for constants, values are subtracted, derivatives stay the same
184  data_[ valuepos_ ] -= other;
185 
186  return *this;
187  }
188 
189  // multiply values and apply chain rule to derivatives: (u*v)' = (v'u + u'v)
190  Evaluation& operator*=(const Evaluation& other)
191  {
192  // while the values are multiplied, the derivatives follow the product rule,
193  // i.e., (u*v)' = (v'u + u'v).
194  const ValueType u = this->value();
195  const ValueType v = other.value();
196 
197  // value
198  data_[valuepos_] *= v ;
199 
200  // derivatives
201  data_[1] = data_[1] * v + other.data_[1] * u;
202 
203  return *this;
204  }
205 
206  // m(c*u)' = c*u'
207  template <class RhsValueType>
208  Evaluation& operator*=(const RhsValueType& other)
209  {
210  data_[0] *= other;
211  data_[1] *= other;
212 
213  return *this;
214  }
215 
216  // m(u*v)' = (vu' - uv')/v^2
217  Evaluation& operator/=(const Evaluation& other)
218  {
219  // values are divided, derivatives follow the rule for division, i.e., (u/v)' = (v'u -
220  // u'v)/v^2.
221  ValueType& u = data_[ valuepos_ ];
222  const ValueType& v = other.value();
223  data_[1] = (v*data_[1] - u*other.data_[1])/(v*v);
224  u /= v;
225 
226  return *this;
227  }
228 
229  // divide value and derivatives by value of other
230  template <class RhsValueType>
231  Evaluation& operator/=(const RhsValueType& other)
232  {
233  const ValueType tmp = 1.0/other;
234 
235  data_[0] *= tmp;
236  data_[1] *= tmp;
237 
238  return *this;
239  }
240 
241  // add two evaluation objects
242  Evaluation operator+(const Evaluation& other) const
243  {
244  Evaluation result(*this);
245 
246  result += other;
247 
248  return result;
249  }
250 
251  // add constant to this object
252  template <class RhsValueType>
253  Evaluation operator+(const RhsValueType& other) const
254  {
255  Evaluation result(*this);
256 
257  result += other;
258 
259  return result;
260  }
261 
262  // subtract two evaluation objects
263  Evaluation operator-(const Evaluation& other) const
264  {
265  Evaluation result(*this);
266 
267  result -= other;
268 
269  return result;
270  }
271 
272  // subtract constant from evaluation object
273  template <class RhsValueType>
274  Evaluation operator-(const RhsValueType& other) const
275  {
276  Evaluation result(*this);
277 
278  result -= other;
279 
280  return result;
281  }
282 
283  // negation (unary minus) operator
284  Evaluation operator-() const
285  {
286  Evaluation result;
287 
288  // set value and derivatives to negative
289  result.data_[0] = - data_[0];
290  result.data_[1] = - data_[1];
291 
292  return result;
293  }
294 
295  Evaluation operator*(const Evaluation& other) const
296  {
297  Evaluation result(*this);
298 
299  result *= other;
300 
301  return result;
302  }
303 
304  template <class RhsValueType>
305  Evaluation operator*(const RhsValueType& other) const
306  {
307  Evaluation result(*this);
308 
309  result *= other;
310 
311  return result;
312  }
313 
314  Evaluation operator/(const Evaluation& other) const
315  {
316  Evaluation result(*this);
317 
318  result /= other;
319 
320  return result;
321  }
322 
323  template <class RhsValueType>
324  Evaluation operator/(const RhsValueType& other) const
325  {
326  Evaluation result(*this);
327 
328  result /= other;
329 
330  return result;
331  }
332 
333  template <class RhsValueType>
334  Evaluation& operator=(const RhsValueType& other)
335  {
336  setValue( other );
337  clearDerivatives();
338 
339  return *this;
340  }
341 
342  // copy assignment from evaluation
343  Evaluation& operator=(const Evaluation& other) = default;
344 
345  template <class RhsValueType>
346  bool operator==(const RhsValueType& other) const
347  { return value() == other; }
348 
349  bool operator==(const Evaluation& other) const
350  {
351  for (int idx = 0; idx < length_; ++idx) {
352  if (data_[idx] != other.data_[idx]) {
353  return false;
354  }
355  }
356  return true;
357  }
358 
359  bool operator!=(const Evaluation& other) const
360  { return !operator==(other); }
361 
362  template <class RhsValueType>
363  bool operator>(RhsValueType other) const
364  { return value() > other; }
365 
366  bool operator>(const Evaluation& other) const
367  { return value() > other.value(); }
368 
369  template <class RhsValueType>
370  bool operator<(RhsValueType other) const
371  { return value() < other; }
372 
373  bool operator<(const Evaluation& other) const
374  { return value() < other.value(); }
375 
376  template <class RhsValueType>
377  bool operator>=(RhsValueType other) const
378  { return value() >= other; }
379 
380  bool operator>=(const Evaluation& other) const
381  { return value() >= other.value(); }
382 
383  template <class RhsValueType>
384  bool operator<=(RhsValueType other) const
385  { return value() <= other; }
386 
387  bool operator<=(const Evaluation& other) const
388  { return value() <= other.value(); }
389 
390  // return value of variable
391  const ValueType& value() const
392  { return data_[valuepos_]; }
393 
394  // set value of variable
395  template <class RhsValueType>
396  void setValue(const RhsValueType& val)
397  { data_[valuepos_] = val; }
398 
399  // return varIdx'th derivative
400  const ValueType& derivative(int varIdx) const
401  {
402  assert(0 <= varIdx && varIdx < size);
403 
404  return data_[dstart_ + varIdx];
405  }
406 
407  // set derivative at position varIdx
408  void setDerivative(int varIdx, const ValueType& derVal)
409  {
410  assert(0 <= varIdx && varIdx < size);
411 
412  data_[dstart_ + varIdx] = derVal;
413  }
414 
415 private:
416  std::array<ValueT, length_> data_;
417 };
418 
419 } } // namespace DenseAd, Opm
420 
421 #endif // OPM_DENSEAD_EVALUATION1_HPP
Evaluation()
default constructor
Definition: Evaluation.hpp:79
ValueT ValueType
field type
Definition: Evaluation1.hpp:56
A number of commonly used algebraic functions for the localized OPM automatic differentiation (AD) fr...
ValueT ValueType
field type
Definition: Evaluation.hpp:61
Evaluation()
default constructor
Definition: Evaluation1.hpp:74
static constexpr int valuepos_
position index for value
Definition: Evaluation.hpp:71
static constexpr int length_
length of internal data vector
Definition: Evaluation.hpp:68
static constexpr int dstart_
start index for derivatives
Definition: Evaluation.hpp:73
Represents a function evaluation and its derivatives w.r.t.
Definition: Evaluation.hpp:57
static constexpr int size
number of derivatives
Definition: Evaluation.hpp:64
static constexpr int dend_
end+1 index for derivatives
Definition: Evaluation.hpp:75
Representation of an evaluation of a function and its derivatives w.r.t.