CppAD: A C++ Algorithmic Differentiation Package 20110419
|
00001 /* $Id$ */ 00002 /* -------------------------------------------------------------------------- 00003 CppAD: C++ Algorithmic Differentiation: Copyright (C) 2003-10 Bradley M. Bell 00004 00005 CppAD is distributed under multiple licenses. This distribution is under 00006 the terms of the 00007 Common 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 "vec_fun_pattern.hpp" 00014 CPPAD_BEGIN_NAMESPACE 00015 /*! 00016 \file vec_fun_pattern.cpp 00017 \brief Determine a sparsity pattern for a vector of AD function objects. 00018 */ 00019 00020 /*! 00021 Determine a sparsity patterns for each function in a vector of functions. 00022 00023 \param K 00024 is the number of functions that we are computing the sparsity pattern for. 00025 00026 \param p 00027 is a vector with size \c K. 00028 For <tt>k = 0 , ... , K-1, p[k]</tt> 00029 is dimension of the range space for \f$ r_k (u) \f$; i.e., 00030 \f$ r_k (u) \in {\bf R}^{p(k)} \f$. 00031 00032 \param q 00033 is a vector with size \c K. 00034 For <tt>k = 0 , ... , K-1, q[k]</tt> 00035 is dimension of the domain space for \f$ r_k (u) \f$; i.e., 00036 \f$ u \in {\bf R}^{q(k)} \f$. 00037 00038 \param retape 00039 is a vector with size \c K. 00040 For <tt>k = 0 , ... , K-1</tt>, 00041 if <tt>retape[k]</tt> is true, 00042 the function object <tt>r[k]</tt> is a valid representation 00043 for \f$ r_k (u) \f$ for all \f$ u \in {\bf R}^{q(k)} \f$. 00044 Otherwise, the function object must be retaped for each 00045 value of \f$ u \f$. 00046 00047 \param r_fun 00048 is the vector of AD function objects which has size size \c K. 00049 For <tt>k = 0 , ... , K-1</tt>, 00050 if <tt>retape[k]</tt> is true, <tt>r_fun[k]</tt> is not used. 00051 If <tt>retape[k]</tt> is false, <tt>r_fun[k]</tt> is not used. 00052 is a CppAD function object correspopnding to the function 00053 \f$ r_k : {\bf R}^{q[k]} \rightarrow {\bf R}^{p[k]} \f$. 00054 The following non-constant member functions will be called: 00055 \verbatim 00056 r_fun[k].ForSparseJac(q[k], pattern_domain) 00057 r_fun[k].RevSparseHes(p[k], pattern_range) 00058 \endverbatim 00059 The following \c const member functions <tt>r_fun[k].Range()</tt> 00060 and <tt>r_fun[k].Domain()</tt> may also be called. 00061 00062 \param pattern_jac_r 00063 is a vector with size \c K. 00064 On input, For <tt>k = 0 , ... , K-1, pattern_jac_r[k]</tt> 00065 is a vector of length p[k] * q[k] 00066 and the value of its elements does not matter. 00067 On output it is a CppAD sparsity pattern for the Jacobian of 00068 \f$ r_k (u) \f$. 00069 00070 \param pattern_hes_r 00071 is a vector with size \c K. 00072 On input, For <tt>k = 0 , ... , K-1, pattern_hes_r[k]</tt> 00073 is a vector of length q[k] * q[k] 00074 and the value of its elements does not matter. 00075 On output it is a CppAD sparsity pattern for the Hessian of 00076 \f$ R : {\bf R}^{q[k]} \rightarrow {\bf R} \f$ which is defined by 00077 \f[ 00078 R(u) = \sum_{i=0}^{p[k]-1} r_k (u)_i 00079 \f] 00080 */ 00081 void vec_fun_pattern( 00082 size_t K , 00083 const CppAD::vector<size_t>& p , 00084 const CppAD::vector<size_t>& q , 00085 const CppAD::vectorBool& retape , 00086 CppAD::vector< CppAD::ADFun<Ipopt::Number> >& r_fun , 00087 CppAD::vector<CppAD::vectorBool>& pattern_jac_r , 00088 CppAD::vector<CppAD::vectorBool>& pattern_hes_r ) 00089 { // check some assumptions 00090 CPPAD_ASSERT_UNKNOWN( K == p.size() ); 00091 CPPAD_ASSERT_UNKNOWN( K == q.size() ); 00092 CPPAD_ASSERT_UNKNOWN( K == retape.size() ); 00093 CPPAD_ASSERT_UNKNOWN( K == r_fun.size() ); 00094 CPPAD_ASSERT_UNKNOWN( K == pattern_jac_r.size() ); 00095 CPPAD_ASSERT_UNKNOWN( K == pattern_hes_r.size() ); 00096 00097 using CppAD::vectorBool; 00098 size_t i, j, k; 00099 00100 for(k = 0; k < K; k++) 00101 { // check some k specific assumptions 00102 CPPAD_ASSERT_UNKNOWN( pattern_jac_r[k].size() == p[k] * q[k] ); 00103 CPPAD_ASSERT_UNKNOWN( pattern_hes_r[k].size() == q[k] * q[k] ); 00104 00105 if( retape[k] ) 00106 { for(i = 0; i < p[k]; i++) 00107 { for(j = 0; j < q[k]; j++) 00108 pattern_jac_r[k][i*q[k] + j] = true; 00109 } 00110 for(i = 0; i < q[k]; i++) 00111 { for(j = 0; j < q[k]; j++) 00112 pattern_hes_r[k][i*q[k] + j] = true; 00113 } 00114 } 00115 else 00116 { // check assumptions about r_k 00117 CPPAD_ASSERT_UNKNOWN( r_fun[k].Range() == p[k] ); 00118 CPPAD_ASSERT_UNKNOWN( r_fun[k].Domain() == q[k] ); 00119 00120 // pattern for the identity matrix 00121 CppAD::vectorBool pattern_domain(q[k] * q[k]); 00122 for(i = 0; i < q[k]; i++) 00123 { for(j = 0; j < q[k]; j++) 00124 pattern_domain[i*q[k] + j] = (i == j); 00125 } 00126 // use forward mode to compute Jacobian sparsity 00127 pattern_jac_r[k] = 00128 r_fun[k].ForSparseJac(q[k], pattern_domain); 00129 // user reverse mode to compute Hessian sparsity 00130 CppAD::vectorBool pattern_ones(p[k]); 00131 for(i = 0; i < p[k]; i++) 00132 pattern_ones[i] = true; 00133 pattern_hes_r[k] = 00134 r_fun[k].RevSparseHes(q[k], pattern_ones); 00135 } 00136 } 00137 } 00138 CPPAD_END_NAMESPACE