NGSolve 5.3
autodiffdiff.hpp
1#ifndef FILE_AUTODIFFDIFF
2#define FILE_AUTODIFFDIFF
3
4/**************************************************************************/
5/* File: autodiffdiff.hpp */
6/* Author: Joachim Schoeberl */
7/* Date: 13. June. 05 */
8/**************************************************************************/
9
10namespace ngstd
11{
12 using ngcore::IfPos;
13
14// Automatic second differentiation datatype
15
16
22template <int D, typename SCAL = double>
24{
25 SCAL val;
26 SCAL dval[D?D:1];
27 SCAL ddval[D?D*D:1];
28public:
29
31
32
35
37 AutoDiffDiff (const AutoDiffDiff & ad2) throw()
38 {
39 val = ad2.val;
40 for (int i = 0; i < D; i++)
41 dval[i] = ad2.dval[i];
42 for (int i = 0; i < D*D; i++)
43 ddval[i] = ad2.ddval[i];
44 }
45
48 {
49 val = aval;
50 for (int i = 0; i < D; i++)
51 dval[i] = 0;
52 for (int i = 0; i < D*D; i++)
53 ddval[i] = 0;
54 }
55
58 {
59 val = ad2.Value();
60 for (int i = 0; i < D; i++)
61 dval[i] = ad2.DValue(i);
62 for (int i = 0; i < D*D; i++)
63 ddval[i] = 0;
64 }
65
68 {
69 val = aval;
70 for (int i = 0; i < D; i++)
71 dval[i] = 0;
72 for (int i = 0; i < D*D; i++)
73 ddval[i] = 0;
74 dval[diffindex] = 1;
75 }
76
78 {
79 val = aval;
80 LoadGradient (grad);
81 for (int i = 0; i < D*D; i++)
82 ddval[i] = 0;
83 }
84
85 INLINE AutoDiffDiff (SCAL aval, const SCAL * grad, const SCAL * hesse)
86 {
87 val = aval;
88 LoadGradient (grad);
89 LoadHessian (hesse);
90 }
91
94 {
95 val = aval;
96 for (int i = 0; i < D; i++)
97 dval[i] = 0;
98 for (int i = 0; i < D*D; i++)
99 ddval[i] = 0;
100 return *this;
101 }
102
103 INLINE void StoreGradient (SCAL * p) const
104 {
105 for (int i = 0; i < D; i++)
106 p[i] = dval[i];
107 }
108
109 INLINE void LoadGradient (const SCAL * p)
110 {
111 for (int i = 0; i < D; i++)
112 dval[i] = p[i];
113 }
114
115 INLINE void StoreHessian (SCAL * p) const
116 {
117 for (int i = 0; i < D*D; i++)
118 p[i] = ddval[i];
119 }
120
121 INLINE void LoadHessian (const SCAL * p)
122 {
123 for (int i = 0; i < D*D; i++)
124 ddval[i] = p[i];
125 }
126
128 SCAL Value() const throw() { return val; }
129
131 SCAL DValue (int i) const throw() { return dval[i]; }
132
133 AutoDiff<D,SCAL> DValueAD (int i) const
134 {
135 AutoDiff<D,SCAL> r(dval[i]);
136 for (int j = 0; j < D; j++)
137 r.DValue(j) = ddval[i*D+j];
138 return r;
139 }
140
142 SCAL DDValue (int i) const throw() { return ddval[i]; }
143
145 SCAL DDValue (int i, int j) const throw() { return ddval[i*D+j]; }
146
148 SCAL & Value() throw() { return val; }
149
151 SCAL & DValue (int i) throw() { return dval[i]; }
152
154 SCAL & DDValue (int i) throw() { return ddval[i]; }
155
157 SCAL & DDValue (int i, int j) throw() { return ddval[i*D+j]; }
158
159 explicit operator AutoDiff<D,SCAL> () const
160 { return AutoDiff<D,SCAL> (val, &dval[0]); }
161
164 {
165 val += y.val;
166 for (int i = 0; i < D; i++)
167 dval[i] += y.dval[i];
168 for (int i = 0; i < D*D; i++)
169 ddval[i] += y.ddval[i];
170 return *this;
171 }
172
175 {
176 val -= y.val;
177 for (int i = 0; i < D; i++)
178 dval[i] -= y.dval[i];
179 for (int i = 0; i < D*D; i++)
180 ddval[i] -= y.ddval[i];
181 return *this;
182 }
183
186 {
187 for (int i = 0; i < D*D; i++)
188 ddval[i] = val * y.ddval[i] + y.val * ddval[i];
189
190 for (int i = 0; i < D; i++)
191 for (int j = 0; j < D; j++)
192 ddval[i*D+j] += dval[i] * y.dval[j] + dval[j] * y.dval[i];
193
194 for (int i = 0; i < D; i++)
195 {
196 dval[i] *= y.val;
197 dval[i] += val * y.dval[i];
198 }
199 val *= y.val;
200 return *this;
201 }
202
205 {
206 for ( int i = 0; i < D*D; i++ )
207 ddval[i] *= y;
208 for (int i = 0; i < D; i++)
209 dval[i] *= y;
210 val *= y;
211 return *this;
212 }
213
216 {
217 SCAL iy = 1.0 / y;
218 for ( int i = 0; i < D*D; i++ )
219 ddval[i] *= iy;
220 for (int i = 0; i < D; i++)
221 dval[i] *= iy;
222 val *= iy;
223 return *this;
224 }
225
227 bool operator== (SCAL val2) throw()
228 {
229 return val == val2;
230 }
231
233 bool operator!= (SCAL val2) throw()
234 {
235 return val != val2;
236 }
237
239 bool operator< (SCAL val2) throw()
240 {
241 return val < val2;
242 }
243
245 bool operator> (SCAL val2) throw()
246 {
247 return val > val2;
248 }
249};
250
251
253
255template<int D, typename SCAL>
257{
258 ost << x.Value() << ", D = ";
259 for (int i = 0; i < D; i++)
260 ost << x.DValue(i) << " ";
261 ost << ", DD = ";
262 for (int i = 0; i < D*D; i++)
263 ost << x.DDValue(i) << " ";
264 return ost;
265}
266
268template<int D, typename SCAL>
269inline AutoDiffDiff<D, SCAL> operator+ (const AutoDiffDiff<D, SCAL> & x, const AutoDiffDiff<D, SCAL> & y) throw()
270{
271 AutoDiffDiff<D, SCAL> res;
272 res.Value () = x.Value()+y.Value();
273 for (int i = 0; i < D; i++)
274 res.DValue(i) = x.DValue(i) + y.DValue(i);
275 for (int i = 0; i < D*D; i++)
276 res.DDValue(i) = x.DDValue(i) + y.DDValue(i);
277 return res;
278}
279
280
282template<int D, typename SCAL>
283inline AutoDiffDiff<D, SCAL> operator- (const AutoDiffDiff<D, SCAL> & x, const AutoDiffDiff<D, SCAL> & y) throw()
284{
285 AutoDiffDiff<D, SCAL> res;
286 res.Value() = x.Value()-y.Value();
287 for (int i = 0; i < D; i++)
288 res.DValue(i) = x.DValue(i) - y.DValue(i);
289 for (int i = 0; i < D*D; i++)
290 res.DDValue(i) = x.DDValue(i) - y.DDValue(i);
291 return res;
292}
293
294
296template<int D, typename SCAL, typename SCAL2,
297 typename std::enable_if<std::is_convertible<SCAL2,SCAL>::value, int>::type = 0>
298inline AutoDiffDiff<D, SCAL> operator+ (SCAL2 x, const AutoDiffDiff<D, SCAL> & y) throw()
299{
300 AutoDiffDiff<D, SCAL> res;
301 res.Value() = x+y.Value();
302 for (int i = 0; i < D; i++)
303 res.DValue(i) = y.DValue(i);
304 for (int i = 0; i < D*D; i++)
305 res.DDValue(i) = y.DDValue(i);
306 return res;
307}
308
310template<int D, typename SCAL, typename SCAL2,
311 typename std::enable_if<std::is_convertible<SCAL2,SCAL>::value, int>::type = 0>
312inline AutoDiffDiff<D, SCAL> operator+ (const AutoDiffDiff<D, SCAL> & y, SCAL2 x) throw()
313{
314 AutoDiffDiff<D, SCAL> res;
315 res.Value() = x+y.Value();
316 for (int i = 0; i < D; i++)
317 res.DValue(i) = y.DValue(i);
318 for (int i = 0; i < D*D; i++)
319 res.DDValue(i) = y.DDValue(i);
320 return res;
321}
322
323
325template<int D, typename SCAL>
326inline AutoDiffDiff<D, SCAL> operator- (const AutoDiffDiff<D, SCAL> & x) throw()
327{
328 AutoDiffDiff<D, SCAL> res;
329 res.Value() = -x.Value();
330 for (int i = 0; i < D; i++)
331 res.DValue(i) = -x.DValue(i);
332 for (int i = 0; i < D*D; i++)
333 res.DDValue(i) = -x.DDValue(i);
334 return res;
335}
336
338template<int D, typename SCAL, typename SCAL2,
339 typename std::enable_if<std::is_convertible<SCAL2,SCAL>::value, int>::type = 0>
340inline AutoDiffDiff<D, SCAL> operator- (const AutoDiffDiff<D, SCAL> & x, SCAL2 y) throw()
341{
342 AutoDiffDiff<D, SCAL> res;
343 res.Value() = x.Value()-y;
344 for (int i = 0; i < D; i++)
345 res.DValue(i) = x.DValue(i);
346 for (int i = 0; i < D*D; i++)
347 res.DDValue(i) = x.DDValue(i);
348 return res;
349}
350
352template<int D, typename SCAL, typename SCAL2,
353 typename std::enable_if<std::is_convertible<SCAL2,SCAL>::value, int>::type = 0>
354inline AutoDiffDiff<D, SCAL> operator- (SCAL2 x, const AutoDiffDiff<D, SCAL> & y) throw()
355{
356 AutoDiffDiff<D, SCAL> res;
357 res.Value() = x-y.Value();
358 for (int i = 0; i < D; i++)
359 res.DValue(i) = -y.DValue(i);
360 for (int i = 0; i < D*D; i++)
361 res.DDValue(i) = -y.DDValue(i);
362 return res;
363}
364
365
367template<int D, typename SCAL, typename SCAL2,
368 typename std::enable_if<std::is_convertible<SCAL2,SCAL>::value, int>::type = 0>
369inline AutoDiffDiff<D, SCAL> operator* (SCAL2 x, const AutoDiffDiff<D, SCAL> & y) throw()
370{
371 AutoDiffDiff<D, SCAL> res;
372 res.Value() = x*y.Value();
373 for (int i = 0; i < D; i++)
374 res.DValue(i) = x*y.DValue(i);
375 for (int i = 0; i < D*D; i++)
376 res.DDValue(i) = x*y.DDValue(i);
377 return res;
378}
379
381template<int D, typename SCAL, typename SCAL2,
382 typename std::enable_if<std::is_convertible<SCAL2,SCAL>::value, int>::type = 0>
383inline AutoDiffDiff<D, SCAL> operator* (const AutoDiffDiff<D, SCAL> & y, SCAL2 x) throw()
384{
385 AutoDiffDiff<D, SCAL> res;
386 res.Value() = x*y.Value();
387 for (int i = 0; i < D; i++)
388 res.DValue(i) = x*y.DValue(i);
389 for (int i = 0; i < D*D; i++)
390 res.DDValue(i) = x*y.DDValue(i);
391 return res;
392}
393
395template<int D, typename SCAL>
396inline AutoDiffDiff<D, SCAL> operator* (const AutoDiffDiff<D, SCAL> & x, const AutoDiffDiff<D, SCAL> & y) throw()
397{
398 AutoDiffDiff<D, SCAL> res;
399 SCAL hx = x.Value();
400 SCAL hy = y.Value();
401
402 res.Value() = hx*hy;
403 for (int i = 0; i < D; i++)
404 res.DValue(i) = hx*y.DValue(i) + hy*x.DValue(i);
405
406 for (int i = 0; i < D; i++)
407 for (int j = 0; j < D; j++)
408 res.DDValue(i,j) = hx * y.DDValue(i,j) + hy * x.DDValue(i,j)
409 + x.DValue(i) * y.DValue(j) + x.DValue(j) * y.DValue(i);
410
411 return res;
412}
413
414
415
416template<int D, typename SCAL>
417inline AutoDiffDiff<D, SCAL> Inv (const AutoDiffDiff<D, SCAL> & x)
418{
419 AutoDiffDiff<D, SCAL> res(1.0 / x.Value());
420 for (int i = 0; i < D; i++)
421 res.DValue(i) = -x.DValue(i) / (x.Value() * x.Value());
422
423 SCAL fac1 = 2/(x.Value()*x.Value()*x.Value());
424 SCAL fac2 = 1/sqr(x.Value());
425 for (int i = 0; i < D; i++)
426 for (int j = 0; j < D; j++)
427 res.DDValue(i,j) = fac1*x.DValue(i)*x.DValue(j) - fac2*x.DDValue(i,j);
428 return res;
429}
430
431
432template<int D, typename SCAL>
433inline AutoDiffDiff<D, SCAL> operator/ (const AutoDiffDiff<D, SCAL> & x, const AutoDiffDiff<D, SCAL> & y)
434{
435 return x * Inv (y);
436}
437
438template<int D, typename SCAL, typename SCAL2,
439 typename std::enable_if<std::is_convertible<SCAL2,SCAL>::value, int>::type = 0>
440inline AutoDiffDiff<D, SCAL> operator/ (const AutoDiffDiff<D, SCAL> & x, SCAL2 y)
441{
442 return (1/y) * x;
443}
444
445template<int D, typename SCAL, typename SCAL2,
446 typename std::enable_if<std::is_convertible<SCAL2,SCAL>::value, int>::type = 0>
447inline AutoDiffDiff<D, SCAL> operator/ (SCAL2 x, const AutoDiffDiff<D, SCAL> & y)
448{
449 return x * Inv(y);
450}
451
452
453template<int D, typename SCAL>
454inline AutoDiffDiff<D, SCAL> sqrt (const AutoDiffDiff<D, SCAL> & x)
455{
456 AutoDiffDiff<D, SCAL> res;
457 res.Value() = sqrt(x.Value());
458 for (int j = 0; j < D; j++)
459 res.DValue(j) = IfZero(x.DValue(j),SCAL{0.},0.5 / res.Value() * x.DValue(j));
460
461
462 for (int i = 0; i < D; i++)
463 for (int j = 0; j < D; j++)
464 res.DDValue(i,j) = IfZero(x.DDValue(i,j)+x.DValue(i) * x.DValue(j),SCAL{0.},0.5/res.Value() * x.DDValue(i,j) - 0.25 / (x.Value()*res.Value()) * x.DValue(i) * x.DValue(j));
465
466 return res;
467}
468
469// df(u)/dx = exp(x) * du/dx
470// d^2 f(u) / dx^2 = exp(x) * (du/dx)^2 + exp(x) * d^2u /dx^2
471template <int D, typename SCAL>
472INLINE AutoDiffDiff<D, SCAL> exp (AutoDiffDiff<D, SCAL> x)
473{
474 AutoDiffDiff<D, SCAL> res;
475 res.Value() = exp(x.Value());
476 for (int k = 0; k < D; k++)
477 res.DValue(k) = x.DValue(k) * res.Value();
478 for (int k = 0; k < D; k++)
479 for (int l = 0; l < D; l++)
480 res.DDValue(k,l) = (x.DValue(k) * x.DValue(l)+x.DDValue(k,l)) * res.Value();
481 return res;
482}
483
484using std::pow;
485template <int D, typename SCAL>
486INLINE AutoDiffDiff<D,SCAL> pow (AutoDiffDiff<D,SCAL> x, AutoDiffDiff<D,SCAL> y )
487{
488 return exp(log(x)*y);
489}
490
491template <int D, typename SCAL>
492INLINE AutoDiffDiff<D, SCAL> log (AutoDiffDiff<D, SCAL> x)
493{
494 AutoDiffDiff<D, SCAL> res;
495 res.Value() = log(x.Value());
496 SCAL xinv = 1.0/x.Value();
497 for (int k = 0; k < D; k++)
498 res.DValue(k) = x.DValue(k) * xinv;
499 for (int k = 0; k < D; k++)
500 for (int l = 0; l < D; l++)
501 res.DDValue(k,l) = -xinv*xinv*x.DValue(k) * x.DValue(l) + xinv * x.DDValue(k,l);
502 return res;
503}
504
505
506
507template <int D, typename SCAL>
508INLINE AutoDiffDiff<D, SCAL> sin (AutoDiffDiff<D, SCAL> x)
509{
510 AutoDiffDiff<D, SCAL> res;
511 SCAL s = sin(x.Value());
512 SCAL c = cos(x.Value());
513
514 res.Value() = s;
515 for (int k = 0; k < D; k++)
516 res.DValue(k) = x.DValue(k) * c;
517 for (int k = 0; k < D; k++)
518 for (int l = 0; l < D; l++)
519 res.DDValue(k,l) = -s * x.DValue(k) * x.DValue(l) + c * x.DDValue(k,l);
520 return res;
521}
522
523
524template <int D, typename SCAL>
525INLINE AutoDiffDiff<D, SCAL> cos (AutoDiffDiff<D, SCAL> x)
526{
527 AutoDiffDiff<D, SCAL> res;
528 SCAL s = sin(x.Value());
529 SCAL c = cos(x.Value());
530
531 res.Value() = c;
532 for (int k = 0; k < D; k++)
533 res.DValue(k) = -s * x.DValue(k);
534 for (int k = 0; k < D; k++)
535 for (int l = 0; l < D; l++)
536 res.DDValue(k,l) = -c * x.DValue(k) * x.DValue(l) - s * x.DDValue(k,l);
537 return res;
538}
539
540template <int D, typename SCAL>
541INLINE AutoDiffDiff<D, SCAL> tan (AutoDiffDiff<D, SCAL> x)
542{ return sin(x) / cos(x); }
543
544
545template <int D, typename SCAL>
546INLINE AutoDiffDiff<D, SCAL> atan (AutoDiffDiff<D, SCAL> x)
547{
548 AutoDiffDiff<D, SCAL> res;
549 SCAL a = atan(x.Value());
550 res.Value() = a;
551 for (int k = 0; k < D; k++)
552 res.DValue(k) = x.DValue(k)/(1+x.Value()*x.Value()) ;
553 for (int k = 0; k < D; k++)
554 for (int l = 0; l < D; l++)
555 res.DDValue(k,l) = -2*x.Value()/((1+x.Value()*x.Value())*(1+x.Value()*x.Value())) * x.DValue(k) * x.DValue(l) + x.DDValue(k,l)/(1+x.Value()*x.Value());
556 return res;
557}
558
559template <int D, typename SCAL>
560INLINE AutoDiffDiff<D, SCAL> atan2 (AutoDiffDiff<D, SCAL> x,AutoDiffDiff<D, SCAL> y)
561{
562 AutoDiffDiff<D, SCAL> res;
563 SCAL a = atan2(x.Value(), y.Value());
564 res.Value() = a;
565 for (int k = 0; k < D; k++)
566 res.DValue(k) = (x.Value()*y.DValue(k)-y.Value()*x.DValue(k))/(y.Value()*y.Value()+x.Value()*x.Value());
567
568 for (int k = 0; k < D; k++)
569 for (int l = 0; l < D; l++)
570 res.DDValue(k,l) = (x.DValue(k)*y.DValue(l)+x.Value()*y.DDValue(l,k) - y.DValue(k)*x.DValue(l) - y.Value()*x.DDValue(l,k))/(y.Value()*y.Value()+x.Value()*x.Value()) - 2 * (x.Value()*y.DValue(k)-y.Value()*x.DValue(k)) * (x.Value()*x.DValue(k) + y.Value()*y.DValue(k))/( (y.Value()*y.Value()+x.Value()*x.Value()) * (y.Value()*y.Value()+x.Value()*x.Value()) );
571 return res;
572}
573
574
575
576using std::acos;
577template <int D, typename SCAL>
578INLINE AutoDiffDiff<D,SCAL> acos (AutoDiffDiff<D,SCAL> x)
579{
580 AutoDiffDiff<D,SCAL> res;
581 SCAL a = acos(x.Value());
582 res.Value() = a;
583 auto omaa = 1-x.Value()*x.Value();
584 auto s = sqrt(omaa);
585 SCAL da = -1 / s;
586 SCAL dda = -x.Value() / (s*omaa);
587 for (int k = 0; k < D; k++)
588 res.DValue(k) = x.DValue(k)*da;
589 for (int k = 0; k < D; k++)
590 for (int l = 0; l < D; l++)
591 res.DDValue(k,l) = dda * x.DValue(k) * x.DValue(l) + da * x.DDValue(k,l);
592
593 return res;
594}
595
596
597using std::acos;
598template <int D, typename SCAL>
599INLINE AutoDiffDiff<D,SCAL> asin (AutoDiffDiff<D,SCAL> x)
600{
601 AutoDiffDiff<D,SCAL> res;
602 SCAL a = asin(x.Value());
603 res.Value() = a;
604 auto omaa = 1-x.Value()*x.Value();
605 auto s = sqrt(omaa);
606 SCAL da = 1 / s;
607 SCAL dda = x.Value() / (s*omaa);
608 for (int k = 0; k < D; k++)
609 res.DValue(k) = x.DValue(k)*da;
610 for (int k = 0; k < D; k++)
611 for (int l = 0; l < D; l++)
612 res.DDValue(k,l) = dda * x.DValue(k) * x.DValue(l) + da * x.DDValue(k,l);
613
614 return res;
615}
616
617
618template <int D, typename SCAL>
619INLINE AutoDiffDiff<D, SCAL> sinh (AutoDiffDiff<D, SCAL> x)
620{
621 AutoDiffDiff<D, SCAL> res;
622 SCAL sh = sinh(x.Value());
623 SCAL ch = cosh(x.Value());
624
625 res.Value() = sh;
626 for (int k = 0; k < D; k++)
627 res.DValue(k) = x.DValue(k) * ch;
628 for (int k = 0; k < D; k++)
629 for (int l = 0; l < D; l++)
630 res.DDValue(k,l) = sh * x.DValue(k) * x.DValue(l) + ch * x.DDValue(k,l);
631 return res;
632}
633
634
635template <int D, typename SCAL>
636INLINE AutoDiffDiff<D, SCAL> cosh (AutoDiffDiff<D, SCAL> x)
637{
638 AutoDiffDiff<D, SCAL> res;
639 SCAL sh = sinh(x.Value());
640 SCAL ch = cosh(x.Value());
641
642 res.Value() = ch;
643 for (int k = 0; k < D; k++)
644 res.DValue(k) = sh * x.DValue(k);
645 for (int k = 0; k < D; k++)
646 for (int l = 0; l < D; l++)
647 res.DDValue(k,l) = ch * x.DValue(k) * x.DValue(l) + sh * x.DDValue(k,l);
648 return res;
649}
650
651template <int D, typename SCAL>
652INLINE AutoDiffDiff<D, SCAL> erf (AutoDiffDiff<D, SCAL> x)
653{
654 AutoDiffDiff<D, SCAL> res;
655 SCAL derf = 2. / sqrt(M_PI) * exp(- x.Value() * x.Value());
656
657 res.Value() = erf(x.Value());
658 for (int k = 0; k < D; k++)
659 res.DValue(k) = - derf * x.DValue(k);
660 for (int k = 0; k < D; k++)
661 for (int l = 0; l < D; l++)
662 res.DDValue(k,l) = derf * (x.DDValue(k, l) - 2 * x.Value() * x.DValue(k) * x.DValue(l));
663 return res;
664}
665
666using std::floor;
667template<int D, typename SCAL>
668INLINE AutoDiffDiff<D,SCAL> floor (const AutoDiffDiff<D,SCAL> & x)
669{
670 return floor(x.Value());
671}
672
673using std::ceil;
674template<int D, typename SCAL>
675INLINE AutoDiffDiff<D,SCAL> ceil (const AutoDiffDiff<D,SCAL> & x)
676{
677 return ceil(x.Value());
678}
679
680
681template <int D, typename SCAL, typename TB, typename TC>
682auto IfPos (AutoDiffDiff<D,SCAL> a, TB b, TC c) -> decltype(IfPos (a.Value(), b, c))
683{
684 return IfPos (a.Value(), b, c);
685}
686
687template <int D, typename SCAL>
688INLINE AutoDiffDiff<D,SCAL> IfPos (SCAL /* SIMD<double> */ a, AutoDiffDiff<D,SCAL> b, AutoDiffDiff<D,SCAL> c)
689{
690 AutoDiffDiff<D,SCAL> res;
691 res.Value() = IfPos (a, b.Value(), c.Value());
692 for (int j = 0; j < D; j++)
693 {
694 res.DValue(j) = IfPos (a, b.DValue(j), c.DValue(j));
695 res.DDValue(j) = IfPos (a, b.DDValue(j), c.DDValue(j));
696 }
697 return res;
698}
699
700template <int D, typename SCAL, typename TC>
701INLINE AutoDiffDiff<D,SCAL> IfPos (SCAL /* SIMD<double> */ a, AutoDiffDiff<D,SCAL> b, TC c)
702{
703 return IfPos (a, b, AutoDiffDiff<D,SCAL> (c));
704}
705
706
707
708
710
711}
712
713
714namespace ngbla
715{
716 template <typename T> struct is_scalar_type;
717 template <int D, typename T>
718 struct is_scalar_type<ngstd::AutoDiffDiff<D,T>> { static constexpr bool value = true; };
719
720
721 // not meaningful for AutoDiff<D,Complex>, since this is
722 // not (complex) differentiable anyway
723 template<int D, typename SCAL>
724 inline auto L2Norm2 (const ngstd::AutoDiffDiff<D,SCAL> & x)
725 {
726 return x*x;
727 }
728
729}
730
731
732
733#endif
Datatype for automatic differentiation.
Definition autodiffdiff.hpp:24
AutoDiffDiff< D, SCAL > & operator/=(const SCAL &y)
divide by scalar
Definition autodiffdiff.hpp:215
AutoDiffDiff< D, SCAL > & operator-=(const AutoDiffDiff< D, SCAL > &y)
subtract autodiffdiff object
Definition autodiffdiff.hpp:174
SCAL DValue(int i) const
returns partial derivative
Definition autodiffdiff.hpp:131
SCAL DDValue(int i) const
returns partial derivative
Definition autodiffdiff.hpp:142
bool operator!=(SCAL val2)
different values
Definition autodiffdiff.hpp:233
AutoDiffDiff(const AutoDiffDiff &ad2)
copy constructor
Definition autodiffdiff.hpp:37
AutoDiffDiff(const AutoDiff< D, SCAL > &ad2)
initial object with value and derivative
Definition autodiffdiff.hpp:57
bool operator<(SCAL val2)
less
Definition autodiffdiff.hpp:239
SCAL & DValue(int i)
accesses partial derivative
Definition autodiffdiff.hpp:151
AutoDiffDiff(SCAL aval)
initial object with constant value
Definition autodiffdiff.hpp:47
SCAL Value() const
returns value
Definition autodiffdiff.hpp:128
SCAL & DDValue(int i)
accesses partial derivative
Definition autodiffdiff.hpp:154
AutoDiffDiff()
elements are undefined
Definition autodiffdiff.hpp:34
AutoDiffDiff< D, SCAL > & operator*=(const AutoDiffDiff< D, SCAL > &y)
multiply with autodiffdiff object
Definition autodiffdiff.hpp:185
SCAL & Value()
access value
Definition autodiffdiff.hpp:148
SCAL DDValue(int i, int j) const
returns partial derivative
Definition autodiffdiff.hpp:145
AutoDiffDiff< D, SCAL > & operator+=(const AutoDiffDiff< D, SCAL > &y)
add autodiffdiff object
Definition autodiffdiff.hpp:163
bool operator>(SCAL val2)
greater
Definition autodiffdiff.hpp:245
AutoDiffDiff & operator=(SCAL aval)
assign constant value
Definition autodiffdiff.hpp:93
SCAL & DDValue(int i, int j)
accesses partial derivative
Definition autodiffdiff.hpp:157
bool operator==(SCAL val2)
same value
Definition autodiffdiff.hpp:227
AutoDiffDiff(SCAL aval, int diffindex)
init object with (val, e_diffindex)
Definition autodiffdiff.hpp:67
Datatype for automatic differentiation.
Definition autodiff.hpp:26
INLINE SCAL DValue(int i) const
returns partial derivative
Definition autodiff.hpp:88
INLINE SCAL Value() const
returns value
Definition autodiff.hpp:85
namespace for standard data types and algorithms.
Definition ngstd.hpp:42
ostream & operator<<(ostream &ost, const AutoDiffDiff< D, SCAL > &x)
Prints AudoDiffDiff.
Definition autodiffdiff.hpp:256
Definition autodiffdiff.hpp:716