CppAD: A C++ Algorithmic Differentiation Package
20130102
|
00001 // $Id$ 00002 # ifndef CPPAD_SPARSE_PACK_INCLUDED 00003 # define CPPAD_SPARSE_PACK_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 # include <cppad/local/cppad_assert.hpp> 00016 # include <cppad/local/pod_vector.hpp> 00017 00018 CPPAD_BEGIN_NAMESPACE 00019 /*! 00020 \defgroup sparse_pack_hpp sparse_pack.hpp 00021 \{ 00022 \file sparse_pack.hpp 00023 Vector of sets of positive integers stored as a packed array of bools. 00024 */ 00025 00026 /*! 00027 Vector of sets of postivie integers, each set stored as a packed boolean array. 00028 */ 00029 00030 class sparse_pack { 00031 private: 00032 /// Type used to pack elements (should be the same as corresponding 00033 /// typedef in multiple_n_bit() in test_more/sparse_hacobian.cpp) 00034 typedef size_t Pack; 00035 /// Number of bits per Pack value 00036 static const size_t n_bit_ = std::numeric_limits<Pack>::digits; 00037 /// Number of sets that we are representing 00038 /// (set by constructor and resize). 00039 size_t n_set_; 00040 /// Possible elements in each set are 0, 1, ..., end_ - 1 00041 /// (set by constructor and resize). 00042 size_t end_; 00043 /// Number of \c Pack values necessary to represent \c end_ bits. 00044 /// (set by constructor and resize). 00045 size_t n_pack_; 00046 /// Data for all the sets. 00047 pod_vector<Pack> data_; 00048 /// index for which we were retrieving next_element 00049 /// (use n_set_ if no such index exists). 00050 size_t next_index_; 00051 /// Next element to start search at 00052 /// (use end_ for no such element exists; i.e., past end of the set). 00053 size_t next_element_; 00054 public: 00055 // ----------------------------------------------------------------- 00056 /*! Default constructor (no sets) 00057 */ 00058 sparse_pack(void) : 00059 n_set_(0) , 00060 end_(0) , 00061 n_pack_(0) , 00062 next_index_(0) , 00063 next_element_(0) 00064 { } 00065 // ----------------------------------------------------------------- 00066 /*! Make use of copy constructor an error 00067 00068 \param v 00069 vector that we are attempting to make a copy of. 00070 */ 00071 sparse_pack(const sparse_pack& v) 00072 { // Error: 00073 // Probably a sparse_pack argument has been passed by value 00074 CPPAD_ASSERT_UNKNOWN(0); 00075 } 00076 // ----------------------------------------------------------------- 00077 /*! Destructor 00078 */ 00079 ~sparse_pack(void) 00080 { } 00081 // ----------------------------------------------------------------- 00082 /*! Change number of sets, set end, and initialize all sets as empty 00083 00084 If \c n_set_in is zero, any memory currently allocated for this object 00085 is freed. Otherwise, new memory may be allocated for the sets (if needed). 00086 00087 \param n_set_in 00088 is the number of sets in this vector of sets. 00089 00090 \param end_in 00091 is the maximum element plus one (the minimum element is 0). 00092 */ 00093 void resize(size_t n_set_in, size_t end_in) 00094 { 00095 n_set_ = n_set_in; 00096 end_ = end_in; 00097 if( n_set_ == 0 ) 00098 { data_.free(); 00099 return; 00100 } 00101 // now start a new vector with empty sets 00102 Pack zero(0); 00103 data_.erase(); 00104 00105 n_pack_ = ( 1 + (end_ - 1) / n_bit_ ); 00106 size_t i = n_set_ * n_pack_; 00107 00108 if( i > 0 ) 00109 { data_.extend(i); 00110 while(i--) 00111 data_[i] = zero; 00112 } 00113 00114 // values that signify past end of list 00115 next_index_ = n_set_; 00116 next_element_ = end_; 00117 } 00118 // ----------------------------------------------------------------- 00119 /*! Add one element to a set. 00120 00121 \param index 00122 is the index for this set in the vector of sets. 00123 00124 \param element 00125 is the element we are adding to the set. 00126 00127 \par Checked Assertions 00128 \li index < n_set_ 00129 \li element < end_ 00130 */ 00131 void add_element(size_t index, size_t element) 00132 { static Pack one(1); 00133 CPPAD_ASSERT_UNKNOWN( index < n_set_ ); 00134 CPPAD_ASSERT_UNKNOWN( element < end_ ); 00135 size_t j = element / n_bit_; 00136 size_t k = element - j * n_bit_; 00137 Pack mask = one << k; 00138 data_[ index * n_pack_ + j] |= mask; 00139 } 00140 // ----------------------------------------------------------------- 00141 /*! Is an element of a set. 00142 00143 \param index 00144 is the index for this set in the vector of sets. 00145 00146 \param element 00147 is the element we are checking to see if it is in the set. 00148 00149 \par Checked Assertions 00150 \li index < n_set_ 00151 \li element < end_ 00152 */ 00153 bool is_element(size_t index, size_t element) 00154 { static Pack one(1); 00155 static Pack zero(0); 00156 CPPAD_ASSERT_UNKNOWN( index < n_set_ ); 00157 CPPAD_ASSERT_UNKNOWN( element < end_ ); 00158 size_t j = element / n_bit_; 00159 size_t k = element - j * n_bit_; 00160 Pack mask = one << k; 00161 return (data_[ index * n_pack_ + j] & mask) != zero; 00162 } 00163 // ----------------------------------------------------------------- 00164 /*! Begin retrieving elements from one of the sets. 00165 00166 \param index 00167 is the index for the set that is going to be retrieved. 00168 The elements of the set are retrieved in increasing order. 00169 00170 \par Checked Assertions 00171 \li index < n_set_ 00172 */ 00173 void begin(size_t index) 00174 { // initialize element to search for in this set 00175 CPPAD_ASSERT_UNKNOWN( index < n_set_ ); 00176 next_index_ = index; 00177 next_element_ = 0; 00178 } 00179 /*! Get the next element from the current retrieval set. 00180 00181 \return 00182 is the next element in the set with index 00183 specified by the previous call to \c begin. 00184 If no such element exists, \c this->end() is returned. 00185 */ 00186 size_t next_element(void) 00187 { static Pack one(1); 00188 CPPAD_ASSERT_UNKNOWN( next_index_ < n_set_ ); 00189 CPPAD_ASSERT_UNKNOWN( next_element_ <= end_ ); 00190 00191 if( next_element_ == end_ ) 00192 return end_; 00193 00194 // initialize packed data index 00195 size_t j = next_element_ / n_bit_; 00196 00197 // initialize bit index 00198 size_t k = next_element_ - j * n_bit_; 00199 00200 // start search at this packed value 00201 Pack check = data_[ next_index_ * n_pack_ + j ]; 00202 while( true ) 00203 { // increment next element before checking this one 00204 next_element_++; 00205 // check if this element is in the set 00206 if( check & (one << k) ) 00207 return next_element_ - 1; 00208 // check if no more elements in the set 00209 if( next_element_ == end_ ) 00210 return end_; 00211 // increment bit index in Pack value so corresponds to 00212 // next element 00213 k++; 00214 CPPAD_ASSERT_UNKNOWN( k <= n_bit_ ); 00215 if( k == n_bit_ ) 00216 { // get next packed value 00217 k = 0; 00218 j++; 00219 CPPAD_ASSERT_UNKNOWN( j < n_pack_ ); 00220 check = data_[ next_index_ * n_pack_ + j ]; 00221 } 00222 } 00223 // should never get here 00224 CPPAD_ASSERT_UNKNOWN(false); 00225 return end_; 00226 } 00227 // ----------------------------------------------------------------- 00228 /*! Assign the empty set to one of the sets. 00229 00230 \param target 00231 is the index of the set we are setting to the empty set. 00232 00233 \par Checked Assertions 00234 \li target < n_set_ 00235 */ 00236 void clear(size_t target) 00237 { // value with all its bits set to false 00238 static Pack zero(0); 00239 CPPAD_ASSERT_UNKNOWN( target < n_set_ ); 00240 size_t t = target * n_pack_; 00241 00242 size_t j = n_pack_; 00243 while(j--) 00244 data_[t++] = zero; 00245 } 00246 // ----------------------------------------------------------------- 00247 /*! Assign one set equal to another set. 00248 00249 \param this_target 00250 is the index (in this \c sparse_pack object) of the set being assinged. 00251 00252 \param other_value 00253 is the index (in the other \c sparse_pack object) of the 00254 that we are using as the value to assign to the target set. 00255 00256 \param other 00257 is the other \c sparse_pack object (which may be the same as this 00258 \c sparse_pack object). 00259 00260 \par Checked Assertions 00261 \li this_target < n_set_ 00262 \li other_value < other.n_set_ 00263 \li n_pack_ == other.n_pack_ 00264 */ 00265 void assignment( 00266 size_t this_target , 00267 size_t other_value , 00268 const sparse_pack& other ) 00269 { CPPAD_ASSERT_UNKNOWN( this_target < n_set_ ); 00270 CPPAD_ASSERT_UNKNOWN( other_value < other.n_set_ ); 00271 CPPAD_ASSERT_UNKNOWN( n_pack_ == other.n_pack_ ); 00272 size_t t = this_target * n_pack_; 00273 size_t v = other_value * n_pack_; 00274 00275 size_t j = n_pack_; 00276 while(j--) 00277 data_[t++] = other.data_[v++]; 00278 } 00279 00280 // ----------------------------------------------------------------- 00281 /*! Assing a set equal to the union of two other sets. 00282 00283 \param this_target 00284 is the index (in this \c sparse_pack object) of the set being assinged. 00285 00286 \param this_left 00287 is the index (in this \c sparse_pack object) of the 00288 left operand for the union operation. 00289 It is OK for \a this_target and \a this_left to be the same value. 00290 00291 \param other_right 00292 is the index (in the other \c sparse_pack object) of the 00293 right operand for the union operation. 00294 It is OK for \a this_target and \a other_right to be the same value. 00295 00296 \param other 00297 is the other \c sparse_pack object (which may be the same as this 00298 \c sparse_pack object). 00299 00300 \par Checked Assertions 00301 \li this_target < n_set_ 00302 \li this_left < n_set_ 00303 \li other_right < other.n_set_ 00304 \li n_pack_ == other.n_pack_ 00305 */ 00306 void binary_union( 00307 size_t this_target , 00308 size_t this_left , 00309 size_t other_right , 00310 const sparse_pack& other ) 00311 { CPPAD_ASSERT_UNKNOWN( this_target < n_set_ ); 00312 CPPAD_ASSERT_UNKNOWN( this_left < n_set_ ); 00313 CPPAD_ASSERT_UNKNOWN( other_right < other.n_set_ ); 00314 CPPAD_ASSERT_UNKNOWN( n_pack_ == other.n_pack_ ); 00315 00316 size_t t = this_target * n_pack_; 00317 size_t l = this_left * n_pack_; 00318 size_t r = other_right * n_pack_; 00319 00320 size_t j = n_pack_; 00321 while(j--) 00322 data_[t++] = ( data_[l++] | other.data_[r++] ); 00323 } 00324 // ----------------------------------------------------------------- 00325 /*! Amount of memory used by this vector of sets 00326 00327 \return 00328 The amount of memory in units of type unsigned char memory. 00329 */ 00330 size_t memory(void) const 00331 { return n_set_ * n_pack_ * sizeof(Pack); 00332 } 00333 // ----------------------------------------------------------------- 00334 /*! Fetch n_set for vector of sets object. 00335 00336 \return 00337 Number of from sets for this vector of sets object 00338 */ 00339 size_t n_set(void) const 00340 { return n_set_; } 00341 // ----------------------------------------------------------------- 00342 /*! Fetch end for this vector of sets object. 00343 00344 \return 00345 is the maximum element value plus one (the minimum element value is 0). 00346 */ 00347 size_t end(void) const 00348 { return end_; } 00349 }; 00350 00351 /*! 00352 Copy a user vector of bools sparsity pattern to an internal sparse_pack object. 00353 00354 \tparam VectorBool 00355 is a simple vector with elements of type bool. 00356 00357 \param internal 00358 The input value of sparisty does not matter. 00359 Upon return it contains the same sparsity pattern as \c user 00360 (or the transposed sparsity pattern). 00361 00362 \param user 00363 sparsity pattern that we are placing \c internal. 00364 00365 \param n_row 00366 number of rows in the sparsity pattern in \c user 00367 (rand dimension). 00368 00369 \param n_col 00370 number of columns in the sparsity pattern in \c user 00371 (domain dimension). 00372 00373 \param transpose 00374 if true, the sparsity pattern in \c internal is the transpose 00375 of the one in \c user. 00376 Otherwise it is the same sparsity pattern. 00377 */ 00378 template<class VectorBool> 00379 void sparsity_user2internal( 00380 sparse_pack& internal , 00381 const VectorBool& user , 00382 size_t n_row , 00383 size_t n_col , 00384 bool transpose ) 00385 { CPPAD_ASSERT_UNKNOWN( n_row * n_col == size_t(user.size()) ); 00386 size_t i, j; 00387 00388 CPPAD_ASSERT_KNOWN( 00389 size_t( user.size() ) == n_row * n_col, 00390 "Size of this vector of bools sparsity pattern is not equal product " 00391 "of the domain and range dimensions for corresponding function." 00392 ); 00393 00394 // transposed pattern case 00395 if( transpose ) 00396 { internal.resize(n_col, n_row); 00397 for(j = 0; j < n_col; j++) 00398 { for(i = 0; i < n_row; i++) 00399 { if( user[ i * n_col + j ] ) 00400 internal.add_element(j, i); 00401 } 00402 } 00403 return; 00404 } 00405 00406 // same pattern case 00407 internal.resize(n_row, n_col); 00408 for(i = 0; i < n_row; i++) 00409 { for(j = 0; j < n_col; j++) 00410 { if( user[ i * n_col + j ] ) 00411 internal.add_element(i, j); 00412 } 00413 } 00414 return; 00415 } 00416 00417 /*! \} */ 00418 CPPAD_END_NAMESPACE 00419 # endif