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