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