Changeset 3727 for trunk/test_more
 Timestamp:
 Sep 22, 2015 1:39:26 PM (4 years ago)
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

trunk/test_more/mul_cond.cpp
r3724 r3727 20 20 bool ok = true; 21 21 using CppAD::vector; 22 using CppAD::NearEqual; 22 23 double eps = 10. * std::numeric_limits<double>::epsilon(); 23 24 // … … 25 26 typedef CppAD::AD<a1double> a2double; 26 27 // 27 a2double a2zero = a2double(0.0); 28 a2double a2one = a2double(1.0); 28 a1double a1zero = 0.0; 29 a2double a2zero = a1zero; 30 a1double a1one = 1.0; 31 a2double a2one = a1one; 29 32 // 30 33 //  … … 32 35 size_t n = 1; 33 36 size_t m = 2; 37 size_t i; 34 38 vector<a2double> a2x(n), a2y(m); 35 39 a2x[0] = a2double( 5.0 ); 36 40 Independent(a2x); 37 41 // 42 i = 0; 38 43 // div 39 a2y[0] = CondExpGt(a2x[0], a2zero, a2one / a2x[0], a2zero); 40 41 // abs 42 a2y[1] = CondExpGt(a2x[0], a2zero, abs( a2one / a2x[0] ), a2zero); 43 44 a2y[i++] = CondExpGt(a2x[0], a2zero, a2one / a2x[0], a2zero); 45 // acosh 46 a2y[i++] = CondExpGt(a2x[0], a2zero, acosh( a2x[0] ), a2zero); 47 // 44 48 CppAD::ADFun<a1double> a1f; 45 49 a1f.Dependent(a2x, a2y); 46 50 //  47 // create g = sum_i f_i '(x) 48 vector<a1double> a1x(n), a1dy(n), a1w(m); 51 // create h = f(x) 52 vector<a1double> a1x(n), a1y(m); 53 a1x[0] = 5.0; 54 // 55 Independent(a1x); 56 i = 0; 57 // div 58 a1y[i++] = CondExpGt(a1x[0], a1zero, a1one / a1x[0], a1zero); 59 // acosh 60 a1y[i++] = CondExpGt(a1x[0], a1zero, acosh( a1x[0] ), a1zero); 61 // 62 CppAD::ADFun<double> h; 63 h.Dependent(a1x, a1y); 64 //  65 // create g = f'(x) 66 vector<a1double> a1dy(m), a1w(m); 49 67 a1x[0] = 2.0; 50 for( size_ti = 0; i < m; i++)51 a1w[i] = 1.0;68 for(i = 0; i < m; i++) 69 a1w[i] = 0.0; 52 70 // 53 71 Independent(a1x); 54 72 a1f.Forward(0, a1x); 55 a1dy = a1f.Reverse(1, a1w); 73 // 74 for(i = 0; i < m; i++) 75 { a1w[i] = 1.0; 76 vector<a1double> dyi_dx = a1f.Reverse(1, a1w); 77 a1dy[i] = dyi_dx[0]; 78 a1w[i] = 0.0; 79 } 56 80 CppAD::ADFun<double> g; 57 81 g.Dependent(a1x, a1dy); 58 82 //  59 // check g where x[0] > 0 60 vector<double> x(1), y(1); 61 x[0] = 2.0; 62 y = g.Forward(0, x); 63 double check = 0.0; 64 check += 1.0 / (x[0] * x[0]); // 1 / x[0] 65 check += 1.0 / (x[0] * x[0]); // abs( 1 / x[0] ) 66 ok &= CppAD::NearEqual(y[0], check, eps, eps); 83 // check case where x[0] > 0 84 vector<double> x(1), dx(1), dg(m), dh(m); 85 x[0] = 2.0; 86 dx[0] = 1.0; 87 h.Forward(0, x); 88 dh = h.Forward(1, dx); 89 dg = g.Forward(0, x); 90 for(i = 0; i < m; i++) 91 ok &= NearEqual(dg[i], dh[i], eps, eps); 67 92 //  68 // check gwhere x[0] = 093 // check case where x[0] = 0 69 94 x[0] = 0.0; 70 y = g.Forward(0, x); 71 ok &= CppAD::NearEqual(y[0], 0.0, eps, eps); 95 dg = g.Forward(0, x); 96 h.Forward(0, x); 97 dh = h.Forward(1, dx); 98 for(i = 0; i < m; i++) 99 { ok &= dg[i] == 0.0; 100 ok &= dh[i] == 0.0; 101 } 72 102 //  73 103 return ok;
Note: See TracChangeset
for help on using the changeset viewer.