CppAD: A C++ Algorithmic Differentiation Package  20130918
div_op.hpp
Go to the documentation of this file.
00001 /* $Id$ */
00002 # ifndef CPPAD_DIV_OP_INCLUDED
00003 # define CPPAD_DIV_OP_INCLUDED
00004 
00005 /* --------------------------------------------------------------------------
00006 CppAD: C++ Algorithmic Differentiation: Copyright (C) 2003-14 Bradley M. Bell
00007 
00008 CppAD is distributed under multiple licenses. This distribution is under
00009 the terms of the 
00010                     Eclipse 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 namespace CppAD { // BEGIN_CPPAD_NAMESPACE
00017 /*!
00018 \file div_op.hpp
00019 Forward and reverse mode calculations for z = x / y.
00020 */
00021 
00022 // --------------------------- Divvv -----------------------------------------
00023 /*!
00024 Compute forward mode Taylor coefficients for result of op = DivvvOp.
00025 
00026 The C++ source code corresponding to this operation is
00027 \verbatim
00028      z = x / y
00029 \endverbatim
00030 In the documentation below,
00031 this operations is for the case where both x and y are variables
00032 and the argument \a parameter is not used.
00033 
00034 \copydetails forward_binary_op
00035 */
00036 
00037 template <class Base>
00038 inline void forward_divvv_op(
00039      size_t        p           , 
00040      size_t        q           , 
00041      size_t        i_z         ,
00042      const addr_t* arg         ,
00043      const Base*   parameter   ,
00044      size_t        nc_taylor   ,
00045      Base*         taylor      )
00046 {
00047      // check assumptions
00048      CPPAD_ASSERT_UNKNOWN( NumArg(DivvvOp) == 2 );
00049      CPPAD_ASSERT_UNKNOWN( NumRes(DivvvOp) == 1 );
00050      CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < i_z );
00051      CPPAD_ASSERT_UNKNOWN( size_t(arg[1]) < i_z );
00052      CPPAD_ASSERT_UNKNOWN( q < nc_taylor );
00053      CPPAD_ASSERT_UNKNOWN( p <= q );
00054 
00055      // Taylor coefficients corresponding to arguments and result
00056      Base* x = taylor + arg[0] * nc_taylor;
00057      Base* y = taylor + arg[1] * nc_taylor;
00058      Base* z = taylor + i_z    * nc_taylor;
00059 
00060 
00061      // Using CondExp, it can make sense to divide by zero,
00062      // so do not make it an error.
00063      size_t k;
00064      for(size_t d = p; d <= q; d++)
00065      {    z[d] = x[d];
00066           for(k = 1; k <= d; k++)
00067                z[d] -= z[d-k] * y[k];
00068           z[d] /= y[0];
00069      }
00070 }
00071 
00072 
00073 /*!
00074 Compute zero order forward mode Taylor coefficients for result of op = DivvvOp.
00075 
00076 The C++ source code corresponding to this operation is
00077 \verbatim
00078      z = x / y
00079 \endverbatim
00080 In the documentation below,
00081 this operations is for the case where both x and y are variables
00082 and the argument \a parameter is not used.
00083 
00084 \copydetails forward_binary_op_0
00085 */
00086 
00087 template <class Base>
00088 inline void forward_divvv_op_0(
00089      size_t        i_z         ,
00090      const addr_t* arg         ,
00091      const Base*   parameter   ,
00092      size_t        nc_taylor   ,
00093      Base*         taylor      )
00094 {
00095      // check assumptions
00096      CPPAD_ASSERT_UNKNOWN( NumArg(DivvvOp) == 2 );
00097      CPPAD_ASSERT_UNKNOWN( NumRes(DivvvOp) == 1 );
00098      CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < i_z );
00099      CPPAD_ASSERT_UNKNOWN( size_t(arg[1]) < i_z );
00100 
00101      // Taylor coefficients corresponding to arguments and result
00102      Base* x = taylor + arg[0] * nc_taylor;
00103      Base* y = taylor + arg[1] * nc_taylor;
00104      Base* z = taylor + i_z    * nc_taylor;
00105 
00106      z[0] = x[0] / y[0];
00107 }
00108 
00109 /*!
00110 Compute reverse mode partial derivatives for result of op = DivvvOp.
00111 
00112 The C++ source code corresponding to this operation is
00113 \verbatim
00114      z = x / y
00115 \endverbatim
00116 In the documentation below,
00117 this operations is for the case where both x and y are variables
00118 and the argument \a parameter is not used.
00119 
00120 \copydetails reverse_binary_op
00121 */
00122 
00123 template <class Base>
00124 inline void reverse_divvv_op(
00125      size_t        d           , 
00126      size_t        i_z         ,
00127      const addr_t* arg         ,
00128      const Base*   parameter   ,
00129      size_t        nc_taylor   ,
00130      const Base*   taylor      ,
00131      size_t        nc_partial  ,
00132      Base*         partial     )
00133 {
00134      // check assumptions
00135      CPPAD_ASSERT_UNKNOWN( NumArg(DivvvOp) == 2 );
00136      CPPAD_ASSERT_UNKNOWN( NumRes(DivvvOp) == 1 );
00137      CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < i_z );
00138      CPPAD_ASSERT_UNKNOWN( size_t(arg[1]) < i_z );
00139      CPPAD_ASSERT_UNKNOWN( d < nc_taylor );
00140      CPPAD_ASSERT_UNKNOWN( d < nc_partial );
00141 
00142      // Arguments
00143      const Base* y  = taylor + arg[1] * nc_taylor;
00144      const Base* z  = taylor + i_z    * nc_taylor;
00145 
00146      // Partial derivatives corresponding to arguments and result
00147      Base* px = partial + arg[0] * nc_partial;
00148      Base* py = partial + arg[1] * nc_partial;
00149      Base* pz = partial + i_z    * nc_partial;
00150 
00151      // Using CondExp, it can make sense to divide by zero
00152      // so do not make it an error.
00153 
00154      size_t k;
00155      // number of indices to access
00156      size_t j = d + 1;
00157      while(j)
00158      {    --j;
00159           // scale partial w.r.t. z[j]
00160           pz[j] /= y[0];
00161 
00162           px[j] += pz[j];
00163           for(k = 1; k <= j; k++)
00164           {    pz[j-k] -= pz[j] * y[k];
00165                py[k]   -= pz[j] * z[j-k];
00166           }    
00167           py[0] -= pz[j] * z[j];
00168      }
00169 }
00170 
00171 // --------------------------- Divpv -----------------------------------------
00172 /*!
00173 Compute forward mode Taylor coefficients for result of op = DivpvOp.
00174 
00175 The C++ source code corresponding to this operation is
00176 \verbatim
00177      z = x / y
00178 \endverbatim
00179 In the documentation below,
00180 this operations is for the case where x is a parameter and y is a variable.
00181 
00182 \copydetails forward_binary_op
00183 */
00184 
00185 template <class Base>
00186 inline void forward_divpv_op(
00187      size_t        p           , 
00188      size_t        q           , 
00189      size_t        i_z         ,
00190      const addr_t* arg         ,
00191      const Base*   parameter   ,
00192      size_t        nc_taylor   ,
00193      Base*         taylor      )
00194 {
00195      // check assumptions
00196      CPPAD_ASSERT_UNKNOWN( NumArg(DivpvOp) == 2 );
00197      CPPAD_ASSERT_UNKNOWN( NumRes(DivpvOp) == 1 );
00198      CPPAD_ASSERT_UNKNOWN( size_t(arg[1]) < i_z );
00199      CPPAD_ASSERT_UNKNOWN( q < nc_taylor );
00200      CPPAD_ASSERT_UNKNOWN( p <= q );
00201 
00202      // Taylor coefficients corresponding to arguments and result
00203      Base* y = taylor + arg[1] * nc_taylor;
00204      Base* z = taylor + i_z    * nc_taylor;
00205 
00206      // Paraemter value
00207      Base x = parameter[ arg[0] ];
00208 
00209      // Using CondExp, it can make sense to divide by zero,
00210      // so do not make it an error.
00211      size_t k;
00212      if( p == 0 )
00213      {    z[0] = x / y[0];
00214           p++;
00215      }
00216      for(size_t d = p; d <= q; d++)
00217      {    z[d] = Base(0);
00218           for(k = 1; k <= d; k++)
00219                z[d] -= z[d-k] * y[k];
00220           z[d] /= y[0];
00221      }
00222 }
00223 
00224 /*!
00225 Compute zero order forward mode Taylor coefficient for result of op = DivpvOp.
00226 
00227 The C++ source code corresponding to this operation is
00228 \verbatim
00229      z = x / y
00230 \endverbatim
00231 In the documentation below,
00232 this operations is for the case where x is a parameter and y is a variable.
00233 
00234 \copydetails forward_binary_op_0
00235 */
00236 
00237 template <class Base>
00238 inline void forward_divpv_op_0(
00239      size_t        i_z         ,
00240      const addr_t* arg         ,
00241      const Base*   parameter   ,
00242      size_t        nc_taylor   ,
00243      Base*         taylor      )
00244 {
00245      // check assumptions
00246      CPPAD_ASSERT_UNKNOWN( NumArg(DivpvOp) == 2 );
00247      CPPAD_ASSERT_UNKNOWN( NumRes(DivpvOp) == 1 );
00248      CPPAD_ASSERT_UNKNOWN( size_t(arg[1]) < i_z );
00249 
00250      // Paraemter value
00251      Base x = parameter[ arg[0] ];
00252 
00253      // Taylor coefficients corresponding to arguments and result
00254      Base* y = taylor + arg[1] * nc_taylor;
00255      Base* z = taylor + i_z    * nc_taylor;
00256 
00257      z[0] = x / y[0];
00258 }
00259 
00260 /*!
00261 Compute reverse mode partial derivative for result of op = DivpvOp.
00262 
00263 The C++ source code corresponding to this operation is
00264 \verbatim
00265      z = x / y
00266 \endverbatim
00267 In the documentation below,
00268 this operations is for the case where x is a parameter and y is a variable.
00269 
00270 \copydetails reverse_binary_op
00271 */
00272 
00273 template <class Base>
00274 inline void reverse_divpv_op(
00275      size_t        d           , 
00276      size_t        i_z         ,
00277      const addr_t* arg         ,
00278      const Base*   parameter   ,
00279      size_t        nc_taylor   ,
00280      const Base*   taylor      ,
00281      size_t        nc_partial  ,
00282      Base*         partial     )
00283 {
00284      // check assumptions
00285      CPPAD_ASSERT_UNKNOWN( NumArg(DivvvOp) == 2 );
00286      CPPAD_ASSERT_UNKNOWN( NumRes(DivvvOp) == 1 );
00287      CPPAD_ASSERT_UNKNOWN( size_t(arg[1]) < i_z );
00288      CPPAD_ASSERT_UNKNOWN( d < nc_taylor );
00289      CPPAD_ASSERT_UNKNOWN( d < nc_partial );
00290 
00291      // Arguments
00292      const Base* y = taylor + arg[1] * nc_taylor;
00293      const Base* z = taylor + i_z    * nc_taylor;
00294 
00295      // Partial derivatives corresponding to arguments and result
00296      Base* py = partial + arg[1] * nc_partial;
00297      Base* pz = partial + i_z    * nc_partial;
00298 
00299      // Using CondExp, it can make sense to divide by zero so do not
00300      // make it an error.
00301 
00302      size_t k;
00303      // number of indices to access
00304      size_t j = d + 1;
00305      while(j)
00306      {    --j;
00307           // scale partial w.r.t z[j]
00308           pz[j] /= y[0];
00309 
00310           for(k = 1; k <= j; k++)
00311           {    pz[j-k] -= pz[j] * y[k];
00312                py[k]   -= pz[j] * z[j-k];
00313           }    
00314           py[0] -= pz[j] * z[j];
00315      }
00316 }
00317 
00318 
00319 // --------------------------- Divvp -----------------------------------------
00320 /*!
00321 Compute forward mode Taylor coefficients for result of op = DivvvOp.
00322 
00323 The C++ source code corresponding to this operation is
00324 \verbatim
00325      z = x / y
00326 \endverbatim
00327 In the documentation below,
00328 this operations is for the case where x is a variable and y is a parameter.
00329 
00330 \copydetails forward_binary_op
00331 */
00332 
00333 template <class Base>
00334 inline void forward_divvp_op(
00335      size_t        p           , 
00336      size_t        q           , 
00337      size_t        i_z         ,
00338      const addr_t* arg         ,
00339      const Base*   parameter   ,
00340      size_t        nc_taylor   ,
00341      Base*         taylor      )
00342 {
00343      // check assumptions
00344      CPPAD_ASSERT_UNKNOWN( NumArg(DivvpOp) == 2 );
00345      CPPAD_ASSERT_UNKNOWN( NumRes(DivvpOp) == 1 );
00346      CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < i_z );
00347      CPPAD_ASSERT_UNKNOWN( q < nc_taylor );
00348      CPPAD_ASSERT_UNKNOWN( p <= q );
00349 
00350      // Taylor coefficients corresponding to arguments and result
00351      Base* x = taylor + arg[0] * nc_taylor;
00352      Base* z = taylor + i_z    * nc_taylor;
00353 
00354      // Parameter value
00355      Base y = parameter[ arg[1] ];
00356 
00357      // Using CondExp and multiple levels of AD, it can make sense 
00358      // to divide by zero so do not make it an error.
00359      for(size_t d = p; d <= q; d++)
00360           z[d] = x[d] / y;
00361 }
00362 
00363 
00364 /*!
00365 Compute zero order forward mode Taylor coefficients for result of op = DivvvOp.
00366 
00367 The C++ source code corresponding to this operation is
00368 \verbatim
00369      z = x / y
00370 \endverbatim
00371 In the documentation below,
00372 this operations is for the case where x is a variable and y is a parameter.
00373 
00374 \copydetails forward_binary_op_0
00375 */
00376 
00377 template <class Base>
00378 inline void forward_divvp_op_0(
00379      size_t        i_z         ,
00380      const addr_t* arg         ,
00381      const Base*   parameter   ,
00382      size_t        nc_taylor   ,
00383      Base*         taylor      )
00384 {
00385      // check assumptions
00386      CPPAD_ASSERT_UNKNOWN( NumArg(DivvpOp) == 2 );
00387      CPPAD_ASSERT_UNKNOWN( NumRes(DivvpOp) == 1 );
00388      CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < i_z );
00389 
00390      // Parameter value
00391      Base y = parameter[ arg[1] ];
00392 
00393      // Taylor coefficients corresponding to arguments and result
00394      Base* x = taylor + arg[0] * nc_taylor;
00395      Base* z = taylor + i_z    * nc_taylor;
00396 
00397      z[0] = x[0] / y;
00398 }
00399 
00400 /*!
00401 Compute reverse mode partial derivative for result of op = DivvpOp.
00402 
00403 The C++ source code corresponding to this operation is
00404 \verbatim
00405      z = x / y
00406 \endverbatim
00407 In the documentation below,
00408 this operations is for the case where x is a variable and y is a parameter.
00409 
00410 \copydetails reverse_binary_op
00411 */
00412 
00413 template <class Base>
00414 inline void reverse_divvp_op(
00415      size_t        d           , 
00416      size_t        i_z         ,
00417      const addr_t* arg         ,
00418      const Base*   parameter   ,
00419      size_t        nc_taylor   ,
00420      const Base*   taylor      ,
00421      size_t        nc_partial  ,
00422      Base*         partial     )
00423 {
00424      // check assumptions
00425      CPPAD_ASSERT_UNKNOWN( NumArg(DivvpOp) == 2 );
00426      CPPAD_ASSERT_UNKNOWN( NumRes(DivvpOp) == 1 );
00427      CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < i_z );
00428      CPPAD_ASSERT_UNKNOWN( d < nc_taylor );
00429      CPPAD_ASSERT_UNKNOWN( d < nc_partial );
00430 
00431      // Argument values
00432      Base  y = parameter[ arg[1] ];
00433 
00434      // Partial derivatives corresponding to arguments and result
00435      Base* px = partial + arg[0] * nc_partial;
00436      Base* pz = partial + i_z    * nc_partial;
00437 
00438      // Using CondExp, it can make sense to divide by zero
00439      // so do not make it an error.
00440 
00441      // number of indices to access
00442      size_t j = d + 1;
00443      while(j)
00444      {    --j;
00445           px[j] += pz[j] / y;
00446      }
00447 }
00448 
00449 } // END_CPPAD_NAMESPACE
00450 # endif
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines