CppAD: A C++ Algorithmic Differentiation Package  20130102
ad_fun.hpp
Go to the documentation of this file.
00001 /* $Id$ */
00002 # ifndef CPPAD_AD_FUN_INCLUDED
00003 # define CPPAD_AD_FUN_INCLUDED
00004 
00005 /* --------------------------------------------------------------------------
00006 CppAD: C++ Algorithmic Differentiation: Copyright (C) 2003-12 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 $begin ADFun$$
00017 $spell 
00018      xk
00019      Ind
00020      bool
00021      taylor_
00022      sizeof
00023      const
00024      std
00025      ind_taddr_
00026      dep_taddr_
00027 $$
00028 
00029 $spell
00030 $$
00031 
00032 $section ADFun Objects$$
00033 
00034 $index ADFun, object$$
00035 $index object, ADFun$$
00036 
00037 $head Purpose$$
00038 An AD of $icode Base$$
00039 $cref/operation sequence/glossary/Operation/Sequence/$$
00040 is stored in an $code ADFun$$ object by its $cref FunConstruct$$.
00041 The $code ADFun$$ object can then be used to calculate function values,
00042 derivative values, and other values related to the corresponding function.
00043 
00044 $childtable%
00045      cppad/local/independent.hpp%
00046      cppad/local/fun_construct.hpp%
00047      cppad/local/dependent.hpp%
00048      cppad/local/abort_recording.hpp%
00049      omh/seq_property.omh%
00050      cppad/local/fun_eval.hpp%
00051      cppad/local/drivers.hpp%
00052      cppad/local/fun_check.hpp%
00053      cppad/local/optimize.hpp
00054 %$$
00055 
00056 $end
00057 */
00058 
00059 CPPAD_BEGIN_NAMESPACE
00060 /*!
00061 \defgroup ad_fun_hpp ad_fun.hpp
00062 \{
00063 \file ad_fun.hpp
00064 File used to define the ADFun<Base> class.
00065 */
00066 
00067 /*!
00068 Class used to hold function objects
00069 
00070 \tparam Base
00071 A function object has a recording of <tt>AD<Base></tt> operations.
00072 It does it calculations using \c Base operations.
00073 */
00074 
00075 template <class Base>
00076 class ADFun {
00077 // ------------------------------------------------------------
00078 // Private member variables
00079 private:
00080      /// debug checking number of comparision operations that changed
00081      size_t compare_change_;
00082 
00083      /// number of taylor_ coefficieint per variable (currently stored)
00084      size_t taylor_per_var_;
00085 
00086      /// number of columns currently allocated for taylor_ array
00087      size_t taylor_col_dim_;
00088 
00089      /// number of rows (variables) in the recording (play_)
00090      size_t total_num_var_;
00091 
00092      /// tape address for the independent variables
00093      CppAD::vector<size_t> ind_taddr_;
00094 
00095      /// tape address and parameter flag for the dependent variables
00096      CppAD::vector<size_t> dep_taddr_;
00097 
00098      /// which dependent variables are actually parameters
00099      CppAD::vector<bool>   dep_parameter_;
00100 
00101      /// the operation sequence corresponding to this object
00102      player<Base> play_;
00103 
00104      /// results of the forward mode calculations
00105      pod_vector<Base> taylor_;
00106 
00107      /// Packed results of the forward mode Jacobian sparsity calculations.
00108      /// for_jac_sparse_pack_.n_set() != 0  implies other sparsity results
00109      /// are empty
00110      sparse_pack      for_jac_sparse_pack_;
00111 
00112      /// Set results of the forward mode Jacobian sparsity calculations
00113      /// for_jac_sparse_set_.n_set() != 0  implies for_sparse_pack_ is empty.
00114      CPPAD_INTERNAL_SPARSE_SET  for_jac_sparse_set_;
00115 
00116 // ------------------------------------------------------------
00117 // Private member functions
00118 
00119      /// change the operation sequence corresponding to this object
00120      template <typename ADvector>
00121      void Dependent(ADTape<Base> *tape, const ADvector &y);
00122 
00123      // ------------------------------------------------------------
00124      // vector of bool version of ForSparseJac
00125      // (see doxygen in for_sparse_jac.hpp)
00126      template <class VectorSet>
00127      void ForSparseJacCase(
00128           bool               set_type  ,
00129           size_t             q         ,
00130           const VectorSet&   r         ,  
00131           VectorSet&         s
00132      );
00133      // vector of std::set<size_t> version of ForSparseJac
00134      // (see doxygen in for_sparse_jac.hpp)
00135      template <class VectorSet>
00136      void ForSparseJacCase(
00137           const std::set<size_t>&  set_type  ,
00138           size_t                   q         ,
00139           const VectorSet&         r         ,  
00140           VectorSet&               s
00141      );
00142      // ------------------------------------------------------------
00143      // vector of bool version of RevSparseJac
00144      // (see doxygen in rev_sparse_jac.hpp)
00145      template <class VectorSet>
00146      void RevSparseJacCase(
00147           bool               set_type  ,
00148           size_t             p         ,
00149           const VectorSet&   s         ,  
00150           VectorSet&         r
00151      );
00152      // vector of std::set<size_t> version of RevSparseJac
00153      // (see doxygen in rev_sparse_jac.hpp)
00154      template <class VectorSet>
00155      void RevSparseJacCase(
00156           const std::set<size_t>&  set_type  ,
00157           size_t                   p         ,
00158           const VectorSet&         s         ,  
00159           VectorSet&               r
00160      );
00161      // ------------------------------------------------------------
00162      // vector of bool version of RevSparseHes
00163      // (see doxygen in rev_sparse_hes.hpp)
00164      template <class VectorSet>
00165      void RevSparseHesCase(
00166           bool               set_type  ,
00167           size_t             q         ,
00168           const VectorSet&   s         ,  
00169           VectorSet&         h
00170      );
00171      // vector of std::set<size_t> version of RevSparseHes
00172      // (see doxygen in rev_sparse_hes.hpp)
00173      template <class VectorSet>
00174      void RevSparseHesCase(
00175           const std::set<size_t>&  set_type  ,
00176           size_t                   q         ,
00177           const VectorSet&         s         ,  
00178           VectorSet&               h
00179      );
00180      // ------------------------------------------------------------
00181      // Forward mode version of SparseJacobian
00182      // (see doxygen in sparse_jacobian.hpp)
00183      template <class VectorBase, class VectorSet>
00184      size_t SparseJacobianFor(
00185           const VectorBase&     x               ,
00186           VectorSet&            p_transpose     ,
00187           VectorBase&           jac             ,
00188           sparse_jacobian_work& work
00189      );
00190      // Reverse mode version of SparseJacobian
00191      // (see doxygen in sparse_jacobian.hpp)
00192      template <class VectorBase, class VectorSet>
00193      size_t SparseJacobianRev(
00194           const VectorBase&     x               ,
00195           VectorSet&            p               ,
00196           VectorBase&           jac             ,
00197           sparse_jacobian_work& work
00198      );
00199      // ------------------------------------------------------------
00200      // vector of bool version of SparseJacobian
00201      // (see doxygen in sparse_jacobian.hpp)
00202      template <class VectorBase, class VectorSet>
00203      size_t SparseJacobianCase(
00204           bool                     set_type    ,
00205           const VectorBase&        x           ,
00206           const VectorSet&         p           ,
00207           VectorBase&              jac         ,
00208           sparse_jacobian_work&    work
00209      );
00210      // vector of std::set<size_t> version of SparseJacobian
00211      // (see doxygen in sparse_jacobian.hpp)
00212      template <class VectorBase, class VectorSet>
00213      size_t SparseJacobianCase(
00214           const std::set<size_t>&  set_type    ,
00215           const VectorBase&        x           ,
00216           const VectorSet&         p           ,
00217           VectorBase&              jac         ,
00218           sparse_jacobian_work&    work
00219      );
00220      // vector of bool version of SparseJacobian
00221      // (see doxygen in sparse_jacobian.hpp)
00222      template <class VectorBase, class VectorSet>
00223      void SparseJacobianCase(
00224           bool                     set_type    ,
00225           const VectorBase&        x           ,
00226           const VectorSet&         p           ,
00227           VectorBase&              jac
00228      );
00229      // vector of std::set<size_t> version of SparseJacobian
00230      // (see doxygen in sparse_jacobian.hpp)
00231      template <class VectorBase, class VectorSet>
00232      void SparseJacobianCase(
00233           const std::set<size_t>&  set_type    ,
00234           const VectorBase&        x           ,
00235           const VectorSet&         p           ,
00236           VectorBase&              jac
00237      );
00238      // ------------------------------------------------------------
00239      // combined sparse_set, sparse_list and sparse_pack version of 
00240      // SparseHessian (see doxygen in sparse_hessian.hpp)
00241      template <class VectorBase, class VectorSet>
00242      size_t SparseHessianCompute(
00243           const VectorBase&        x           ,
00244           const VectorBase&        w           ,
00245           VectorSet&               sparsity    ,
00246           VectorBase&              hes         ,
00247           sparse_hessian_work&     work
00248      );
00249      // vector of bool version of SparseHessian
00250      // (see doxygen in sparse_hessian.hpp)
00251      template <class VectorBase, class VectorSet>
00252      size_t SparseHessianCase(
00253           bool                     set_type    ,
00254           const VectorBase&        x           ,
00255           const VectorBase&        w           ,
00256           const VectorSet&         p           ,
00257           VectorBase&              hes         ,
00258           sparse_hessian_work&     work
00259      );
00260      // vector of std::set<size_t> version of SparseHessian
00261      // (see doxygen in sparse_hessian.hpp)
00262      template <class VectorBase, class VectorSet>
00263      size_t SparseHessianCase(
00264           const std::set<size_t>&  set_type    ,
00265           const VectorBase&        x           ,
00266           const VectorBase&        w           ,
00267           const VectorSet&         p           ,
00268           VectorBase&              hes         ,
00269           sparse_hessian_work&     work
00270      );
00271      // vector of bool version of SparseHessian
00272      // (see doxygen in sparse_hessian.hpp)
00273      template <class VectorBase, class VectorSet>
00274      void SparseHessianCase(
00275           bool                     set_type    ,
00276           const VectorBase&        x           ,
00277           const VectorBase&        w           ,
00278           const VectorSet&         p           ,
00279           VectorBase&              hes
00280      );
00281      // vector of std::set<size_t> version of SparseHessian
00282      // (see doxygen in sparse_hessian.hpp)
00283      template <class VectorBase, class VectorSet>
00284      void SparseHessianCase(
00285           const std::set<size_t>&  set_type    ,
00286           const VectorBase&        x           ,
00287           const VectorBase&        w           ,
00288           const VectorSet&         p           ,
00289           VectorBase&              hes
00290      );
00291 // ------------------------------------------------------------
00292 public:
00293      /// copy constructor
00294      ADFun(const ADFun& g) 
00295      : total_num_var_(0)
00296      {    CppAD::ErrorHandler::Call(
00297           true,
00298           __LINE__,
00299           __FILE__,
00300           "ADFun(const ADFun& g)",
00301           "Attempting to use the ADFun<Base> copy constructor.\n"
00302           "Perhaps you are passing an ADFun<Base> object "
00303           "by value instead of by reference."
00304           );
00305       }
00306 
00307      /// default constructor
00308      ADFun(void); 
00309 
00310      // assignment operator
00311      // (see doxygen in fun_construct.hpp)
00312      void operator=(const ADFun& f);
00313 
00314      /// sequence constructor
00315      template <typename ADvector>
00316      ADFun(const ADvector &x, const ADvector &y);
00317 
00318      /// destructor
00319      ~ADFun(void)
00320      { }
00321 
00322      /// deprecated: assign a new operation sequence
00323      template <typename ADvector>
00324      void Dependent(const ADvector &y);
00325 
00326      /// assign a new operation sequence
00327      template <typename ADvector>
00328      void Dependent(const ADvector &x, const ADvector &y);
00329 
00330      /// forward mode sweep
00331      template <typename VectorBase>
00332      VectorBase Forward(
00333           size_t p, const VectorBase& x, std::ostream& s = std::cout);
00334 
00335      /// reverse mode sweep
00336      template <typename VectorBase>
00337      VectorBase Reverse(size_t p, const VectorBase &v);
00338 
00339      // forward mode Jacobian sparsity 
00340      // (see doxygen documentation in for_sparse_jac.hpp)
00341      template <typename VectorSet>
00342      VectorSet ForSparseJac(
00343           size_t q, const VectorSet &r
00344      );
00345      // reverse mode Jacobian sparsity 
00346      // (see doxygen documentation in rev_sparse_jac.hpp)
00347      template <typename VectorSet>
00348      VectorSet RevSparseJac(
00349           size_t q, const VectorSet &s
00350      );
00351      // reverse mode Hessian sparsity 
00352      // (see doxygen documentation in rev_sparse_hes.hpp)
00353      template <typename VectorSet>
00354      VectorSet RevSparseHes(
00355           size_t q, const VectorSet &s
00356      );
00357 
00358      /// amount of memeory used for Jacobain sparsity pattern
00359      size_t size_forward_bool(void) const
00360      {    return for_jac_sparse_pack_.memory(); }
00361 
00362      /// free memeory used for Jacobain sparsity pattern
00363      void size_forward_bool(size_t zero) 
00364      {    CPPAD_ASSERT_KNOWN(
00365                zero == 0,
00366                "size_forward_bool: argument not equal to zero"
00367           );
00368           for_jac_sparse_pack_.resize(0, 0); 
00369      }
00370 
00371      /// total number of elements used for Jacobian sparsity pattern
00372      size_t size_forward_set(void) const
00373      {    return for_jac_sparse_set_.number_elements(); }
00374 
00375      /// free memeory used for Jacobain sparsity pattern
00376      void size_forward_set(size_t zero)
00377      {    CPPAD_ASSERT_KNOWN(
00378                zero == 0,
00379                "size_forward_bool: argument not equal to zero"
00380           );
00381           for_jac_sparse_set_.resize(0, 0); 
00382      }
00383 
00384      /// number of operators in the operation sequence
00385      size_t size_op(void) const
00386      {    return play_.num_rec_op(); }
00387 
00388      /// number of operator arguments in the operation sequence
00389      size_t size_op_arg(void) const
00390      {    return play_.num_rec_op_arg(); }
00391 
00392      /// amount of memory required for the operation sequence
00393      size_t size_op_seq(void) const
00394      {    return play_.Memory(); }
00395 
00396      /// number of parameters in the operation sequence
00397      size_t size_par(void) const
00398      {    return play_.num_rec_par(); }
00399 
00400      /// number of taylor_ coefficients currently calculated (per variable)
00401      size_t size_taylor(void) const
00402      {    return taylor_per_var_; } 
00403 
00404      /// number of characters in the operation sequence
00405      size_t size_text(void) const
00406      {    return play_.num_rec_text(); }
00407 
00408      /// number of variables in opertion sequence
00409      size_t size_var(void) const
00410      {    return total_num_var_; }
00411 
00412      /// number of VecAD indices in the operation sequence
00413      size_t size_VecAD(void) const
00414      {    return play_.num_rec_vecad_ind(); }
00415 
00416      /// set number of coefficients currently allocated (per variable)
00417      void capacity_taylor(size_t per_var);   
00418 
00419      /// number of independent variables
00420      size_t Domain(void) const
00421      {    return ind_taddr_.size(); }
00422 
00423      /// number of dependent variables
00424      size_t Range(void) const
00425      {    return dep_taddr_.size(); }
00426 
00427      /// is variable a parameter
00428      bool Parameter(size_t i)
00429      {    CPPAD_ASSERT_KNOWN(
00430                i < dep_taddr_.size(),
00431                "Argument to Parameter is >= dimension of range space"
00432           );
00433           return dep_parameter_[i]; 
00434      }
00435 
00436 # ifndef NDEBUG
00437      /// in not NDEBUG case, number of comparison operations that change
00438      size_t CompareChange(void) const
00439      {    return compare_change_; }
00440 # endif
00441 
00442      /// calculate entire Jacobian
00443      template <typename VectorBase>
00444      VectorBase Jacobian(const VectorBase &x); 
00445 
00446      /// calculate Hessian for one component of f
00447      template <typename VectorBase>
00448      VectorBase Hessian(const VectorBase &x, const VectorBase &w); 
00449      template <typename VectorBase>
00450      VectorBase Hessian(const VectorBase &x, size_t i); 
00451 
00452      /// forward mode calculation of partial w.r.t one domain component
00453      template <typename VectorBase>
00454      VectorBase ForOne(
00455           const VectorBase   &x ,
00456           size_t              j );
00457 
00458      /// reverse mode calculation of derivative of one range component
00459      template <typename VectorBase>
00460      VectorBase RevOne(
00461           const VectorBase   &x ,
00462           size_t              i );
00463 
00464      /// forward mode calculation of a subset of second order partials
00465      template <typename VectorBase, typename VectorSize_t>
00466      VectorBase ForTwo(
00467           const VectorBase   &x ,
00468           const VectorSize_t &J ,
00469           const VectorSize_t &K );
00470 
00471      /// reverse mode calculation of a subset of second order partials
00472      template <typename VectorBase, typename VectorSize_t>
00473      VectorBase RevTwo(
00474           const VectorBase   &x ,
00475           const VectorSize_t &I ,
00476           const VectorSize_t &J );
00477 
00478      /// calculate sparse Jacobians 
00479      template <typename VectorBase>
00480      VectorBase SparseJacobian(
00481           const VectorBase &x
00482      ); 
00483      template <typename VectorBase, typename VectorSet>
00484      VectorBase SparseJacobian(
00485           const VectorBase &x , 
00486           const VectorSet  &p
00487      ); 
00488      template <class VectorBase, class VectorSet, class VectorSize>
00489      size_t SparseJacobianForward(
00490           const VectorBase&     x     ,
00491           const VectorSet&      p     ,
00492           const VectorSize&     r     ,
00493           const VectorSize&     c     ,
00494           VectorBase&           jac   ,
00495           sparse_jacobian_work& work
00496      );
00497      template <class VectorBase, class VectorSet, class VectorSize>
00498      size_t SparseJacobianReverse(
00499           const VectorBase&     x    ,
00500           const VectorSet&      p    ,
00501           const VectorSize&     r    ,
00502           const VectorSize&     c    ,
00503           VectorBase&           jac  ,
00504           sparse_jacobian_work& work
00505      );
00506 
00507      /// calculate sparse Hessians 
00508      template <typename VectorBase>
00509      VectorBase SparseHessian(
00510           const VectorBase&    x  , 
00511           const VectorBase&    w
00512      ); 
00513      template <typename VectorBase, typename VectorBool>
00514      VectorBase SparseHessian(
00515           const VectorBase&    x  ,
00516           const VectorBase&    w  ,
00517           const VectorBool&    p
00518      ); 
00519      template <class VectorBase, class VectorSet, class VectorSize>
00520      size_t SparseHessian(
00521           const VectorBase&    x   ,
00522           const VectorBase&    w   ,
00523           const VectorSet&     p   ,
00524           const VectorSize&    r   ,
00525           const VectorSize&    c   ,
00526           VectorBase&          hes ,
00527           sparse_hessian_work& work
00528      ); 
00529 
00530      // Optimize the tape
00531      // (see doxygen documentation in optimize.hpp)
00532      void optimize(void);
00533      // ------------------- Deprecated -----------------------------
00534 
00535      /// deprecated: number of variables in opertion sequence
00536      size_t Size(void) const
00537      {    return total_num_var_; }
00538 
00539      /// deprecated: # taylor_ coefficients currently stored (per variable)
00540      size_t Order(void) const
00541      {    return taylor_per_var_ - 1; }
00542 
00543      /// Deprecated: amount of memory for this object 
00544      /// Note that an approximation is used for the std::set<size_t> memory
00545      size_t Memory(void) const
00546      {    size_t pervar  = taylor_col_dim_ * sizeof(Base)
00547           + for_jac_sparse_pack_.memory()
00548           + 3 * sizeof(size_t) * for_jac_sparse_set_.number_elements();
00549           size_t total   = total_num_var_ * pervar + play_.Memory();
00550           return total;
00551      }
00552 
00553      /// deprecated: # taylor_ coefficients stored (per variable)
00554      size_t taylor_size(void) const
00555      {    return taylor_per_var_; } 
00556 
00557      /// deprecated: Does this AD operation sequence use 
00558      //VecAD<Base>::reference operands
00559      bool use_VecAD(void) const
00560      {    return play_.num_rec_vecad_ind() > 0; }
00561 };
00562 // ---------------------------------------------------------------------------
00563 
00564 /*! \} */
00565 CPPAD_END_NAMESPACE
00566 
00567 // non-user interfaces
00568 # include <cppad/local/forward0sweep.hpp>
00569 # include <cppad/local/forward_sweep.hpp>
00570 # include <cppad/local/reverse_sweep.hpp>
00571 # include <cppad/local/for_jac_sweep.hpp>
00572 # include <cppad/local/rev_jac_sweep.hpp>
00573 # include <cppad/local/rev_hes_sweep.hpp>
00574 
00575 // user interfaces
00576 # include <cppad/local/parallel_ad.hpp>
00577 # include <cppad/local/independent.hpp>
00578 # include <cppad/local/dependent.hpp>
00579 # include <cppad/local/fun_construct.hpp>
00580 # include <cppad/local/abort_recording.hpp>
00581 # include <cppad/local/fun_eval.hpp>
00582 # include <cppad/local/drivers.hpp>
00583 # include <cppad/local/fun_check.hpp>
00584 # include <cppad/local/omp_max_thread.hpp>
00585 # include <cppad/local/optimize.hpp> 
00586 
00587 # endif
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines