1 /* --------------------------------------------------------------------------
2 CppAD: C++ Algorithmic Differentiation: Copyright (C) 2003-17 Bradley M. Bell
3 
4 CppAD is distributed under the terms of the
5              Eclipse Public License Version 2.0.
6 
7 This Source Code may also be made available under the following
8 Secondary License when the conditions for such availability set forth
9 in the Eclipse Public License, Version 2.0 are satisfied:
10       GNU General Public License, Version 2.0 or later.
11 ---------------------------------------------------------------------------- */
12 
13 
14 /*
15 Two old log examples now used just for validation testing
16 */
17 
18 # include <cppad/cppad.hpp>
19 
20 namespace { // BEGIN empty namespace
21 
LogTestOne(void)22 bool LogTestOne(void)
23 {   bool ok = true;
24     using CppAD::log;
25     using namespace CppAD;
26     double eps99 = 99.0 * std::numeric_limits<double>::epsilon();
27 
28     // independent variable vector, indices, values, and declaration
29     CPPAD_TESTVECTOR(AD<double>) U(1);
30     size_t s = 0;
31     U[s]     = 2.;
32     Independent(U);
33 
34     // dependent variable vector, indices, and values
35     CPPAD_TESTVECTOR(AD<double>) Z(2);
36     size_t x = 0;
37     size_t y = 1;
38     Z[x]     = log(U[s]);
39     Z[y]     = log(Z[x]);
40 
41     // define f : U -> Z and vectors for derivative calculations
42     ADFun<double> f(U, Z);
43     CPPAD_TESTVECTOR(double) v( f.Domain() );
44     CPPAD_TESTVECTOR(double) w( f.Range() );
45 
46     // check values
47     ok &= NearEqual(Z[x] , log(2.),  eps99 , eps99);
48     ok &= NearEqual(Z[y] , log( log(2.) ),  eps99 , eps99);
49 
50     // forward computation of partials w.r.t. s
51     v[s] = 1.;
52     w    = f.Forward(1, v);
53     ok &= NearEqual(w[x], 1. / U[s],          eps99 , eps99); // dx/ds
54     ok &= NearEqual(w[y], 1. / (U[s] * Z[x]), eps99 , eps99); // dy/ds
55 
56     // reverse computation of partials of y
57     w[x] = 0.;
58     w[y] = 1.;
59     v    = f.Reverse(1,w);
60     ok &= NearEqual(v[s], 1. / (U[s] * Z[x]), eps99 , eps99); // dy/ds
61 
62     // forward computation of second partials w.r.t. s
63     v[s] = 1.;
64     w    = f.Forward(1, v);
65     v[s] = 0.;
66     w    = f.Forward(2, v);
67     ok &= NearEqual(
68         2. * w[y] ,
69         - 1. / (Z[x]*Z[x]*U[s]*U[s]) - 1. / (Z[x]*U[s]*U[s]),
70         eps99 ,
71         eps99
72     );
73 
74     // reverse computation of second partials of y
75     CPPAD_TESTVECTOR(double) r( f.Domain() * 2 );
76     w[x] = 0.;
77     w[y] = 1.;
78     r    = f.Reverse(2, w);
79     ok &= NearEqual(
80         r[2 * s + 1] ,
81         - 1. / (Z[x]*Z[x]*U[s]*U[s]) - 1. / (Z[x]*U[s]*U[s]),
82         eps99 ,
83         eps99
84     );
85 
86     return ok;
87 }
LogTestTwo(void)88 bool LogTestTwo(void)
89 {   bool ok = true;
90     using CppAD::log;
91     using namespace CppAD;
92     double eps99 = 99.0 * std::numeric_limits<double>::epsilon();
93 
94     // independent variable vector
95     CPPAD_TESTVECTOR(AD<double>) U(1);
96     U[0]     = 1.;
97     Independent(U);
98 
99     // a temporary values
100     AD<double> x = exp(U[0]);
101 
102     // dependent variable vector
103     CPPAD_TESTVECTOR(AD<double>) Z(1);
104     Z[0] = log(x);       // log( exp(u) )
105 
106     // create f: U -> Z and vectors used for derivative calculations
107     ADFun<double> f(U, Z);
108     CPPAD_TESTVECTOR(double) v(1);
109     CPPAD_TESTVECTOR(double) w(1);
110 
111     // check value
112     ok &= NearEqual(U[0] , Z[0],  eps99 , eps99);
113 
114     // forward computation of partials w.r.t. u
115     size_t j;
116     size_t p     = 5;
117     double jfac  = 1.;
118     double value = 1.;
119     v[0]         = 1.;
120     for(j = 1; j < p; j++)
121     {   jfac *= double(j);
122         w     = f.Forward(j, v);
123         ok &= NearEqual(w[0], value/jfac, eps99, eps99); // d^jz/du^j
124         v[0]  = 0.;
125         value = 0.;
126     }
127 
128     // reverse computation of partials of Taylor coefficients
129     CPPAD_TESTVECTOR(double) r(p);
130     w[0]  = 1.;
131     r     = f.Reverse(p, w);
132     jfac  = 1.;
133     value = 1.;
134     for(j = 0; j < p; j++)
135     {   ok &= NearEqual(r[j], value/jfac, eps99, eps99); // d^jz/du^j
136         jfac *= double(j + 1);
137         value = 0.;
138     }
139 
140     return ok;
141 }
142 
143 } // END empty namespace
144 
log(void)145 bool log(void)
146 {   bool ok = true;
147     ok &= LogTestOne();
148     ok &= LogTestTwo();
149     return ok;
150 }
151