Changeset 1659 for trunk/test_more
 Timestamp:
 Mar 10, 2010 10:48:31 PM (10 years ago)
 Location:
 trunk/test_more
 Files:

 2 edited
Legend:
 Unmodified
 Added
 Removed

trunk/test_more/reverse.cpp
r1370 r1659 1 1 /* $Id$ */ 2 2 /*  3 CppAD: C++ Algorithmic Differentiation: Copyright (C) 2003 07Bradley M. Bell3 CppAD: C++ Algorithmic Differentiation: Copyright (C) 200310 Bradley M. Bell 4 4 5 5 CppAD is distributed under multiple licenses. This distribution is under … … 16 16 17 17 # include <cppad/cppad.hpp> 18 namespace { //  18 19 19 20 bool Reverse(void) … … 119 120 return ok; 120 121 } 122 123 // define the template function reverse_any_cases<Vector> in empty namespace 124 template <typename Vector> 125 bool reverse_any_cases(void) 126 { bool ok = true; 127 using CppAD::AD; 128 using CppAD::NearEqual; 129 130 // domain space vector 131 size_t n = 3; 132 CPPAD_TEST_VECTOR< AD<double> > X(n); 133 X[0] = 0.; 134 X[1] = 1.; 135 X[2] = 2.; 136 137 // declare independent variables and start recording 138 CppAD::Independent(X); 139 140 // range space vector 141 size_t m = 1; 142 CPPAD_TEST_VECTOR< AD<double> > Y(m); 143 Y[0] = X[0] * X[1] * X[2]; 144 145 // create f : X > Y and stop recording 146 CppAD::ADFun<double> f(X, Y); 147 148 // define W(t, u) = (u_0 + dx_0*t)*(u_1 + dx_1*t)*(u_2 + dx_2*t) 149 // use zero order forward to evaluate W0(u) = W(0, u) 150 Vector u(n), W0(m); 151 u[0] = 2.; 152 u[1] = 3.; 153 u[2] = 4.; 154 W0 = f.Forward(0, u); 155 double check; 156 check = u[0]*u[1]*u[2]; 157 ok &= NearEqual(W0[0] , check, 1e10, 1e10); 158 159 // define W_t(t, u) = partial W(t, u) w.r.t t 160 // W_t(t, u) = (u_0 + dx_0*t)*(u_1 + dx_1*t)*dx_2 161 // + (u_0 + dx_0*t)*(u_2 + dx_2*t)*dx_1 162 // + (u_1 + dx_1*t)*(u_2 + dx_2*t)*dx_0 163 // use first order forward mode to evaluate W1(u) = W_t(0, u) 164 Vector dx(n), W1(m); 165 dx[0] = .2; 166 dx[1] = .3; 167 dx[2] = .4; 168 W1 = f.Forward(1, dx); 169 check = u[0]*u[1]*dx[2] + u[0]*u[2]*dx[1] + u[1]*u[2]*dx[0]; 170 ok &= NearEqual(W1[0], check, 1e10, 1e10); 171 172 // define W_tt (t, u) = partial W_t(t, u) w.r.t t 173 // W_tt(t, u) = 2*(u_0 + dx_0*t)*dx_1*dx_2 174 // + 2*(u_1 + dx_1*t)*dx_0*dx_2 175 // + 2*(u_3 + dx_3*t)*dx_0*dx_1 176 // use second order forward to evaluate W2(u) = 1/2 * W_tt(0, u) 177 Vector ddx(n), W2(m); 178 ddx[0] = ddx[1] = ddx[2] = 0.; 179 W2 = f.Forward(2, ddx); 180 check = u[0]*dx[1]*dx[2] + u[1]*dx[0]*dx[2] + u[2]*dx[0]*dx[1]; 181 ok &= NearEqual(W2[0], check, 1e10, 1e10); 182 183 // use third order reverse mode to evaluate derivatives 184 size_t p = 3; 185 Vector w(m), dw(n * p); 186 w[0] = 1.; 187 dw = f.Reverse(p, w); 188 189 // check derivative of W0(u) w.r.t. u 190 ok &= NearEqual(dw[0*p+0], u[1]*u[2], 1e10, 1e10); 191 ok &= NearEqual(dw[1*p+0], u[0]*u[2], 1e10, 1e10); 192 ok &= NearEqual(dw[2*p+0], u[0]*u[1], 1e10, 1e10); 193 194 // check derivative of W1(u) w.r.t. u 195 ok &= NearEqual(dw[0*p+1], u[1]*dx[2] + u[2]*dx[1], 1e10, 1e10); 196 ok &= NearEqual(dw[1*p+1], u[0]*dx[2] + u[2]*dx[0], 1e10, 1e10); 197 ok &= NearEqual(dw[2*p+1], u[0]*dx[1] + u[1]*dx[0], 1e10, 1e10); 198 199 // check derivative of W2(u) w.r.t u 200 ok &= NearEqual(dw[0*p+2], dx[1]*dx[2], 1e10, 1e10); 201 ok &= NearEqual(dw[1*p+2], dx[0]*dx[2], 1e10, 1e10); 202 ok &= NearEqual(dw[2*p+2], dx[0]*dx[1], 1e10, 1e10); 203 204 return ok; 205 } 206 } // End empty namespace 207 208 # include <vector> 209 # include <valarray> 210 bool reverse(void) 211 { bool ok = true; 212 ok &= Reverse(); 213 214 ok &= reverse_any_cases< CppAD::vector <double> >(); 215 ok &= reverse_any_cases< std::vector <double> >(); 216 ok &= reverse_any_cases< std::valarray <double> >(); 217 return ok; 218 } 
trunk/test_more/test_more.cpp
r1633 r1659 60 60 extern bool Pow(void); 61 61 extern bool PowInt(void); 62 extern bool Reverse(void);62 extern bool reverse(void); 63 63 extern bool rev_sparse_hes(void); 64 64 extern bool rev_sparse_jac(void); … … 157 157 ok &= Run( Pow, "Pow" ); 158 158 ok &= Run( PowInt, "PowInt" ); 159 ok &= Run( Reverse, "Reverse" );159 ok &= Run( reverse, "reverse" ); 160 160 ok &= Run( rev_sparse_hes, "rev_sparse_hes" ); 161 161 ok &= Run( rev_sparse_jac, "rev_sparse_jac" );
Note: See TracChangeset
for help on using the changeset viewer.