Theoretica
A C++ numerical and automatic mathematical library
Loading...
Searching...
No Matches
multidual.h
Go to the documentation of this file.
1#ifndef THEORETICA_MULTIDUAL_H
2#define THEORETICA_MULTIDUAL_H
3
4
8
9#ifndef THEORETICA_NO_PRINT
10#include <sstream>
11#include <ostream>
12#endif
13
14#include "../algebra/vec.h"
15#include "../algebra/mat.h"
16
17
18namespace theoretica {
19
25 template<unsigned int N = 0>
26 class multidual {
27
28 public:
29
32
36
37 // The template argument of the vector type used
38 static constexpr size_t vector_argument = N;
39
40
43 multidual() : a(0) {}
44
45
49
50
53 multidual(real r) : a(r), v(vec<real, N>()) {}
54
55 ~multidual() = default;
56
57
60 a = x;
61 v = vec<real, N>();
62 return *this;
63 }
64
65
67 inline real Re() const {
68 return a;
69 }
70
71
73 inline real& Re() {
74 return a;
75 }
76
77
79 inline vec<real, N> Dual() const {
80 return v;
81 }
82
83
85 inline vec<real, N>& Dual() {
86 return v;
87 }
88
89
91 inline multidual conjugate() const {
92 return multidual(a, -v);
93 }
94
95
97 inline multidual inverse() const {
98
99 if(a == 0) {
100 TH_MATH_ERROR("multidual::inverse", 0, DIV_BY_ZERO);
101 return multidual(nan(), vec<real, N>(nan(), N));
102 }
103
104 return multidual(1.0 / a, v * (-1 / (a * a)));
105 }
106
107
109 inline multidual operator+() const {
110 return multidual(a, v);
111 }
112
113
115 inline multidual operator+(const multidual& other) const {
116 return multidual(a + other.a, v + other.v);
117 }
118
119
121 inline multidual operator+(real r) const {
122 return multidual(a + r, v);
123 }
124
125
127 inline multidual operator-() const {
128 return multidual(-a, -v);
129 }
130
131
133 inline multidual operator-(const multidual& other) const {
134 return multidual(a - other.a, v - other.v);
135 }
136
137
139 inline multidual operator-(real r) const {
140 return multidual(a - r, v);
141 }
142
143
145 inline multidual operator*(const multidual& other) const {
146 return multidual(a * other.a, other.v * a + v * other.a);
147 }
148
149
151 inline multidual operator*(real r) const {
152 return multidual(a * r, v * r);
153 }
154
155
157 inline multidual operator/(const multidual& other) const {
158
159 if(a == 0) {
160 TH_MATH_ERROR("multidual::operator/", 0, DIV_BY_ZERO);
161 return multidual(nan(), vec<real, N>(nan(), N));
162 }
163
164 return multidual(a / other.a,
165 (v * other.a - other.v * a) / (other.a * other.a));
166 }
167
169 inline multidual operator/(real r) const {
170 return multidual(a / r, v / r);
171 }
172
173
176
177 a += other.a;
178 v += other.v;
179 return *this;
180 }
181
182
185
186 a += r;
187 return *this;
188 }
189
192
193 a -= other.a;
194 v -= other.v;
195 return *this;
196 }
197
200
201 a -= r;
202 return *this;
203 }
204
207
208 a = (a * other.a);
209 v = (other.v * a) + (v * other.a);
210 return *this;
211 }
212
215
216 a *= r;
217 v *= r;
218 return *this;
219 }
220
223
224 a = (a / other.a);
225 v = (v * other.a - other.v * a) / (other.a * other.a);
226 return *this;
227 }
228
231
232 if(r == 0) {
233 TH_MATH_ERROR("multidual::operator/=", 0, DIV_BY_ZERO);
234 a = nan();
235 v = vec<real, N>(nan(), N);
236 return *this;
237 }
238
239 a /= r;
240 v /= r;
241
242 return *this;
243 }
244
245
248 inline bool operator==(const multidual& other) {
249 return (a == other.a) && (v == other.v);
250 }
251
252
257 const vec<real, N>& x) {
258
259 vec<multidual<N>, N> arg;
260 arg.resize(x.size());
261
262 for (unsigned int i = 0; i < x.size(); ++i)
263 arg[i] = multidual<N>(
265 );
266
267 return arg;
268 }
269
270
274 const vec<multidual<N>, N>& v) {
275
276 vec<real, N> x;
277 x.resize(v.size());
278
279 for (unsigned int i = 0; i < N; ++i)
280 x[i] = v.get(i).Re();
281
282 return x;
283 }
284
285
290 const vec<multidual<N>, N>& v) {
291
293 J.resize(v.size(), v.size());
294
295 for (unsigned int i = 0; i < N; ++i)
296 for (unsigned int j = 0; j < N; ++j)
297 J.at(j, i) = v.get(j).Dual().get(i);
298
299 return J;
300 }
301
302
306 inline static void extract(
307 const vec<multidual<N>, N>& v,
308 vec<real, N>& x,
310
311 for (unsigned int i = 0; i < N; ++i) {
312
313 for (unsigned int j = 0; j < N; ++j)
314 J.at(j, i) = v.get(j).Dual().get(i);
315
316 x[i] = v.get(i).Re();
317 }
318 }
319
320
323 unsigned int size() const {
324 return v.size();
325 }
326
327
330 inline void resize(unsigned int size) {
331 v.resize(size);
332 }
333
334
335 // Friend operators to enable equations of the form
336 // (real) op. (multidual)
337
338 inline friend multidual operator+(real a, const multidual& d) {
339 return d + a;
340 }
341
342 inline friend multidual operator-(real a, const multidual& d) {
343 return -d + a;
344 }
345
346 inline friend multidual operator*(real a, const multidual& d) {
347 return d * a;
348 }
349
350 inline friend multidual operator/(real a, const multidual& d) {
351 return a * d.inverse();
352 }
353
354
355#ifndef THEORETICA_NO_PRINT
356
359 inline std::string to_string(const std::string& epsilon = "e") const {
360
361 std::stringstream res;
362 res << a << " + " << v << epsilon;
363
364 return res.str();
365 }
366
367
369 inline operator std::string() {
370 return to_string();
371 }
372
373
376 inline friend std::ostream& operator<<(std::ostream& out, const multidual& obj) {
377 return out << obj.to_string();
378 }
379
380#endif
381
382 };
383
384}
385
386#endif
Multidual number algebra for functions of the form .
Definition multidual.h:26
multidual operator+(const multidual &other) const
Sum two multidual numbers.
Definition multidual.h:115
multidual operator+() const
Identity (for consistency)
Definition multidual.h:109
multidual operator/(const multidual &other) const
Dual division.
Definition multidual.h:157
bool operator==(const multidual &other)
Check whether two multidual numbers have the same real and multidual parts.
Definition multidual.h:248
multidual & operator*=(real r)
Multiply this multidual number by a real number.
Definition multidual.h:214
multidual & operator-=(real r)
Subtract a real number from this multidual number.
Definition multidual.h:199
multidual & operator*=(const multidual &other)
Multiply this multidual number by another one.
Definition multidual.h:206
multidual & operator/=(const multidual &other)
Divide this multidual number by another one.
Definition multidual.h:222
multidual & operator+=(const multidual &other)
Sum a real number to this one.
Definition multidual.h:175
vec< real, N > v
The dual part of the multidimensional dual number as a real vector.
Definition multidual.h:35
multidual & operator-=(const multidual &other)
Subtract a real number from this one.
Definition multidual.h:191
friend std::ostream & operator<<(std::ostream &out, const multidual &obj)
Stream the multidual number in string representation to an output stream (std::ostream)
Definition multidual.h:376
unsigned int size() const
Get the number of independent variables associated with the multidual number.
Definition multidual.h:323
multidual operator-(real r) const
Subtract a real number from a multidual number.
Definition multidual.h:139
multidual operator*(const multidual &other) const
Multiply two multidual numbers.
Definition multidual.h:145
multidual & operator+=(real r)
Sum a real number to this multidual number.
Definition multidual.h:184
multidual operator*(real r) const
Multiply a multidual number by a real number.
Definition multidual.h:151
multidual(real r, vec< real, N > u)
Construct a multidual number from a real number and an N dimensional vector.
Definition multidual.h:48
real a
The real part of the multidimensional dual number.
Definition multidual.h:31
multidual operator/(real r) const
Divide a multidual number by a real number.
Definition multidual.h:169
multidual & operator/=(real r)
Divide a multidual number by a real number.
Definition multidual.h:230
vec< real, N > & Dual()
Access the multidual part.
Definition multidual.h:85
static vec< real, N > extract_real(const vec< multidual< N >, N > &v)
Extract the real vector from a vector of multidual numbers as a vec<real, N>
Definition multidual.h:273
static vec< multidual< N >, N > make_argument(const vec< real, N > &x)
Construct an N-dimensional vector of multidual numbers to be passed as argument to a multidual functi...
Definition multidual.h:256
vec< real, N > Dual() const
Get the multidual part.
Definition multidual.h:79
multidual(real r)
Construct a multidual number from a real number.
Definition multidual.h:53
real Re() const
Get the real part.
Definition multidual.h:67
std::string to_string(const std::string &epsilon="e") const
Convert the multidual number to string representation.
Definition multidual.h:359
static mat< real, N, N > extract_dual(const vec< multidual< N >, N > &v)
Extract the dual matrix (Jacobian) from a vector of multidual numbers as a mat<N, N>
Definition multidual.h:289
multidual inverse() const
Get the inverse of a multidual number.
Definition multidual.h:97
multidual operator+(real r) const
Sum a real number to a multidual number.
Definition multidual.h:121
multidual operator-() const
Get the opposite of a multidual number.
Definition multidual.h:127
static void extract(const vec< multidual< N >, N > &v, vec< real, N > &x, mat< real, N, N > &J)
Extract the real vector and dual matrix from a vector of multidual numbers as a vec<real,...
Definition multidual.h:306
void resize(unsigned int size)
Change the size of the dual part of the number (only for dynamically allocated vectors)
Definition multidual.h:330
multidual conjugate() const
Get the multidual conjugate.
Definition multidual.h:91
real & Re()
Access the real part.
Definition multidual.h:73
multidual()
Construct a multidual number as .
Definition multidual.h:43
multidual & operator=(real x)
Initialize a multidual number from a real number.
Definition multidual.h:59
multidual operator-(const multidual &other) const
Subtract two multidual numbers.
Definition multidual.h:133
A statically allocated N-dimensional vector with elements of the given type.
Definition vec.h:92
void resize(size_t n) const
Compatibility function to allow for allocation or resizing of dynamic vectors.
Definition vec.h:416
Type get(unsigned int i) const
Getters and setters.
Definition vec.h:326
TH_CONSTEXPR unsigned int size() const
Returns the size of the vector (N)
Definition vec.h:406
#define TH_MATH_ERROR(F_NAME, VALUE, EXCEPTION)
TH_MATH_ERROR is a macro which throws exceptions or modifies errno (depending on which compiling opti...
Definition error.h:219
Main namespace of the library which contains all functions and objects.
Definition algebra.h:27
double real
A real number, defined as a floating point type.
Definition constants.h:198
std::remove_reference_t< decltype(std::declval< Structure >()[0])> vector_element_t
Extract the type of a vector (or any indexable container) from its operator[].
Definition core_traits.h:134
real nan()
Return a quiet NaN number in floating point representation.
Definition error.h:54