source: trunk/cppad/local/cos_op.hpp @ 3301

Last change on this file since 3301 was 3301, checked in by bradbell, 6 years ago

merge in multiple forward direcitons from branches/forward_dir

  • Property svn:keywords set to Id
File size: 6.3 KB
Line 
1/* $Id: cos_op.hpp 3301 2014-05-24 05:20:21Z bradbell $ */
2# ifndef CPPAD_COS_OP_INCLUDED
3# define CPPAD_COS_OP_INCLUDED
4
5/* --------------------------------------------------------------------------
6CppAD: C++ Algorithmic Differentiation: Copyright (C) 2003-14 Bradley M. Bell
7
8CppAD is distributed under multiple licenses. This distribution is under
9the terms of the
10                    Eclipse Public License Version 1.0.
11
12A copy of this license is included in the COPYING file of this distribution.
13Please visit http://www.coin-or.org/CppAD/ for information on other licenses.
14-------------------------------------------------------------------------- */
15
16
17namespace CppAD { // BEGIN_CPPAD_NAMESPACE
18/*!
19\file cos_op.hpp
20Forward and reverse mode calculations for z = cos(x).
21*/
22
23/*!
24Compute forward mode Taylor coefficient for result of op = CosOp.
25
26The C++ source code corresponding to this operation is
27\verbatim
28        z = cos(x)
29\endverbatim
30The auxillary result is
31\verbatim
32        y = sin(x)
33\endverbatim
34The value of y, and its derivatives, are computed along with the value
35and derivatives of z.
36
37\copydetails forward_unary2_op
38*/
39template <class Base>
40inline void forward_cos_op(
41        size_t p           ,
42        size_t q           ,
43        size_t i_z         ,
44        size_t i_x         ,
45        size_t cap_order   , 
46        Base*  taylor      )
47{       
48        // check assumptions
49        CPPAD_ASSERT_UNKNOWN( NumArg(CosOp) == 1 );
50        CPPAD_ASSERT_UNKNOWN( NumRes(CosOp) == 2 );
51        CPPAD_ASSERT_UNKNOWN( i_x + 1 < i_z );
52        CPPAD_ASSERT_UNKNOWN( q < cap_order );
53        CPPAD_ASSERT_UNKNOWN( p <= q );
54
55        // Taylor coefficients corresponding to argument and result
56        Base* x = taylor + i_x * cap_order;
57        Base* c = taylor + i_z * cap_order;
58        Base* s = c      -       cap_order;
59
60
61        // rest of this routine is identical for the following cases:
62        // forward_sin_op, forward_cos_op, forward_sinh_op, forward_cosh_op.
63        // (except that there is a sign difference for the hyperbolic case).
64        size_t k;
65        if( p == 0 )
66        {       s[0] = sin( x[0] );
67                c[0] = cos( x[0] );
68                p++;
69        }
70        for(size_t j = p; j <= q; j++)
71        {
72                s[j] = Base(0);
73                c[j] = Base(0);
74                for(k = 1; k <= j; k++)
75                {       s[j] += Base(k) * x[k] * c[j-k];
76                        c[j] -= Base(k) * x[k] * s[j-k];
77                }
78                s[j] /= Base(j);
79                c[j] /= Base(j);
80        }
81}
82/*!
83Compute forward mode Taylor coefficient for result of op = CosOp.
84
85The C++ source code corresponding to this operation is
86\verbatim
87        z = cos(x)
88\endverbatim
89The auxillary result is
90\verbatim
91        y = sin(x)
92\endverbatim
93The value of y, and its derivatives, are computed along with the value
94and derivatives of z.
95
96\copydetails forward_unary2_op_dir
97*/
98template <class Base>
99inline void forward_cos_op_dir(
100        size_t q           ,
101        size_t r           ,
102        size_t i_z         ,
103        size_t i_x         ,
104        size_t cap_order   , 
105        Base*  taylor      )
106{       
107        // check assumptions
108        CPPAD_ASSERT_UNKNOWN( NumArg(CosOp) == 1 );
109        CPPAD_ASSERT_UNKNOWN( NumRes(CosOp) == 2 );
110        CPPAD_ASSERT_UNKNOWN( i_x + 1 < i_z );
111        CPPAD_ASSERT_UNKNOWN( 0 < q );
112        CPPAD_ASSERT_UNKNOWN( q < cap_order );
113
114        // Taylor coefficients corresponding to argument and result
115        size_t num_taylor_per_var = (cap_order-1) * r + 1;
116        Base* x = taylor + i_x * num_taylor_per_var;
117        Base* c = taylor + i_z * num_taylor_per_var;
118        Base* s = c      -       num_taylor_per_var;
119
120
121        // rest of this routine is identical for the following cases:
122        // forward_sin_op, forward_cos_op, forward_sinh_op, forward_cosh_op
123        // (except that there is a sign difference for the hyperbolic case).
124        size_t m = (q-1) * r + 1;
125        for(size_t ell = 0; ell < r; ell++)
126        {       s[m+ell] =   Base(q) * x[m + ell] * c[0];
127                c[m+ell] = - Base(q) * x[m + ell] * s[0];
128                for(size_t k = 1; k < q; k++)
129                {       s[m+ell] += Base(k) * x[(k-1)*r+1+ell] * c[(q-k-1)*r+1+ell];
130                        c[m+ell] -= Base(k) * x[(k-1)*r+1+ell] * s[(q-k-1)*r+1+ell];
131                }
132                s[m+ell] /= Base(q);
133                c[m+ell] /= Base(q);
134        }
135}
136
137/*!
138Compute zero order forward mode Taylor coefficient for result of op = CosOp.
139
140The C++ source code corresponding to this operation is
141\verbatim
142        z = cos(x)
143\endverbatim
144The auxillary result is
145\verbatim
146        y = sin(x)
147\endverbatim
148The value of y is computed along with the value of z.
149
150\copydetails forward_unary2_op_0
151*/
152template <class Base>
153inline void forward_cos_op_0(
154        size_t i_z         ,
155        size_t i_x         ,
156        size_t cap_order   , 
157        Base*  taylor      )
158{
159        // check assumptions
160        CPPAD_ASSERT_UNKNOWN( NumArg(CosOp) == 1 );
161        CPPAD_ASSERT_UNKNOWN( NumRes(CosOp) == 2 );
162        CPPAD_ASSERT_UNKNOWN( i_x + 1 < i_z );
163        CPPAD_ASSERT_UNKNOWN( 0 < cap_order );
164
165        // Taylor coefficients corresponding to argument and result
166        Base* x = taylor + i_x * cap_order;
167        Base* c = taylor + i_z * cap_order;  // called z in documentation
168        Base* s = c      -       cap_order;  // called y in documentation
169
170        c[0] = cos( x[0] );
171        s[0] = sin( x[0] );
172}
173/*!
174Compute reverse mode partial derivatives for result of op = CosOp.
175
176The C++ source code corresponding to this operation is
177\verbatim
178        z = cos(x)
179\endverbatim
180The auxillary result is
181\verbatim
182        y = sin(x)
183\endverbatim
184The value of y is computed along with the value of z.
185
186\copydetails reverse_unary2_op
187*/
188
189template <class Base>
190inline void reverse_cos_op(
191        size_t      d            ,
192        size_t      i_z          ,
193        size_t      i_x          ,
194        size_t      cap_order    , 
195        const Base* taylor       ,
196        size_t      nc_partial   ,
197        Base*       partial      )
198{
199        // check assumptions
200        CPPAD_ASSERT_UNKNOWN( NumArg(CosOp) == 1 );
201        CPPAD_ASSERT_UNKNOWN( NumRes(CosOp) == 2 );
202        CPPAD_ASSERT_UNKNOWN( i_x + 1 < i_z );
203        CPPAD_ASSERT_UNKNOWN( d < cap_order );
204        CPPAD_ASSERT_UNKNOWN( d < nc_partial );
205
206        // Taylor coefficients and partials corresponding to argument
207        const Base* x  = taylor  + i_x * cap_order;
208        Base* px       = partial + i_x * nc_partial;
209
210        // Taylor coefficients and partials corresponding to first result
211        const Base* c  = taylor  + i_z * cap_order; // called z in doc
212        Base* pc       = partial + i_z * nc_partial;
213
214        // Taylor coefficients and partials corresponding to auxillary result
215        const Base* s  = c  - cap_order; // called y in documentation
216        Base* ps       = pc - nc_partial;
217
218
219        // rest of this routine is identical for the following cases:
220        // reverse_sin_op, reverse_cos_op, reverse_sinh_op, reverse_cosh_op.
221        size_t j = d;
222        size_t k;
223        while(j)
224        {
225                ps[j]   /= Base(j);
226                pc[j]   /= Base(j);
227                for(k = 1; k <= j; k++)
228                {
229                        px[k]   += ps[j] * Base(k) * c[j-k];
230                        px[k]   -= pc[j] * Base(k) * s[j-k];
231       
232                        ps[j-k] -= pc[j] * Base(k) * x[k];
233                        pc[j-k] += ps[j] * Base(k) * x[k];
234
235                }
236                --j;
237        }
238        px[0] += ps[0] * c[0];
239        px[0] -= pc[0] * s[0];
240}
241
242} // END_CPPAD_NAMESPACE
243# endif
Note: See TracBrowser for help on using the repository browser.