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