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
46namespace Opm {
47namespace DenseAd {
48
49template <class ValueT>
50class Evaluation<ValueT, 7>
51{
52public:
55 static const int numVars = 7;
56
58 typedef ValueT ValueType;
59
61 constexpr int size() const
62 { return 7; };
63
64protected:
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
90public:
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
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
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
589private:
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