CppAD: A C++ Algorithmic Differentiation Package  20130918
store_op.hpp
Go to the documentation of this file.
00001 /* $Id$ */
00002 # ifndef CPPAD_STORE_OP_INCLUDED
00003 # define CPPAD_STORE_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 store_op.hpp
00019 Changing the current value of a VecAD element.
00020 */
00021 /*!
00022 Shared documentation for zero order forward implementation of 
00023 op = StppOp, StpvOp, StvpOp, or StvvOp (not called).
00024 
00025 The C++ source code corresponding to this operation is
00026 \verbatim
00027      v[x] = y
00028 \endverbatim
00029 where v is a VecAD<Base> vector, x is an AD<Base> object,
00030 and y is AD<Base> or Base objects. 
00031 We define the index corresponding to v[x] by
00032 \verbatim
00033      i_v_x = index_by_ind[ arg[0] + i_vec ]
00034 \endverbatim
00035 where i_vec is defined under the heading arg[1] below:
00036 
00037 \tparam Base
00038 base type for the operator; i.e., this operation was recorded
00039 using AD<Base> and computations by this routine are done using type Base.
00040 
00041 \param i_z
00042 is the index corresponding to the previous variable on the tape
00043 (only used for error checking).
00044 
00045 \param arg
00046 \n
00047 arg[0]
00048 \n
00049 is the offset of this VecAD vector relative to the beginning 
00050 of the isvar_by_ind and index_by_ind arrays.
00051 \n
00052 \n 
00053 arg[1] 
00054 \n
00055 If this is a StppOp or StpvOp operation (if x is a parameter), 
00056 i_vec is defined by
00057 \verbatim
00058      i_vec = arg[1]
00059 \endverbatim
00060 If this is a StvpOp or StvvOp operation (if x is a variable),
00061 i_vec is defined by
00062 \verbatim
00063      i_vec = floor( taylor[ arg[1] * nc_taylor + 0 ] )
00064 \endverbatim
00065 where floor(c) is the greatest integer less that or equal c.
00066 \n
00067 \n
00068 arg[2]
00069 \n
00070 index corresponding to the third operand for this operator;
00071 i.e. the index corresponding to y.
00072 
00073 \param num_par
00074 is the total number of parameters on the tape
00075 (only used for error checking).
00076 
00077 \param nc_taylor
00078 number of columns in the matrix containing the Taylor coefficients.
00079 
00080 \param taylor
00081 In StvpOp and StvvOp cases, <code><taylor[ arg[1] * nc_taylor + 0 ]</code>
00082 is used to compute the index in the definition of i_vec above.
00083 
00084 \param isvar_by_ind
00085 If y is a varable (StpvOp and StvvOp cases), 
00086 <code>isvar_by_ind[ arg[0] + i_vec ] </code> is set to true.
00087 Otherwise y is a paraemter (StppOp and StvpOp cases) and 
00088 <code>isvar_by_ind[ arg[0] + i_vec ] </code> is set to false.
00089 
00090 \param index_by_ind
00091 <code>index_by_ind[ arg[0] - 1 ]</code>
00092 is the number of elements in the user vector containing this element.
00093 The value <code>index_by_ind[ arg[0] + i_vec]</code>
00094 is set equal to arg[2].
00095 
00096 \par Check User Errors
00097 \li Check that the index is with in range; i.e.
00098 <code>i_vec < index_by_ind[ arg[0] - 1 ]</code>
00099 Note that, if x is a parameter, 
00100 the corresponding vector index and it does not change.
00101 In this case, the error above should be detected during tape recording.
00102 
00103 \par Checked Assertions 
00104 \li NumArg(op) == 3
00105 \li NumRes(op) == 0
00106 \li 0 <  arg[0]
00107 \li if y is a parameter, arg[2] < num_par
00108 \li if x is a variable, arg[1] <= i_z
00109 \li if y is a variable, arg[2] <= i_z
00110 */
00111 template <class Base>
00112 inline void forward_store_op_0(
00113      size_t         i_z         ,
00114      const addr_t*  arg         , 
00115      size_t         num_par     ,
00116      size_t         nc_taylor   ,
00117      Base*          taylor      ,
00118      bool*          isvar_by_ind   ,
00119      size_t*        index_by_ind   )
00120 {
00121      // This routine is only for documentaiton, it should not be used
00122      CPPAD_ASSERT_UNKNOWN( false );
00123 }
00124 /*!
00125 Shared documnetation for sparsity operations corresponding to 
00126 op = StpvOp or StvvOp (not called).
00127 
00128 <!-- define sparse_store_op -->
00129 The C++ source code corresponding to this operation is
00130 \verbatim
00131      v[x] = y
00132 \endverbatim
00133 where v is a VecAD<Base> vector, x is an AD<Base> object,
00134 and y is AD<Base> or Base objects. 
00135 We define the index corresponding to v[x] by
00136 \verbatim
00137      i_v_x = combined[ arg[0] + i_vec ]
00138 \endverbatim
00139 where i_vec is defined under the heading \a arg[1] below:
00140 
00141 \tparam Vector_set
00142 is the type used for vectors of sets. It can be either
00143 \c sparse_pack, \c sparse_set, or \c sparse_list.
00144 
00145 \param op
00146 is the code corresponding to this operator; i.e., StpvOp or StvvOp
00147 (only used for error checking).
00148 
00149 \param arg
00150 \n
00151 \a arg[0]
00152 is the offset corresponding to this VecAD vector in the combined array.
00153 \n
00154 \n 
00155 \a arg[2]
00156 \n
00157 The set with index \a arg[2] in \a var_sparsity 
00158 is the sparsity pattern corresponding to y.
00159 (Note that \a arg[2] > 0 because y is a variable.) 
00160 
00161 \param num_combined
00162 is the total number of elements in the VecAD address array.
00163 
00164 \param combined
00165 \a combined [ arg[0] - 1 ]
00166 is the index of the set in \a vecad_sparsity corresponding
00167 to the sparsity pattern for the vector v.
00168 We use the notation i_v below which is defined by
00169 \verbatim
00170      i_v = combined[ \a arg[0] - 1 ]
00171 \endverbatim
00172 
00173 \param var_sparsity
00174 The set  with index \a arg[2] in \a var_sparsity 
00175 is the sparsity pattern for y.
00176 This is an input for forward mode operations.
00177 For reverse mode operations:
00178 The sparsity pattern for v is added to the spartisy pattern for y.
00179 
00180 \param vecad_sparsity
00181 The set with index \a i_v in \a vecad_sparsity 
00182 is the sparsity pattern for v.
00183 This is an input for reverse mode operations.
00184 For forward mode operations, the sparsity pattern for y is added
00185 to the sparsity pattern for the vector v.
00186 
00187 \par Checked Assertions 
00188 \li NumArg(op) == 3
00189 \li NumRes(op) == 0
00190 \li 0 <  \a arg[0]
00191 \li \a arg[0] < \a num_combined
00192 \li \a arg[2] < \a var_sparsity.n_set()
00193 \li i_v       < \a vecad_sparsity.n_set()
00194 <!-- end sparse_store_op -->
00195 */
00196 template <class Vector_set>
00197 inline void sparse_store_op(
00198      OpCode         op             ,
00199      const addr_t*  arg            , 
00200      size_t         num_combined   ,
00201      const size_t*  combined       ,
00202      Vector_set&    var_sparsity   ,
00203      Vector_set&    vecad_sparsity )
00204 {
00205      // This routine is only for documentaiton, it should not be used
00206      CPPAD_ASSERT_UNKNOWN( false );
00207 }
00208 
00209 
00210 /*!
00211 Zero order forward mode implementation of op = StppOp.
00212 
00213 \copydetails forward_store_op_0
00214 */
00215 template <class Base>
00216 inline void forward_store_pp_op_0(
00217      size_t         i_z         ,
00218      const addr_t*  arg         , 
00219      size_t         num_par     ,
00220      size_t         nc_taylor   ,
00221      Base*          taylor      ,
00222      bool*          isvar_by_ind   ,
00223      size_t*        index_by_ind   )
00224 {    size_t i_vec = arg[1];
00225 
00226      // Because the index is a parameter, this indexing error should be 
00227      // caught and reported to the user when the tape is recording.
00228      CPPAD_ASSERT_UNKNOWN( i_vec < index_by_ind[ arg[0] - 1 ] );
00229 
00230      CPPAD_ASSERT_UNKNOWN( NumArg(StppOp) == 3 );
00231      CPPAD_ASSERT_UNKNOWN( NumRes(StppOp) == 0 );
00232      CPPAD_ASSERT_UNKNOWN( 0 < arg[0] );
00233      CPPAD_ASSERT_UNKNOWN( size_t(arg[2]) < num_par );
00234 
00235      isvar_by_ind[ arg[0] + i_vec ]  = false;
00236      index_by_ind[ arg[0] + i_vec ]  = arg[2];
00237 }
00238 
00239 /*!
00240 Zero order forward mode implementation of op = StpvOp.
00241 
00242 \copydetails forward_store_op_0
00243 */
00244 template <class Base>
00245 inline void forward_store_pv_op_0(
00246      size_t         i_z         ,
00247      const addr_t*  arg         , 
00248      size_t         num_par     ,
00249      size_t         nc_taylor   ,
00250      Base*          taylor      ,
00251      bool*          isvar_by_ind   ,
00252      size_t*        index_by_ind   )
00253 {    size_t i_vec = arg[1];
00254 
00255      // Because the index is a parameter, this indexing error should be 
00256      // caught and reported to the user when the tape is recording.
00257      CPPAD_ASSERT_UNKNOWN( i_vec < index_by_ind[ arg[0] - 1 ] );
00258 
00259      CPPAD_ASSERT_UNKNOWN( NumArg(StpvOp) == 3 );
00260      CPPAD_ASSERT_UNKNOWN( NumRes(StpvOp) == 0 );
00261      CPPAD_ASSERT_UNKNOWN( 0 < arg[0] );
00262      CPPAD_ASSERT_UNKNOWN( size_t(arg[2]) <= i_z );
00263 
00264      isvar_by_ind[ arg[0] + i_vec ]  = true;
00265      index_by_ind[ arg[0] + i_vec ]  = arg[2];
00266 }
00267 
00268 /*!
00269 Zero order forward mode implementation of op = StvpOp.
00270 
00271 \copydetails forward_store_op_0
00272 */
00273 template <class Base>
00274 inline void forward_store_vp_op_0(
00275      size_t         i_z         ,
00276      const addr_t*  arg         , 
00277      size_t         num_par     ,
00278      size_t         nc_taylor   ,
00279      Base*          taylor      ,
00280      bool*          isvar_by_ind   ,
00281      size_t*        index_by_ind   )
00282 {    
00283      CPPAD_ASSERT_UNKNOWN( size_t(arg[1]) <= i_z );
00284      size_t i_vec = Integer( taylor[ arg[1] * nc_taylor + 0 ] );
00285      CPPAD_ASSERT_KNOWN( 
00286           i_vec < index_by_ind[ arg[0] - 1 ] ,
00287           "VecAD: index during zero order forward sweep is out of range"
00288      );
00289 
00290      CPPAD_ASSERT_UNKNOWN( NumArg(StvpOp) == 3 );
00291      CPPAD_ASSERT_UNKNOWN( NumRes(StvpOp) == 0 );
00292      CPPAD_ASSERT_UNKNOWN( 0 < arg[0] );
00293      CPPAD_ASSERT_UNKNOWN( size_t(arg[2]) < num_par );
00294 
00295      isvar_by_ind[ arg[0] + i_vec ]  = false;
00296      index_by_ind[ arg[0] + i_vec ]  = arg[2];
00297 }
00298 
00299 /*!
00300 Zero order forward mode implementation of op = StvvOp.
00301 
00302 \copydetails forward_store_op_0
00303 */
00304 template <class Base>
00305 inline void forward_store_vv_op_0(
00306      size_t         i_z         ,
00307      const addr_t*  arg         , 
00308      size_t         num_par     ,
00309      size_t         nc_taylor   ,
00310      Base*          taylor      ,
00311      bool*          isvar_by_ind   ,
00312      size_t*        index_by_ind   )
00313 {    
00314      CPPAD_ASSERT_UNKNOWN( size_t(arg[1]) <= i_z );
00315      size_t i_vec = Integer( taylor[ arg[1] * nc_taylor + 0 ] );
00316      CPPAD_ASSERT_KNOWN( 
00317           i_vec < index_by_ind[ arg[0] - 1 ] ,
00318           "VecAD: index during zero order forward sweep is out of range"
00319      );
00320 
00321      CPPAD_ASSERT_UNKNOWN( NumArg(StvpOp) == 3 );
00322      CPPAD_ASSERT_UNKNOWN( NumRes(StvpOp) == 0 );
00323      CPPAD_ASSERT_UNKNOWN( 0 < arg[0] );
00324      CPPAD_ASSERT_UNKNOWN( size_t(arg[2]) <= i_z );
00325 
00326      isvar_by_ind[ arg[0] + i_vec ]  = true;
00327      index_by_ind[ arg[0] + i_vec ]  = arg[2];
00328 }
00329 
00330 /*!
00331 Forward mode sparsity operations for StpvOp and StvvOp
00332 
00333 \copydetails sparse_store_op
00334 */
00335 template <class Vector_set>
00336 inline void forward_sparse_store_op(
00337      OpCode              op             ,
00338      const addr_t*       arg            , 
00339      size_t              num_combined   ,
00340      const size_t*       combined       ,
00341      Vector_set&         var_sparsity   ,
00342      Vector_set&         vecad_sparsity )
00343 {
00344      CPPAD_ASSERT_UNKNOWN( NumArg(op) == 3 );
00345      CPPAD_ASSERT_UNKNOWN( NumRes(op) == 0 );
00346      CPPAD_ASSERT_UNKNOWN( 0 < arg[0] );
00347      CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < num_combined );
00348      size_t i_v = combined[ arg[0] - 1 ];
00349      CPPAD_ASSERT_UNKNOWN( i_v < vecad_sparsity.n_set() );
00350      CPPAD_ASSERT_UNKNOWN( size_t(arg[2]) < var_sparsity.n_set() );
00351 
00352      vecad_sparsity.binary_union(i_v, i_v, arg[2], var_sparsity);
00353 
00354      return;
00355 }
00356 
00357 /*!
00358 Reverse mode sparsity operations for StpvOp and StvvOp
00359 
00360 This routine is given the sparsity patterns for
00361 G(v[x], y , w , u ... ) and it uses them to compute the 
00362 sparsity patterns for  
00363 \verbatim
00364      H(y , w , u , ... ) = G[ v[x], y , w , u , ... ]
00365 \endverbatim
00366 
00367 <!-- replace sparse_store_op -->
00368 The C++ source code corresponding to this operation is
00369 \verbatim
00370      v[x] = y
00371 \endverbatim
00372 where v is a VecAD<Base> vector, x is an AD<Base> object,
00373 and y is AD<Base> or Base objects. 
00374 We define the index corresponding to v[x] by
00375 \verbatim
00376      i_v_x = combined[ arg[0] + i_vec ]
00377 \endverbatim
00378 where i_vec is defined under the heading \a arg[1] below:
00379 
00380 \tparam Vector_set
00381 is the type used for vectors of sets. It can be either
00382 \c sparse_pack, \c sparse_set, or \c sparse_list.
00383 
00384 \param op
00385 is the code corresponding to this operator; i.e., StpvOp or StvvOp
00386 (only used for error checking).
00387 
00388 \param arg
00389 \n
00390 \a arg[0]
00391 is the offset corresponding to this VecAD vector in the combined array.
00392 \n
00393 \n 
00394 \a arg[2]
00395 \n
00396 The set with index \a arg[2] in \a var_sparsity 
00397 is the sparsity pattern corresponding to y.
00398 (Note that \a arg[2] > 0 because y is a variable.) 
00399 
00400 \param num_combined
00401 is the total number of elements in the VecAD address array.
00402 
00403 \param combined
00404 \a combined [ arg[0] - 1 ]
00405 is the index of the set in \a vecad_sparsity corresponding
00406 to the sparsity pattern for the vector v.
00407 We use the notation i_v below which is defined by
00408 \verbatim
00409      i_v = combined[ \a arg[0] - 1 ]
00410 \endverbatim
00411 
00412 \param var_sparsity
00413 The set  with index \a arg[2] in \a var_sparsity 
00414 is the sparsity pattern for y.
00415 This is an input for forward mode operations.
00416 For reverse mode operations:
00417 The sparsity pattern for v is added to the spartisy pattern for y.
00418 
00419 \param vecad_sparsity
00420 The set with index \a i_v in \a vecad_sparsity 
00421 is the sparsity pattern for v.
00422 This is an input for reverse mode operations.
00423 For forward mode operations, the sparsity pattern for y is added
00424 to the sparsity pattern for the vector v.
00425 
00426 \par Checked Assertions 
00427 \li NumArg(op) == 3
00428 \li NumRes(op) == 0
00429 \li 0 <  \a arg[0]
00430 \li \a arg[0] < \a num_combined
00431 \li \a arg[2] < \a var_sparsity.n_set()
00432 \li i_v       < \a vecad_sparsity.n_set()
00433 <!-- end sparse_store_op -->
00434 */
00435 template <class Vector_set>
00436 inline void reverse_sparse_jacobian_store_op(
00437      OpCode             op              ,
00438      const addr_t*      arg             , 
00439      size_t             num_combined    ,
00440      const size_t*      combined        ,
00441      Vector_set&        var_sparsity    ,
00442      Vector_set&        vecad_sparsity  )
00443 {
00444      CPPAD_ASSERT_UNKNOWN( NumArg(op) == 3 );
00445      CPPAD_ASSERT_UNKNOWN( NumRes(op) == 0 );
00446      CPPAD_ASSERT_UNKNOWN( 0 < arg[0] );
00447      CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < num_combined );
00448      size_t i_v = combined[ arg[0] - 1 ];
00449      CPPAD_ASSERT_UNKNOWN( i_v < vecad_sparsity.n_set() );
00450      CPPAD_ASSERT_UNKNOWN( size_t(arg[2]) < var_sparsity.n_set() );
00451 
00452      var_sparsity.binary_union(arg[2], arg[2], i_v, vecad_sparsity);
00453 
00454      return;
00455 }
00456 
00457 /*!
00458 Reverse mode sparsity operations for StpvOp and StvvOp
00459 
00460 This routine is given the sparsity patterns for
00461 G(v[x], y , w , u ... )
00462 and it uses them to compute the sparsity patterns for 
00463 \verbatim
00464      H(y , w , u , ... ) = G[ v[x], y , w , u , ... ]
00465 \endverbatim
00466 
00467 <!-- replace sparse_store_op -->
00468 The C++ source code corresponding to this operation is
00469 \verbatim
00470      v[x] = y
00471 \endverbatim
00472 where v is a VecAD<Base> vector, x is an AD<Base> object,
00473 and y is AD<Base> or Base objects. 
00474 We define the index corresponding to v[x] by
00475 \verbatim
00476      i_v_x = combined[ arg[0] + i_vec ]
00477 \endverbatim
00478 where i_vec is defined under the heading \a arg[1] below:
00479 
00480 \tparam Vector_set
00481 is the type used for vectors of sets. It can be either
00482 \c sparse_pack, \c sparse_set, or \c sparse_list.
00483 
00484 \param op
00485 is the code corresponding to this operator; i.e., StpvOp or StvvOp
00486 (only used for error checking).
00487 
00488 \param arg
00489 \n
00490 \a arg[0]
00491 is the offset corresponding to this VecAD vector in the combined array.
00492 \n
00493 \n 
00494 \a arg[2]
00495 \n
00496 The set with index \a arg[2] in \a var_sparsity 
00497 is the sparsity pattern corresponding to y.
00498 (Note that \a arg[2] > 0 because y is a variable.) 
00499 
00500 \param num_combined
00501 is the total number of elements in the VecAD address array.
00502 
00503 \param combined
00504 \a combined [ arg[0] - 1 ]
00505 is the index of the set in \a vecad_sparsity corresponding
00506 to the sparsity pattern for the vector v.
00507 We use the notation i_v below which is defined by
00508 \verbatim
00509      i_v = combined[ \a arg[0] - 1 ]
00510 \endverbatim
00511 
00512 \param var_sparsity
00513 The set  with index \a arg[2] in \a var_sparsity 
00514 is the sparsity pattern for y.
00515 This is an input for forward mode operations.
00516 For reverse mode operations:
00517 The sparsity pattern for v is added to the spartisy pattern for y.
00518 
00519 \param vecad_sparsity
00520 The set with index \a i_v in \a vecad_sparsity 
00521 is the sparsity pattern for v.
00522 This is an input for reverse mode operations.
00523 For forward mode operations, the sparsity pattern for y is added
00524 to the sparsity pattern for the vector v.
00525 
00526 \par Checked Assertions 
00527 \li NumArg(op) == 3
00528 \li NumRes(op) == 0
00529 \li 0 <  \a arg[0]
00530 \li \a arg[0] < \a num_combined
00531 \li \a arg[2] < \a var_sparsity.n_set()
00532 \li i_v       < \a vecad_sparsity.n_set()
00533 <!-- end sparse_store_op -->
00534 
00535 \param var_jacobian
00536 \a var_jacobian[ \a arg[2] ] 
00537 is false (true) if the Jacobian of G with respect to y is always zero 
00538 (may be non-zero).
00539 
00540 \param vecad_jacobian
00541 \a vecad_jacobian[i_v] 
00542 is false (true) if the Jacobian with respect to x is always zero 
00543 (may be non-zero).
00544 On input, it corresponds to the function G,
00545 and on output it corresponds to the function H.
00546 */
00547 template <class Vector_set>
00548 inline void reverse_sparse_hessian_store_op(
00549      OpCode             op           ,
00550      const addr_t*      arg          , 
00551      size_t             num_combined ,
00552      const size_t*      combined     ,
00553      Vector_set&        var_sparsity ,
00554      Vector_set&        vecad_sparsity ,
00555      bool*              var_jacobian   ,
00556      bool*              vecad_jacobian )
00557 {
00558      CPPAD_ASSERT_UNKNOWN( NumArg(op) == 3 );
00559      CPPAD_ASSERT_UNKNOWN( NumRes(op) == 0 );
00560      CPPAD_ASSERT_UNKNOWN( 0 < arg[0] );
00561      CPPAD_ASSERT_UNKNOWN( size_t(arg[0]) < num_combined );
00562      size_t i_v = combined[ arg[0] - 1 ];
00563      CPPAD_ASSERT_UNKNOWN( i_v < vecad_sparsity.n_set() );
00564      CPPAD_ASSERT_UNKNOWN( size_t(arg[2]) < var_sparsity.n_set() );
00565 
00566      var_sparsity.binary_union(arg[2], arg[2], i_v, vecad_sparsity);
00567 
00568      var_jacobian[ arg[2] ] |= vecad_jacobian[i_v];
00569 
00570      return;
00571 }
00572 
00573 } // END_CPPAD_NAMESPACE
00574 # endif
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines