CppAD: A C++ Algorithmic Differentiation Package  20130918
rev_sparse_hes.hpp
Go to the documentation of this file.
00001 /* $Id$ */
00002 # ifndef CPPAD_REV_SPARSE_HES_INCLUDED
00003 # define CPPAD_REV_SPARSE_HES_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 /*
00017 $begin RevSparseHes$$
00018 $spell
00019      std
00020      VecAD
00021      Jacobian
00022      Jac
00023      Hessian
00024      Hes
00025      const
00026      Bool
00027      Dep
00028      proportional
00029      var
00030 $$
00031 
00032 $section Hessian Sparsity Pattern: Reverse Mode$$ 
00033 
00034 $index RevSparseHes$$
00035 $index reverse, sparse Hessian$$
00036 $index sparse, reverse Hessian$$
00037 $index pattern, reverse Hessian$$
00038 
00039 $head Syntax$$
00040 $icode%h% = %f%.RevSparseHes(%q%, %s%)
00041 %$$
00042 $icode%h% = %f%.RevSparseHes(%q%, %s%, %transpose%)%$$
00043 
00044 
00045 $head Purpose$$
00046 We use $latex F : \B{R}^n \rightarrow \B{R}^m$$ to denote the
00047 $cref/AD function/glossary/AD Function/$$ corresponding to $icode f$$.
00048 For a fixed matrix $latex R \in \B{R}^{n \times q}$$ 
00049 and a fixed vector $latex S \in \B{R}^{1 \times m}$$,
00050 we define 
00051 $latex \[
00052 \begin{array}{rcl}
00053 H(x) 
00054 & = & \partial_x \left[ \partial_u S * F[ x + R * u ] \right]_{u=0}
00055 \\
00056 & = & R^\R{T} * (S * F)^{(2)} ( x )
00057 \\
00058 H(x)^\R{T} 
00059 & = & (S * F)^{(2)} ( x ) * R
00060 \end{array}
00061 \] $$
00062 Given a
00063 $cref/sparsity pattern/glossary/Sparsity Pattern/$$ 
00064 for the matrix $latex R$$ and the vector $latex S$$,
00065 $code RevSparseHes$$ returns a sparsity pattern for the $latex H(x)$$.
00066 
00067 $head f$$
00068 The object $icode f$$ has prototype
00069 $codei%
00070      const ADFun<%Base%> %f%
00071 %$$
00072 
00073 $head x$$
00074 the sparsity pattern is valid for all values of the independent
00075 variables in $latex x \in \B{R}^n$$
00076 (even if it has $cref CondExp$$ or $cref VecAD$$ operations).
00077 
00078 $head q$$
00079 The argument $icode q$$ has prototype
00080 $codei%
00081      size_t %q%
00082 %$$
00083 It specifies the number of columns in $latex R \in \B{R}^{n \times q}$$
00084 and the number of rows in $latex H(x) \in \B{R}^{q \times n}$$.
00085 It must be the same value as in the previous $cref ForSparseJac$$ call 
00086 $codei%
00087      %f%.ForSparseJac(%q%, %r%, %r_transpose%)
00088 %$$
00089 Note that if $icode r_transpose$$ is true, $icode r$$ in the call above
00090 corresponding to $latex R^\R{T} \in \B{R}^{q \times n}$$
00091 
00092 $head transpose$$
00093 The argument $icode transpose$$ has prototype
00094 $codei%
00095      bool %transpose%
00096 %$$
00097 The default value $code false$$ is used when $icode transpose$$ is not present.
00098 
00099 
00100 $head r$$
00101 The matrix $latex R$$ is specified by the previous call
00102 $codei%
00103      %f%.ForSparseJac(%q%, %r%, %transpose%)
00104 %$$
00105 see $cref/r/ForSparseJac/r/$$.
00106 The type of the elements of
00107 $cref/VectorSet/RevSparseHes/VectorSet/$$ must be the 
00108 same as the type of the elements of $icode r$$.
00109 
00110 $head s$$
00111 The argument $icode s$$ has prototype
00112 $codei%
00113      const %VectorSet%& %s%
00114 %$$
00115 (see $cref/VectorSet/RevSparseHes/VectorSet/$$ below)
00116 If it has elements of type $code bool$$,
00117 its size is $latex m$$.
00118 If it has elements of type $code std::set<size_t>$$,
00119 its size is one and all the elements of $icode%s%[0]%$$
00120 are between zero and $latex m - 1$$.
00121 It specifies a 
00122 $cref/sparsity pattern/glossary/Sparsity Pattern/$$ 
00123 for the vector $icode S$$.
00124 
00125 $head h$$
00126 The result $icode h$$ has prototype
00127 $codei%
00128      %VectorSet%& %h%
00129 %$$
00130 (see $cref/VectorSet/RevSparseHes/VectorSet/$$ below).
00131 
00132 $subhead transpose false$$
00133 If $icode h$$ has elements of type $code bool$$,
00134 its size is $latex q * n$$.
00135 If it has elements of type $code std::set<size_t>$$,
00136 its size is $latex q$$ and all the set elements are between
00137 zero and $icode%n%-1%$$ inclusive.
00138 It specifies a 
00139 $cref/sparsity pattern/glossary/Sparsity Pattern/$$ 
00140 for the matrix $latex H(x)$$.
00141 
00142 $subhead transpose true$$
00143 If $icode h$$ has elements of type $code bool$$,
00144 its size is $latex n * q$$.
00145 If it has elements of type $code std::set<size_t>$$,
00146 its size is $latex n$$ and all the set elements are between
00147 zero and $icode%q%-1%$$ inclusive.
00148 It specifies a 
00149 $cref/sparsity pattern/glossary/Sparsity Pattern/$$ 
00150 for the matrix $latex H(x)^\R{T}$$.
00151 
00152 $head VectorSet$$
00153 The type $icode VectorSet$$ must be a $cref SimpleVector$$ class with
00154 $cref/elements of type/SimpleVector/Elements of Specified Type/$$
00155 $code bool$$ or $code std::set<size_t>$$;
00156 see $cref/sparsity pattern/glossary/Sparsity Pattern/$$ for a discussion
00157 of the difference.
00158 The type of the elements of
00159 $cref/VectorSet/RevSparseHes/VectorSet/$$ must be the 
00160 same as the type of the elements of $icode r$$.
00161 
00162 $head Entire Sparsity Pattern$$
00163 Suppose that $latex q = n$$ and
00164 $latex R \in \B{R}^{n \times n}$$ is the $latex n \times n$$ identity matrix.
00165 Further suppose that the $latex S$$ is the $th k$$
00166 $cref/elementary vector/glossary/Elementary Vector/$$; i.e.
00167 $latex \[
00168 S_j = \left\{ \begin{array}{ll}
00169      1  & {\rm if} \; j = k 
00170      \\
00171      0  & {\rm otherwise}
00172 \end{array} \right. 
00173 \] $$
00174 In this case,
00175 the corresponding value $icode h$$ is a 
00176 sparsity pattern for the Hessian matrix
00177 $latex F_k^{(2)} (x) \in \B{R}^{n \times n}$$.
00178 
00179 $head Example$$
00180 $children%
00181      example/rev_sparse_hes.cpp
00182 %$$
00183 The file
00184 $cref rev_sparse_hes.cpp$$
00185 contains an example and test of this operation.
00186 It returns true if it succeeds and false otherwise.
00187 
00188 $end
00189 -----------------------------------------------------------------------------
00190 */
00191 # include <algorithm>
00192 # include <cppad/local/pod_vector.hpp>
00193 # include <cppad/local/std_set.hpp>
00194 
00195 namespace CppAD { // BEGIN_CPPAD_NAMESPACE
00196 /*!
00197 \file rev_sparse_hes.hpp
00198 Reverse mode Hessian sparsity patterns.
00199 */
00200 
00201 /*!
00202 Calculate Hessian sparsity patterns using reverse mode.
00203 
00204 The C++ source code corresponding to this operation is
00205 \verbatim
00206      h = f.RevSparseHes(q, s)
00207 \endverbatim
00208 
00209 \tparam Base
00210 is the base type for this recording.
00211 
00212 \tparam VectorSet
00213 is a simple vector with elements of type \c bool.
00214 
00215 \tparam Sparsity
00216 is either \c sparse_pack, \c sparse_set, or \c sparse_list. 
00217 
00218 \param transpose
00219 is true (false) if \c is is equal to \f$ H(x) \f$ (\f$ H(x)^T \f$)
00220 where
00221 \f[
00222      H(x) = R^T (S * F)^{(2)} (x)
00223 \f]
00224 where \f$ F \f$ is the function corresponding to the operation sequence
00225 and \a x is any argument value.
00226 
00227 \param q
00228 is the value of \a q in the 
00229 by the previous call of the form 
00230 \verbatim
00231      f.ForSparseJac(q, r)
00232 \endverbatim
00233 The value \c r in this call is a sparsity pattern for the matrix \f$ R \f$.
00234 
00235 \param s
00236 is a vector with size \c m that specifies the sparsity pattern
00237 for the vector \f$ S \f$,
00238 where \c m is the number of dependent variables
00239 corresponding to the operation sequence stored in \a play. 
00240 
00241 \param h
00242 the input value of \a h must be a vector with size \c q*n.
00243 The input value of its elements does not matter.
00244 On output, \a h is the sparsity pattern for the matrix \f$ H(x) \f$
00245 or \f$ H(x)^T \f$ depending on \c transpose.
00246 
00247 \param total_num_var
00248 is the total number of variables in this recording.
00249 
00250 \param dep_taddr
00251 maps dependendent variable index
00252 to the corresponding variable in the tape.
00253 
00254 \param ind_taddr
00255 maps independent variable index
00256 to the corresponding variable in the tape.
00257 
00258 \param play
00259 is the recording that defines the function we are computing the sparsity 
00260 pattern for.
00261 
00262 \param for_jac_sparsity
00263 is a vector of sets containing the 
00264 the forward Jacobian sparsity pattern corresponding to 
00265 $latex R$$ for all of the variables on the tape. 
00266 */
00267 
00268 template <class Base, class VectorSet, class Sparsity>
00269 void RevSparseHesBool(
00270      bool                      transpose         ,
00271      size_t                    q                 ,
00272      const VectorSet&          s                 ,
00273      VectorSet&                h                 ,
00274      size_t                    total_num_var     ,
00275      CppAD::vector<size_t>&    dep_taddr         ,
00276      CppAD::vector<size_t>&    ind_taddr         ,
00277      CppAD::player<Base>&      play              ,
00278      Sparsity&                 for_jac_sparsity  )
00279 {
00280      // temporary indices
00281      size_t i, j;
00282 
00283      // check Vector is Simple VectorSet class with bool elements
00284      CheckSimpleVector<bool, VectorSet>();
00285 
00286      // range and domain dimensions for F
00287      size_t m = dep_taddr.size();
00288      size_t n = ind_taddr.size();
00289 
00290      CPPAD_ASSERT_KNOWN(
00291           q == for_jac_sparsity.end(),
00292           "RevSparseHes: q is not equal to its value\n"
00293           "in the previous call to ForSparseJac with this ADFun object."
00294      );
00295      CPPAD_ASSERT_KNOWN(
00296           size_t(s.size()) == m,
00297           "RevSparseHes: size of s is not equal to\n"
00298           "range dimension for ADFun object."
00299      );
00300 
00301      // Array that will hold reverse Jacobian dependency flag.
00302      // Initialize as true for the dependent variables.
00303      pod_vector<bool> RevJac;
00304      RevJac.extend(total_num_var); 
00305      for(i = 0; i < total_num_var; i++)
00306           RevJac[i] = false;
00307      for(i = 0; i < m; i++)
00308      {    CPPAD_ASSERT_UNKNOWN( dep_taddr[i] < total_num_var );
00309           RevJac[ dep_taddr[i] ] = s[i];
00310      }
00311 
00312      // vector of sets that will hold reverse Hessain values
00313      Sparsity       rev_hes_sparsity;
00314      rev_hes_sparsity.resize(total_num_var, q);
00315 
00316      // compute the Hessian sparsity patterns
00317      RevHesSweep(
00318           n,
00319           total_num_var,
00320           &play,
00321           for_jac_sparsity, 
00322           RevJac.data(),
00323           rev_hes_sparsity
00324      );
00325 
00326      // return values corresponding to independent variables
00327      CPPAD_ASSERT_UNKNOWN( size_t(h.size()) == n * q );
00328      for(j = 0; j < n; j++)
00329      {    for(i = 0; i < q; i++) 
00330           {    if( transpose )
00331                     h[ j * q + i ] = false;
00332                else h[ i * n + j ] = false;
00333           }
00334      }
00335 
00336      // j is index corresponding to reverse mode partial
00337      for(j = 0; j < n; j++)
00338      {    CPPAD_ASSERT_UNKNOWN( ind_taddr[j] < total_num_var );
00339 
00340           // ind_taddr[j] is operator taddr for j-th independent variable
00341           CPPAD_ASSERT_UNKNOWN( ind_taddr[j] == j + 1 );
00342           CPPAD_ASSERT_UNKNOWN( play.GetOp( ind_taddr[j] ) == InvOp );
00343 
00344           // extract the result from rev_hes_sparsity
00345           CPPAD_ASSERT_UNKNOWN( rev_hes_sparsity.end() == q );
00346           rev_hes_sparsity.begin(j + 1);
00347           i = rev_hes_sparsity.next_element();
00348           while( i < q )
00349           {    if( transpose )
00350                     h[ j * q + i ] = true;
00351                else h[ i * n + j ] = true;
00352                i = rev_hes_sparsity.next_element();
00353           }
00354      }
00355 
00356      return;
00357 }
00358 
00359 /*!
00360 Calculate Hessian sparsity patterns using reverse mode.
00361 
00362 The C++ source code corresponding to this operation is
00363 \verbatim
00364      h = f.RevSparseHes(q, s)
00365 \endverbatim
00366 
00367 \tparam Base
00368 is the base type for this recording.
00369 
00370 \tparam VectorSet
00371 is a simple vector with elements of type \c std::set<size_t>.
00372 
00373 \tparam Sparsity
00374 is either \c sparse_pack, \c sparse_set, or \c sparse_list. 
00375 
00376 \param transpose
00377 is true (false) if \c is is equal to \f$ H(x) \f$ (\f$ H(x)^T \f$)
00378 where
00379 \f[
00380      H(x) = R^T (S * F)^{(2)} (x)
00381 \f]
00382 where \f$ F \f$ is the function corresponding to the operation sequence
00383 and \a x is any argument value.
00384 
00385 \param q
00386 is the value of \a q in the 
00387 by the previous call of the form 
00388 \verbatim
00389      f.ForSparseJac(q, r)
00390 \endverbatim
00391 The value \c r in this call is a sparsity pattern for the matrix \f$ R \f$.
00392 
00393 \param s
00394 is a vector with size \c m that specifies the sparsity pattern
00395 for the vector \f$ S \f$,
00396 where \c m is the number of dependent variables
00397 corresponding to the operation sequence stored in \a play. 
00398 
00399 \param h
00400 If \c transpose, the input value of \a h must be a vector with size \a q.
00401 Otherwise, its input value must have size \c n;
00402 On input, each element of \a h must be an empty set.
00403 On output, \a h is the sparsity pattern for the matrix \f$ H(x) \f$
00404 or \f$ H(x)^T \f$ depending on \c transpose.
00405 
00406 \param total_num_var
00407 is the total number of variables in this recording.
00408 
00409 \param dep_taddr
00410 maps dependendent variable index
00411 to the corresponding variable in the tape.
00412 
00413 \param ind_taddr
00414 maps independent variable index
00415 to the corresponding variable in the tape.
00416 
00417 \param play
00418 is the recording that defines the function we are computing the sparsity 
00419 pattern for.
00420 
00421 \param for_jac_sparsity
00422 is a vector of sets containing the 
00423 the forward Jacobian sparsity pattern corresponding to 
00424 $latex R$$ for all of the variables on the tape. 
00425 */
00426 
00427 template <class Base, class VectorSet, class Sparsity>
00428 void RevSparseHesSet(
00429      bool                      transpose         ,
00430      size_t                    q                 ,
00431      const VectorSet&          s                 ,
00432      VectorSet&                h                 ,
00433      size_t                    total_num_var     ,
00434      CppAD::vector<size_t>&    dep_taddr         ,
00435      CppAD::vector<size_t>&    ind_taddr         ,
00436      CppAD::player<Base>&      play              ,
00437      Sparsity&                 for_jac_sparsity  )
00438 {
00439      // temporary indices
00440      size_t i, j;
00441      std::set<size_t>::const_iterator itr;
00442 
00443      // check VectorSet is Simple Vector class with sets for elements
00444      CheckSimpleVector<std::set<size_t>, VectorSet>(
00445           one_element_std_set<size_t>(), two_element_std_set<size_t>()
00446      );
00447 
00448      // range and domain dimensions for F
00449 # ifndef NDEBUG
00450      size_t m = dep_taddr.size();
00451 # endif
00452      size_t n = ind_taddr.size();
00453 
00454      CPPAD_ASSERT_KNOWN(
00455           q == for_jac_sparsity.end(),
00456           "RevSparseHes: q is not equal to its value\n"
00457           "in the previous call to ForSparseJac with this ADFun object."
00458      );
00459      CPPAD_ASSERT_KNOWN(
00460           s.size() == 1,
00461           "RevSparseHes: size of s is not equal to one."
00462      );
00463 
00464      // Array that will hold reverse Jacobian dependency flag.
00465      // Initialize as true for the dependent variables.
00466      pod_vector<bool> RevJac;
00467      RevJac.extend(total_num_var); 
00468      for(i = 0; i < total_num_var; i++)
00469           RevJac[i] = false;
00470      itr = s[0].begin();
00471      while( itr != s[0].end() )
00472      {    i = *itr++;
00473           CPPAD_ASSERT_KNOWN(
00474                i < m,
00475                "RevSparseHes: an element of the set s[0] has value "
00476                "greater than or equal m"
00477           );
00478           CPPAD_ASSERT_UNKNOWN( dep_taddr[i] < total_num_var );
00479           RevJac[ dep_taddr[i] ] = true;
00480      }
00481 
00482 
00483      // vector of sets that will hold reverse Hessain values
00484      Sparsity       rev_hes_sparsity;
00485      rev_hes_sparsity.resize(total_num_var, q);
00486 
00487      // compute the Hessian sparsity patterns
00488      RevHesSweep(
00489           n,
00490           total_num_var,
00491           &play,
00492           for_jac_sparsity, 
00493           RevJac.data(),
00494           rev_hes_sparsity
00495      );
00496 
00497      // return values corresponding to independent variables
00498      // j is index corresponding to reverse mode partial
00499      CPPAD_ASSERT_UNKNOWN( size_t(h.size()) == q || transpose );
00500      CPPAD_ASSERT_UNKNOWN( size_t(h.size()) == n || ! transpose );
00501      for(j = 0; j < n; j++)
00502      {    CPPAD_ASSERT_UNKNOWN( ind_taddr[j] < total_num_var );
00503           CPPAD_ASSERT_UNKNOWN( ind_taddr[j] == j + 1 );
00504           CPPAD_ASSERT_UNKNOWN( play.GetOp( ind_taddr[j] ) == InvOp );
00505 
00506           // extract the result from rev_hes_sparsity
00507           // and add corresponding elements to result sets in h
00508           CPPAD_ASSERT_UNKNOWN( rev_hes_sparsity.end() == q );
00509           rev_hes_sparsity.begin(j+1);
00510           i = rev_hes_sparsity.next_element();
00511           while( i < q )
00512           {    if( transpose )
00513                     h[j].insert(i);
00514                else h[i].insert(j);
00515                i = rev_hes_sparsity.next_element();
00516           }
00517      }
00518 
00519      return;
00520 }
00521 
00522 
00523 
00524 /*!
00525 User API for Hessian sparsity patterns using reverse mode.
00526 
00527 The C++ source code corresponding to this operation is
00528 \verbatim
00529      h = f.RevSparseHes(q, r)
00530 \endverbatim
00531 
00532 \tparam Base
00533 is the base type for this recording.
00534 
00535 \tparam VectorSet
00536 is a simple vector with elements of type \c bool
00537 or \c std::set<size_t>.
00538 
00539 \param transpose
00540 is true (false) if \c is is equal to \f$ H(x) \f$ (\f$ H(x)^T \f$)
00541 where
00542 \f[
00543      H(x) = R^T (S * F)^{(2)} (x)
00544 \f]
00545 where \f$ F \f$ is the function corresponding to the operation sequence
00546 and \a x is any argument value.
00547 
00548 \param q
00549 is the value of \a q in the 
00550 by the previous call of the form 
00551 \verbatim
00552      f.ForSparseJac(q, r, packed)
00553 \endverbatim
00554 The value \c r in this call is a sparsity pattern for the matrix \f$ R \f$.
00555 The type of the element of \c r for the previous call to \c ForSparseJac
00556 must be the same as the type of the elements of \c s.
00557 
00558 \param s
00559 is a vector with size \c m that specifies the sparsity pattern
00560 for the vector \f$ S \f$,
00561 where \c m is the number of dependent variables
00562 corresponding to the operation sequence stored in \a play. 
00563 
00564 \return
00565 If \c transpose is false (true), 
00566 the return vector is a sparsity pattern for \f$ H(x) \f$ (\f$ H(x)^T \f$).
00567 \f[
00568      H(x) = R^T ( S * F)^{(2)} (x)
00569 \f]
00570 where \f$ F \f$ is the function corresponding to the operation sequence
00571 and \a x is any argument value.
00572 */
00573 
00574 template <class Base>
00575 template <class VectorSet>
00576 VectorSet ADFun<Base>::RevSparseHes(
00577      size_t q,  const VectorSet& s, bool transpose
00578 )
00579 {    VectorSet h;
00580      typedef typename VectorSet::value_type Set_type;
00581 
00582      // Should check to make sure q is same as in previous call to
00583      // forward sparse Jacobian.
00584      RevSparseHesCase(
00585           Set_type()    ,
00586           transpose     ,
00587           q             ,
00588           s             ,
00589           h
00590      );
00591 
00592      return h;
00593 }
00594 
00595 /*!
00596 Private helper function for RevSparseHes(q, s).
00597 
00598 All of the description in the public member function RevSparseHes(q, s)
00599 applies.
00600 
00601 \param set_type
00602 is a \c bool value. This argument is used to dispatch to the proper source
00603 code depending on the vlaue of \c VectorSet::value_type.
00604 
00605 \param transpose
00606 See \c RevSparseHes(q, s).
00607 
00608 \param q
00609 See \c RevSparseHes(q, s).
00610 
00611 \param s
00612 See \c RevSparseHes(q, s).
00613 
00614 \param h
00615 is the return value for the corresponging call to \c RevSparseJac(q, s).
00616 */
00617 template <class Base>
00618 template <class VectorSet>
00619 void ADFun<Base>::RevSparseHesCase(
00620      bool              set_type         ,
00621      bool              transpose        ,  
00622      size_t            q                ,  
00623      const VectorSet&  s                ,
00624      VectorSet&        h                )
00625 {    size_t n = Domain();     
00626      h.resize(q * n );
00627 
00628      CPPAD_ASSERT_KNOWN( 
00629           for_jac_sparse_pack_.n_set() > 0,
00630           "RevSparseHes: previous stored call to ForSparseJac did not "
00631           "use bool for the elements of r."
00632      );
00633      CPPAD_ASSERT_UNKNOWN( for_jac_sparse_set_.n_set() == 0 );
00634      CPPAD_ASSERT_UNKNOWN( for_jac_sparse_pack_.n_set() == num_var_tape_  );
00635      
00636      // use sparse_pack for the calculation
00637      CppAD::RevSparseHesBool( 
00638           transpose                ,
00639           q                        ,
00640           s                        ,
00641           h                        ,
00642           num_var_tape_            ,
00643           dep_taddr_               ,
00644           ind_taddr_               ,
00645           play_                    ,
00646           for_jac_sparse_pack_ 
00647      );
00648 }
00649 
00650 /*!
00651 Private helper function for RevSparseHes(q, s).
00652 
00653 All of the description in the public member function RevSparseHes(q, s)
00654 applies.
00655 
00656 \param set_type
00657 is a \c std::set<size_t> value. 
00658 This argument is used to dispatch to the proper source
00659 code depending on the vlaue of \c VectorSet::value_type.
00660 
00661 \param transpose
00662 See \c RevSparseHes(q, s).
00663 
00664 \param q
00665 See \c RevSparseHes(q, s).
00666 
00667 \param s
00668 See \c RevSparseHes(q, s).
00669 
00670 \param h
00671 is the return value for the corresponging call to \c RevSparseJac(q, s).
00672 */
00673 template <class Base>
00674 template <class VectorSet>
00675 void ADFun<Base>::RevSparseHesCase(
00676      const std::set<size_t>&   set_type         ,
00677      bool                      transpose        ,  
00678      size_t                    q                ,  
00679      const VectorSet&          s                ,
00680      VectorSet&                h                )
00681 {    size_t n = Domain();
00682      if( transpose )
00683           h.resize(n);
00684      else h.resize(q);
00685 
00686      CPPAD_ASSERT_KNOWN( 
00687           for_jac_sparse_set_.n_set() > 0,
00688           "RevSparseHes: previous stored call to ForSparseJac did not "
00689           "use std::set<size_t> for the elements of r."
00690      );
00691      CPPAD_ASSERT_UNKNOWN( for_jac_sparse_pack_.n_set() == 0 );
00692      CPPAD_ASSERT_UNKNOWN( for_jac_sparse_set_.n_set() == num_var_tape_  );
00693      
00694      // use sparse_pack for the calculation
00695      CppAD::RevSparseHesSet( 
00696           transpose                ,
00697           q                        ,
00698           s                        ,
00699           h                        ,
00700           num_var_tape_            ,
00701           dep_taddr_               ,
00702           ind_taddr_               ,
00703           play_                    ,
00704           for_jac_sparse_set_ 
00705      );
00706 }
00707 
00708 } // END_CPPAD_NAMESPACE
00709 # endif
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines