source: trunk/test_more/sin_cos.cpp @ 1613

Last change on this file since 1613 was 1370, checked in by bradbell, 11 years ago

trunk: Fix svn_add_id.sh and use it set Id property for some missed files.

  • Property svn:keywords set to Id
File size: 8.7 KB
Line 
1/* $Id: sin_cos.cpp 1370 2009-05-31 05:31:50Z bradbell $ */
2/* --------------------------------------------------------------------------
3CppAD: C++ Algorithmic Differentiation: Copyright (C) 2003-07 Bradley M. Bell
4
5CppAD is distributed under multiple licenses. This distribution is under
6the terms of the
7                    Common Public License Version 1.0.
8
9A copy of this license is included in the COPYING file of this distribution.
10Please visit http://www.coin-or.org/CppAD/ for information on other licenses.
11-------------------------------------------------------------------------- */
12
13/*
14Comprehensive test of Trigonometric and Hyperbolic Sine and Cosine
15*/
16
17# include <cppad/cppad.hpp>
18# include <cmath>
19
20namespace { // Begin empty namespace
21
22bool Sin(void)
23{       bool ok = true;
24
25        using CppAD::sin;
26        using CppAD::cos;
27        using namespace CppAD;
28
29        // independent variable vector
30        double x = .5;
31        double y = .8;
32        CPPAD_TEST_VECTOR< AD<double> > X(2);
33        X[0]     = x;
34        X[1]     = y;
35        Independent(X);
36
37        // dependent variable vector
38        CPPAD_TEST_VECTOR< AD<double> > Z(1);
39        AD<double> U = X[0] * X[1];
40        Z[0] = sin( U ); 
41
42        // create f: X -> Z and vectors used for derivative calculations
43        // f(x, y) = sin(x, y)
44        ADFun<double> f(X, Z); 
45        CPPAD_TEST_VECTOR<double> v( 2 );
46        CPPAD_TEST_VECTOR<double> w( 1 );
47
48        // check value
49        double sin_u = sin( Value(U) );
50        double cos_u = cos( Value(U) );
51
52        ok &= NearEqual(sin_u, Value(Z[0]),  1e-10 , 1e-10);
53
54        // forward computation of partials w.r.t. u
55        size_t j;
56        size_t p     = 5;
57        double jfac  = 1.;
58        v[0]         = 1.;  // differential w.r.t. x
59        v[1]         = 0;   // differential w.r.t. y
60        double yj    = 1;   // y^j
61        for(j = 1; j < p; j++)
62        {       w      = f.Forward(j, v);       
63
64                // compute j-th power of y
65                yj *= y ;
66
67                // compute j-th derivartive of sin function
68                double sinj;
69                if( j % 4 == 1 )
70                        sinj = cos_u;
71                else if( j % 4 == 2 )
72                        sinj = -sin_u;
73                else if( j % 4 == 3 )
74                        sinj = -cos_u;
75                else    sinj = sin_u;
76
77                jfac *= j;
78
79                // check j-th derivative of z w.r.t x
80                ok &= NearEqual(jfac*w[0], sinj * yj, 1e-10 , 1e-10); 
81
82                v[0]  = 0.;
83        }
84
85        // reverse computation of partials of Taylor coefficients
86        CPPAD_TEST_VECTOR<double> r( 2 * p); 
87        w[0]  = 1.;
88        r     = f.Reverse(p, w);
89        jfac  = 1.;
90        yj    = 1.;
91        double sinjp = 0.;
92        for(j = 0; j < p; j++)
93        {
94                double sinj = sinjp;
95
96                // compute j+1 derivative of sin funciton
97                if( j % 4 == 0 )
98                        sinjp = cos_u;
99                else if( j % 4 == 1 )
100                        sinjp = -sin_u;
101                else if( j % 4 == 2 )
102                        sinjp = -cos_u;
103                else    sinjp = sin_u;
104
105                // derivative w.r.t x of sin^{(j)} (x * y) * y^j
106                ok &= NearEqual(jfac*r[0+j], sinjp * yj * y, 1e-10 , 1e-10);
107
108                // derivative w.r.t y of sin^{(j)} (x * y) * y^j
109                double value = sinjp * yj * x + j * sinj * yj / y;
110                ok &= NearEqual(jfac*r[p+j], value , 1e-10 , 1e-10);
111
112                jfac  *= (j + 1);
113                yj    *= y;
114        }
115
116        return ok;
117}
118
119bool Cos(void)
120{       bool ok = true;
121
122        using CppAD::sin;
123        using CppAD::cos;
124        using namespace CppAD;
125
126        // independent variable vector
127        double x = .5;
128        double y = .8;
129        CPPAD_TEST_VECTOR< AD<double> > X(2);
130        X[0]     = x;
131        X[1]     = y;
132        Independent(X);
133
134        // dependent variable vector
135        CPPAD_TEST_VECTOR< AD<double> > Z(1);
136        AD<double> U = X[0] * X[1];
137        Z[0] = cos( U ); 
138
139        // create f: X -> Z and vectors used for derivative calculations
140        // f(x, y) = cos(x, y)
141        ADFun<double> f(X, Z); 
142        CPPAD_TEST_VECTOR<double> v( 2 );
143        CPPAD_TEST_VECTOR<double> w( 1 );
144
145        // check value
146        double sin_u = sin( Value(U) );
147        double cos_u = cos( Value(U) );
148
149        ok &= NearEqual(cos_u, Value(Z[0]),  1e-10 , 1e-10);
150
151        // forward computation of partials w.r.t. u
152        size_t j;
153        size_t p     = 5;
154        double jfac  = 1.;
155        v[0]         = 1.;  // differential w.r.t. x
156        v[1]         = 0;   // differential w.r.t. y
157        double yj    = 1;   // y^j
158        for(j = 1; j < p; j++)
159        {       w      = f.Forward(j, v);       
160
161                // compute j-th power of y
162                yj *= y ;
163
164                // compute j-th derivartive of cos function
165                double cosj;
166                if( j % 4 == 1 )
167                        cosj = -sin_u;
168                else if( j % 4 == 2 )
169                        cosj = -cos_u;
170                else if( j % 4 == 3 )
171                        cosj = sin_u;
172                else    cosj = cos_u;
173
174                jfac *= j;
175
176                // check j-th derivative of z w.r.t x
177                ok &= NearEqual(jfac*w[0], cosj * yj, 1e-10 , 1e-10); 
178
179                v[0]  = 0.;
180        }
181
182        // reverse computation of partials of Taylor coefficients
183        CPPAD_TEST_VECTOR<double> r( 2 * p); 
184        w[0]  = 1.;
185        r     = f.Reverse(p, w);
186        jfac  = 1.;
187        yj    = 1.;
188        double cosjp = 0.;
189        for(j = 0; j < p; j++)
190        {
191                double cosj = cosjp;
192
193                // compute j+1 derivative of cos funciton
194                if( j % 4 == 0 )
195                        cosjp = -sin_u;
196                else if( j % 4 == 1 )
197                        cosjp = -cos_u;
198                else if( j % 4 == 2 )
199                        cosjp = sin_u;
200                else    cosjp = cos_u;
201
202                // derivative w.r.t x of cos^{(j)} (x * y) * y^j
203                ok &= NearEqual(jfac*r[0+j], cosjp * yj * y, 1e-10 , 1e-10);
204
205                // derivative w.r.t y of cos^{(j)} (x * y) * y^j
206                double value = cosjp * yj * x + j * cosj * yj / y;
207                ok &= NearEqual(jfac*r[p+j], value , 1e-10 , 1e-10);
208
209                jfac  *= (j + 1);
210                yj    *= y;
211        }
212
213        return ok;
214}
215
216bool Cosh(void)
217{       bool ok = true;
218
219        using CppAD::sinh;
220        using CppAD::cosh;
221        using namespace CppAD;
222
223        // independent variable vector
224        double x = .5;
225        double y = .8;
226        CPPAD_TEST_VECTOR< AD<double> > X(2);
227        X[0]     = x;
228        X[1]     = y;
229        Independent(X);
230
231        // dependent variable vector
232        CPPAD_TEST_VECTOR< AD<double> > Z(1);
233        AD<double> U = X[0] * X[1];
234        Z[0] = cosh( U ); 
235
236        // create f: X -> Z and vectors used for derivative calculations
237        // f(x, y) = cosh(x, y)
238        ADFun<double> f(X, Z); 
239        CPPAD_TEST_VECTOR<double> v( 2 );
240        CPPAD_TEST_VECTOR<double> w( 1 );
241
242        // check value
243        double sinh_u = sinh( Value(U) );
244        double cosh_u = cosh( Value(U) );
245
246        ok &= NearEqual(cosh_u, Value(Z[0]),  1e-10 , 1e-10);
247
248        // forward computation of partials w.r.t. u
249        size_t j;
250        size_t p     = 5;
251        double jfac  = 1.;
252        v[0]         = 1.;  // differential w.r.t. x
253        v[1]         = 0;   // differential w.r.t. y
254        double yj    = 1;   // y^j
255        for(j = 1; j < p; j++)
256        {       w      = f.Forward(j, v);       
257
258                // compute j-th power of y
259                yj *= y ;
260
261                // compute j-th derivartive of cosh function
262                double coshj;
263                if( j % 2 == 1 )
264                        coshj = sinh_u;
265                else    coshj = cosh_u;
266
267                jfac *= j;
268
269                // check j-th derivative of z w.r.t x
270                ok &= NearEqual(jfac*w[0], coshj * yj, 1e-10 , 1e-10); 
271
272                v[0]  = 0.;
273        }
274
275        // reverse computation of partials of Taylor coefficients
276        CPPAD_TEST_VECTOR<double> r( 2 * p); 
277        w[0]  = 1.;
278        r     = f.Reverse(p, w);
279        jfac  = 1.;
280        yj    = 1.;
281        double coshjp = 0.;
282        for(j = 0; j < p; j++)
283        {
284                double coshj = coshjp;
285
286                // compute j+1 derivative of cosh funciton
287                if( j % 2 == 0 )
288                        coshjp = sinh_u;
289                else    coshjp = cosh_u;
290
291                // derivative w.r.t x of cosh^{(j)} (x * y) * y^j
292                ok &= NearEqual(jfac*r[0+j], coshjp * yj * y, 1e-10 , 1e-10);
293
294                // derivative w.r.t y of cosh^{(j)} (x * y) * y^j
295                double value = coshjp * yj * x + j * coshj * yj / y;
296                ok &= NearEqual(jfac*r[p+j], value , 1e-10 , 1e-10);
297
298                jfac  *= (j + 1);
299                yj    *= y;
300        }
301
302        return ok;
303}
304
305bool Sinh(void)
306{       bool ok = true;
307
308        using CppAD::sinh;
309        using CppAD::cosh;
310        using namespace CppAD;
311
312        // independent variable vector
313        double x = .5;
314        double y = .8;
315        CPPAD_TEST_VECTOR< AD<double> > X(2);
316        X[0]     = x;
317        X[1]     = y;
318        Independent(X);
319
320        // dependent variable vector
321        CPPAD_TEST_VECTOR< AD<double> > Z(1);
322        AD<double> U = X[0] * X[1];
323        Z[0] = sinh( U ); 
324
325        // create f: X -> Z and vectors used for derivative calculations
326        // f(x, y) = sinh(x, y)
327        ADFun<double> f(X, Z); 
328        CPPAD_TEST_VECTOR<double> v( 2 );
329        CPPAD_TEST_VECTOR<double> w( 1 );
330
331        // check value
332        double sinh_u = sinh( Value(U) );
333        double cosh_u = cosh( Value(U) );
334
335        ok &= NearEqual(sinh_u, Value(Z[0]),  1e-10 , 1e-10);
336
337        // forward computation of partials w.r.t. u
338        size_t j;
339        size_t p     = 5;
340        double jfac  = 1.;
341        v[0]         = 1.;  // differential w.r.t. x
342        v[1]         = 0;   // differential w.r.t. y
343        double yj    = 1;   // y^j
344        for(j = 1; j < p; j++)
345        {       w      = f.Forward(j, v);       
346
347                // compute j-th power of y
348                yj *= y ;
349
350                // compute j-th derivartive of sinh function
351                double sinhj;
352                if( j % 2 == 1 )
353                        sinhj = cosh_u;
354                else    sinhj = sinh_u;
355
356                jfac *= j;
357
358                // check j-th derivative of z w.r.t x
359                ok &= NearEqual(jfac*w[0], sinhj * yj, 1e-10 , 1e-10); 
360
361                v[0]  = 0.;
362        }
363
364        // reverse computation of partials of Taylor coefficients
365        CPPAD_TEST_VECTOR<double> r( 2 * p); 
366        w[0]  = 1.;
367        r     = f.Reverse(p, w);
368        jfac  = 1.;
369        yj    = 1.;
370        double sinhjp = 0.;
371        for(j = 0; j < p; j++)
372        {
373                double sinhj = sinhjp;
374
375                // compute j+1 derivative of sinh funciton
376                if( j % 2 == 0 )
377                        sinhjp = cosh_u;
378                else    sinhjp = sinh_u;
379
380                // derivative w.r.t x of sinh^{(j)} (x * y) * y^j
381                ok &= NearEqual(jfac*r[0+j], sinhjp * yj * y, 1e-10 , 1e-10);
382
383                // derivative w.r.t y of sinh^{(j)} (x * y) * y^j
384                double value = sinhjp * yj * x + j * sinhj * yj / y;
385                ok &= NearEqual(jfac*r[p+j], value , 1e-10 , 1e-10);
386
387                jfac  *= (j + 1);
388                yj    *= y;
389        }
390
391        return ok;
392}
393
394} // End empty namespace
395
396bool SinCos(void)
397{       bool ok = Sin() & Cos() & Cosh() & Sinh();
398        return ok;
399}
400       
Note: See TracBrowser for help on using the repository browser.