GCC Code Coverage Report

Directory: ../../../builds/dumux-repositories/
File: /builds/dumux-repositories/dumux/dumux/common/numericdifferentiation.hh
Date: 2024-09-21 20:52:54
Exec Total Coverage
Lines: 13 15 86.7%
Functions: 712 743 95.8%
Branches: 7 16 43.8%
Line Branch Exec Source
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 // SPDX-FileCopyrightInfo: Copyright © DuMux Project contributors, see AUTHORS.md in root folder
5 // SPDX-License-Identifier: GPL-3.0-or-later
6 //
7 /*!
8 * \file
9 * \ingroup Core
10 * \brief A class for numeric differentiation
11 *
12 */
13 #ifndef DUMUX_NUMERIC_DIFFERENTIATION_HH
14 #define DUMUX_NUMERIC_DIFFERENTIATION_HH
15
16 #include <cmath>
17 #include <limits>
18
19 namespace Dumux {
20
21 /*!
22 * \ingroup Core
23 * \brief A class for numeric differentiation with respect to a scalar parameter
24 */
25 class NumericDifferentiation
26 {
27 public:
28
29 /*!
30 * \brief Computes the epsilon used for numeric differentiation
31 * \param value The value of the variable with respect to which we are differentiating
32 * \param baseEps The step width which we are using for differentiation
33 */
34 template<class Scalar>
35 static Scalar epsilon(const Scalar value, const Scalar baseEps = 1e-10)
36 {
37
2/4
✗ Branch 0 not taken.
✓ Branch 1 taken 177927155 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 9328012 times.
 187255167 assert(std::numeric_limits<Scalar>::epsilon()*1e4 < baseEps);
38 // the epsilon value used for the numeric differentiation is
39 // now scaled by the absolute value of the primary variable...
40 using std::abs;
41 374510336 return baseEps*(abs(value) + 1.0);
42 }
43
44 /*!
45 * \brief Computes the derivative of a function with respect to a function parameter
46 * \note Overload using default epsilon computation
47 */
48 template<class Function, class Scalar, class FunctionEvalType>
49 static void partialDerivative(const Function& function, Scalar x0,
50 FunctionEvalType& derivative,
51 const FunctionEvalType& fx0,
52 const int numericDifferenceMethod = 1)
53
2/8
✓ Branch 1 taken 70288 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 70288 times.
✗ Branch 5 not taken.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
 140576 { partialDerivative(function, x0, derivative, fx0, epsilon(x0), numericDifferenceMethod); }
54
55 /*!
56 * \brief Computes the derivative of a function with respect to a function parameter
57 * \param function The function to derive
58 * \param x0 The parameter at which the derivative is ought to be evaluated
59 * \param derivative The partial derivative (output)
60 * \param fx0 The result of the function evaluated at x0
61 * \param eps The numeric epsilon used in the differentiation
62 * \param numericDifferenceMethod The numeric difference method
63 * (1: forward differences (default), 0: central differences, -1: backward differences)
64 */
65 template<class Function, class Scalar, class FunctionEvalType>
66 586567457 static void partialDerivative(const Function& function, Scalar x0,
67 FunctionEvalType& derivative,
68 const FunctionEvalType& fx0,
69 const Scalar eps,
70 const int numericDifferenceMethod = 1)
71 {
72 586567457 Scalar delta = 0.0;
73 // we are using forward or central differences, i.e. we need to calculate f(x + \epsilon)
74
1/2
✓ Branch 0 taken 350770782 times.
✗ Branch 1 not taken.
 586567457 if (numericDifferenceMethod >= 0)
75 {
76 586567457 delta += eps;
77 // calculate the function evaluated with the deflected variable
78 586567457 derivative = function(x0 + eps);
79 }
80
81 // we are using backward differences,
82 // i.e. we don't need to calculate f(x + \epsilon)
83 // we can recycle the (possibly cached) f(x)
84 ✗ else derivative = fx0;
85
86 // we are using backward or central differences,
87 // i.e. we need to calculate f(x - \epsilon)
88
2/2
✓ Branch 0 taken 292669869 times.
✓ Branch 1 taken 58100913 times.
 586567457 if (numericDifferenceMethod <= 0)
89 {
90 158358500 delta += eps;
91 // subtract the function evaluated with the deflected variable
92 194349552 derivative -= function(x0 - eps);
93 }
94
95 // we are using forward differences,
96 // i.e. we don't need to calculate f(x - \epsilon)
97 // we can recycle the (possibly cached) f(x)
98 ✗ else derivative -= fx0;
99
100 // divide difference in residuals by the magnitude of the
101 // deflections between the two function evaluation
102 787603599 derivative /= delta;
103 586567457 }
104 };
105
106 } // end namespace Dumux
107
108 #endif
109