1 | /* $Id: poly.cpp 3311 2014-05-28 16:21:08Z bradbell $ */ |
---|
2 | /* -------------------------------------------------------------------------- |
---|
3 | CppAD: C++ Algorithmic Differentiation: Copyright (C) 2003-14 Bradley M. Bell |
---|
4 | |
---|
5 | CppAD is distributed under multiple licenses. This distribution is under |
---|
6 | the terms of the |
---|
7 | Eclipse Public License Version 1.0. |
---|
8 | |
---|
9 | A copy of this license is included in the COPYING file of this distribution. |
---|
10 | Please visit http://www.coin-or.org/CppAD/ for information on other licenses. |
---|
11 | -------------------------------------------------------------------------- */ |
---|
12 | /* |
---|
13 | $begin cppad_poly.cpp$$ |
---|
14 | $spell |
---|
15 | onetape |
---|
16 | coef |
---|
17 | ddp |
---|
18 | ADScalar |
---|
19 | dz |
---|
20 | ddz |
---|
21 | Taylor |
---|
22 | vector Vector |
---|
23 | typedef |
---|
24 | cppad |
---|
25 | Lu |
---|
26 | CppAD |
---|
27 | det |
---|
28 | hpp |
---|
29 | const |
---|
30 | CPPAD_TESTVECTOR |
---|
31 | bool |
---|
32 | var |
---|
33 | std |
---|
34 | cout |
---|
35 | endl |
---|
36 | $$ |
---|
37 | |
---|
38 | $section CppAD Speed: Second Derivative of a Polynomial$$ |
---|
39 | |
---|
40 | $index link_poly, cppad$$ |
---|
41 | $index cppad, link_poly$$ |
---|
42 | $index speed, cppad$$ |
---|
43 | $index cppad, speed$$ |
---|
44 | $index polynomial, speed cppad$$ |
---|
45 | |
---|
46 | $head Specifications$$ |
---|
47 | See $cref link_poly$$. |
---|
48 | |
---|
49 | $head Implementation$$ |
---|
50 | |
---|
51 | $codep */ |
---|
52 | # include <cppad/cppad.hpp> |
---|
53 | # include <cppad/speed/uniform_01.hpp> |
---|
54 | |
---|
55 | // Note that CppAD uses global_memory at the main program level |
---|
56 | extern bool |
---|
57 | global_onetape, global_atomic, global_optimize; |
---|
58 | |
---|
59 | bool link_poly( |
---|
60 | size_t size , |
---|
61 | size_t repeat , |
---|
62 | CppAD::vector<double> &a , // coefficients of polynomial |
---|
63 | CppAD::vector<double> &z , // polynomial argument value |
---|
64 | CppAD::vector<double> &ddp ) // second derivative w.r.t z |
---|
65 | { |
---|
66 | // speed test global option values |
---|
67 | if( global_atomic ) |
---|
68 | return false; |
---|
69 | |
---|
70 | // ----------------------------------------------------- |
---|
71 | // setup |
---|
72 | typedef CppAD::AD<double> ADScalar; |
---|
73 | typedef CppAD::vector<ADScalar> ADVector; |
---|
74 | |
---|
75 | size_t i; // temporary index |
---|
76 | size_t m = 1; // number of dependent variables |
---|
77 | size_t n = 1; // number of independent variables |
---|
78 | ADVector Z(n); // AD domain space vector |
---|
79 | ADVector P(m); // AD range space vector |
---|
80 | |
---|
81 | // choose the polynomial coefficients |
---|
82 | CppAD::uniform_01(size, a); |
---|
83 | |
---|
84 | // AD copy of the polynomial coefficients |
---|
85 | ADVector A(size); |
---|
86 | for(i = 0; i < size; i++) |
---|
87 | A[i] = a[i]; |
---|
88 | |
---|
89 | // forward mode first and second differentials |
---|
90 | CppAD::vector<double> p(1), dp(1), dz(1), ddz(1); |
---|
91 | dz[0] = 1.; |
---|
92 | ddz[0] = 0.; |
---|
93 | |
---|
94 | // AD function object |
---|
95 | CppAD::ADFun<double> f; |
---|
96 | |
---|
97 | // -------------------------------------------------------------------- |
---|
98 | if( ! global_onetape ) while(repeat--) |
---|
99 | { |
---|
100 | // choose an argument value |
---|
101 | CppAD::uniform_01(1, z); |
---|
102 | Z[0] = z[0]; |
---|
103 | |
---|
104 | // declare independent variables |
---|
105 | Independent(Z); |
---|
106 | |
---|
107 | // AD computation of the function value |
---|
108 | P[0] = CppAD::Poly(0, A, Z[0]); |
---|
109 | |
---|
110 | // create function object f : A -> detA |
---|
111 | f.Dependent(Z, P); |
---|
112 | |
---|
113 | if( global_optimize ) |
---|
114 | f.optimize(); |
---|
115 | |
---|
116 | // pre-allocate memory for three forward mode calculations |
---|
117 | f.capacity_order(3); |
---|
118 | |
---|
119 | // evaluate the polynomial |
---|
120 | p = f.Forward(0, z); |
---|
121 | |
---|
122 | // evaluate first order Taylor coefficient |
---|
123 | dp = f.Forward(1, dz); |
---|
124 | |
---|
125 | // second derivative is twice second order Taylor coef |
---|
126 | ddp = f.Forward(2, ddz); |
---|
127 | ddp[0] *= 2.; |
---|
128 | } |
---|
129 | else |
---|
130 | { |
---|
131 | // choose an argument value |
---|
132 | CppAD::uniform_01(1, z); |
---|
133 | Z[0] = z[0]; |
---|
134 | |
---|
135 | // declare independent variables |
---|
136 | Independent(Z); |
---|
137 | |
---|
138 | // AD computation of the function value |
---|
139 | P[0] = CppAD::Poly(0, A, Z[0]); |
---|
140 | |
---|
141 | // create function object f : A -> detA |
---|
142 | f.Dependent(Z, P); |
---|
143 | |
---|
144 | if( global_optimize ) |
---|
145 | f.optimize(); |
---|
146 | |
---|
147 | while(repeat--) |
---|
148 | { // sufficient memory is allocated by second repetition |
---|
149 | |
---|
150 | // get the next argument value |
---|
151 | CppAD::uniform_01(1, z); |
---|
152 | |
---|
153 | // evaluate the polynomial at the new argument value |
---|
154 | p = f.Forward(0, z); |
---|
155 | |
---|
156 | // evaluate first order Taylor coefficient |
---|
157 | dp = f.Forward(1, dz); |
---|
158 | |
---|
159 | // second derivative is twice second order Taylor coef |
---|
160 | ddp = f.Forward(2, ddz); |
---|
161 | ddp[0] *= 2.; |
---|
162 | } |
---|
163 | } |
---|
164 | return true; |
---|
165 | } |
---|
166 | /* $$ |
---|
167 | $end |
---|
168 | */ |
---|