ADOL-C: Automatic Differentiation in C/C++ 2.7.3
Loading...
Searching...
No Matches
adoublecuda.h
Go to the documentation of this file.
1/* ---------------------------------------------------------------------------
2 ADOL-C -- Automatic Differentiation by Overloading in C++
3
4 Revision: $Id$
5 Contents: adoublecuda.h contains the class of adouble specifically
6 suited to be used within CUDA environment
7
8 Copyright (c) Alina Koniaeva, Andrea Walther
9
10 This file is part of ADOL-C. This software is provided as open source.
11 Any use, reproduction, or distribution of the software constitutes
12 recipient's acceptance of the terms of the accompanying license file.
13
14---------------------------------------------------------------------------*/
15
16#if !defined(ADOLC_ADOUBLECUDA_H)
17#define ADOLC_ADOUBLECUDA_H 1
18
19#include <cmath>
20#include <cstdlib>
21#include <iostream>
22using std::cout;
23using std::istream;
24using std::ostream;
25
26#include <cuda_runtime.h>
27#include <math_constants.h>
28
29namespace adtlc {
30
31#if defined(NUMBER_DIRECTIONS)
32__managed__ size_t ADOLC_numDir = NUMBER_DIRECTIONS;
33#if defined(DYNAMIC_DIRECTIONS)
34#define ADVAL_DECL *adval
35#define ADVAL_TYPE_ADV const double *adv
36#else
37#define ADVAL_DECL adval[NUMBER_DIRECTIONS]
38#define ADVAL_TYPE_ADV const double adv[NUMBER_DIRECTIONS]
39#endif
40#define ADVAL_TYPE double *
41#define FOR_I_EQ_0_LT_NUMDIR for (size_t i = 0; i < ADOLC_numDir; i++)
42#define ADVAL_I adval[i]
43#define ADV_I adv[i]
44#define V_I v[i]
45#else
46#define ADVAL_DECL adval
47#define ADVAL adval
48#define ADVAL_TYPE_ADV double adv
49#define ADVAL_TYPE double
50#define FOR_I_EQ_0_LT_NUMDIR
51#define ADVAL_I adval
52#define ADV_I adv
53#define V_I v
54#endif
55
56#define ADOLC_MATH_NSP std
57
58inline __device__ double makeNaN() { return CUDART_NAN; }
59
60inline __device__ double makeInf() { return CUDART_INF; }
61
62#define CUDADEV __device__ inline
63#define CUDAHOST __host__ inline
64#define CUDAHOSTDEV __host__ __device__ inline
65
72class adouble {
73public:
74 // ctors
75 CUDADEV adouble();
76 CUDADEV adouble(const double v);
77 CUDADEV adouble(const double v, ADVAL_TYPE_ADV);
78 CUDADEV adouble(const adouble &a);
79#if defined(DYNAMIC_DIRECTIONS)
80 CUDADEV ~adouble();
81#endif
82 /******************* temporary results ******************************/
83 // sign
84 CUDADEV adouble operator-() const;
85 CUDADEV adouble operator+() const;
86
87 // addition
88 CUDADEV adouble operator+(const double v) const;
89 CUDADEV adouble operator+(const adouble &a) const;
90 CUDADEV friend adouble operator+(const double v, const adouble &a);
91
92 // subtraction
93 CUDADEV adouble operator-(const double v) const;
94 CUDADEV adouble operator-(const adouble &a) const;
95 CUDADEV friend adouble operator-(const double v, const adouble &a);
96
97 // multiplication
98 CUDADEV adouble operator*(const double v) const;
99 CUDADEV adouble operator*(const adouble &a) const;
100 CUDADEV friend adouble operator*(const double v, const adouble &a);
101
102 // division
103 CUDADEV adouble operator/(const double v) const;
104 CUDADEV adouble operator/(const adouble &a) const;
105 CUDADEV friend adouble operator/(const double v, const adouble &a);
106
107 // inc/dec
108 CUDADEV adouble operator++();
109 CUDADEV adouble operator++(int);
110 CUDADEV adouble operator--();
111 CUDADEV adouble operator--(int);
112
113 // functions
114 CUDADEV friend adouble tan(const adouble &a);
115 CUDADEV friend adouble exp(const adouble &a);
116 CUDADEV friend adouble log(const adouble &a);
117 CUDADEV friend adouble sqrt(const adouble &a);
118 CUDADEV friend adouble sin(const adouble &a);
119 CUDADEV friend adouble cos(const adouble &a);
120 CUDADEV friend adouble asin(const adouble &a);
121 CUDADEV friend adouble acos(const adouble &a);
122 CUDADEV friend adouble atan(const adouble &a);
123
124 CUDADEV friend adouble atan2(const adouble &a, const adouble &b);
125 CUDADEV friend adouble pow(const adouble &a, double v);
126 CUDADEV friend adouble pow(const adouble &a, const adouble &b);
127 CUDADEV friend adouble pow(double v, const adouble &a);
128 CUDADEV friend adouble log10(const adouble &a);
129
130 CUDADEV friend adouble sinh(const adouble &a);
131 CUDADEV friend adouble cosh(const adouble &a);
132 CUDADEV friend adouble tanh(const adouble &a);
133 CUDADEV friend adouble asinh(const adouble &a);
134 CUDADEV friend adouble acosh(const adouble &a);
135 CUDADEV friend adouble atanh(const adouble &a);
136 CUDADEV friend adouble fabs(const adouble &a);
137 CUDADEV friend adouble ceil(const adouble &a);
138 CUDADEV friend adouble floor(const adouble &a);
139 CUDADEV friend adouble fmax(const adouble &a, const adouble &b);
140 CUDADEV friend adouble fmax(double v, const adouble &a);
141 CUDADEV friend adouble fmax(const adouble &a, double v);
142 CUDADEV friend adouble fmin(const adouble &a, const adouble &b);
143 CUDADEV friend adouble fmin(double v, const adouble &a);
144 CUDADEV friend adouble fmin(const adouble &a, double v);
145 CUDADEV friend adouble ldexp(const adouble &a, const adouble &b);
146 CUDADEV friend adouble ldexp(const adouble &a, const double v);
147 CUDADEV friend adouble ldexp(const double v, const adouble &a);
148 CUDADEV friend double frexp(const adouble &a, int *v);
149 CUDADEV friend adouble erf(const adouble &a);
150 CUDADEV friend adouble erfc(const adouble &a);
151
152 /******************* nontemporary results ***************************/
153 // assignment
154 CUDADEV void operator=(const double v);
155 CUDADEV void operator=(const adouble &a);
156
157 // addition
158 CUDADEV void operator+=(const double v);
159 CUDADEV void operator+=(const adouble &a);
160
161 // subtraction
162 CUDADEV void operator-=(const double v);
163 CUDADEV void operator-=(const adouble &a);
164
165 // multiplication
166 CUDADEV void operator*=(const double v);
167 CUDADEV void operator*=(const adouble &a);
168
169 // division
170 CUDADEV void operator/=(const double v);
171 CUDADEV void operator/=(const adouble &a);
172
173 // not
174 CUDADEV int operator!() const;
175
176 // comparison
177 CUDADEV int operator!=(const adouble &) const;
178 CUDADEV int operator!=(const double) const;
179 CUDADEV friend int operator!=(const double, const adouble &);
180
181 CUDADEV int operator==(const adouble &) const;
182 CUDADEV int operator==(const double) const;
183 CUDADEV friend int operator==(const double, const adouble &);
184
185 CUDADEV int operator<=(const adouble &) const;
186 CUDADEV int operator<=(const double) const;
187 CUDADEV friend int operator<=(const double, const adouble &);
188
189 CUDADEV int operator>=(const adouble &) const;
190 CUDADEV int operator>=(const double) const;
191 CUDADEV friend int operator>=(const double, const adouble &);
192
193 CUDADEV int operator>(const adouble &) const;
194 CUDADEV int operator>(const double) const;
195 CUDADEV friend int operator>(const double, const adouble &);
196
197 CUDADEV int operator<(const adouble &) const;
198 CUDADEV int operator<(const double) const;
199 CUDADEV friend int operator<(const double, const adouble &);
200
201 /******************* getter / setter ********************************/
202 CUDAHOSTDEV double getValue() const;
203 CUDAHOSTDEV void setValue(const double v);
204 CUDAHOSTDEV ADVAL_TYPE getADValue() const;
205 CUDAHOSTDEV void setADValue(ADVAL_TYPE v);
206#if defined(NUMBER_DIRECTIONS)
207 CUDAHOSTDEV double getADValue(const unsigned int p) const;
208 CUDAHOSTDEV void setADValue(const unsigned int p, const double v);
209#endif
210
211 /******************* i/o operations *********************************/
212 CUDAHOST friend ostream &operator<<(ostream &, const adouble &);
213 CUDAHOST friend istream &operator>>(istream &, adouble &);
214
215private:
216 // internal variables
217 double val;
218 double ADVAL_DECL;
219};
220
221/******************************* ctors ************************************/
222CUDADEV adouble::adouble() {
223#if defined(DYNAMIC_DIRECTIONS)
224 adval = new double[ADOLC_numDir];
225#endif
226}
227
228CUDADEV adouble::adouble(const double v) : val(v) {
229#if defined(DYNAMIC_DIRECTIONS)
230 adval = new double[ADOLC_numDir];
231#endif
232 FOR_I_EQ_0_LT_NUMDIR
233 ADVAL_I = 0.0;
234}
235
236CUDADEV adouble::adouble(const double v, ADVAL_TYPE_ADV) : val(v) {
237#if defined(DYNAMIC_DIRECTIONS)
238 adval = new double[ADOLC_numDir];
239#endif
240 FOR_I_EQ_0_LT_NUMDIR
241 ADVAL_I = ADV_I;
242}
243
244CUDADEV adouble::adouble(const adouble &a) : val(a.val) {
245#if defined(DYNAMIC_DIRECTIONS)
246 adval = new double[ADOLC_numDir];
247#endif
248 FOR_I_EQ_0_LT_NUMDIR
249 ADVAL_I = a.ADVAL_I;
250}
251
252/******************************* dtors ************************************/
253#if defined(DYNAMIC_DIRECTIONS)
254CUDADEV adouble::~adouble() { delete[] adval; }
255#endif
256
257/************************* temporary results ******************************/
258// sign
259CUDADEV adouble adouble::operator-() const {
260 adouble tmp;
261 tmp.val = -val;
262 FOR_I_EQ_0_LT_NUMDIR
263 tmp.ADVAL_I = -ADVAL_I;
264 return tmp;
265}
266
267CUDADEV adouble adouble::operator+() const { return *this; }
268
269// addition
270CUDADEV adouble adouble::operator+(const double v) const {
271 return adouble(val + v, adval);
272}
273
274CUDADEV adouble adouble::operator+(const adouble &a) const {
275 adouble tmp;
276 tmp.val = val + a.val;
277 FOR_I_EQ_0_LT_NUMDIR
278 tmp.ADVAL_I = ADVAL_I + a.ADVAL_I;
279 return tmp;
280}
281
282CUDADEV adouble operator+(const double v, const adouble &a) {
283 return adouble(v + a.val, a.adval);
284}
285
286// subtraction
287CUDADEV adouble adouble::operator-(const double v) const {
288 return adouble(val - v, adval);
289}
290
291CUDADEV adouble adouble::operator-(const adouble &a) const {
292 adouble tmp;
293 tmp.val = val - a.val;
294 FOR_I_EQ_0_LT_NUMDIR
295 tmp.ADVAL_I = ADVAL_I - a.ADVAL_I;
296 return tmp;
297}
298
299CUDADEV adouble operator-(const double v, const adouble &a) {
300 adouble tmp;
301 tmp.val = v - a.val;
302 FOR_I_EQ_0_LT_NUMDIR
303 tmp.ADVAL_I = -a.ADVAL_I;
304 return tmp;
305}
306
307// multiplication
308CUDADEV adouble adouble::operator*(const double v) const {
309 adouble tmp;
310 tmp.val = val * v;
311 FOR_I_EQ_0_LT_NUMDIR
312 tmp.ADVAL_I = ADVAL_I * v;
313 return tmp;
314}
315
316CUDADEV adouble adouble::operator*(const adouble &a) const {
317 adouble tmp;
318 tmp.val = val * a.val;
319 FOR_I_EQ_0_LT_NUMDIR
320 tmp.ADVAL_I = ADVAL_I * a.val + val * a.ADVAL_I;
321 return tmp;
322}
323
324CUDADEV adouble operator*(const double v, const adouble &a) {
325 adouble tmp;
326 tmp.val = v * a.val;
327 FOR_I_EQ_0_LT_NUMDIR
328 tmp.ADVAL_I = v * a.ADVAL_I;
329 return tmp;
330}
331
332// division
333CUDADEV adouble adouble::operator/(const double v) const {
334 adouble tmp;
335 tmp.val = val / v;
336 FOR_I_EQ_0_LT_NUMDIR
337 tmp.ADVAL_I = ADVAL_I / v;
338 return tmp;
339}
340
341CUDADEV adouble adouble::operator/(const adouble &a) const {
342 adouble tmp;
343 tmp.val = val / a.val;
344 FOR_I_EQ_0_LT_NUMDIR
345 tmp.ADVAL_I = (ADVAL_I * a.val - val * a.ADVAL_I) / (a.val * a.val);
346 return tmp;
347}
348
349CUDADEV adouble operator/(const double v, const adouble &a) {
350 adouble tmp;
351 tmp.val = v / a.val;
352 FOR_I_EQ_0_LT_NUMDIR
353 tmp.ADVAL_I = (-v * a.ADVAL_I) / (a.val * a.val);
354 return tmp;
355}
356
357// inc/dec
358CUDADEV adouble adouble::operator++() {
359 ++val;
360 return *this;
361}
362
363CUDADEV adouble adouble::operator++(int) {
364 adouble tmp;
365 tmp.val = val++;
366 FOR_I_EQ_0_LT_NUMDIR
367 tmp.ADVAL_I = ADVAL_I;
368 return tmp;
369}
370
371CUDADEV adouble adouble::operator--() {
372 --val;
373 return *this;
374}
375
376CUDADEV adouble adouble::operator--(int) {
377 adouble tmp;
378 tmp.val = val--;
379 FOR_I_EQ_0_LT_NUMDIR
380 tmp.ADVAL_I = ADVAL_I;
381 return tmp;
382}
383
384// functions
385CUDADEV adouble tan(const adouble &a) {
386 adouble tmp;
387 double tmp2;
388 tmp.val = ADOLC_MATH_NSP::tan(a.val);
389 tmp2 = ADOLC_MATH_NSP::cos(a.val);
390 tmp2 *= tmp2;
391 FOR_I_EQ_0_LT_NUMDIR
392 tmp.ADVAL_I = a.ADVAL_I / tmp2;
393 return tmp;
394}
395
396CUDADEV adouble exp(const adouble &a) {
397 adouble tmp;
398 tmp.val = ADOLC_MATH_NSP::exp(a.val);
399 FOR_I_EQ_0_LT_NUMDIR
400 tmp.ADVAL_I = tmp.val * a.ADVAL_I;
401 return tmp;
402}
403
404CUDADEV adouble log(const adouble &a) {
405 adouble tmp;
406 tmp.val = ADOLC_MATH_NSP::log(a.val);
407 FOR_I_EQ_0_LT_NUMDIR
408 if (a.val > 0)
409 tmp.ADVAL_I = a.ADVAL_I / a.val;
410 else if (a.val == 0 && a.ADVAL_I != 0.0) {
411 int sign = (a.ADVAL_I < 0) ? -1 : 1;
412 tmp.ADVAL_I = sign * makeInf();
413 } else
414 tmp.ADVAL_I = makeNaN();
415 return tmp;
416}
417
418CUDADEV adouble sqrt(const adouble &a) {
419 adouble tmp;
420 tmp.val = ADOLC_MATH_NSP::sqrt(a.val);
421 FOR_I_EQ_0_LT_NUMDIR
422 if (a.val > 0)
423 tmp.ADVAL_I = a.ADVAL_I / (tmp.val * 2);
424 else if (a.val == 0.0 && a.ADVAL_I != 0.0) {
425 int sign = (a.ADVAL_I < 0) ? -1 : 1;
426 tmp.ADVAL_I = sign * makeInf();
427 } else
428 tmp.ADVAL_I = makeNaN();
429 return tmp;
430}
431
432CUDADEV adouble sin(const adouble &a) {
433 adouble tmp;
434 double tmp2;
435 tmp.val = ADOLC_MATH_NSP::sin(a.val);
436 tmp2 = ADOLC_MATH_NSP::cos(a.val);
437 FOR_I_EQ_0_LT_NUMDIR
438 tmp.ADVAL_I = tmp2 * a.ADVAL_I;
439 return tmp;
440}
441
442CUDADEV adouble cos(const adouble &a) {
443 adouble tmp;
444 double tmp2;
445 tmp.val = ADOLC_MATH_NSP::cos(a.val);
446 tmp2 = -ADOLC_MATH_NSP::sin(a.val);
447 FOR_I_EQ_0_LT_NUMDIR
448 tmp.ADVAL_I = tmp2 * a.ADVAL_I;
449 return tmp;
450}
451
452CUDADEV adouble asin(const adouble &a) {
453 adouble tmp;
454 tmp.val = ADOLC_MATH_NSP::asin(a.val);
455 double tmp2 = ADOLC_MATH_NSP::sqrt(1 - a.val * a.val);
456 FOR_I_EQ_0_LT_NUMDIR
457 tmp.ADVAL_I = a.ADVAL_I / tmp2;
458 return tmp;
459}
460
461CUDADEV adouble acos(const adouble &a) {
462 adouble tmp;
463 tmp.val = ADOLC_MATH_NSP::acos(a.val);
464 double tmp2 = -ADOLC_MATH_NSP::sqrt(1 - a.val * a.val);
465 FOR_I_EQ_0_LT_NUMDIR
466 tmp.ADVAL_I = a.ADVAL_I / tmp2;
467 return tmp;
468}
469
470CUDADEV adouble atan(const adouble &a) {
471 adouble tmp;
472 tmp.val = ADOLC_MATH_NSP::atan(a.val);
473 double tmp2 = 1 + a.val * a.val;
474 tmp2 = 1 / tmp2;
475 if (tmp2 != 0)
476 FOR_I_EQ_0_LT_NUMDIR
477 tmp.ADVAL_I = a.ADVAL_I * tmp2;
478 else FOR_I_EQ_0_LT_NUMDIR tmp.ADVAL_I = 0.0;
479 return tmp;
480}
481
482CUDADEV adouble atan2(const adouble &a, const adouble &b) {
483 adouble tmp;
484 tmp.val = ADOLC_MATH_NSP::atan2(a.val, b.val);
485 double tmp2 = a.val * a.val;
486 double tmp3 = b.val * b.val;
487 double tmp4 = tmp3 / (tmp2 + tmp3);
488 if (tmp4 != 0)
489 FOR_I_EQ_0_LT_NUMDIR
490 tmp.ADVAL_I = (a.ADVAL_I * b.val - a.val * b.ADVAL_I) / tmp3 * tmp4;
491 else FOR_I_EQ_0_LT_NUMDIR tmp.ADVAL_I = 0.0;
492 return tmp;
493}
494
495CUDADEV adouble pow(const adouble &a, double v) {
496 adouble tmp;
497 tmp.val = ADOLC_MATH_NSP::pow(a.val, v);
498 double tmp2 = v * ADOLC_MATH_NSP::pow(a.val, v - 1);
499 FOR_I_EQ_0_LT_NUMDIR
500 tmp.ADVAL_I = tmp2 * a.ADVAL_I;
501 return tmp;
502}
503
504CUDADEV adouble pow(const adouble &a, const adouble &b) {
505 adouble tmp;
506 tmp.val = ADOLC_MATH_NSP::pow(a.val, b.val);
507 double tmp2 = b.val * ADOLC_MATH_NSP::pow(a.val, b.val - 1);
508 double tmp3 = ADOLC_MATH_NSP::log(a.val) * tmp.val;
509 FOR_I_EQ_0_LT_NUMDIR
510 tmp.ADVAL_I = tmp2 * a.ADVAL_I + tmp3 * b.ADVAL_I;
511 return tmp;
512}
513
514CUDADEV adouble pow(double v, const adouble &a) {
515 adouble tmp;
516 tmp.val = ADOLC_MATH_NSP::pow(v, a.val);
517 double tmp2 = tmp.val * ADOLC_MATH_NSP::log(v);
518 FOR_I_EQ_0_LT_NUMDIR
519 tmp.ADVAL_I = tmp2 * a.ADVAL_I;
520 return tmp;
521}
522
523CUDADEV adouble log10(const adouble &a) {
524 adouble tmp;
525 tmp.val = ADOLC_MATH_NSP::log10(a.val);
526 double tmp2 = ADOLC_MATH_NSP::log((double)10) * a.val;
527 FOR_I_EQ_0_LT_NUMDIR
528 tmp.ADVAL_I = a.ADVAL_I / tmp2;
529 return tmp;
530}
531
532CUDADEV adouble sinh(const adouble &a) {
533 adouble tmp;
534 tmp.val = ADOLC_MATH_NSP::sinh(a.val);
535 double tmp2 = ADOLC_MATH_NSP::cosh(a.val);
536 FOR_I_EQ_0_LT_NUMDIR
537 tmp.ADVAL_I = a.ADVAL_I * tmp2;
538 return tmp;
539}
540
541CUDADEV adouble cosh(const adouble &a) {
542 adouble tmp;
543 tmp.val = ADOLC_MATH_NSP::cosh(a.val);
544 double tmp2 = ADOLC_MATH_NSP::sinh(a.val);
545 FOR_I_EQ_0_LT_NUMDIR
546 tmp.ADVAL_I = a.ADVAL_I * tmp2;
547 return tmp;
548}
549
550CUDADEV adouble tanh(const adouble &a) {
551 adouble tmp;
552 tmp.val = ADOLC_MATH_NSP::tanh(a.val);
553 double tmp2 = ADOLC_MATH_NSP::cosh(a.val);
554 tmp2 *= tmp2;
555 FOR_I_EQ_0_LT_NUMDIR
556 tmp.ADVAL_I = a.ADVAL_I / tmp2;
557 return tmp;
558}
559
560CUDADEV adouble asinh(const adouble &a) {
561 adouble tmp;
562 tmp.val = ADOLC_MATH_NSP_ERF::asinh(a.val);
563 double tmp2 = ADOLC_MATH_NSP::sqrt(a.val * a.val + 1);
564 FOR_I_EQ_0_LT_NUMDIR
565 tmp.ADVAL_I = a.ADVAL_I / tmp2;
566 return tmp;
567}
568
569CUDADEV adouble acosh(const adouble &a) {
570 adouble tmp;
571 tmp.val = ADOLC_MATH_NSP_ERF::acosh(a.val);
572 double tmp2 = ADOLC_MATH_NSP::sqrt(a.val * a.val - 1);
573 FOR_I_EQ_0_LT_NUMDIR
574 tmp.ADVAL_I = a.ADVAL_I / tmp2;
575 return tmp;
576}
577
578CUDADEV adouble atanh(const adouble &a) {
579 adouble tmp;
580 tmp.val = ADOLC_MATH_NSP_ERF::atanh(a.val);
581 double tmp2 = 1 - a.val * a.val;
582 FOR_I_EQ_0_LT_NUMDIR
583 tmp.ADVAL_I = a.ADVAL_I / tmp2;
584 return tmp;
585}
586
587CUDADEV adouble fabs(const adouble &a) {
588 adouble tmp;
589 tmp.val = ADOLC_MATH_NSP::fabs(a.val);
590 int as = 0;
591 if (a.val > 0)
592 as = 1;
593 if (a.val < 0)
594 as = -1;
595 if (as != 0)
596 FOR_I_EQ_0_LT_NUMDIR
597 tmp.ADVAL_I = a.ADVAL_I * as;
598 else FOR_I_EQ_0_LT_NUMDIR {
599 as = 0;
600 if (a.ADVAL_I > 0)
601 as = 1;
602 if (a.ADVAL_I < 0)
603 as = -1;
604 tmp.ADVAL_I = a.ADVAL_I * as;
605 }
606 return tmp;
607}
608
609CUDADEV adouble ceil(const adouble &a) {
610 adouble tmp;
611 tmp.val = ADOLC_MATH_NSP::ceil(a.val);
612 FOR_I_EQ_0_LT_NUMDIR
613 tmp.ADVAL_I = 0.0;
614 return tmp;
615}
616
617CUDADEV adouble floor(const adouble &a) {
618 adouble tmp;
619 tmp.val = ADOLC_MATH_NSP::floor(a.val);
620 FOR_I_EQ_0_LT_NUMDIR
621 tmp.ADVAL_I = 0.0;
622 return tmp;
623}
624
625CUDADEV adouble fmax(const adouble &a, const adouble &b) {
626 adouble tmp;
627 double tmp2 = a.val - b.val;
628 if (tmp2 < 0) {
629 tmp.val = b.val;
630 FOR_I_EQ_0_LT_NUMDIR
631 tmp.ADVAL_I = b.ADVAL_I;
632 } else {
633 tmp.val = a.val;
634 if (tmp2 > 0) {
635 FOR_I_EQ_0_LT_NUMDIR
636 tmp.ADVAL_I = a.ADVAL_I;
637 } else {
638 FOR_I_EQ_0_LT_NUMDIR {
639 if (a.ADVAL_I < b.ADVAL_I)
640 tmp.ADVAL_I = b.ADVAL_I;
641 else
642 tmp.ADVAL_I = a.ADVAL_I;
643 }
644 }
645 }
646 return tmp;
647}
648
649CUDADEV adouble fmax(double v, const adouble &a) {
650 adouble tmp;
651 double tmp2 = v - a.val;
652 if (tmp2 < 0) {
653 tmp.val = a.val;
654 FOR_I_EQ_0_LT_NUMDIR
655 tmp.ADVAL_I = a.ADVAL_I;
656 } else {
657 tmp.val = v;
658 if (tmp2 > 0) {
659 FOR_I_EQ_0_LT_NUMDIR
660 tmp.ADVAL_I = 0.0;
661 } else {
662 FOR_I_EQ_0_LT_NUMDIR {
663 if (a.ADVAL_I > 0)
664 tmp.ADVAL_I = a.ADVAL_I;
665 else
666 tmp.ADVAL_I = 0.0;
667 }
668 }
669 }
670 return tmp;
671}
672
673CUDADEV adouble fmax(const adouble &a, double v) {
674 adouble tmp;
675 double tmp2 = a.val - v;
676 if (tmp2 < 0) {
677 tmp.val = v;
678 FOR_I_EQ_0_LT_NUMDIR
679 tmp.ADVAL_I = 0.0;
680 } else {
681 tmp.val = a.val;
682 if (tmp2 > 0) {
683 FOR_I_EQ_0_LT_NUMDIR
684 tmp.ADVAL_I = a.ADVAL_I;
685 } else {
686 FOR_I_EQ_0_LT_NUMDIR {
687 if (a.ADVAL_I > 0)
688 tmp.ADVAL_I = a.ADVAL_I;
689 else
690 tmp.ADVAL_I = 0.0;
691 }
692 }
693 }
694 return tmp;
695}
696
697CUDADEV adouble fmin(const adouble &a, const adouble &b) {
698 adouble tmp;
699 double tmp2 = a.val - b.val;
700 if (tmp2 < 0) {
701 tmp.val = a.val;
702 FOR_I_EQ_0_LT_NUMDIR
703 tmp.ADVAL_I = a.ADVAL_I;
704 } else {
705 tmp.val = b.val;
706 if (tmp2 > 0) {
707 FOR_I_EQ_0_LT_NUMDIR
708 tmp.ADVAL_I = b.ADVAL_I;
709 } else {
710 FOR_I_EQ_0_LT_NUMDIR {
711 if (a.ADVAL_I < b.ADVAL_I)
712 tmp.ADVAL_I = a.ADVAL_I;
713 else
714 tmp.ADVAL_I = b.ADVAL_I;
715 }
716 }
717 }
718 return tmp;
719}
720
721CUDADEV adouble fmin(double v, const adouble &a) {
722 adouble tmp;
723 double tmp2 = v - a.val;
724 if (tmp2 < 0) {
725 tmp.val = v;
726 FOR_I_EQ_0_LT_NUMDIR
727 tmp.ADVAL_I = 0.0;
728 } else {
729 tmp.val = a.val;
730 if (tmp2 > 0) {
731 FOR_I_EQ_0_LT_NUMDIR
732 tmp.ADVAL_I = a.ADVAL_I;
733 } else {
734 FOR_I_EQ_0_LT_NUMDIR {
735 if (a.ADVAL_I < 0)
736 tmp.ADVAL_I = a.ADVAL_I;
737 else
738 tmp.ADVAL_I = 0.0;
739 }
740 }
741 }
742 return tmp;
743}
744
745CUDADEV adouble fmin(const adouble &a, double v) {
746 adouble tmp;
747 double tmp2 = a.val - v;
748 if (tmp2 < 0) {
749 tmp.val = a.val;
750 FOR_I_EQ_0_LT_NUMDIR
751 tmp.ADVAL_I = a.ADVAL_I;
752 } else {
753 tmp.val = v;
754 if (tmp2 > 0) {
755 FOR_I_EQ_0_LT_NUMDIR
756 tmp.ADVAL_I = 0.0;
757 } else {
758 FOR_I_EQ_0_LT_NUMDIR {
759 if (a.ADVAL_I < 0)
760 tmp.ADVAL_I = a.ADVAL_I;
761 else
762 tmp.ADVAL_I = 0.0;
763 }
764 }
765 }
766 return tmp;
767}
768
769CUDADEV adouble ldexp(const adouble &a, const adouble &b) {
770 return a * pow(2., b);
771}
772
773CUDADEV adouble ldexp(const adouble &a, const double v) {
774 return a * ADOLC_MATH_NSP::pow(2., v);
775}
776
777CUDADEV adouble ldexp(const double v, const adouble &a) {
778 return v * pow(2., a);
779}
780
781CUDADEV double frexp(const adouble &a, int *v) {
782 return ADOLC_MATH_NSP::frexp(a.val, v);
783}
784
785CUDADEV adouble erf(const adouble &a) {
786 adouble tmp;
787 tmp.val = ADOLC_MATH_NSP_ERF::erf(a.val);
788 double tmp2 = 2.0 / ADOLC_MATH_NSP_ERF::sqrt(ADOLC_MATH_NSP::acos(-1.0)) *
789 ADOLC_MATH_NSP_ERF::exp(-a.val * a.val);
790 FOR_I_EQ_0_LT_NUMDIR
791 tmp.ADVAL_I = tmp2 * a.ADVAL_I;
792 return tmp;
793}
794
795CUDADEV adouble erfc(const adouble &a) {
796 adouble tmp;
797 tmp.val = ADOLC_MATH_NSP_ERF::erfc(a.val);
798 double tmp2 = -2.0 / ADOLC_MATH_NSP_ERF::sqrt(ADOLC_MATH_NSP::acos(-1.0)) *
799 ADOLC_MATH_NSP_ERF::exp(-a.val * a.val);
800 FOR_I_EQ_0_LT_NUMDIR
801 tmp.ADVAL_I = tmp2 * a.ADVAL_I;
802 return tmp;
803}
804
805/******************* nontemporary results *********************************/
806CUDADEV void adouble::operator=(const double v) {
807 val = v;
808 FOR_I_EQ_0_LT_NUMDIR
809 ADVAL_I = 0.0;
810}
811
812CUDADEV void adouble::operator=(const adouble &a) {
813 val = a.val;
814 FOR_I_EQ_0_LT_NUMDIR
815 ADVAL_I = a.ADVAL_I;
816}
817
818CUDADEV void adouble::operator+=(const double v) { val += v; }
819
820CUDADEV void adouble::operator+=(const adouble &a) {
821 val = val + a.val;
822 FOR_I_EQ_0_LT_NUMDIR
823 ADVAL_I += a.ADVAL_I;
824}
825
826CUDADEV void adouble::operator-=(const double v) { val -= v; }
827
828CUDADEV void adouble::operator-=(const adouble &a) {
829 val = val - a.val;
830 FOR_I_EQ_0_LT_NUMDIR
831 ADVAL_I -= a.ADVAL_I;
832}
833
834CUDADEV void adouble::operator*=(const double v) {
835 val = val * v;
836 FOR_I_EQ_0_LT_NUMDIR
837 ADVAL_I *= v;
838}
839
840CUDADEV void adouble::operator*=(const adouble &a) {
841 FOR_I_EQ_0_LT_NUMDIR
842 ADVAL_I = ADVAL_I * a.val + val * a.ADVAL_I;
843 val *= a.val;
844}
845
846CUDADEV void adouble::operator/=(const double v) {
847 val /= v;
848 FOR_I_EQ_0_LT_NUMDIR
849 ADVAL_I /= v;
850}
851
852CUDADEV void adouble::operator/=(const adouble &a) {
853 FOR_I_EQ_0_LT_NUMDIR
854 ADVAL_I = (ADVAL_I * a.val - val * a.ADVAL_I) / (a.val * a.val);
855 val = val / a.val;
856}
857
858// not
859CUDADEV int adouble::operator!() const { return val == 0.0; }
860
861// comparison
862CUDADEV int adouble::operator!=(const adouble &a) const { return val != a.val; }
863
864CUDADEV int adouble::operator!=(const double v) const { return val != v; }
865
866CUDADEV int operator!=(const double v, const adouble &a) { return v != a.val; }
867
868CUDADEV int adouble::operator==(const adouble &a) const { return val == a.val; }
869
870CUDADEV int adouble::operator==(const double v) const { return val == v; }
871
872CUDADEV int operator==(const double v, const adouble &a) { return v == a.val; }
873
874CUDADEV int adouble::operator<=(const adouble &a) const { return val <= a.val; }
875
876CUDADEV int adouble::operator<=(const double v) const { return val <= v; }
877
878CUDADEV int operator<=(const double v, const adouble &a) { return v <= a.val; }
879
880CUDADEV int adouble::operator>=(const adouble &a) const { return val >= a.val; }
881
882CUDADEV int adouble::operator>=(const double v) const { return val >= v; }
883
884CUDADEV int operator>=(const double v, const adouble &a) { return v >= a.val; }
885
886CUDADEV int adouble::operator>(const adouble &a) const { return val > a.val; }
887
888CUDADEV int adouble::operator>(const double v) const { return val > v; }
889
890CUDADEV int operator>(const double v, const adouble &a) { return v > a.val; }
891
892CUDADEV int adouble::operator<(const adouble &a) const { return val < a.val; }
893
894CUDADEV int adouble::operator<(const double v) const { return val < v; }
895
896CUDADEV int operator<(const double v, const adouble &a) { return v < a.val; }
897
898/******************* getter / setter **************************************/
899CUDAHOSTDEV double adouble::getValue() const { return val; }
900
901CUDAHOSTDEV void adouble::setValue(const double v) { val = v; }
902
903CUDAHOSTDEV ADVAL_TYPE adouble::getADValue() const { return (ADVAL_TYPE)adval; }
904
909
910#if defined(NUMBER_DIRECTIONS)
911CUDAHOSTDEV double adouble::getADValue(const unsigned int p) const {
912 unsigned int locp = p;
913 if (locp >= ADOLC_numDir) {
914 locp = ADOLC_numDir - 1;
915 }
916 return adval[locp];
917}
918
919CUDAHOSTDEV void adouble::setADValue(const unsigned int p, const double v) {
920 unsigned int locp = p;
921 if (locp >= ADOLC_numDir) {
922 locp = ADOLC_numDir - 1;
923 }
924 adval[locp] = v;
925}
926#endif
927
928#if defined(NUMBER_DIRECTIONS)
929void setNumDir(const unsigned int p) {
930#if !defined(DYNAMIC_DIRECTIONS)
931 if (p > NUMBER_DIRECTIONS)
932 ADOLC_numDir = NUMBER_DIRECTIONS;
933 else
934 ADOLC_numDir = p;
935#else
936 ADOLC_numDir = p;
937#endif
938}
939#endif
940
941/******************* i/o operations ***************************************/
942CUDAHOST ostream &operator<<(ostream &out, const adouble &a) {
943 out << "Value: " << a.val;
944#if !defined(NUMBER_DIRECTIONS)
945 out << " ADValue: ";
946#else
947 out << " ADValues (" << ADOLC_numDir << "): ";
948#endif
949 FOR_I_EQ_0_LT_NUMDIR
950 out << a.ADVAL_I << " ";
951 out << "(a)";
952 return out;
953}
954
955CUDAHOST istream &operator>>(istream &in, adouble &a) {
956 char c;
957 do {
958 in >> c;
959 } while (c != ':' && !in.eof());
960 in >> a.val;
961#if !defined(NUMBER_DIRECTIONS)
962 do
963 in >> c;
964 while (c != ':' && !in.eof());
965#else
966 unsigned int num;
967 do
968 in >> c;
969 while (c != '(' && !in.eof());
970 in >> num;
971 if (num > NUMBER_DIRECTIONS) {
972 cout << "ADOL-C error: to many directions in input\n";
973 exit(-1);
974 }
975 do
976 in >> c;
977 while (c != ')' && !in.eof());
978#endif
979 FOR_I_EQ_0_LT_NUMDIR
980 in >> a.ADVAL_I;
981 do
982 in >> c;
983 while (c != ')' && !in.eof());
984 return in;
985}
986} // namespace adtlc
987/****************************************************************************/
988/* end traceless gpu implementation first order derivatives */
989/****************************************************************************/
990#endif
FOR_I_EQ_0_LT_NUMDIR
#define FOR_I_EQ_0_LT_NUMDIR
Definition adoublecuda.h:50
ADVAL_TYPE_ADV
#define ADVAL_TYPE_ADV
Definition adoublecuda.h:48
CUDADEV
#define CUDADEV
Definition adoublecuda.h:62
ADVAL_I
#define ADVAL_I
Definition adoublecuda.h:51
CUDAHOSTDEV
#define CUDAHOSTDEV
Definition adoublecuda.h:64
CUDAHOST
#define CUDAHOST
Definition adoublecuda.h:63
ADV_I
#define ADV_I
Definition adoublecuda.h:52
ADVAL_TYPE
#define ADVAL_TYPE
Definition adoublecuda.h:49
V_I
#define V_I
Definition adoublecuda.h:53
adouble::~adouble
~adouble()=default
Default destructor.
adtlc::adouble
CUDA-compatible tapeless active scalar for forward derivatives.
Definition adoublecuda.h:72
adtlc::adouble::acos
CUDADEV friend adouble acos(const adouble &a)
Definition adoublecuda.h:461
adtlc::adouble::operator<=
CUDADEV friend int operator<=(const double, const adouble &)
Definition adoublecuda.h:878
adtlc::adouble::operator/=
CUDADEV void operator/=(const double v)
Definition adoublecuda.h:846
adtlc::adouble::atan2
CUDADEV friend adouble atan2(const adouble &a, const adouble &b)
Definition adoublecuda.h:482
adtlc::adouble::getADValue
CUDAHOSTDEV ADVAL_TYPE getADValue() const
Definition adoublecuda.h:903
adtlc::adouble::asinh
CUDADEV friend adouble asinh(const adouble &a)
Definition adoublecuda.h:560
adtlc::adouble::operator/
CUDADEV friend adouble operator/(const double v, const adouble &a)
Definition adoublecuda.h:349
adtlc::adouble::operator-
CUDADEV adouble operator-() const
Definition adoublecuda.h:259
adtlc::adouble::operator>
CUDADEV friend int operator>(const double, const adouble &)
Definition adoublecuda.h:890
adtlc::adouble::ldexp
CUDADEV friend adouble ldexp(const adouble &a, const adouble &b)
Definition adoublecuda.h:769
adtlc::adouble::exp
CUDADEV friend adouble exp(const adouble &a)
Definition adoublecuda.h:396
adtlc::adouble::operator>>
CUDAHOST friend istream & operator>>(istream &, adouble &)
Definition adoublecuda.h:955
adtlc::adouble::fabs
CUDADEV friend adouble fabs(const adouble &a)
Definition adoublecuda.h:587
adtlc::adouble::log10
CUDADEV friend adouble log10(const adouble &a)
Definition adoublecuda.h:523
adtlc::adouble::floor
CUDADEV friend adouble floor(const adouble &a)
Definition adoublecuda.h:617
adtlc::adouble::setADValue
CUDAHOSTDEV void setADValue(ADVAL_TYPE v)
Definition adoublecuda.h:905
adtlc::adouble::sqrt
CUDADEV friend adouble sqrt(const adouble &a)
Definition adoublecuda.h:418
adtlc::adouble::atan
CUDADEV friend adouble atan(const adouble &a)
Definition adoublecuda.h:470
adtlc::adouble::pow
CUDADEV friend adouble pow(const adouble &a, double v)
Definition adoublecuda.h:495
adtlc::adouble::atanh
CUDADEV friend adouble atanh(const adouble &a)
Definition adoublecuda.h:578
adtlc::adouble::operator--
CUDADEV adouble operator--()
Definition adoublecuda.h:371
adtlc::adouble::sin
CUDADEV friend adouble sin(const adouble &a)
Definition adoublecuda.h:432
adtlc::adouble::cosh
CUDADEV friend adouble cosh(const adouble &a)
Definition adoublecuda.h:541
adtlc::adouble::getValue
CUDAHOSTDEV double getValue() const
Definition adoublecuda.h:899
adtlc::adouble::asin
CUDADEV friend adouble asin(const adouble &a)
Definition adoublecuda.h:452
adtlc::adouble::adouble
CUDADEV adouble()
Definition adoublecuda.h:222
adtlc::adouble::operator+=
CUDADEV void operator+=(const double v)
Definition adoublecuda.h:818
adtlc::adouble::acosh
CUDADEV friend adouble acosh(const adouble &a)
Definition adoublecuda.h:569
adtlc::adouble::setValue
CUDAHOSTDEV void setValue(const double v)
Definition adoublecuda.h:901
adtlc::adouble::ceil
CUDADEV friend adouble ceil(const adouble &a)
Definition adoublecuda.h:609
adtlc::adouble::operator=
CUDADEV void operator=(const double v)
Definition adoublecuda.h:806
adtlc::adouble::operator<
CUDADEV friend int operator<(const double, const adouble &)
Definition adoublecuda.h:896
adtlc::adouble::cos
CUDADEV friend adouble cos(const adouble &a)
Definition adoublecuda.h:442
adtlc::adouble::operator*=
CUDADEV void operator*=(const double v)
Definition adoublecuda.h:834
adtlc::adouble::operator==
CUDADEV friend int operator==(const double, const adouble &)
Definition adoublecuda.h:872
adtlc::adouble::fmax
CUDADEV friend adouble fmax(const adouble &a, const adouble &b)
Definition adoublecuda.h:625
adtlc::adouble::fmin
CUDADEV friend adouble fmin(const adouble &a, const adouble &b)
Definition adoublecuda.h:697
adtlc::adouble::log
CUDADEV friend adouble log(const adouble &a)
Definition adoublecuda.h:404
adtlc::adouble::operator!=
CUDADEV friend int operator!=(const double, const adouble &)
Definition adoublecuda.h:866
adtlc::adouble::operator+
CUDADEV adouble operator+() const
Definition adoublecuda.h:267
adtlc::adouble::operator<<
CUDAHOST friend ostream & operator<<(ostream &, const adouble &)
Definition adoublecuda.h:942
adtlc::adouble::operator*
CUDADEV friend adouble operator*(const double v, const adouble &a)
Definition adoublecuda.h:324
adtlc::adouble::tan
CUDADEV friend adouble tan(const adouble &a)
Definition adoublecuda.h:385
adtlc::adouble::frexp
CUDADEV friend double frexp(const adouble &a, int *v)
Definition adoublecuda.h:781
adtlc::adouble::erf
CUDADEV friend adouble erf(const adouble &a)
Definition adoublecuda.h:785
adtlc::adouble::erfc
CUDADEV friend adouble erfc(const adouble &a)
Definition adoublecuda.h:795
adtlc::adouble::tanh
CUDADEV friend adouble tanh(const adouble &a)
Definition adoublecuda.h:550
adtlc::adouble::sinh
CUDADEV friend adouble sinh(const adouble &a)
Definition adoublecuda.h:532
adtlc::adouble::operator-=
CUDADEV void operator-=(const double v)
Definition adoublecuda.h:826
adtlc::adouble::operator!
CUDADEV int operator!() const
Definition adoublecuda.h:859
adtlc::adouble::operator>=
CUDADEV friend int operator>=(const double, const adouble &)
Definition adoublecuda.h:884
adtlc::adouble::operator++
CUDADEV adouble operator++()
Definition adoublecuda.h:358
adtlc
Definition adoublecuda.h:29
adtlc::sinh
CUDADEV adouble sinh(const adouble &a)
Definition adoublecuda.h:532
adtlc::sin
CUDADEV adouble sin(const adouble &a)
Definition adoublecuda.h:432
adtlc::operator-
CUDADEV adouble operator-(const double v, const adouble &a)
Definition adoublecuda.h:299
adtlc::operator==
CUDADEV int operator==(const double v, const adouble &a)
Definition adoublecuda.h:872
adtlc::atan2
CUDADEV adouble atan2(const adouble &a, const adouble &b)
Definition adoublecuda.h:482
adtlc::acosh
CUDADEV adouble acosh(const adouble &a)
Definition adoublecuda.h:569
adtlc::asinh
CUDADEV adouble asinh(const adouble &a)
Definition adoublecuda.h:560
adtlc::operator<
CUDADEV int operator<(const double v, const adouble &a)
Definition adoublecuda.h:896
adtlc::operator>>
CUDAHOST istream & operator>>(istream &in, adouble &a)
Definition adoublecuda.h:955
adtlc::makeInf
__device__ double makeInf()
Definition adoublecuda.h:60
adtlc::frexp
CUDADEV double frexp(const adouble &a, int *v)
Definition adoublecuda.h:781
adtlc::operator>=
CUDADEV int operator>=(const double v, const adouble &a)
Definition adoublecuda.h:884
adtlc::log10
CUDADEV adouble log10(const adouble &a)
Definition adoublecuda.h:523
adtlc::operator+
CUDADEV adouble operator+(const double v, const adouble &a)
Definition adoublecuda.h:282
adtlc::cos
CUDADEV adouble cos(const adouble &a)
Definition adoublecuda.h:442
adtlc::atanh
CUDADEV adouble atanh(const adouble &a)
Definition adoublecuda.h:578
adtlc::asin
CUDADEV adouble asin(const adouble &a)
Definition adoublecuda.h:452
adtlc::operator/
CUDADEV adouble operator/(const double v, const adouble &a)
Definition adoublecuda.h:349
adtlc::erfc
CUDADEV adouble erfc(const adouble &a)
Definition adoublecuda.h:795
adtlc::fabs
CUDADEV adouble fabs(const adouble &a)
Definition adoublecuda.h:587
adtlc::atan
CUDADEV adouble atan(const adouble &a)
Definition adoublecuda.h:470
adtlc::log
CUDADEV adouble log(const adouble &a)
Definition adoublecuda.h:404
adtlc::sqrt
CUDADEV adouble sqrt(const adouble &a)
Definition adoublecuda.h:418
adtlc::operator!=
CUDADEV int operator!=(const double v, const adouble &a)
Definition adoublecuda.h:866
adtlc::acos
CUDADEV adouble acos(const adouble &a)
Definition adoublecuda.h:461
adtlc::operator<<
CUDAHOST ostream & operator<<(ostream &out, const adouble &a)
Definition adoublecuda.h:942
adtlc::fmin
CUDADEV adouble fmin(const adouble &a, const adouble &b)
Definition adoublecuda.h:697
adtlc::pow
CUDADEV adouble pow(const adouble &a, double v)
Definition adoublecuda.h:495
adtlc::erf
CUDADEV adouble erf(const adouble &a)
Definition adoublecuda.h:785
adtlc::cosh
CUDADEV adouble cosh(const adouble &a)
Definition adoublecuda.h:541
adtlc::tanh
CUDADEV adouble tanh(const adouble &a)
Definition adoublecuda.h:550
adtlc::floor
CUDADEV adouble floor(const adouble &a)
Definition adoublecuda.h:617
adtlc::operator<=
CUDADEV int operator<=(const double v, const adouble &a)
Definition adoublecuda.h:878
adtlc::tan
CUDADEV adouble tan(const adouble &a)
Definition adoublecuda.h:385
adtlc::operator>
CUDADEV int operator>(const double v, const adouble &a)
Definition adoublecuda.h:890
adtlc::fmax
CUDADEV adouble fmax(const adouble &a, const adouble &b)
Definition adoublecuda.h:625
adtlc::operator*
CUDADEV adouble operator*(const double v, const adouble &a)
Definition adoublecuda.h:324
adtlc::makeNaN
__device__ double makeNaN()
Definition adoublecuda.h:58
adtlc::exp
CUDADEV adouble exp(const adouble &a)
Definition adoublecuda.h:396
adtlc::ceil
CUDADEV adouble ceil(const adouble &a)
Definition adoublecuda.h:609
adtlc::ldexp
CUDADEV adouble ldexp(const adouble &a, const adouble &b)
Definition adoublecuda.h:769
Generated by doxygen 1.9.8