SHOGUN
v2.0.0
|
00001 /* 00002 * This program is free software; you can redistribute it and/or modify 00003 * it under the terms of the GNU General Public License as published by 00004 * the Free Software Foundation; either version 3 of the License, or 00005 * (at your option) any later version. 00006 * 00007 * Written (W) 2012 Fernando José Iglesias García 00008 * Written (W) 2010,2012 Soeren Sonnenburg 00009 * Copyright (C) 2010 Berlin Institute of Technology 00010 * Copyright (C) 2012 Soeren Sonnenburg 00011 */ 00012 #include <shogun/lib/config.h> 00013 #include <shogun/mathematics/Math.h> 00014 #include <shogun/lib/SGVector.h> 00015 #include <shogun/lib/SGReferencedData.h> 00016 #include <shogun/mathematics/lapack.h> 00017 00018 #ifdef HAVE_EIGEN3 00019 #include <shogun/mathematics/eigen3.h> 00020 #endif 00021 00022 namespace shogun { 00023 00024 template<class T> SGVector<T>::SGVector() : SGReferencedData() 00025 { 00026 init_data(); 00027 } 00028 00029 template<class T> SGVector<T>::SGVector(T* v, index_t len, bool ref_counting) 00030 : SGReferencedData(ref_counting), vector(v), vlen(len) 00031 { 00032 } 00033 00034 template<class T> SGVector<T>::SGVector(index_t len, bool ref_counting) 00035 : SGReferencedData(ref_counting), vlen(len) 00036 { 00037 vector=SG_MALLOC(T, len); 00038 } 00039 00040 template<class T> SGVector<T>::SGVector(const SGVector &orig) : SGReferencedData(orig) 00041 { 00042 copy_data(orig); 00043 } 00044 00045 template<class T> SGVector<T>::~SGVector() 00046 { 00047 unref(); 00048 } 00049 00050 template<class T> void SGVector<T>::zero() 00051 { 00052 if (vector && vlen) 00053 set_const(0); 00054 } 00055 00056 template<class T> void SGVector<T>::set_const(T const_elem) 00057 { 00058 for (index_t i=0; i<vlen; i++) 00059 vector[i]=const_elem ; 00060 } 00061 00062 template<class T> void SGVector<T>::range_fill(T start) 00063 { 00064 range_fill_vector(vector, vlen, start); 00065 } 00066 00067 template<class T> void SGVector<T>::random(T min_value, T max_value) 00068 { 00069 random_vector(vector, vlen, min_value, max_value); 00070 } 00071 00072 template<class T> void SGVector<T>::randperm() 00073 { 00074 randperm(vector, vlen); 00075 } 00076 00077 template<class T> SGVector<T> SGVector<T>::clone() const 00078 { 00079 return SGVector<T>(clone_vector(vector, vlen), vlen); 00080 } 00081 00082 template<class T> const T& SGVector<T>::get_element(index_t index) 00083 { 00084 ASSERT(vector && (index>=0) && (index<vlen)); 00085 return vector[index]; 00086 } 00087 00088 template<class T> void SGVector<T>::set_element(const T& p_element, index_t index) 00089 { 00090 ASSERT(vector && (index>=0) && (index<vlen)); 00091 vector[index]=p_element; 00092 } 00093 00094 template<class T> void SGVector<T>::resize_vector(int32_t n) 00095 { 00096 vector=SG_REALLOC(T, vector, n); 00097 00098 if (n > vlen) 00099 memset(&vector[vlen], 0, (n-vlen)*sizeof(T)); 00100 vlen=n; 00101 } 00102 00103 template<class T> void SGVector<T>::add(const SGVector<T> x) 00104 { 00105 ASSERT(x.vector && vector); 00106 ASSERT(x.vlen == vlen); 00107 00108 for (int32_t i=0; i<vlen; i++) 00109 vector[i]+=x.vector[i]; 00110 } 00111 00112 template<class T> void SGVector<T>::add(const T x) 00113 { 00114 ASSERT(vector); 00115 00116 for (int32_t i=0; i<vlen; i++) 00117 vector[i]+=x; 00118 } 00119 00120 template<class T> void SGVector<T>::display_size() const 00121 { 00122 SG_SPRINT("SGVector '%p' of size: %d\n", vector, vlen); 00123 } 00124 00125 template<class T> void SGVector<T>::copy_data(const SGReferencedData &orig) 00126 { 00127 vector=((SGVector*)(&orig))->vector; 00128 vlen=((SGVector*)(&orig))->vlen; 00129 } 00130 00131 template<class T> void SGVector<T>::init_data() 00132 { 00133 vector=NULL; 00134 vlen=0; 00135 } 00136 00137 template<class T> void SGVector<T>::free_data() 00138 { 00139 SG_FREE(vector); 00140 vector=NULL; 00141 vlen=0; 00142 } 00143 00144 template<class T> void SGVector<T>::display_vector(const char* name, 00145 const char* prefix) const 00146 { 00147 display_vector(vector, vlen, name, prefix); 00148 } 00149 00150 template <class T> 00151 void SGVector<T>::display_vector(const SGVector<T> vector, const char* name, 00152 const char* prefix) 00153 { 00154 vector.display_vector(prefix); 00155 } 00156 00157 template <> 00158 void SGVector<bool>::display_vector(const bool* vector, int32_t n, const char* name, 00159 const char* prefix) 00160 { 00161 ASSERT(n>=0); 00162 SG_SPRINT("%s%s=[", prefix, name); 00163 for (int32_t i=0; i<n; i++) 00164 SG_SPRINT("%s%d%s", prefix, vector[i] ? 1 : 0, i==n-1? "" : ","); 00165 SG_SPRINT("%s]\n", prefix); 00166 } 00167 00168 template <> 00169 void SGVector<char>::display_vector(const char* vector, int32_t n, const char* name, 00170 const char* prefix) 00171 { 00172 ASSERT(n>=0); 00173 SG_SPRINT("%s%s=[", prefix, name); 00174 for (int32_t i=0; i<n; i++) 00175 SG_SPRINT("%s%c%s", prefix, vector[i], i==n-1? "" : ","); 00176 SG_SPRINT("%s]\n", prefix); 00177 } 00178 00179 template <> 00180 void SGVector<uint8_t>::display_vector(const uint8_t* vector, int32_t n, const char* name, 00181 const char* prefix) 00182 { 00183 ASSERT(n>=0); 00184 SG_SPRINT("%s%s=[", prefix, name); 00185 for (int32_t i=0; i<n; i++) 00186 SG_SPRINT("%s%d%s", prefix, vector[i], i==n-1? "" : ","); 00187 SG_SPRINT("%s]\n", prefix); 00188 } 00189 00190 template <> 00191 void SGVector<int8_t>::display_vector(const int8_t* vector, int32_t n, const char* name, 00192 const char* prefix) 00193 { 00194 ASSERT(n>=0); 00195 SG_SPRINT("%s%s=[", prefix, name); 00196 for (int32_t i=0; i<n; i++) 00197 SG_SPRINT("%s%d%s", prefix, vector[i], i==n-1? "" : ","); 00198 SG_SPRINT("%s]\n", prefix); 00199 } 00200 00201 template <> 00202 void SGVector<uint16_t>::display_vector(const uint16_t* vector, int32_t n, const char* name, 00203 const char* prefix) 00204 { 00205 ASSERT(n>=0); 00206 SG_SPRINT("%s%s=[", prefix, name); 00207 for (int32_t i=0; i<n; i++) 00208 SG_SPRINT("%s%d%s", prefix, vector[i], i==n-1? "" : ","); 00209 SG_SPRINT("%s]\n", prefix); 00210 } 00211 00212 template <> 00213 void SGVector<int16_t>::display_vector(const int16_t* vector, int32_t n, const char* name, 00214 const char* prefix) 00215 { 00216 ASSERT(n>=0); 00217 SG_SPRINT("%s%s=[", prefix, name); 00218 for (int32_t i=0; i<n; i++) 00219 SG_SPRINT("%s%d%s", prefix, vector[i], i==n-1? "" : ","); 00220 SG_SPRINT("%s]\n", prefix); 00221 } 00222 00223 template <> 00224 void SGVector<int32_t>::display_vector(const int32_t* vector, int32_t n, const char* name, 00225 const char* prefix) 00226 { 00227 ASSERT(n>=0); 00228 SG_SPRINT("%s%s=[", prefix, name); 00229 for (int32_t i=0; i<n; i++) 00230 SG_SPRINT("%s%d%s", prefix, vector[i], i==n-1? "" : ","); 00231 SG_SPRINT("%s]\n", prefix); 00232 } 00233 00234 template <> 00235 void SGVector<uint32_t>::display_vector(const uint32_t* vector, int32_t n, const char* name, 00236 const char* prefix) 00237 { 00238 ASSERT(n>=0); 00239 SG_SPRINT("%s%s=[", prefix, name); 00240 for (int32_t i=0; i<n; i++) 00241 SG_SPRINT("%s%d%s", prefix, vector[i], i==n-1? "" : ","); 00242 SG_SPRINT("%s]\n", prefix); 00243 } 00244 00245 00246 template <> 00247 void SGVector<int64_t>::display_vector(const int64_t* vector, int32_t n, const char* name, 00248 const char* prefix) 00249 { 00250 ASSERT(n>=0); 00251 SG_SPRINT("%s%s=[", prefix, name); 00252 for (int32_t i=0; i<n; i++) 00253 SG_SPRINT("%s%lld%s", prefix, vector[i], i==n-1? "" : ","); 00254 SG_SPRINT("%s]\n", prefix); 00255 } 00256 00257 template <> 00258 void SGVector<uint64_t>::display_vector(const uint64_t* vector, int32_t n, const char* name, 00259 const char* prefix) 00260 { 00261 ASSERT(n>=0); 00262 SG_SPRINT("%s%s=[", prefix, name); 00263 for (int32_t i=0; i<n; i++) 00264 SG_SPRINT("%s%llu%s", prefix, vector[i], i==n-1? "" : ","); 00265 SG_SPRINT("%s]\n", prefix); 00266 } 00267 00268 template <> 00269 void SGVector<float32_t>::display_vector(const float32_t* vector, int32_t n, const char* name, 00270 const char* prefix) 00271 { 00272 ASSERT(n>=0); 00273 SG_SPRINT("%s%s=[", prefix, name); 00274 for (int32_t i=0; i<n; i++) 00275 SG_SPRINT("%s%g%s", prefix, vector[i], i==n-1? "" : ","); 00276 SG_SPRINT("%s]\n", prefix); 00277 } 00278 00279 template <> 00280 void SGVector<float64_t>::display_vector(const float64_t* vector, int32_t n, const char* name, 00281 const char* prefix) 00282 { 00283 ASSERT(n>=0); 00284 SG_SPRINT("%s%s=[", prefix, name); 00285 for (int32_t i=0; i<n; i++) 00286 SG_SPRINT("%s%.18g%s", prefix, vector[i], i==n-1? "" : ","); 00287 SG_SPRINT("%s]\n", prefix); 00288 } 00289 00290 template <> 00291 void SGVector<floatmax_t>::display_vector(const floatmax_t* vector, int32_t n, 00292 const char* name, const char* prefix) 00293 { 00294 ASSERT(n>=0); 00295 SG_SPRINT("%s%s=[", prefix, name); 00296 for (int32_t i=0; i<n; i++) 00297 { 00298 SG_SPRINT("%s%.36Lg%s", prefix, (long double) vector[i], 00299 i==n-1? "" : ","); 00300 } 00301 SG_SPRINT("%s]\n", prefix); 00302 } 00303 00304 template <class T> 00305 void SGVector<T>::vec1_plus_scalar_times_vec2(float64_t* vec1, 00306 const float64_t scalar, const float64_t* vec2, int32_t n) 00307 { 00308 #ifdef HAVE_LAPACK 00309 int32_t skip=1; 00310 cblas_daxpy(n, scalar, vec2, skip, vec1, skip); 00311 #else 00312 for (int32_t i=0; i<n; i++) 00313 vec1[i]+=scalar*vec2[i]; 00314 #endif 00315 } 00316 00317 template <class T> 00318 void SGVector<T>::vec1_plus_scalar_times_vec2(float32_t* vec1, 00319 const float32_t scalar, const float32_t* vec2, int32_t n) 00320 { 00321 #ifdef HAVE_LAPACK 00322 int32_t skip=1; 00323 cblas_saxpy(n, scalar, vec2, skip, vec1, skip); 00324 #else 00325 for (int32_t i=0; i<n; i++) 00326 vec1[i]+=scalar*vec2[i]; 00327 #endif 00328 } 00329 00330 template <class T> 00331 float64_t SGVector<T>::dot(const float64_t* v1, const float64_t* v2, int32_t n) 00332 { 00333 float64_t r=0; 00334 #ifdef HAVE_EIGEN3 00335 Eigen::Map<const Eigen::VectorXd> ev1(v1,n); 00336 Eigen::Map<const Eigen::VectorXd> ev2(v2,n); 00337 r = ev1.dot(ev2); 00338 #else 00339 #ifdef HAVE_LAPACK 00340 int32_t skip=1; 00341 r = cblas_ddot(n, v1, skip, v2, skip); 00342 #else 00343 for (int32_t i=0; i<n; i++) 00344 r+=v1[i]*v2[i]; 00345 #endif /* HAVE_LAPACK */ 00346 #endif /* HAVE_EIGEN3 */ 00347 return r; 00348 } 00349 00350 template <class T> 00351 float32_t SGVector<T>::dot(const float32_t* v1, const float32_t* v2, int32_t n) 00352 { 00353 float32_t r=0; 00354 #ifdef HAVE_EIGEN3 00355 Eigen::Map<const Eigen::VectorXf> ev1(v1,n); 00356 Eigen::Map<const Eigen::VectorXf> ev2(v2,n); 00357 r = ev1.dot(ev2); 00358 #else 00359 #ifdef HAVE_LAPACK 00360 int32_t skip=1; 00361 r = cblas_sdot(n, v1, skip, v2, skip); 00362 #else 00363 for (int32_t i=0; i<n; i++) 00364 r+=v1[i]*v2[i]; 00365 #endif /* HAVE_LAPACK */ 00366 #endif /* HAVE_EIGEN3 */ 00367 return r; 00368 } 00369 00371 template <class T> 00372 void SGVector<T>::random_vector(T* vec, int32_t len, T min_value, T max_value) 00373 { 00374 for (int32_t i=0; i<len; i++) 00375 vec[i]=CMath::random(min_value, max_value); 00376 } 00377 00379 template <class T> 00380 void SGVector<T>::randperm(T* perm, int32_t n) 00381 { 00382 for (int32_t i = 0; i < n; i++) 00383 perm[i] = i; 00384 permute(perm,n); 00385 } 00386 00388 template <class T> 00389 void SGVector<T>::permute(T* vec, int32_t n) 00390 { 00391 for (int32_t i = 0; i < n; i++) 00392 CMath::swap(vec[i], vec[CMath::random(i, n-1)]); 00393 } 00394 00395 template<class T> void SGVector<T>::permute() 00396 { 00397 SGVector<T>::permute(vector, vlen); 00398 } 00399 00400 template <class T> 00401 void SGVector<T>::permute_vector(SGVector<T> vec) 00402 { 00403 for (index_t i=0; i<vec.vlen; ++i) 00404 { 00405 CMath::swap(vec.vector[i], 00406 vec.vector[CMath::random(i, vec.vlen-1)]); 00407 } 00408 } 00409 00410 template <> 00411 bool SGVector<bool>::twonorm(const bool* x, int32_t len) 00412 { 00413 SG_SNOTIMPLEMENTED; 00414 return false; 00415 } 00416 00417 template <> 00418 char SGVector<char>::twonorm(const char* x, int32_t len) 00419 { 00420 SG_SNOTIMPLEMENTED; 00421 return '\0'; 00422 } 00423 00424 template <> 00425 int8_t SGVector<int8_t>::twonorm(const int8_t* x, int32_t len) 00426 { 00427 float64_t result=0; 00428 for (int32_t i=0; i<len; i++) 00429 result+=x[i]*x[i]; 00430 00431 return CMath::sqrt(result); 00432 } 00433 00434 template <> 00435 uint8_t SGVector<uint8_t>::twonorm(const uint8_t* x, int32_t len) 00436 { 00437 float64_t result=0; 00438 for (int32_t i=0; i<len; i++) 00439 result+=x[i]*x[i]; 00440 00441 return CMath::sqrt(result); 00442 } 00443 00444 template <> 00445 int16_t SGVector<int16_t>::twonorm(const int16_t* x, int32_t len) 00446 { 00447 float64_t result=0; 00448 for (int32_t i=0; i<len; i++) 00449 result+=x[i]*x[i]; 00450 00451 return CMath::sqrt(result); 00452 } 00453 00454 template <> 00455 uint16_t SGVector<uint16_t>::twonorm(const uint16_t* x, int32_t len) 00456 { 00457 float64_t result=0; 00458 for (int32_t i=0; i<len; i++) 00459 result+=x[i]*x[i]; 00460 00461 return CMath::sqrt(result); 00462 } 00463 00464 template <> 00465 int32_t SGVector<int32_t>::twonorm(const int32_t* x, int32_t len) 00466 { 00467 float64_t result=0; 00468 for (int32_t i=0; i<len; i++) 00469 result+=x[i]*x[i]; 00470 00471 return CMath::sqrt(result); 00472 } 00473 00474 template <> 00475 uint32_t SGVector<uint32_t>::twonorm(const uint32_t* x, int32_t len) 00476 { 00477 float64_t result=0; 00478 for (int32_t i=0; i<len; i++) 00479 result+=x[i]*x[i]; 00480 00481 return CMath::sqrt(result); 00482 } 00483 00484 template <> 00485 int64_t SGVector<int64_t>::twonorm(const int64_t* x, int32_t len) 00486 { 00487 float64_t result=0; 00488 for (int32_t i=0; i<len; i++) 00489 result+=x[i]*x[i]; 00490 00491 return CMath::sqrt(result); 00492 } 00493 00494 template <> 00495 uint64_t SGVector<uint64_t>::twonorm(const uint64_t* x, int32_t len) 00496 { 00497 float64_t result=0; 00498 for (int32_t i=0; i<len; i++) 00499 result+=x[i]*x[i]; 00500 00501 return CMath::sqrt(result); 00502 } 00503 00504 template <> 00505 float32_t SGVector<float32_t>::twonorm(const float32_t* x, int32_t len) 00506 { 00507 float64_t result=0; 00508 for (int32_t i=0; i<len; i++) 00509 result+=x[i]*x[i]; 00510 00511 return CMath::sqrt(result); 00512 } 00513 00514 template <> 00515 float64_t SGVector<float64_t>::twonorm(const float64_t* v, int32_t n) 00516 { 00517 float64_t norm = 0.0; 00518 #ifdef HAVE_LAPACK 00519 norm = cblas_dnrm2(n, v, 1); 00520 #else 00521 norm = CMath::sqrt(SGVector::dot(v, v, n)); 00522 #endif 00523 return norm; 00524 } 00525 00526 template <> 00527 floatmax_t SGVector<floatmax_t>::twonorm(const floatmax_t* x, int32_t len) 00528 { 00529 float64_t result=0; 00530 for (int32_t i=0; i<len; i++) 00531 result+=x[i]*x[i]; 00532 00533 return CMath::sqrt(result); 00534 } 00535 00536 00537 template <class T> 00538 float64_t SGVector<T>::onenorm(T* x, int32_t len) 00539 { 00540 float64_t result=0; 00541 for (int32_t i=0;i<len; ++i) 00542 result+=CMath::abs(x[i]); 00543 00544 return result; 00545 } 00546 00548 template <class T> 00549 T SGVector<T>::qsq(T* x, int32_t len, float64_t q) 00550 { 00551 float64_t result=0; 00552 for (int32_t i=0; i<len; i++) 00553 result+=CMath::pow(fabs(x[i]), q); 00554 00555 return result; 00556 } 00557 00559 template <class T> 00560 T SGVector<T>::qnorm(T* x, int32_t len, float64_t q) 00561 { 00562 ASSERT(q!=0); 00563 return CMath::pow((float64_t) qsq(x, len, q), 1.0/q); 00564 } 00565 00567 template <class T> 00568 T SGVector<T>::min(T* vec, int32_t len) 00569 { 00570 ASSERT(len>0); 00571 T minv=vec[0]; 00572 00573 for (int32_t i=1; i<len; i++) 00574 minv=CMath::min(vec[i], minv); 00575 00576 return minv; 00577 } 00578 00580 template <class T> 00581 T SGVector<T>::max(T* vec, int32_t len) 00582 { 00583 ASSERT(len>0); 00584 T maxv=vec[0]; 00585 00586 for (int32_t i=1; i<len; i++) 00587 maxv=CMath::max(vec[i], maxv); 00588 00589 return maxv; 00590 } 00591 00593 template <class T> 00594 T SGVector<T>::sum_abs(T* vec, int32_t len) 00595 { 00596 T result=0; 00597 for (int32_t i=0; i<len; i++) 00598 result+=CMath::abs(vec[i]); 00599 00600 return result; 00601 } 00602 00604 template <class T> 00605 bool SGVector<T>::fequal(T x, T y, float64_t precision) 00606 { 00607 return CMath::abs(x-y)<precision; 00608 } 00609 00610 template <class T> 00611 int32_t SGVector<T>::unique(T* output, int32_t size) 00612 { 00613 CMath::qsort(output, size); 00614 int32_t j=0; 00615 00616 for (int32_t i=0; i<size; i++) 00617 { 00618 if (i==0 || output[i]!=output[i-1]) 00619 output[j++]=output[i]; 00620 } 00621 return j; 00622 } 00623 00624 template <class T> 00625 SGVector<index_t> SGVector<T>::find(T elem) 00626 { 00627 SGVector<index_t> idx(vlen); 00628 index_t k=0; 00629 00630 for (index_t i=0; i < vlen; ++i) 00631 if (vector[i] == elem) 00632 idx[k++] = i; 00633 idx.vlen = k; 00634 return idx; 00635 } 00636 00637 template<class T> void SGVector<T>::scale(T alpha) 00638 { 00639 scale_vector(alpha, vector, vlen); 00640 } 00641 00642 template<class T> float64_t SGVector<T>::mean() const 00643 { 00644 float64_t cum = 0; 00645 00646 for ( index_t i = 0 ; i < vlen ; ++i ) 00647 cum += vector[i]; 00648 00649 return cum/vlen; 00650 } 00651 00652 template class SGVector<bool>; 00653 template class SGVector<char>; 00654 template class SGVector<int8_t>; 00655 template class SGVector<uint8_t>; 00656 template class SGVector<int16_t>; 00657 template class SGVector<uint16_t>; 00658 template class SGVector<int32_t>; 00659 template class SGVector<uint32_t>; 00660 template class SGVector<int64_t>; 00661 template class SGVector<uint64_t>; 00662 template class SGVector<float32_t>; 00663 template class SGVector<float64_t>; 00664 template class SGVector<floatmax_t>; 00665 }