00001
00029 #ifndef VEC_H
00030 #define VEC_H
00031
00032 #ifndef _MSC_VER
00033 # include <itpp/config.h>
00034 #else
00035 # include <itpp/config_msvc.h>
00036 #endif
00037
00038 #include <itpp/base/itassert.h>
00039 #include <itpp/base/math/misc.h>
00040 #include <itpp/base/copy_vector.h>
00041 #include <itpp/base/factory.h>
00042 #include <vector>
00043
00044
00045 namespace itpp
00046 {
00047
00048
00049 template<class Num_T> class Vec;
00050
00051 template<class Num_T> class Mat;
00052
00053 class bin;
00054
00055
00056
00057
00058
00060 template<class Num_T>
00061 Vec<Num_T> operator+(const Vec<Num_T> &v1, const Vec<Num_T> &v2);
00063 template<class Num_T>
00064 Vec<Num_T> operator+(const Vec<Num_T> &v, Num_T t);
00066 template<class Num_T>
00067 Vec<Num_T> operator+(Num_T t, const Vec<Num_T> &v);
00068
00070 template<class Num_T>
00071 Vec<Num_T> operator-(const Vec<Num_T> &v1, const Vec<Num_T> &v2);
00073 template<class Num_T>
00074 Vec<Num_T> operator-(const Vec<Num_T> &v, Num_T t);
00076 template<class Num_T>
00077 Vec<Num_T> operator-(Num_T t, const Vec<Num_T> &v);
00079 template<class Num_T>
00080 Vec<Num_T> operator-(const Vec<Num_T> &v);
00081
00083 template<class Num_T>
00084 Num_T dot(const Vec<Num_T> &v1, const Vec<Num_T> &v2);
00086 template<class Num_T>
00087 Num_T operator*(const Vec<Num_T> &v1, const Vec<Num_T> &v2);
00096 template<class Num_T>
00097 Mat<Num_T> outer_product(const Vec<Num_T> &v1, const Vec<Num_T> &v2,
00098 bool hermitian = false);
00100 template<class Num_T>
00101 Vec<Num_T> operator*(const Vec<Num_T> &v, Num_T t);
00103 template<class Num_T>
00104 Vec<Num_T> operator*(Num_T t, const Vec<Num_T> &v);
00105
00107 template<class Num_T>
00108 Vec<Num_T> elem_mult(const Vec<Num_T> &a, const Vec<Num_T> &b);
00110 template<class Num_T>
00111 Vec<Num_T> elem_mult(const Vec<Num_T> &a, const Vec<Num_T> &b,
00112 const Vec<Num_T> &c);
00114 template<class Num_T>
00115 Vec<Num_T> elem_mult(const Vec<Num_T> &a, const Vec<Num_T> &b,
00116 const Vec<Num_T> &c, const Vec<Num_T> &d);
00117
00119 template<class Num_T>
00120 void elem_mult_out(const Vec<Num_T> &a, const Vec<Num_T> &b,
00121 Vec<Num_T> &out);
00123 template<class Num_T>
00124 void elem_mult_out(const Vec<Num_T> &a, const Vec<Num_T> &b,
00125 const Vec<Num_T> &c, Vec<Num_T> &out);
00127 template<class Num_T>
00128 void elem_mult_out(const Vec<Num_T> &a, const Vec<Num_T> &b,
00129 const Vec<Num_T> &c, const Vec<Num_T> &d,
00130 Vec<Num_T> &out);
00131
00133 template<class Num_T>
00134 void elem_mult_inplace(const Vec<Num_T> &a, Vec<Num_T> &b);
00136 template<class Num_T>
00137 Num_T elem_mult_sum(const Vec<Num_T> &a, const Vec<Num_T> &b);
00138
00140 template<class Num_T>
00141 Vec<Num_T> operator/(const Vec<Num_T> &v, Num_T t);
00143 template<class Num_T>
00144 Vec<Num_T> operator/(Num_T t, const Vec<Num_T> &v);
00145
00147 template<class Num_T>
00148 Vec<Num_T> elem_div(const Vec<Num_T> &a, const Vec<Num_T> &b);
00150 template<class Num_T>
00151 Vec<Num_T> elem_div(Num_T t, const Vec<Num_T> &v);
00153 template<class Num_T>
00154 void elem_div_out(const Vec<Num_T> &a, const Vec<Num_T> &b, Vec<Num_T> &out);
00156 template<class Num_T>
00157 Num_T elem_div_sum(const Vec<Num_T> &a, const Vec<Num_T> &b);
00158
00160 template<class Num_T>
00161 Vec<Num_T> concat(const Vec<Num_T> &v, Num_T a);
00163 template<class Num_T>
00164 Vec<Num_T> concat(Num_T a, const Vec<Num_T> &v);
00166 template<class Num_T>
00167 Vec<Num_T> concat(const Vec<Num_T> &v1, const Vec<Num_T> &v2);
00169 template<class Num_T>
00170 Vec<Num_T> concat(const Vec<Num_T> &v1, const Vec<Num_T> &v2,
00171 const Vec<Num_T> &v3);
00173 template<class Num_T>
00174 Vec<Num_T> concat(const Vec<Num_T> &v1, const Vec<Num_T> &v2,
00175 const Vec<Num_T> &v3, const Vec<Num_T> &v4);
00177 template<class Num_T>
00178 Vec<Num_T> concat(const Vec<Num_T> &v1, const Vec<Num_T> &v2,
00179 const Vec<Num_T> &v3, const Vec<Num_T> &v4,
00180 const Vec<Num_T> &v5);
00181
00182
00183
00184
00185
00249 template<class Num_T>
00250 class Vec
00251 {
00252 public:
00254 typedef Num_T value_type;
00255
00257 explicit Vec(const Factory &f = DEFAULT_FACTORY);
00259 explicit Vec(int size, const Factory &f = DEFAULT_FACTORY);
00261 Vec(const Vec<Num_T> &v);
00263 Vec(const Vec<Num_T> &v, const Factory &f);
00265 Vec(const char *str, const Factory &f = DEFAULT_FACTORY);
00267 Vec(const std::string &str, const Factory &f = DEFAULT_FACTORY);
00269 Vec(const Num_T *c_array, int size, const Factory &f = DEFAULT_FACTORY);
00270
00272 ~Vec();
00273
00275 int length() const { return datasize; }
00277 int size() const { return datasize; }
00278
00280 void set_size(int size, bool copy = false);
00282 void set_length(int size, bool copy = false) { set_size(size, copy); }
00284 void zeros();
00286 void clear() { zeros(); }
00288 void ones();
00290 void set(const char *str);
00292 void set(const std::string &str);
00293
00295 const Num_T &operator[](int i) const;
00297 const Num_T &operator()(int i) const;
00299 Num_T &operator[](int i);
00301 Num_T &operator()(int i);
00303 Vec<Num_T> operator()(int i1, int i2) const;
00305 Vec<Num_T> operator()(const Vec<int> &indexlist) const;
00307 Vec<Num_T> operator()(const Vec<bin> &binlist) const;
00308
00310 const Num_T &get(int i) const;
00312 Vec<Num_T> get(int i1, int i2) const;
00314 Vec<Num_T> get(const Vec<int> &indexlist) const;
00316 Vec<Num_T> get(const Vec<bin> &binlist) const;
00317
00319 void set(int i, Num_T t);
00320
00322 Mat<Num_T> transpose() const;
00324 Mat<Num_T> T() const { return this->transpose(); }
00326 Mat<Num_T> hermitian_transpose() const;
00328 Mat<Num_T> H() const { return this->hermitian_transpose(); }
00329
00331 Vec<Num_T>& operator+=(const Vec<Num_T> &v);
00333 Vec<Num_T>& operator+=(Num_T t);
00335 friend Vec<Num_T> operator+<>(const Vec<Num_T> &v1, const Vec<Num_T> &v2);
00337 friend Vec<Num_T> operator+<>(const Vec<Num_T> &v, Num_T t);
00339 friend Vec<Num_T> operator+<>(Num_T t, const Vec<Num_T> &v);
00340
00342 Vec<Num_T>& operator-=(const Vec<Num_T> &v);
00344 Vec<Num_T>& operator-=(Num_T t);
00346 friend Vec<Num_T> operator-<>(const Vec<Num_T> &v1, const Vec<Num_T> &v2);
00348 friend Vec<Num_T> operator-<>(const Vec<Num_T> &v, Num_T t);
00350 friend Vec<Num_T> operator-<>(Num_T t, const Vec<Num_T> &v);
00352 friend Vec<Num_T> operator-<>(const Vec<Num_T> &v);
00353
00355 Vec<Num_T>& operator*=(Num_T t);
00357 friend Num_T operator*<>(const Vec<Num_T> &v1, const Vec<Num_T> &v2);
00359 friend Num_T dot<>(const Vec<Num_T> &v1, const Vec<Num_T> &v2);
00361 friend Mat<Num_T> outer_product<>(const Vec<Num_T> &v1,
00362 const Vec<Num_T> &v2, bool hermitian);
00364 friend Vec<Num_T> operator*<>(const Vec<Num_T> &v, Num_T t);
00366 friend Vec<Num_T> operator*<>(Num_T t, const Vec<Num_T> &v);
00367
00369 friend Vec<Num_T> elem_mult<>(const Vec<Num_T> &a, const Vec<Num_T> &b);
00371 friend Vec<Num_T> elem_mult<>(const Vec<Num_T> &a, const Vec<Num_T> &b,
00372 const Vec<Num_T> &c);
00374 friend Vec<Num_T> elem_mult<>(const Vec<Num_T> &a, const Vec<Num_T> &b,
00375 const Vec<Num_T> &c, const Vec<Num_T> &d);
00376
00378 friend void elem_mult_out<>(const Vec<Num_T> &a, const Vec<Num_T> &b,
00379 Vec<Num_T> &out);
00381 friend void elem_mult_out<>(const Vec<Num_T> &a, const Vec<Num_T> &b,
00382 const Vec<Num_T> &c, Vec<Num_T> &out);
00384 friend void elem_mult_out<>(const Vec<Num_T> &a, const Vec<Num_T> &b,
00385 const Vec<Num_T> &c, const Vec<Num_T> &d,
00386 Vec<Num_T> &out);
00387
00389 friend void elem_mult_inplace<>(const Vec<Num_T> &a, Vec<Num_T> &b);
00391 friend Num_T elem_mult_sum<>(const Vec<Num_T> &a, const Vec<Num_T> &b);
00392
00394 Vec<Num_T>& operator/=(Num_T t);
00396 Vec<Num_T>& operator/=(const Vec<Num_T> &v);
00397
00399 friend Vec<Num_T> operator/<>(const Vec<Num_T> &v, Num_T t);
00401 friend Vec<Num_T> operator/<>(Num_T t, const Vec<Num_T> &v);
00402
00404 friend Vec<Num_T> elem_div<>(const Vec<Num_T> &v1, const Vec<Num_T> &v2);
00406 friend Vec<Num_T> elem_div<>(Num_T t, const Vec<Num_T> &v);
00408 friend void elem_div_out<>(const Vec<Num_T> &v1, const Vec<Num_T> &v2,
00409 Vec<Num_T> &out);
00411 friend Num_T elem_div_sum<>(const Vec<Num_T> &a, const Vec<Num_T> &b);
00412
00414 Vec<Num_T> right(int nr) const;
00416 Vec<Num_T> left(int nr) const;
00418 Vec<Num_T> mid(int start, int nr) const;
00426 Vec<Num_T> split(int pos);
00428 void shift_right(Num_T t, int n = 1);
00430 void shift_right(const Vec<Num_T> &v);
00432 void shift_left(Num_T t, int n = 1);
00434 void shift_left(const Vec<Num_T> &v);
00435
00437 friend Vec<Num_T> concat<>(const Vec<Num_T> &v, Num_T t);
00439 friend Vec<Num_T> concat<>(Num_T t, const Vec<Num_T> &v);
00441 friend Vec<Num_T> concat<>(const Vec<Num_T> &v1, const Vec<Num_T> &v2);
00443 friend Vec<Num_T> concat<>(const Vec<Num_T> &v1, const Vec<Num_T> &v2,
00444 const Vec<Num_T> &v3);
00446 friend Vec<Num_T> concat<>(const Vec<Num_T> &v1, const Vec<Num_T> &v2,
00447 const Vec<Num_T> &v3, const Vec<Num_T> &v4);
00449 friend Vec<Num_T> concat<>(const Vec<Num_T> &v1, const Vec<Num_T> &v2,
00450 const Vec<Num_T> &v3, const Vec<Num_T> &v4,
00451 const Vec<Num_T> &v5);
00452
00454 void set_subvector(int i1, int i2, const Vec<Num_T> &v);
00456 void set_subvector(int i, const Vec<Num_T> &v);
00458 void set_subvector(int i1, int i2, Num_T t);
00460 void replace_mid(int i, const Vec<Num_T> &v);
00462 void del(int i);
00464 void del(int i1, int i2);
00466 void ins(int i, Num_T t);
00468 void ins(int i, const Vec<Num_T> &v);
00469
00471 Vec<Num_T>& operator=(Num_T t);
00473 Vec<Num_T>& operator=(const Vec<Num_T> &v);
00475 Vec<Num_T>& operator=(const Mat<Num_T> &m);
00477 Vec<Num_T>& operator=(const char *str);
00479 Vec<Num_T>& operator=(const std::string &str);
00480
00482 Vec<bin> operator==(Num_T t) const;
00484 Vec<bin> operator!=(Num_T t) const;
00486 Vec<bin> operator<(Num_T t) const;
00488 Vec<bin> operator<=(Num_T t) const;
00490 Vec<bin> operator>(Num_T t) const;
00492 Vec<bin> operator>=(Num_T t) const;
00493
00495 bool operator==(const Vec<Num_T> &v) const;
00497 bool operator!=(const Vec<Num_T> &v) const;
00498
00500 Num_T &_elem(int i) { return data[i]; }
00502 const Num_T &_elem(int i) const { return data[i]; }
00503
00505 Num_T *_data() { return data; }
00507 const Num_T *_data() const { return data; }
00508
00509 protected:
00511 void alloc(int size);
00513 void free();
00514
00516 int datasize;
00518 Num_T *data;
00520 const Factory &factory;
00521 private:
00522
00523 std::vector<std::string> tokenize(const std::string &str,
00524 bool &abc_format) const;
00525
00526 Num_T parse_token(const std::string &s) const;
00527
00528 void parse_abc_token(const std::string &s, Num_T &a, Num_T &b,
00529 Num_T &c) const;
00531 bool in_range(int i) const { return ((i < datasize) && (i >= 0)); }
00532 };
00533
00534
00535
00536
00537
00542 typedef Vec<double> vec;
00543
00548 typedef Vec<std::complex<double> > cvec;
00549
00554 typedef Vec<int> ivec;
00555
00560 typedef Vec<short int> svec;
00561
00566 typedef Vec<bin> bvec;
00567
00568 }
00569
00570
00571 #include <itpp/base/mat.h>
00572
00573 namespace itpp
00574 {
00575
00576
00577
00578
00579
00584 template<class Num_T>
00585 std::ostream &operator<<(std::ostream &os, const Vec<Num_T> &v);
00586
00598 template<class Num_T>
00599 std::istream &operator>>(std::istream &is, Vec<Num_T> &v);
00600
00601
00602
00603
00604
00605 template<class Num_T> inline
00606 void Vec<Num_T>::alloc(int size)
00607 {
00608 if (size > 0) {
00609 create_elements(data, size, factory);
00610 datasize = size;
00611 }
00612 else {
00613 data = 0;
00614 datasize = 0;
00615 }
00616 }
00617
00618 template<class Num_T> inline
00619 void Vec<Num_T>::free()
00620 {
00621 destroy_elements(data, datasize);
00622 datasize = 0;
00623 }
00624
00625
00626 template<class Num_T> inline
00627 Vec<Num_T>::Vec(const Factory &f) : datasize(0), data(0), factory(f) {}
00628
00629 template<class Num_T> inline
00630 Vec<Num_T>::Vec(int size, const Factory &f) : datasize(0), data(0), factory(f)
00631 {
00632 it_assert_debug(size >= 0, "Negative size in Vec::Vec(int)");
00633 alloc(size);
00634 }
00635
00636 template<class Num_T> inline
00637 Vec<Num_T>::Vec(const Vec<Num_T> &v) : datasize(0), data(0), factory(v.factory)
00638 {
00639 alloc(v.datasize);
00640 copy_vector(datasize, v.data, data);
00641 }
00642
00643 template<class Num_T> inline
00644 Vec<Num_T>::Vec(const Vec<Num_T> &v, const Factory &f) : datasize(0), data(0), factory(f)
00645 {
00646 alloc(v.datasize);
00647 copy_vector(datasize, v.data, data);
00648 }
00649
00650 template<class Num_T> inline
00651 Vec<Num_T>::Vec(const char *str, const Factory &f) : datasize(0), data(0), factory(f)
00652 {
00653 set(std::string(str));
00654 }
00655
00656 template<class Num_T> inline
00657 Vec<Num_T>::Vec(const std::string &str, const Factory &f) : datasize(0), data(0), factory(f)
00658 {
00659 set(str);
00660 }
00661
00662 template<class Num_T> inline
00663 Vec<Num_T>::Vec(const Num_T *c_array, int size, const Factory &f) : datasize(0), data(0), factory(f)
00664 {
00665 alloc(size);
00666 copy_vector(size, c_array, data);
00667 }
00668
00669 template<class Num_T> inline
00670 Vec<Num_T>::~Vec()
00671 {
00672 free();
00673 }
00674
00675 template<class Num_T>
00676 void Vec<Num_T>::set_size(int size, bool copy)
00677 {
00678 it_assert_debug(size >= 0, "Vec::set_size(): New size must not be negative");
00679 if (datasize == size)
00680 return;
00681 if (copy) {
00682
00683 Num_T* tmp = data;
00684
00685 int old_datasize = datasize;
00686
00687 int min = datasize < size ? datasize : size;
00688
00689 alloc(size);
00690
00691 copy_vector(min, tmp, data);
00692
00693 for (int i = min; i < size; ++i)
00694 data[i] = Num_T(0);
00695
00696 destroy_elements(tmp, old_datasize);
00697 }
00698 else {
00699 free();
00700 alloc(size);
00701 }
00702 }
00703
00704 template<class Num_T> inline
00705 const Num_T& Vec<Num_T>::operator[](int i) const
00706 {
00707 it_assert_debug(in_range(i), "Vec<>::operator[]: Index out of range");
00708 return data[i];
00709 }
00710
00711 template<class Num_T> inline
00712 const Num_T& Vec<Num_T>::operator()(int i) const
00713 {
00714 return (*this)[i];
00715 }
00716
00717 template<class Num_T> inline
00718 Num_T& Vec<Num_T>::operator[](int i)
00719 {
00720 it_assert_debug(in_range(i), "Vec<>::operator[]: Index out of range");
00721 return data[i];
00722 }
00723
00724 template<class Num_T> inline
00725 Num_T& Vec<Num_T>::operator()(int i)
00726 {
00727 return (*this)[i];
00728 }
00729
00730 template<class Num_T> inline
00731 Vec<Num_T> Vec<Num_T>::operator()(int i1, int i2) const
00732 {
00733 if (i1 == -1) i1 = datasize - 1;
00734 if (i2 == -1) i2 = datasize - 1;
00735
00736 it_assert_debug((i1 >= 0) && (i1 <= i2) && (i2 < datasize),
00737 "Vec<>::operator()(i1, i2): Indexing out of range");
00738
00739 Vec<Num_T> s(i2 - i1 + 1);
00740 copy_vector(s.datasize, data + i1, s.data);
00741
00742 return s;
00743 }
00744
00745 template<class Num_T>
00746 Vec<Num_T> Vec<Num_T>::operator()(const Vec<int> &indexlist) const
00747 {
00748 int size = indexlist.size();
00749 Vec<Num_T> temp(size);
00750 for (int i = 0; i < size; ++i) {
00751 it_assert_debug(in_range(indexlist(i)), "Vec<>::operator()(ivec &): "
00752 "Index i=" << i << " out of range");
00753 temp(i) = data[indexlist(i)];
00754 }
00755 return temp;
00756 }
00757
00758 template<class Num_T>
00759 Vec<Num_T> Vec<Num_T>::operator()(const Vec<bin> &binlist) const
00760 {
00761 int size = binlist.size();
00762 it_assert_debug(datasize == size, "Vec<>::operator()(bvec &): "
00763 "Wrong size of binlist vector");
00764 Vec<Num_T> temp(size);
00765 int j = 0;
00766 for (int i = 0; i < size; ++i)
00767 if (binlist(i) == bin(1))
00768 temp(j++) = data[i];
00769 temp.set_size(j, true);
00770 return temp;
00771 }
00772
00773
00774 template<class Num_T> inline
00775 const Num_T& Vec<Num_T>::get(int i) const
00776 {
00777 return (*this)[i];
00778 }
00779
00780 template<class Num_T> inline
00781 Vec<Num_T> Vec<Num_T>::get(int i1, int i2) const
00782 {
00783 return (*this)(i1, i2);
00784 }
00785
00786 template<class Num_T> inline
00787 Vec<Num_T> Vec<Num_T>::get(const Vec<int> &indexlist) const
00788 {
00789 return (*this)(indexlist);
00790 }
00791
00792 template<class Num_T> inline
00793 Vec<Num_T> Vec<Num_T>::get(const Vec<bin> &binlist) const
00794 {
00795 return (*this)(binlist);
00796 }
00797
00798 template<class Num_T> inline
00799 void Vec<Num_T>::zeros()
00800 {
00801 for (int i = 0; i < datasize; i++)
00802 data[i] = Num_T(0);
00803 }
00804
00805 template<class Num_T> inline
00806 void Vec<Num_T>::ones()
00807 {
00808 for (int i = 0; i < datasize; i++)
00809 data[i] = Num_T(1);
00810 }
00811
00812 template<class Num_T> inline
00813 void Vec<Num_T>::set(int i, Num_T t)
00814 {
00815 it_assert_debug(in_range(i), "Vec<>::set(i, t): Index out of range");
00816 data[i] = t;
00817 }
00818
00820 template<>
00821 void Vec<double>::set(const std::string &str);
00822 template<>
00823 void Vec<std::complex<double> >::set(const std::string &str);
00824 template<>
00825 void Vec<int>::set(const std::string &str);
00826 template<>
00827 void Vec<short int>::set(const std::string &str);
00828 template<>
00829 void Vec<bin>::set(const std::string &str);
00831
00832 template<class Num_T>
00833 void Vec<Num_T>::set(const std::string &str)
00834 {
00835 it_error("Vec::set(): Only `double', `complex<double>', `int', "
00836 "`short int' and `bin' types supported");
00837 }
00838
00839 template<class Num_T> inline
00840 void Vec<Num_T>::set(const char *str)
00841 {
00842 set(std::string(str));
00843 }
00844
00845
00846 template<class Num_T>
00847 Mat<Num_T> Vec<Num_T>::transpose() const
00848 {
00849 Mat<Num_T> temp(1, datasize);
00850 copy_vector(datasize, data, temp._data());
00851 return temp;
00852 }
00853
00855 template<>
00856 Mat<std::complex<double> > Vec<std::complex<double> >::hermitian_transpose() const;
00858
00859 template<class Num_T>
00860 Mat<Num_T> Vec<Num_T>::hermitian_transpose() const
00861 {
00862 Mat<Num_T> temp(1, datasize);
00863 copy_vector(datasize, data, temp._data());
00864 return temp;
00865 }
00866
00867 template<class Num_T>
00868 Vec<Num_T>& Vec<Num_T>::operator+=(const Vec<Num_T> &v)
00869 {
00870 if (datasize == 0) {
00871 if (this != &v) {
00872 alloc(v.datasize);
00873 copy_vector(datasize, v.data, data);
00874 }
00875 }
00876 else {
00877 it_assert_debug(datasize == v.datasize, "Vec::operator+=: Wrong sizes");
00878 for (int i = 0; i < datasize; i++)
00879 data[i] += v.data[i];
00880 }
00881 return *this;
00882 }
00883
00884 template<class Num_T> inline
00885 Vec<Num_T>& Vec<Num_T>::operator+=(Num_T t)
00886 {
00887 for (int i = 0;i < datasize;i++)
00888 data[i] += t;
00889 return *this;
00890 }
00891
00892 template<class Num_T>
00893 Vec<Num_T> operator+(const Vec<Num_T> &v1, const Vec<Num_T> &v2)
00894 {
00895 int i;
00896 Vec<Num_T> r(v1.datasize);
00897
00898 it_assert_debug(v1.datasize == v2.datasize, "Vec::operator+: wrong sizes");
00899 for (i = 0; i < v1.datasize; i++)
00900 r.data[i] = v1.data[i] + v2.data[i];
00901
00902 return r;
00903 }
00904
00905 template<class Num_T>
00906 Vec<Num_T> operator+(const Vec<Num_T> &v, Num_T t)
00907 {
00908 int i;
00909 Vec<Num_T> r(v.datasize);
00910
00911 for (i = 0; i < v.datasize; i++)
00912 r.data[i] = v.data[i] + t;
00913
00914 return r;
00915 }
00916
00917 template<class Num_T>
00918 Vec<Num_T> operator+(Num_T t, const Vec<Num_T> &v)
00919 {
00920 int i;
00921 Vec<Num_T> r(v.datasize);
00922
00923 for (i = 0; i < v.datasize; i++)
00924 r.data[i] = t + v.data[i];
00925
00926 return r;
00927 }
00928
00929 template<class Num_T>
00930 Vec<Num_T>& Vec<Num_T>::operator-=(const Vec<Num_T> &v)
00931 {
00932 if (datasize == 0) {
00933 if (this != &v) {
00934 alloc(v.datasize);
00935 for (int i = 0; i < v.datasize; i++)
00936 data[i] = -v.data[i];
00937 }
00938 }
00939 else {
00940 it_assert_debug(datasize == v.datasize, "Vec::operator-=: Wrong sizes");
00941 for (int i = 0; i < datasize; i++)
00942 data[i] -= v.data[i];
00943 }
00944 return *this;
00945 }
00946
00947 template<class Num_T> inline
00948 Vec<Num_T>& Vec<Num_T>::operator-=(Num_T t)
00949 {
00950 for (int i = 0;i < datasize;i++)
00951 data[i] -= t;
00952 return *this;
00953 }
00954
00955 template<class Num_T>
00956 Vec<Num_T> operator-(const Vec<Num_T> &v1, const Vec<Num_T> &v2)
00957 {
00958 int i;
00959 Vec<Num_T> r(v1.datasize);
00960
00961 it_assert_debug(v1.datasize == v2.datasize, "Vec::operator-: wrong sizes");
00962 for (i = 0; i < v1.datasize; i++)
00963 r.data[i] = v1.data[i] - v2.data[i];
00964
00965 return r;
00966 }
00967
00968 template<class Num_T>
00969 Vec<Num_T> operator-(const Vec<Num_T> &v, Num_T t)
00970 {
00971 int i;
00972 Vec<Num_T> r(v.datasize);
00973
00974 for (i = 0; i < v.datasize; i++)
00975 r.data[i] = v.data[i] - t;
00976
00977 return r;
00978 }
00979
00980 template<class Num_T>
00981 Vec<Num_T> operator-(Num_T t, const Vec<Num_T> &v)
00982 {
00983 int i;
00984 Vec<Num_T> r(v.datasize);
00985
00986 for (i = 0; i < v.datasize; i++)
00987 r.data[i] = t - v.data[i];
00988
00989 return r;
00990 }
00991
00992 template<class Num_T>
00993 Vec<Num_T> operator-(const Vec<Num_T> &v)
00994 {
00995 int i;
00996 Vec<Num_T> r(v.datasize);
00997
00998 for (i = 0; i < v.datasize; i++)
00999 r.data[i] = -v.data[i];
01000
01001 return r;
01002 }
01003
01004 template<class Num_T> inline
01005 Vec<Num_T>& Vec<Num_T>::operator*=(Num_T t)
01006 {
01007 scal_vector(datasize, t, data);
01008 return *this;
01009 }
01010
01011 #if defined(HAVE_BLAS)
01012 template<> inline
01013 double dot(const vec &v1, const vec &v2)
01014 {
01015 it_assert_debug(v1.datasize == v2.datasize, "vec::dot: wrong sizes");
01016 int incr = 1;
01017 return blas::ddot_(&v1.datasize, v1.data, &incr, v2.data, &incr);
01018 }
01019
01020 #if defined(HAVE_ZDOTUSUB) || defined(HAVE_ZDOTU_VOID)
01021 template<> inline
01022 std::complex<double> dot(const cvec &v1, const cvec &v2)
01023 {
01024 it_assert_debug(v1.datasize == v2.datasize, "cvec::dot: wrong sizes");
01025 int incr = 1;
01026 std::complex<double> output;
01027 blas::zdotusub_(&output, &v1.datasize, v1.data, &incr, v2.data, &incr);
01028 return output;
01029 }
01030 #endif // HAVE_ZDOTUSUB || HAVE_ZDOTU_VOID
01031 #endif // HAVE_BLAS
01032
01033 template<class Num_T>
01034 Num_T dot(const Vec<Num_T> &v1, const Vec<Num_T> &v2)
01035 {
01036 int i;
01037 Num_T r = Num_T(0);
01038
01039 it_assert_debug(v1.datasize == v2.datasize, "Vec::dot: wrong sizes");
01040 for (i = 0; i < v1.datasize; i++)
01041 r += v1.data[i] * v2.data[i];
01042
01043 return r;
01044 }
01045
01046 template<class Num_T> inline
01047 Num_T operator*(const Vec<Num_T> &v1, const Vec<Num_T> &v2)
01048 {
01049 return dot(v1, v2);
01050 }
01051
01052
01053 #if defined(HAVE_BLAS)
01054 template<> inline
01055 mat outer_product(const vec &v1, const vec &v2, bool)
01056 {
01057 it_assert_debug((v1.datasize > 0) && (v2.datasize > 0),
01058 "Vec::outer_product():: Input vector of zero size");
01059
01060 mat out(v1.datasize, v2.datasize);
01061 out.zeros();
01062 double alpha = 1.0;
01063 int incr = 1;
01064 blas::dger_(&v1.datasize, &v2.datasize, &alpha, v1.data, &incr,
01065 v2.data, &incr, out._data(), &v1.datasize);
01066 return out;
01067 }
01068
01069 template<> inline
01070 cmat outer_product(const cvec &v1, const cvec &v2, bool hermitian)
01071 {
01072 it_assert_debug((v1.datasize > 0) && (v2.datasize > 0),
01073 "Vec::outer_product():: Input vector of zero size");
01074
01075 cmat out(v1.datasize, v2.datasize);
01076 out.zeros();
01077 std::complex<double> alpha(1.0);
01078 int incr = 1;
01079 if (hermitian) {
01080 blas::zgerc_(&v1.datasize, &v2.datasize, &alpha, v1.data, &incr,
01081 v2.data, &incr, out._data(), &v1.datasize);
01082 }
01083 else {
01084 blas::zgeru_(&v1.datasize, &v2.datasize, &alpha, v1.data, &incr,
01085 v2.data, &incr, out._data(), &v1.datasize);
01086 }
01087 return out;
01088 }
01089 #else
01091 template<> inline
01092 cmat outer_product(const cvec &v1, const cvec &v2, bool hermitian)
01093 {
01094 it_assert_debug((v1.datasize > 0) && (v2.datasize > 0),
01095 "Vec::outer_product():: Input vector of zero size");
01096
01097 cmat out(v1.datasize, v2.datasize);
01098 if (hermitian) {
01099 for (int i = 0; i < v1.datasize; ++i) {
01100 for (int j = 0; j < v2.datasize; ++j) {
01101 out(i, j) = v1.data[i] * conj(v2.data[j]);
01102 }
01103 }
01104 }
01105 else {
01106 for (int i = 0; i < v1.datasize; ++i) {
01107 for (int j = 0; j < v2.datasize; ++j) {
01108 out(i, j) = v1.data[i] * v2.data[j];
01109 }
01110 }
01111 }
01112 return out;
01113 }
01114 #endif // HAVE_BLAS
01115
01116 template<class Num_T>
01117 Mat<Num_T> outer_product(const Vec<Num_T> &v1, const Vec<Num_T> &v2, bool)
01118 {
01119 int i, j;
01120
01121 it_assert_debug((v1.datasize > 0) && (v2.datasize > 0),
01122 "Vec::outer_product:: Input vector of zero size");
01123
01124 Mat<Num_T> r(v1.datasize, v2.datasize);
01125 for (i = 0; i < v1.datasize; i++) {
01126 for (j = 0; j < v2.datasize; j++) {
01127 r(i, j) = v1.data[i] * v2.data[j];
01128 }
01129 }
01130
01131 return r;
01132 }
01133
01134 template<class Num_T>
01135 Vec<Num_T> operator*(const Vec<Num_T> &v, Num_T t)
01136 {
01137 int i;
01138 Vec<Num_T> r(v.datasize);
01139 for (i = 0; i < v.datasize; i++)
01140 r.data[i] = v.data[i] * t;
01141
01142 return r;
01143 }
01144
01145 template<class Num_T> inline
01146 Vec<Num_T> operator*(Num_T t, const Vec<Num_T> &v)
01147 {
01148 return operator*(v, t);
01149 }
01150
01151 template<class Num_T> inline
01152 Vec<Num_T> elem_mult(const Vec<Num_T> &a, const Vec<Num_T> &b)
01153 {
01154 Vec<Num_T> out;
01155 elem_mult_out(a, b, out);
01156 return out;
01157 }
01158
01159 template<class Num_T> inline
01160 Vec<Num_T> elem_mult(const Vec<Num_T> &a, const Vec<Num_T> &b,
01161 const Vec<Num_T> &c)
01162 {
01163 Vec<Num_T> out;
01164 elem_mult_out(a, b, c, out);
01165 return out;
01166 }
01167
01168 template<class Num_T> inline
01169 Vec<Num_T> elem_mult(const Vec<Num_T> &a, const Vec<Num_T> &b,
01170 const Vec<Num_T> &c, const Vec<Num_T> &d)
01171 {
01172 Vec<Num_T> out;
01173 elem_mult_out(a, b, c, d, out);
01174 return out;
01175 }
01176
01177 template<class Num_T>
01178 void elem_mult_out(const Vec<Num_T> &a, const Vec<Num_T> &b, Vec<Num_T> &out)
01179 {
01180 it_assert_debug(a.datasize == b.datasize,
01181 "Vec<>::elem_mult_out(): Wrong sizes");
01182 out.set_size(a.datasize);
01183 for (int i = 0; i < a.datasize; i++)
01184 out.data[i] = a.data[i] * b.data[i];
01185 }
01186
01187 template<class Num_T>
01188 void elem_mult_out(const Vec<Num_T> &a, const Vec<Num_T> &b,
01189 const Vec<Num_T> &c, Vec<Num_T> &out)
01190 {
01191 it_assert_debug((a.datasize == b.datasize) && (a.datasize == c.datasize),
01192 "Vec<>::elem_mult_out(): Wrong sizes");
01193 out.set_size(a.datasize);
01194 for (int i = 0; i < a.datasize; i++)
01195 out.data[i] = a.data[i] * b.data[i] * c.data[i];
01196 }
01197
01198 template<class Num_T>
01199 void elem_mult_out(const Vec<Num_T> &a, const Vec<Num_T> &b,
01200 const Vec<Num_T> &c, const Vec<Num_T> &d, Vec<Num_T> &out)
01201 {
01202 it_assert_debug((a.datasize == b.datasize) && (a.datasize == c.datasize)
01203 && (a.datasize == d.datasize),
01204 "Vec<>::elem_mult_out(): Wrong sizes");
01205 out.set_size(a.datasize);
01206 for (int i = 0; i < a.datasize; i++)
01207 out.data[i] = a.data[i] * b.data[i] * c.data[i] * d.data[i];
01208 }
01209
01210 template<class Num_T>
01211 #ifndef _MSC_VER
01212 inline
01213 #endif
01214 void elem_mult_inplace(const Vec<Num_T> &a, Vec<Num_T> &b)
01215 {
01216 it_assert_debug(a.datasize == b.datasize,
01217 "Vec<>::elem_mult_inplace(): Wrong sizes");
01218 for (int i = 0; i < a.datasize; i++)
01219 b.data[i] *= a.data[i];
01220 }
01221
01222 template<class Num_T> inline
01223 Num_T elem_mult_sum(const Vec<Num_T> &a, const Vec<Num_T> &b)
01224 {
01225 it_assert_debug(a.datasize == b.datasize,
01226 "Vec<>::elem_mult_sum(): Wrong sizes");
01227 Num_T acc = 0;
01228 for (int i = 0; i < a.datasize; i++)
01229 acc += a.data[i] * b.data[i];
01230 return acc;
01231 }
01232
01233 template<class Num_T>
01234 Vec<Num_T> operator/(const Vec<Num_T> &v, Num_T t)
01235 {
01236 int i;
01237 Vec<Num_T> r(v.datasize);
01238
01239 for (i = 0; i < v.datasize; i++)
01240 r.data[i] = v.data[i] / t;
01241
01242 return r;
01243 }
01244
01245 template<class Num_T>
01246 Vec<Num_T> operator/(Num_T t, const Vec<Num_T> &v)
01247 {
01248 int i;
01249 Vec<Num_T> r(v.datasize);
01250
01251 for (i = 0; i < v.datasize; i++)
01252 r.data[i] = t / v.data[i];
01253
01254 return r;
01255 }
01256
01257 template<class Num_T>
01258 Vec<Num_T> elem_div(Num_T t, const Vec<Num_T> &v)
01259 {
01260 it_warning("Vec<>::elem_div(Num_T, const Vec<Num_T> &): This function is "
01261 "deprecated and might be removed from future IT++ releases. "
01262 "Please use Vec<>::operator/(Num_T, const Vec<Num_T> &) "
01263 "instead.");
01264 return operator/(t, v);
01265 }
01266
01267 template<class Num_T> inline
01268 Vec<Num_T>& Vec<Num_T>::operator/=(Num_T t)
01269 {
01270 for (int i = 0; i < datasize; ++i) {
01271 data[i] /= t;
01272 }
01273 return *this;
01274 }
01275
01276 template<class Num_T> inline
01277 Vec<Num_T>& Vec<Num_T>::operator/=(const Vec<Num_T> &v)
01278 {
01279 it_assert_debug(datasize == v.datasize, "Vec::operator/=(): wrong sizes");
01280 for (int i = 0; i < datasize; ++i) {
01281 data[i] /= v.data[i];
01282 }
01283 return *this;
01284 }
01285
01286 template<class Num_T> inline
01287 Vec<Num_T> elem_div(const Vec<Num_T> &a, const Vec<Num_T> &b)
01288 {
01289 Vec<Num_T> out;
01290 elem_div_out(a, b, out);
01291 return out;
01292 }
01293
01294 template<class Num_T>
01295 void elem_div_out(const Vec<Num_T> &a, const Vec<Num_T> &b, Vec<Num_T> &out)
01296 {
01297 it_assert_debug(a.datasize == b.datasize, "Vecelem_div_out: wrong sizes");
01298
01299 out.set_size(a.size());
01300
01301 for (int i = 0; i < a.datasize; i++)
01302 out.data[i] = a.data[i] / b.data[i];
01303 }
01304
01305 template<class Num_T> inline
01306 Num_T elem_div_sum(const Vec<Num_T> &a, const Vec<Num_T> &b)
01307 {
01308 it_assert_debug(a.datasize == b.datasize, "Vec::elem_div_sum: wrong sizes");
01309
01310 Num_T acc = 0;
01311
01312 for (int i = 0; i < a.datasize; i++)
01313 acc += a.data[i] / b.data[i];
01314
01315 return acc;
01316 }
01317
01318 template<class Num_T>
01319 Vec<Num_T> Vec<Num_T>::right(int nr) const
01320 {
01321 it_assert_debug(nr <= datasize, "Vec::right(): index out of range");
01322 Vec<Num_T> temp(nr);
01323 if (nr > 0) {
01324 copy_vector(nr, &data[datasize-nr], temp.data);
01325 }
01326 return temp;
01327 }
01328
01329 template<class Num_T>
01330 Vec<Num_T> Vec<Num_T>::left(int nr) const
01331 {
01332 it_assert_debug(nr <= datasize, "Vec::left(): index out of range");
01333 Vec<Num_T> temp(nr);
01334 if (nr > 0) {
01335 copy_vector(nr, data, temp.data);
01336 }
01337 return temp;
01338 }
01339
01340 template<class Num_T>
01341 Vec<Num_T> Vec<Num_T>::mid(int start, int nr) const
01342 {
01343 it_assert_debug((start >= 0) && ((start + nr) <= datasize),
01344 "Vec::mid(): indexing out of range");
01345 Vec<Num_T> temp(nr);
01346 if (nr > 0) {
01347 copy_vector(nr, &data[start], temp.data);
01348 }
01349 return temp;
01350 }
01351
01352 template<class Num_T>
01353 Vec<Num_T> Vec<Num_T>::split(int pos)
01354 {
01355 it_assert_debug((pos >= 0) && (pos <= datasize),
01356 "Vec<>::split(): Index out of range");
01357 Vec<Num_T> temp1(pos);
01358 if (pos > 0) {
01359 copy_vector(pos, data, temp1.data);
01360 if (pos < datasize) {
01361 Vec<Num_T> temp2(datasize - pos);
01362 copy_vector(datasize - pos, &data[pos], temp2.data);
01363 (*this) = temp2;
01364 }
01365 else {
01366 set_size(0);
01367 }
01368 }
01369 return temp1;
01370 }
01371
01372 template<class Num_T>
01373 void Vec<Num_T>::shift_right(Num_T t, int n)
01374 {
01375 int i = datasize;
01376
01377 it_assert_debug(n >= 0, "Vec::shift_right: index out of range");
01378 while (--i >= n)
01379 data[i] = data[i-n];
01380 while (i >= 0)
01381 data[i--] = t;
01382 }
01383
01384 template<class Num_T>
01385 void Vec<Num_T>::shift_right(const Vec<Num_T> &v)
01386 {
01387 for (int i = datasize - 1; i >= v.datasize; i--)
01388 data[i] = data[i-v.datasize];
01389 for (int i = 0; i < v.datasize; i++)
01390 data[i] = v[i];
01391 }
01392
01393 template<class Num_T>
01394 void Vec<Num_T>::shift_left(Num_T t, int n)
01395 {
01396 int i;
01397
01398 it_assert_debug(n >= 0, "Vec::shift_left: index out of range");
01399 for (i = 0; i < datasize - n; i++)
01400 data[i] = data[i+n];
01401 while (i < datasize)
01402 data[i++] = t;
01403 }
01404
01405 template<class Num_T>
01406 void Vec<Num_T>::shift_left(const Vec<Num_T> &v)
01407 {
01408 for (int i = 0; i < datasize - v.datasize; i++)
01409 data[i] = data[i+v.datasize];
01410 for (int i = datasize - v.datasize; i < datasize; i++)
01411 data[i] = v[i-datasize+v.datasize];
01412 }
01413
01414 template<class Num_T>
01415 Vec<Num_T> concat(const Vec<Num_T> &v, Num_T t)
01416 {
01417 int size = v.size();
01418 Vec<Num_T> temp(size + 1);
01419 copy_vector(size, v.data, temp.data);
01420 temp(size) = t;
01421 return temp;
01422 }
01423
01424 template<class Num_T>
01425 Vec<Num_T> concat(Num_T t, const Vec<Num_T> &v)
01426 {
01427 int size = v.size();
01428 Vec<Num_T> temp(size + 1);
01429 temp(0) = t;
01430 copy_vector(size, v.data, &temp.data[1]);
01431 return temp;
01432 }
01433
01434 template<class Num_T>
01435 Vec<Num_T> concat(const Vec<Num_T> &v1, const Vec<Num_T> &v2)
01436 {
01437 int size1 = v1.size();
01438 int size2 = v2.size();
01439 Vec<Num_T> temp(size1 + size2);
01440 copy_vector(size1, v1.data, temp.data);
01441 copy_vector(size2, v2.data, &temp.data[size1]);
01442 return temp;
01443 }
01444
01445 template<class Num_T>
01446 Vec<Num_T> concat(const Vec<Num_T> &v1, const Vec<Num_T> &v2,
01447 const Vec<Num_T> &v3)
01448 {
01449 int size1 = v1.size();
01450 int size2 = v2.size();
01451 int size3 = v3.size();
01452 Vec<Num_T> temp(size1 + size2 + size3);
01453 copy_vector(size1, v1.data, temp.data);
01454 copy_vector(size2, v2.data, &temp.data[size1]);
01455 copy_vector(size3, v3.data, &temp.data[size1+size2]);
01456 return temp;
01457 }
01458
01459 template<class Num_T>
01460 Vec<Num_T> concat(const Vec<Num_T> &v1, const Vec<Num_T> &v2,
01461 const Vec<Num_T> &v3, const Vec<Num_T> &v4)
01462 {
01463 int size1 = v1.size();
01464 int size2 = v2.size();
01465 int size3 = v3.size();
01466 int size4 = v4.size();
01467 Vec<Num_T> temp(size1 + size2 + size3 + size4);
01468 copy_vector(size1, v1.data, temp.data);
01469 copy_vector(size2, v2.data, &temp.data[size1]);
01470 copy_vector(size3, v3.data, &temp.data[size1+size2]);
01471 copy_vector(size4, v4.data, &temp.data[size1+size2+size3]);
01472 return temp;
01473 }
01474
01475 template<class Num_T>
01476 Vec<Num_T> concat(const Vec<Num_T> &v1, const Vec<Num_T> &v2,
01477 const Vec<Num_T> &v3, const Vec<Num_T> &v4,
01478 const Vec<Num_T> &v5)
01479 {
01480 int size1 = v1.size();
01481 int size2 = v2.size();
01482 int size3 = v3.size();
01483 int size4 = v4.size();
01484 int size5 = v5.size();
01485 Vec<Num_T> temp(size1 + size2 + size3 + size4 + size5);
01486 copy_vector(size1, v1.data, temp.data);
01487 copy_vector(size2, v2.data, &temp.data[size1]);
01488 copy_vector(size3, v3.data, &temp.data[size1+size2]);
01489 copy_vector(size4, v4.data, &temp.data[size1+size2+size3]);
01490 copy_vector(size5, v5.data, &temp.data[size1+size2+size3+size4]);
01491 return temp;
01492 }
01493
01494 template<class Num_T>
01495 void Vec<Num_T>::set_subvector(int i1, int, const Vec<Num_T> &v)
01496 {
01497 it_warning("Vec<>::set_subvector(int, int, const Vec<> &): This function "
01498 "is deprecated and might be removed from future IT++ releases. "
01499 "Please use Vec<>::set_subvector(int, const Vec<> &) instead.");
01500 set_subvector(i1, v);
01501 }
01502
01503 template<class Num_T> inline
01504 void Vec<Num_T>:: set_subvector(int i, const Vec<Num_T> &v)
01505 {
01506 it_assert_debug((i >= 0) && (i + v.datasize <= datasize),
01507 "Vec<>::set_subvector(int, const Vec<> &): "
01508 "Index out of range or too long input vector");
01509 copy_vector(v.datasize, v.data, data + i);
01510 }
01511
01512 template<class Num_T> inline
01513 void Vec<Num_T>::set_subvector(int i1, int i2, Num_T t)
01514 {
01515 if (i1 == -1) i1 = datasize - 1;
01516 if (i2 == -1) i2 = datasize - 1;
01517 it_assert_debug((i1 >= 0) && (i1 <= i2) && (i2 < datasize),
01518 "Vec<>::set_subvector(int, int, Num_T): Indexing out "
01519 "of range");
01520 for (int i = i1; i <= i2; i++)
01521 data[i] = t;
01522 }
01523
01524 template<class Num_T> inline
01525 void Vec<Num_T>::replace_mid(int i, const Vec<Num_T> &v)
01526 {
01527 set_subvector(i, v);
01528 }
01529
01530 template<class Num_T>
01531 void Vec<Num_T>::del(int index)
01532 {
01533 it_assert_debug(in_range(index), "Vec<>::del(int): Index out of range");
01534 Vec<Num_T> temp(*this);
01535 set_size(datasize - 1, false);
01536 copy_vector(index, temp.data, data);
01537 copy_vector(datasize - index, &temp.data[index+1], &data[index]);
01538 }
01539
01540 template<class Num_T>
01541 void Vec<Num_T>::del(int i1, int i2)
01542 {
01543 if (i1 == -1) i1 = datasize - 1;
01544 if (i2 == -1) i2 = datasize - 1;
01545 it_assert_debug((i1 >= 0) && (i1 <= i2) && (i2 < datasize),
01546 "Vec<>::del(int, int): Indexing out of range");
01547 Vec<Num_T> temp(*this);
01548 int new_size = datasize - (i2 - i1 + 1);
01549 set_size(new_size, false);
01550 copy_vector(i1, temp.data, data);
01551 copy_vector(datasize - i1, &temp.data[i2+1], &data[i1]);
01552 }
01553
01554 template<class Num_T>
01555 void Vec<Num_T>::ins(int index, const Num_T t)
01556 {
01557 it_assert_debug((index >= 0) && (index <= datasize),
01558 "Vec<>::ins(): Index out of range");
01559 Vec<Num_T> Temp(*this);
01560
01561 set_size(datasize + 1, false);
01562 copy_vector(index, Temp.data, data);
01563 data[index] = t;
01564 copy_vector(Temp.datasize - index, Temp.data + index, data + index + 1);
01565 }
01566
01567 template<class Num_T>
01568 void Vec<Num_T>::ins(int index, const Vec<Num_T> &v)
01569 {
01570 it_assert_debug((index >= 0) && (index <= datasize),
01571 "Vec<>::ins(): Index out of range");
01572 Vec<Num_T> Temp(*this);
01573
01574 set_size(datasize + v.length(), false);
01575 copy_vector(index, Temp.data, data);
01576 copy_vector(v.size(), v.data, &data[index]);
01577 copy_vector(Temp.datasize - index, Temp.data + index, data + index + v.size());
01578 }
01579
01580 template<class Num_T> inline
01581 Vec<Num_T>& Vec<Num_T>::operator=(Num_T t)
01582 {
01583 for (int i = 0;i < datasize;i++)
01584 data[i] = t;
01585 return *this;
01586 }
01587
01588 template<class Num_T> inline
01589 Vec<Num_T>& Vec<Num_T>::operator=(const Vec<Num_T> &v)
01590 {
01591 if (this != &v) {
01592 set_size(v.datasize, false);
01593 copy_vector(datasize, v.data, data);
01594 }
01595 return *this;
01596 }
01597
01598 template<class Num_T>
01599 Vec<Num_T>& Vec<Num_T>::operator=(const Mat<Num_T> &m)
01600 {
01601 if (m.cols() == 1) {
01602 set_size(m.rows(), false);
01603 copy_vector(m.rows(), m._data(), data);
01604 }
01605 else if (m.rows() == 1) {
01606 set_size(m.cols(), false);
01607 copy_vector(m.cols(), m._data(), m.rows(), data, 1);
01608 }
01609 else
01610 it_error("Vec<>::operator=(Mat<Num_T> &): Wrong size of input matrix");
01611 return *this;
01612 }
01613
01614 template<class Num_T> inline
01615 Vec<Num_T>& Vec<Num_T>::operator=(const char *str)
01616 {
01617 set(std::string(str));
01618 return *this;
01619 }
01620
01621 template<class Num_T> inline
01622 Vec<Num_T>& Vec<Num_T>::operator=(const std::string &str)
01623 {
01624 set(str);
01625 return *this;
01626 }
01627
01629 template<>
01630 bvec Vec<std::complex<double> >::operator==(std::complex<double>) const;
01632
01633 template<class Num_T>
01634 bvec Vec<Num_T>::operator==(Num_T t) const
01635 {
01636 it_assert_debug(datasize > 0, "Vec<>::operator==(): Wrong size");
01637 bvec temp(datasize);
01638 for (int i = 0; i < datasize; i++)
01639 temp(i) = (data[i] == t);
01640 return temp;
01641 }
01642
01644 template<>
01645 bvec Vec<std::complex<double> >::operator!=(std::complex<double>) const;
01647
01648 template<class Num_T>
01649 bvec Vec<Num_T>::operator!=(Num_T t) const
01650 {
01651 it_assert_debug(datasize > 0, "Vec<>::operator!=(): Wrong size");
01652 bvec temp(datasize);
01653 for (int i = 0; i < datasize; i++)
01654 temp(i) = (data[i] != t);
01655 return temp;
01656 }
01657
01659 template<>
01660 bvec Vec<std::complex<double> >::operator<(std::complex<double>) const;
01662
01663 template<class Num_T>
01664 bvec Vec<Num_T>::operator<(Num_T t) const
01665 {
01666 it_assert_debug(datasize > 0, "Vec<>::operator<(): Wrong size");
01667 bvec temp(datasize);
01668 for (int i = 0; i < datasize; i++)
01669 temp(i) = (data[i] < t);
01670 return temp;
01671 }
01672
01674 template<>
01675 bvec Vec<std::complex<double> >::operator<=(std::complex<double>) const;
01677
01678 template<class Num_T>
01679 bvec Vec<Num_T>::operator<=(Num_T t) const
01680 {
01681 it_assert_debug(datasize > 0, "Vec<>::operator<=(): Wrong size");
01682 bvec temp(datasize);
01683 for (int i = 0; i < datasize; i++)
01684 temp(i) = (data[i] <= t);
01685 return temp;
01686 }
01687
01689 template<>
01690 bvec Vec<std::complex<double> >::operator>(std::complex<double>) const;
01692
01693 template<class Num_T>
01694 bvec Vec<Num_T>::operator>(Num_T t) const
01695 {
01696 it_assert_debug(datasize > 0, "Vec<>::operator>(): Wrong size");
01697 bvec temp(datasize);
01698 for (int i = 0; i < datasize; i++)
01699 temp(i) = (data[i] > t);
01700 return temp;
01701 }
01702
01704 template<>
01705 bvec Vec<std::complex<double> >::operator>=(std::complex<double>) const;
01707
01708 template<class Num_T>
01709 bvec Vec<Num_T>::operator>=(Num_T t) const
01710 {
01711 it_assert_debug(datasize > 0, "Vec<>::operator>=(): Wrong size");
01712 bvec temp(datasize);
01713 for (int i = 0; i < datasize; i++)
01714 temp(i) = (data[i] >= t);
01715 return temp;
01716 }
01717
01718 template<class Num_T>
01719 bool Vec<Num_T>::operator==(const Vec<Num_T> &invector) const
01720 {
01721
01722 if (datasize != invector.datasize) return false;
01723 for (int i = 0;i < datasize;i++) {
01724 if (data[i] != invector.data[i]) return false;
01725 }
01726 return true;
01727 }
01728
01729 template<class Num_T>
01730 bool Vec<Num_T>::operator!=(const Vec<Num_T> &invector) const
01731 {
01732 if (datasize != invector.datasize) return true;
01733 for (int i = 0;i < datasize;i++) {
01734 if (data[i] != invector.data[i]) return true;
01735 }
01736 return false;
01737 }
01738
01740 template<class Num_T>
01741 std::ostream &operator<<(std::ostream &os, const Vec<Num_T> &v)
01742 {
01743 int i, sz = v.length();
01744
01745 os << "[" ;
01746 for (i = 0; i < sz; i++) {
01747 os << v(i) ;
01748 if (i < sz - 1)
01749 os << " ";
01750 }
01751 os << "]" ;
01752
01753 return os;
01754 }
01755
01757 template<class Num_T>
01758 std::istream &operator>>(std::istream &is, Vec<Num_T> &v)
01759 {
01760 std::ostringstream buffer;
01761 bool started = false;
01762 bool finished = false;
01763 bool brackets = false;
01764 char c;
01765
01766 while (!finished) {
01767 if (is.eof()) {
01768 finished = true;
01769 }
01770 else {
01771 is.get(c);
01772
01773 if (is.eof() || (c == '\n')) {
01774 if (brackets) {
01775
01776 is.setstate(std::ios_base::failbit);
01777 finished = true;
01778 }
01779 else if (!((c == '\n') && !started)) {
01780 finished = true;
01781 }
01782 }
01783 else if ((c == ' ') || (c == '\t')) {
01784 if (started) {
01785 buffer << ' ';
01786 }
01787 }
01788 else if (c == '[') {
01789 if (started) {
01790
01791 is.setstate(std::ios_base::failbit);
01792 finished = true;
01793 }
01794 else {
01795 started = true;
01796 brackets = true;
01797 }
01798 }
01799 else if (c == ']') {
01800 if (!started || !brackets) {
01801
01802 is.setstate(std::ios_base::failbit);
01803 finished = true;
01804 }
01805 else {
01806 finished = true;
01807 }
01808 while (!is.eof() && (((c = static_cast<char>(is.peek())) == ' ')
01809 || (c == '\t'))) {
01810 is.get();
01811 }
01812 if (!is.eof() && (c == '\n')) {
01813 is.get();
01814 }
01815 }
01816 else {
01817 started = true;
01818 buffer << c;
01819 }
01820 }
01821 }
01822
01823 if (!started) {
01824 v.set_size(0, false);
01825 }
01826 else {
01827 v.set(buffer.str());
01828 }
01829
01830 return is;
01831 }
01832
01834
01835
01836
01837
01838
01839 template<class Num_T>
01840 void Vec<Num_T>::parse_abc_token(const std::string &s, Num_T &a, Num_T &b,
01841 Num_T &c) const
01842 {
01843 std::string::size_type beg = 0;
01844 std::string::size_type end = s.find(':', 0);
01845 a = parse_token(s.substr(beg, end-beg));
01846 beg = end + 1;
01847 end = s.find(':', beg);
01848 if (end != std::string::npos) {
01849 b = parse_token(s.substr(beg, end-beg));
01850 c = parse_token(s.substr(end+1, s.size()-end));
01851 }
01852 else {
01853 b = Num_T(1);
01854 c = parse_token(s.substr(beg, end-beg-1));
01855 }
01856 }
01857
01858 template<class Num_T>
01859 Num_T Vec<Num_T>::parse_token(const std::string &s) const
01860 {
01861 it_error("Vec::parse_token(): Only `double' and `int' types are supported");
01862 return 0;
01863 }
01864
01865 template<>
01866 double Vec<double>::parse_token(const std::string &s) const;
01867 template<>
01868 int Vec<int>::parse_token(const std::string &s) const;
01869
01870
01871
01872
01873
01874
01875 #ifdef HAVE_EXTERN_TEMPLATE
01876
01877 extern template class Vec<double>;
01878 extern template class Vec<int>;
01879 extern template class Vec<short int>;
01880 extern template class Vec<std::complex<double> >;
01881 extern template class Vec<bin>;
01882
01883
01884
01885 extern template vec operator+(const vec &v1, const vec &v2);
01886 extern template cvec operator+(const cvec &v1, const cvec &v2);
01887 extern template ivec operator+(const ivec &v1, const ivec &v2);
01888 extern template svec operator+(const svec &v1, const svec &v2);
01889 extern template bvec operator+(const bvec &v1, const bvec &v2);
01890
01891 extern template vec operator+(const vec &v1, double t);
01892 extern template cvec operator+(const cvec &v1, std::complex<double> t);
01893 extern template ivec operator+(const ivec &v1, int t);
01894 extern template svec operator+(const svec &v1, short t);
01895 extern template bvec operator+(const bvec &v1, bin t);
01896
01897 extern template vec operator+(double t, const vec &v1);
01898 extern template cvec operator+(std::complex<double> t, const cvec &v1);
01899 extern template ivec operator+(int t, const ivec &v1);
01900 extern template svec operator+(short t, const svec &v1);
01901 extern template bvec operator+(bin t, const bvec &v1);
01902
01903
01904
01905 extern template vec operator-(const vec &v1, const vec &v2);
01906 extern template cvec operator-(const cvec &v1, const cvec &v2);
01907 extern template ivec operator-(const ivec &v1, const ivec &v2);
01908 extern template svec operator-(const svec &v1, const svec &v2);
01909 extern template bvec operator-(const bvec &v1, const bvec &v2);
01910
01911 extern template vec operator-(const vec &v, double t);
01912 extern template cvec operator-(const cvec &v, std::complex<double> t);
01913 extern template ivec operator-(const ivec &v, int t);
01914 extern template svec operator-(const svec &v, short t);
01915 extern template bvec operator-(const bvec &v, bin t);
01916
01917 extern template vec operator-(double t, const vec &v);
01918 extern template cvec operator-(std::complex<double> t, const cvec &v);
01919 extern template ivec operator-(int t, const ivec &v);
01920 extern template svec operator-(short t, const svec &v);
01921 extern template bvec operator-(bin t, const bvec &v);
01922
01923
01924
01925 extern template vec operator-(const vec &v);
01926 extern template cvec operator-(const cvec &v);
01927 extern template ivec operator-(const ivec &v);
01928 extern template svec operator-(const svec &v);
01929 extern template bvec operator-(const bvec &v);
01930
01931
01932
01933 #if !defined(HAVE_BLAS)
01934 extern template double dot(const vec &v1, const vec &v2);
01935 #if !(defined(HAVE_ZDOTUSUB) || defined(HAVE_ZDOTU_VOID))
01936 extern template std::complex<double> dot(const cvec &v1, const cvec &v2);
01937 #endif // !(HAVE_ZDOTUSUB || HAVE_ZDOTU_VOID)
01938 #endif // HAVE_BLAS
01939 extern template int dot(const ivec &v1, const ivec &v2);
01940 extern template short dot(const svec &v1, const svec &v2);
01941 extern template bin dot(const bvec &v1, const bvec &v2);
01942
01943 #if !defined(HAVE_BLAS)
01944 extern template double operator*(const vec &v1, const vec &v2);
01945 extern template std::complex<double> operator*(const cvec &v1,
01946 const cvec &v2);
01947 #endif
01948 extern template int operator*(const ivec &v1, const ivec &v2);
01949 extern template short operator*(const svec &v1, const svec &v2);
01950 extern template bin operator*(const bvec &v1, const bvec &v2);
01951
01952 #if !defined(HAVE_BLAS)
01953 extern template mat outer_product(const vec &v1, const vec &v2,
01954 bool hermitian);
01955 #endif
01956 extern template imat outer_product(const ivec &v1, const ivec &v2,
01957 bool hermitian);
01958 extern template smat outer_product(const svec &v1, const svec &v2,
01959 bool hermitian);
01960 extern template bmat outer_product(const bvec &v1, const bvec &v2,
01961 bool hermitian);
01962
01963 extern template vec operator*(const vec &v, double t);
01964 extern template cvec operator*(const cvec &v, std::complex<double> t);
01965 extern template ivec operator*(const ivec &v, int t);
01966 extern template svec operator*(const svec &v, short t);
01967 extern template bvec operator*(const bvec &v, bin t);
01968
01969 extern template vec operator*(double t, const vec &v);
01970 extern template cvec operator*(std::complex<double> t, const cvec &v);
01971 extern template ivec operator*(int t, const ivec &v);
01972 extern template svec operator*(short t, const svec &v);
01973 extern template bvec operator*(bin t, const bvec &v);
01974
01975
01976
01977 extern template vec elem_mult(const vec &a, const vec &b);
01978 extern template cvec elem_mult(const cvec &a, const cvec &b);
01979 extern template ivec elem_mult(const ivec &a, const ivec &b);
01980 extern template svec elem_mult(const svec &a, const svec &b);
01981 extern template bvec elem_mult(const bvec &a, const bvec &b);
01982
01983 extern template void elem_mult_out(const vec &a, const vec &b, vec &out);
01984 extern template void elem_mult_out(const cvec &a, const cvec &b, cvec &out);
01985 extern template void elem_mult_out(const ivec &a, const ivec &b, ivec &out);
01986 extern template void elem_mult_out(const svec &a, const svec &b, svec &out);
01987 extern template void elem_mult_out(const bvec &a, const bvec &b, bvec &out);
01988
01989 extern template vec elem_mult(const vec &a, const vec &b, const vec &c);
01990 extern template cvec elem_mult(const cvec &a, const cvec &b, const cvec &c);
01991 extern template ivec elem_mult(const ivec &a, const ivec &b, const ivec &c);
01992 extern template svec elem_mult(const svec &a, const svec &b, const svec &c);
01993 extern template bvec elem_mult(const bvec &a, const bvec &b, const bvec &c);
01994
01995 extern template void elem_mult_out(const vec &a, const vec &b,
01996 const vec &c, vec &out);
01997 extern template void elem_mult_out(const cvec &a, const cvec &b,
01998 const cvec &c, cvec &out);
01999 extern template void elem_mult_out(const ivec &a, const ivec &b,
02000 const ivec &c, ivec &out);
02001 extern template void elem_mult_out(const svec &a, const svec &b,
02002 const svec &c, svec &out);
02003 extern template void elem_mult_out(const bvec &a, const bvec &b,
02004 const bvec &c, bvec &out);
02005
02006 extern template vec elem_mult(const vec &a, const vec &b,
02007 const vec &c, const vec &d);
02008 extern template cvec elem_mult(const cvec &a, const cvec &b,
02009 const cvec &c, const cvec &d);
02010 extern template ivec elem_mult(const ivec &a, const ivec &b,
02011 const ivec &c, const ivec &d);
02012 extern template svec elem_mult(const svec &a, const svec &b,
02013 const svec &c, const svec &d);
02014 extern template bvec elem_mult(const bvec &a, const bvec &b,
02015 const bvec &c, const bvec &d);
02016
02017 extern template void elem_mult_out(const vec &a, const vec &b, const vec &c,
02018 const vec &d, vec &out);
02019 extern template void elem_mult_out(const cvec &a, const cvec &b,
02020 const cvec &c, const cvec &d, cvec &out);
02021 extern template void elem_mult_out(const ivec &a, const ivec &b,
02022 const ivec &c, const ivec &d, ivec &out);
02023 extern template void elem_mult_out(const svec &a, const svec &b,
02024 const svec &c, const svec &d, svec &out);
02025 extern template void elem_mult_out(const bvec &a, const bvec &b,
02026 const bvec &c, const bvec &d, bvec &out);
02027
02028
02029
02030 extern template void elem_mult_inplace(const vec &a, vec &b);
02031 extern template void elem_mult_inplace(const cvec &a, cvec &b);
02032 extern template void elem_mult_inplace(const ivec &a, ivec &b);
02033 extern template void elem_mult_inplace(const svec &a, svec &b);
02034 extern template void elem_mult_inplace(const bvec &a, bvec &b);
02035
02036
02037
02038 extern template double elem_mult_sum(const vec &a, const vec &b);
02039 extern template std::complex<double> elem_mult_sum(const cvec &a,
02040 const cvec &b);
02041 extern template int elem_mult_sum(const ivec &a, const ivec &b);
02042 extern template short elem_mult_sum(const svec &a, const svec &b);
02043 extern template bin elem_mult_sum(const bvec &a, const bvec &b);
02044
02045
02046
02047 extern template vec operator/(const vec &v, double t);
02048 extern template cvec operator/(const cvec &v, std::complex<double> t);
02049 extern template ivec operator/(const ivec &v, int t);
02050 extern template svec operator/(const svec &v, short t);
02051 extern template bvec operator/(const bvec &v, bin t);
02052
02053 extern template vec operator/(double t, const vec &v);
02054 extern template cvec operator/(std::complex<double> t, const cvec &v);
02055 extern template ivec operator/(int t, const ivec &v);
02056 extern template svec operator/(short t, const svec &v);
02057 extern template bvec operator/(bin t, const bvec &v);
02058
02059
02060
02061 extern template vec elem_div(const vec &a, const vec &b);
02062 extern template cvec elem_div(const cvec &a, const cvec &b);
02063 extern template ivec elem_div(const ivec &a, const ivec &b);
02064 extern template svec elem_div(const svec &a, const svec &b);
02065 extern template bvec elem_div(const bvec &a, const bvec &b);
02066
02067 extern template vec elem_div(double t, const vec &v);
02068 extern template cvec elem_div(std::complex<double> t, const cvec &v);
02069 extern template ivec elem_div(int t, const ivec &v);
02070 extern template svec elem_div(short t, const svec &v);
02071 extern template bvec elem_div(bin t, const bvec &v);
02072
02073 extern template void elem_div_out(const vec &a, const vec &b, vec &out);
02074 extern template void elem_div_out(const cvec &a, const cvec &b, cvec &out);
02075 extern template void elem_div_out(const ivec &a, const ivec &b, ivec &out);
02076 extern template void elem_div_out(const svec &a, const svec &b, svec &out);
02077 extern template void elem_div_out(const bvec &a, const bvec &b, bvec &out);
02078
02079
02080
02081 extern template double elem_div_sum(const vec &a, const vec &b);
02082 extern template std::complex<double> elem_div_sum(const cvec &a,
02083 const cvec &b);
02084 extern template int elem_div_sum(const ivec &a, const ivec &b);
02085 extern template short elem_div_sum(const svec &a, const svec &b);
02086 extern template bin elem_div_sum(const bvec &a, const bvec &b);
02087
02088
02089
02090 extern template vec concat(const vec &v, double a);
02091 extern template cvec concat(const cvec &v, std::complex<double> a);
02092 extern template ivec concat(const ivec &v, int a);
02093 extern template svec concat(const svec &v, short a);
02094 extern template bvec concat(const bvec &v, bin a);
02095
02096 extern template vec concat(double a, const vec &v);
02097 extern template cvec concat(std::complex<double> a, const cvec &v);
02098 extern template ivec concat(int a, const ivec &v);
02099 extern template svec concat(short a, const svec &v);
02100 extern template bvec concat(bin a, const bvec &v);
02101
02102 extern template vec concat(const vec &v1, const vec &v2);
02103 extern template cvec concat(const cvec &v1, const cvec &v2);
02104 extern template ivec concat(const ivec &v1, const ivec &v2);
02105 extern template svec concat(const svec &v1, const svec &v2);
02106 extern template bvec concat(const bvec &v1, const bvec &v2);
02107
02108 extern template vec concat(const vec &v1, const vec &v2, const vec &v3);
02109 extern template cvec concat(const cvec &v1, const cvec &v2, const cvec &v3);
02110 extern template ivec concat(const ivec &v1, const ivec &v2, const ivec &v3);
02111 extern template svec concat(const svec &v1, const svec &v2, const svec &v3);
02112 extern template bvec concat(const bvec &v1, const bvec &v2, const bvec &v3);
02113
02114 extern template vec concat(const vec &v1, const vec &v2,
02115 const vec &v3, const vec &v4);
02116 extern template cvec concat(const cvec &v1, const cvec &v2,
02117 const cvec &v3, const cvec &v4);
02118 extern template ivec concat(const ivec &v1, const ivec &v2,
02119 const ivec &v3, const ivec &v4);
02120 extern template svec concat(const svec &v1, const svec &v2,
02121 const svec &v3, const svec &v4);
02122 extern template bvec concat(const bvec &v1, const bvec &v2,
02123 const bvec &v3, const bvec &v4);
02124
02125 extern template vec concat(const vec &v1, const vec &v2, const vec &v3,
02126 const vec &v4, const vec &v5);
02127 extern template cvec concat(const cvec &v1, const cvec &v2, const cvec &v3,
02128 const cvec &v4, const cvec &v5);
02129 extern template ivec concat(const ivec &v1, const ivec &v2, const ivec &v3,
02130 const ivec &v4, const ivec &v5);
02131 extern template svec concat(const svec &v1, const svec &v2, const svec &v3,
02132 const svec &v4, const svec &v5);
02133 extern template bvec concat(const bvec &v1, const bvec &v2, const bvec &v3,
02134 const bvec &v4, const bvec &v5);
02135
02136
02137
02138 extern template std::ostream &operator<<(std::ostream& os, const vec &vect);
02139 extern template std::ostream &operator<<(std::ostream& os, const cvec &vect);
02140 extern template std::ostream &operator<<(std::ostream& os, const svec &vect);
02141 extern template std::ostream &operator<<(std::ostream& os, const ivec &vect);
02142 extern template std::ostream &operator<<(std::ostream& os, const bvec &vect);
02143 extern template std::istream &operator>>(std::istream& is, vec &vect);
02144 extern template std::istream &operator>>(std::istream& is, cvec &vect);
02145 extern template std::istream &operator>>(std::istream& is, svec &vect);
02146 extern template std::istream &operator>>(std::istream& is, ivec &vect);
02147 extern template std::istream &operator>>(std::istream& is, bvec &vect);
02148
02149 #endif // HAVE_EXTERN_TEMPLATE
02150
02152
02153 }
02154
02155 #endif // #ifndef VEC_H
