CppAD: A C++ Algorithmic Differentiation Package  20130918
hes_fg_map.cpp
Go to the documentation of this file.
00001 /* $Id$ */
00002 /* --------------------------------------------------------------------------
00003 CppAD: C++ Algorithmic Differentiation: Copyright (C) 2003-14 Bradley M. Bell
00004 
00005 CppAD is distributed under multiple licenses. This distribution is under
00006 the terms of the 
00007                     Eclipse Public License Version 1.0.
00008 
00009 A copy of this license is included in the COPYING file of this distribution.
00010 Please visit http://www.coin-or.org/CppAD/ for information on other licenses.
00011 -------------------------------------------------------------------------- */
00012 # include "cppad_ipopt_nlp.hpp"
00013 # include "hes_fg_map.hpp"
00014 
00015 // ---------------------------------------------------------------------------
00016 namespace cppad_ipopt {
00017 // ---------------------------------------------------------------------------
00018 /*!
00019 \{
00020 \file hes_fg_map.cpp
00021 \brief Creates a mapping between two representations for Hessian of fg. 
00022 */
00023 
00024 
00025 /*!
00026 Create mapping from CppAD to Ipopt sparse representations of Hessian of F(x).
00027 
00028 The functions 
00029 \f$ f : {\bf R}^n \rightarrow {\bf R} \f$ and
00030 \f$ g : {\bf R}^n \rightarrow {\bf R}^m \f$ are defined by the 
00031 \ref Users_Representation.
00032 We define the function 
00033 \f$ F : {\bf R}^n \rightarrow {\bf R} \f$ by
00034 \f[
00035      F(x) = \sum_{i=0}^m fg(x)_i
00036 \f]
00037 
00038 \param fg_info
00039 For <tt>k = 0 , ... , K-1</tt>,
00040 for <tt>ell = 0 , ... , L[k]</tt>,
00041 the function call 
00042 \verbatim
00043      fg_info->index(k, ell, I, J); 
00044 \endverbatim
00045 is made by \c hes_fg_map. 
00046 The values \c k and \c ell are inputs.
00047 The input size of \c I ( \c J ) 
00048 is greater than or equal <tt>p[k] ( q[k] )</tt>
00049 and this size is not changed.
00050 The input values of the elements of \c I and \c J are not specified.
00051 The output value of the elements of \c I define
00052 \f[
00053 I_{k, \ell} = ( {\rm I[0]} , \cdots , {\rm I[p[k]-1]} )
00054 \f]
00055 The output value of the elements of \c J define
00056 \f[
00057 J_{k, \ell} = ( {\rm J[0]} , \cdots , {\rm J[q[k]-1]} )
00058 \f]
00059 
00060 \param m
00061 is the dimension of the range space for \f$ g(x) \f$; i.e.,
00062 \f$ g(x) \in {\bf R}^m \f$.
00063 
00064 \param n
00065 is the dimension of the domain space for \f$ f(x) \f$ and \f$ g(x) \f$;
00066 i.e., \f$ x \in {\bf R}^n \f$.
00067 
00068 \param K
00069 is the number of functions \f$ r_k ( u ) \f$ used for the representation of 
00070 \f$ f(x) \f$ and \f$ g(x) \f$.
00071 
00072 \param L
00073 is a vector with size \c K.
00074 For <tt>k = 0 , ... , K-1, L[k]</tt>
00075 is the number of terms that use \f$ r_k (u) \f$
00076 in the representation of \f$ f(x) \f$ and \f$ g(x) \f$.
00077 
00078 \param p
00079 is a vector with size \c K.
00080 For <tt>k = 0 , ... , K-1, p[k]</tt>
00081 is dimension of the range space for \f$ r_k (u) \f$; i.e.,
00082 \f$ r_k (u) \in {\bf R}^{p(k)} \f$.
00083 
00084 \param q
00085 is a vector with size \c K.
00086 For <tt>k = 0 , ... , K-1, q[k]</tt>
00087 is dimension of the domain space for \f$ r_k (u) \f$; i.e.,
00088 \f$ u \in {\bf R}^{q(k)} \f$.
00089 
00090 \param pattern_hes_r
00091 is a vector with size \c K.
00092 For <tt>k = 0 , ... , K-1, pattern_jac_r[k]</tt>
00093 is a CppAD sparsity pattern for the Hessian of the function 
00094 \f[
00095      R(u) = \sum_{i=0}^{p[k]-1} r_k (u)_i
00096 \f]
00097 As such, <tt>pattern_hes_r[k].size() == q[k] * q[k]</tt>.
00098 
00099 \param I
00100 is a work vector of length greater than or equal <tt>p[k]</tt> for all \c k.
00101 The input and output value of its elements are unspecified.
00102 The size of \c I is not changed.
00103 
00104 \param J
00105 is a work vector of length greater than or equal <tt>q[k]</tt> for all \c k.
00106 The input and output value of its elements are unspecified.
00107 The size of \c J is not changed.
00108 
00109 \param index_hes_fg:
00110 On input, this is empty; i.e., <tt>index_jac_g.size() == 0</tt>.
00111 On output, it is the index mapping from \f$ (i, j) \f$ in the Jacobian of 
00112 \f$ g(x) \f$ to the corresponding index value used by Ipopt to represent
00113 the Jacobian.
00114 Furthermore, if <tt>index_jac_g[i].find(j) == index_jac_g[i].end()</tt>,
00115 then the \f$ (i, j)\f$ entry in the Jacobian of \f$ g(x) \f$ is always zero.
00116 */
00117 void hes_fg_map(
00118      cppad_ipopt_fg_info*  fg_info                                  , 
00119      size_t                                          m              ,
00120      size_t                                          n              ,
00121      size_t                                          K              ,
00122      const CppAD::vector<size_t>&                    L              ,
00123      const CppAD::vector<size_t>&                    p              ,
00124      const CppAD::vector<size_t>&                    q              ,
00125      const CppAD::vector<CppAD::vectorBool>&         pattern_hes_r  ,
00126      CppAD::vector<size_t>&                          I              ,
00127      CppAD::vector<size_t>&                          J              ,
00128      CppAD::vector< std::map<size_t,size_t> >&       index_hes_fg   )
00129 {
00130      using CppAD::vectorBool;
00131      size_t i, j, ij, k, ell;
00132 
00133      CPPAD_ASSERT_UNKNOWN( K == L.size() );
00134      CPPAD_ASSERT_UNKNOWN( K == p.size() );
00135      CPPAD_ASSERT_UNKNOWN( K == q.size() );
00136      CPPAD_ASSERT_UNKNOWN( K == pattern_hes_r.size() );
00137 # ifndef NDEBUG
00138      for(k = 0; k < K; k++)
00139      {     CPPAD_ASSERT_UNKNOWN( p[k] <= I.size() );
00140            CPPAD_ASSERT_UNKNOWN( q[k] <= J.size() );
00141            CPPAD_ASSERT_UNKNOWN( q[k]*q[k] == pattern_hes_r[k].size() );
00142      }
00143 # endif
00144      
00145      // Now compute pattern for fg
00146      // (use standard set representation because can be huge).
00147      CppAD::vector< std::set<size_t> > pattern_hes_fg(n);
00148      for(k = 0; k < K; k++) for(ell = 0; ell < L[k]; ell++)
00149      {    fg_info->index(k, ell, I, J); 
00150           for(i = 0; i < q[k]; i++) 
00151           {    for(j = 0; j < q[k]; j++)
00152                {    ij  = i * q[k] + j;
00153                     if( pattern_hes_r[k][ij] )
00154                          pattern_hes_fg[J[i]].insert(J[j]);
00155                }
00156           }
00157      }
00158 
00159      // Now compute the mapping from (i, j) in the Hessian of fg to the 
00160      // corresponding index value used by Ipopt to represent the Hessian.
00161      CPPAD_ASSERT_UNKNOWN( index_hes_fg.size() == 0 );
00162      index_hes_fg.resize(n);
00163      std::set<size_t>::const_iterator itr;
00164      ell = 0;
00165      for(i = 0; i < n; i++)
00166      {    for( itr = pattern_hes_fg[i].begin(); 
00167                itr != pattern_hes_fg[i].end(); itr++)
00168           {
00169                index_hes_fg[i][*itr] = ell++;
00170           }
00171      }
00172      return;
00173 }
00174 
00175 // ---------------------------------------------------------------------------
00176 } // end namespace cppad_ipopt
00177 // ---------------------------------------------------------------------------
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines