v0.14.0
Loading...
Searching...
No Matches
adol-c_atom.cpp
Go to the documentation of this file.
1/*----------------------------------------------------------------------------
2 ADOL-C -- Automatic Differentiation by Overloading in C++
3 File: powexam.cpp
4 Revision: $Id: powexam.cpp 295 2009-02-25 13:32:25Z awalther $
5 Contents: computation of n-th power, described in the manual
6
7 Copyright (c) Andrea Walther, Andreas Griewank, Andreas Kowarz,
8 Hristo Mitev, Sebastian Schlenkrich, Jean Utke, Olaf Vogel
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/****************************************************************************/
17/* INCLUDES */
18#include <adolc/adolc.h> /* use of ALL ADOL-C interfaces */
19
20#include <iostream>
21using namespace std;
22
23/****************************************************************************/
24/* ADOUBLE ROUTINE */
26 adouble z = 1;
27
28 if (n>0) /* Recursion and branches */
29 { int nh = n/2; /* that do not depend on */
30 z = power(x,nh); /* adoubles are fine !!!! */
31 z *= z;
32 if (2*nh != n)
33 z *= x;
34 return z;
35 } /* end if */
36 else {
37 if (n==0) /* The local adouble z dies */
38 return z; /* as it goes out of scope. */
39 else
40 return 1/power(x,-n);
41 } /* end else */
42} /* end power */
43
44/****************************************************************************/
45/* MAIN PROGRAM */
46int main() {
47 int i,tag = 1;
48 int n = 6;
49
50 //cout << "COMPUTATION OF N-TH POWER (ADOL-C Documented Example)\n\n";
51 //cout << "monomial degree=? \n"; /* input the desired degree */
52 //cin >> n;
53 /* allocations and initializations */
54 double** X;
55 double** Y;
56 X = myalloc2(1,n+4);
57 Y = myalloc2(1,n+4);
58 X[0][0] = 0.5; /* function value = 0. coefficient */
59 X[0][1] = 1.0; /* first derivative = 1. coefficient */
60 for(i=0; i<n+2; i++)
61 X[0][i+2] = 0; /* further coefficients */
62 double** Z; /* used for checking consistency */
63 Z = myalloc2(1,n+2); /* between forward and reverse */
64
65 adouble y,x; /* declare active variables */
66 /* beginning of active section */
67 trace_on(tag); /* tag = 1 and keep = 0 */
68 x <<= X[0][0]; /* only one independent var */
69 y = power(x,n); /* actual function call */
70 y >>= Y[0][0]; /* only one dependent adouble */
71 trace_off(); /* no global adouble has died */
72 /* end of active section */
73 double u[1]; /* weighting vector */
74 u[0]=1; /* for reverse call */
75 for(i=0; i<n+2; i++) /* note that keep = i+1 in call */
76 { forward(tag,1,1,i,i+1,X,Y); /* evaluate the i-the derivative */
77 if (i==0)
78 cout << Y[0][i] << " - " << y.value() << " = " << Y[0][i]-y.value()
79 << " (should be 0)\n";
80 else {
81 Z[0][i] = Z[0][i-1]/i; /* scale derivative to Taylorcoeff. */
82 cout << Y[0][i] << " - " << Z[0][i] << " = " << Y[0][i]-Z[0][i]
83 << " (should be 0)\n";
84 }
85 reverse(tag,1,1,i,u,Z); /* evaluate the (i+1)-st deriv. */
86 } /* end for */
87
88 return 0;
89} /* end main */
adouble power(adouble x, int n)
int main()
FTensor::Index< 'n', SPACE_DIM > n
FTensor::Index< 'i', SPACE_DIM > i