CppAD: A C++ Algorithmic Differentiation Package 20110419
sin_op.hpp
Go to the documentation of this file.
00001 /* $Id$ */
00002 # ifndef CPPAD_SIN_OP_INCLUDED
00003 # define CPPAD_SIN_OP_INCLUDED
00004 
00005 /* --------------------------------------------------------------------------
00006 CppAD: C++ Algorithmic Differentiation: Copyright (C) 2003-10 Bradley M. Bell
00007 
00008 CppAD is distributed under multiple licenses. This distribution is under
00009 the terms of the 
00010                     Common Public License Version 1.0.
00011 
00012 A copy of this license is included in the COPYING file of this distribution.
00013 Please visit http://www.coin-or.org/CppAD/ for information on other licenses.
00014 -------------------------------------------------------------------------- */
00015 
00016 
00017 CPPAD_BEGIN_NAMESPACE
00018 /*!
00019 \file sin_op.hpp
00020 Forward and reverse mode calculations for z = sin(x).
00021 */
00022 
00023 
00024 /*!
00025 Compute forward mode Taylor coefficient for result of op = SinOp.
00026 
00027 The C++ source code corresponding to this operation is
00028 \verbatim
00029         z = sin(x)
00030 \endverbatim
00031 The auxillary result is
00032 \verbatim
00033         y = cos(x)
00034 \endverbatim
00035 The value of y, and its derivatives, are computed along with the value
00036 and derivatives of z.
00037 
00038 \copydetails forward_unary2_op
00039 */
00040 template <class Base>
00041 inline void forward_sin_op(
00042         size_t j           ,
00043         size_t i_z         ,
00044         size_t i_x         ,
00045         size_t nc_taylor   , 
00046         Base*  taylor      )
00047 {       
00048         // check assumptions
00049         CPPAD_ASSERT_UNKNOWN( NumArg(SinOp) == 1 );
00050         CPPAD_ASSERT_UNKNOWN( NumRes(SinOp) == 2 );
00051         CPPAD_ASSERT_UNKNOWN( i_x + 1 < i_z );
00052         CPPAD_ASSERT_UNKNOWN( j < nc_taylor );
00053 
00054         // Taylor coefficients corresponding to argument and result
00055         Base* x = taylor + i_x * nc_taylor;
00056         Base* s = taylor + i_z * nc_taylor;
00057         Base* c = s      -       nc_taylor;
00058 
00059         // rest of this routine is identical for the following cases:
00060         // forward_sin_op, forward_cos_op, forward_sinh_op, forward_cosh_op.
00061         size_t k;
00062         if( j == 0 )
00063         {       s[j] = sin( x[0] );
00064                 c[j] = cos( x[0] );
00065         }
00066         else
00067         {
00068                 s[j] = Base(0);
00069                 c[j] = Base(0);
00070                 for(k = 1; k <= j; k++)
00071                 {       s[j] += Base(k) * x[k] * c[j-k];
00072                         c[j] -= Base(k) * x[k] * s[j-k];
00073                 }
00074                 s[j] /= Base(j);
00075                 c[j] /= Base(j);
00076         }
00077 }
00078 
00079 
00080 /*!
00081 Compute zero order forward mode Taylor coefficient for result of op = SinOp.
00082 
00083 The C++ source code corresponding to this operation is
00084 \verbatim
00085         z = sin(x)
00086 \endverbatim
00087 The auxillary result is
00088 \verbatim
00089         y = cos(x)
00090 \endverbatim
00091 The value of y is computed along with the value of z.
00092 
00093 \copydetails forward_unary2_op_0
00094 */
00095 template <class Base>
00096 inline void forward_sin_op_0(
00097         size_t i_z         ,
00098         size_t i_x         ,
00099         size_t nc_taylor   , 
00100         Base*  taylor      )
00101 {
00102         // check assumptions
00103         CPPAD_ASSERT_UNKNOWN( NumArg(SinOp) == 1 );
00104         CPPAD_ASSERT_UNKNOWN( NumRes(SinOp) == 2 );
00105         CPPAD_ASSERT_UNKNOWN( i_x + 1 < i_z );
00106         CPPAD_ASSERT_UNKNOWN( 0 < nc_taylor );
00107 
00108         // Taylor coefficients corresponding to argument and result
00109         Base* x = taylor + i_x * nc_taylor;
00110         Base* s = taylor + i_z * nc_taylor;  // called z in documentation
00111         Base* c = s      -       nc_taylor;  // called y in documentation
00112 
00113         s[0] = sin( x[0] );
00114         c[0] = cos( x[0] );
00115 }
00116 
00117 /*!
00118 Compute reverse mode partial derivatives for result of op = SinOp.
00119 
00120 The C++ source code corresponding to this operation is
00121 \verbatim
00122         z = sin(x)
00123 \endverbatim
00124 The auxillary result is
00125 \verbatim
00126         y = cos(x)
00127 \endverbatim
00128 The value of y is computed along with the value of z.
00129 
00130 \copydetails reverse_unary2_op
00131 */
00132 
00133 template <class Base>
00134 inline void reverse_sin_op(
00135         size_t      d            ,
00136         size_t      i_z          ,
00137         size_t      i_x          ,
00138         size_t      nc_taylor    , 
00139         const Base* taylor       ,
00140         size_t      nc_partial   ,
00141         Base*       partial      )
00142 {
00143         // check assumptions
00144         CPPAD_ASSERT_UNKNOWN( NumArg(SinOp) == 1 );
00145         CPPAD_ASSERT_UNKNOWN( NumRes(SinOp) == 2 );
00146         CPPAD_ASSERT_UNKNOWN( i_x + 1 < i_z );
00147         CPPAD_ASSERT_UNKNOWN( d < nc_taylor );
00148         CPPAD_ASSERT_UNKNOWN( d < nc_partial );
00149 
00150         // Taylor coefficients and partials corresponding to argument
00151         const Base* x  = taylor  + i_x * nc_taylor;
00152         Base* px       = partial + i_x * nc_partial;
00153 
00154         // Taylor coefficients and partials corresponding to first result
00155         const Base* s  = taylor  + i_z * nc_taylor; // called z in doc
00156         Base* ps       = partial + i_z * nc_partial;
00157 
00158         // Taylor coefficients and partials corresponding to auxillary result
00159         const Base* c  = s  - nc_taylor; // called y in documentation
00160         Base* pc       = ps - nc_partial;
00161 
00162         // rest of this routine is identical for the following cases:
00163         // reverse_sin_op, reverse_cos_op, reverse_sinh_op, reverse_cosh_op.
00164         size_t j = d;
00165         size_t k;
00166         while(j)
00167         {
00168                 ps[j]   /= Base(j);
00169                 pc[j]   /= Base(j);
00170                 for(k = 1; k <= j; k++)
00171                 {
00172                         px[k]   += ps[j] * Base(k) * c[j-k];
00173                         px[k]   -= pc[j] * Base(k) * s[j-k];
00174         
00175                         ps[j-k] -= pc[j] * Base(k) * x[k];
00176                         pc[j-k] += ps[j] * Base(k) * x[k];
00177 
00178                 }
00179                 --j;
00180         }
00181         px[0] += ps[0] * c[0];
00182         px[0] -= pc[0] * s[0];
00183 }
00184 
00185 CPPAD_END_NAMESPACE
00186 # endif