21 #ifndef OPM_AUTODIFF_HPP_HEADER
22 #define OPM_AUTODIFF_HPP_HEADER
35 template <
typename Scalar>
44 return function(x, Scalar(0));
52 return function(x, Scalar(1));
58 function(
const Scalar x,
const Scalar dx)
64 operator +=(
const Scalar& rhs)
70 operator +=(
const AutoDiff& rhs)
77 operator -=(
const Scalar& rhs)
83 operator -=(
const AutoDiff& rhs)
90 operator *=(
const Scalar& rhs)
97 operator *=(
const AutoDiff& rhs)
99 der_ = der_*rhs.val_ + val_*rhs.der_;
104 operator /=(
const Scalar& rhs)
111 operator /=(
const AutoDiff& rhs)
113 der_ = (der_*rhs.val_ - val_*rhs.der_) / (rhs.val_ * rhs.val_);
117 template <
class Ostream>
119 print(Ostream& os)
const
121 os <<
"(x,dx) = (" << val_ <<
',' << der_ <<
")";
126 const Scalar val()
const {
return val_; }
127 const Scalar der()
const {
return der_; }
130 AutoDiff(
const Scalar x,
const Scalar dx)
139 template <
class Ostream,
typename Scalar>
141 operator<<(Ostream& os, const AutoDiff<Scalar>& fw)
146 template <
typename Scalar>
148 operator +(
const AutoDiff<Scalar>& lhs,
149 const AutoDiff<Scalar>& rhs)
151 AutoDiff<Scalar> ret = lhs;
157 template <
typename Scalar,
typename T>
159 operator +(
const T lhs,
160 const AutoDiff<Scalar>& rhs)
162 AutoDiff<Scalar> ret = rhs;
168 template <
typename Scalar,
typename T>
170 operator +(
const AutoDiff<Scalar>& lhs,
173 AutoDiff<Scalar> ret = lhs;
179 template <
typename Scalar>
181 operator -(
const AutoDiff<Scalar>& lhs,
182 const AutoDiff<Scalar>& rhs)
184 AutoDiff<Scalar> ret = lhs;
190 template <
typename Scalar,
typename T>
192 operator -(
const T lhs,
193 const AutoDiff<Scalar>& rhs)
198 template <
typename Scalar,
typename T>
200 operator -(
const AutoDiff<Scalar>& lhs,
203 AutoDiff<Scalar> ret = lhs;
209 template <
typename Scalar>
211 operator *(
const AutoDiff<Scalar>& lhs,
212 const AutoDiff<Scalar>& rhs)
214 AutoDiff<Scalar> ret = lhs;
220 template <
typename Scalar,
typename T>
222 operator *(
const T lhs,
223 const AutoDiff<Scalar>& rhs)
225 AutoDiff<Scalar> ret = rhs;
231 template <
typename Scalar,
typename T>
233 operator *(
const AutoDiff<Scalar>& lhs,
236 AutoDiff<Scalar> ret = lhs;
242 template <
typename Scalar>
244 operator /(
const AutoDiff<Scalar>& lhs,
245 const AutoDiff<Scalar>& rhs)
247 AutoDiff<Scalar> ret = lhs;
253 template <
typename Scalar,
typename T>
255 operator /(
const T lhs,
256 const AutoDiff<Scalar>& rhs)
258 Scalar a = Scalar(lhs) / rhs.val();
259 Scalar b = (-Scalar(lhs) / (rhs.val() * rhs.val())) * rhs.der();
264 template <
typename Scalar,
typename T>
266 operator /(
const AutoDiff<Scalar>& lhs,
269 Scalar a = lhs.val() / Scalar(rhs);
270 Scalar b = lhs.der() / Scalar(rhs);
275 template <
typename Scalar>
277 cos(
const AutoDiff<Scalar>& x)
279 Scalar a = std::cos(x.val());
280 Scalar b = -std::sin(x.val()) * x.der();
285 template <
typename Scalar>
287 sqrt(
const AutoDiff<Scalar>& x)
289 Scalar a = std::sqrt(x.val());
290 Scalar b = (Scalar(1.0) / (Scalar(2.0) * a)) * x.der();
static AutoDiff constant(const Scalar x)
Create an AutoDiff object representing a constant, that is, its derivative is zero.
Definition: AutoDiff.hpp:42
A simple class for forward-mode automatic differentiation.
Definition: AutoDiff.hpp:36
static AutoDiff variable(const Scalar x)
Create an AutoDiff object representing a primary variable, that is, its derivative is one...
Definition: AutoDiff.hpp:50
static AutoDiff function(const Scalar x, const Scalar dx)
Create an AutoDiff object representing a function value and its derivative.
Definition: AutoDiff.hpp:58