SHOGUN  v2.0.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines
SGVector.cpp
Go to the documentation of this file.
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 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines

SHOGUN Machine Learning Toolbox - Documentation