IT++ Logo

mat.h

Go to the documentation of this file.
00001 
00030 #ifndef MAT_H
00031 #define MAT_H
00032 
00033 #ifndef _MSC_VER
00034 #  include <itpp/config.h>
00035 #else
00036 #  include <itpp/config_msvc.h>
00037 #endif
00038 
00039 #include <itpp/base/itassert.h>
00040 #include <itpp/base/math/misc.h>
00041 #include <itpp/base/factory.h>
00042 
00043 
00044 namespace itpp {
00045 
00046   // Declaration of Vec
00047   template<class Num_T> class Vec;
00048   // Declaration of Mat
00049   template<class Num_T> class Mat;
00050   // Declaration of bin
00051   class bin;
00052 
00054   template<class Num_T> const Mat<Num_T> concat_horizontal(const Mat<Num_T> &m1, const Mat<Num_T> &m2);
00056   template<class Num_T> const Mat<Num_T> concat_vertical(const Mat<Num_T> &m1, const Mat<Num_T> &m2);
00057 
00059   template<class Num_T> const Mat<Num_T> operator+(const Mat<Num_T> &m1, const Mat<Num_T> &m2);
00061   template<class Num_T> const Mat<Num_T> operator+(const Mat<Num_T> &m, Num_T t);
00063   template<class Num_T> const Mat<Num_T> operator+(Num_T t, const Mat<Num_T> &m);
00064 
00066   template<class Num_T> const Mat<Num_T> operator-(const Mat<Num_T> &m1, const Mat<Num_T> &m2);
00068   template<class Num_T> const Mat<Num_T> operator-(const Mat<Num_T> &m, Num_T t);
00070   template<class Num_T> const Mat<Num_T> operator-(Num_T t, const Mat<Num_T> &m);
00072   template<class Num_T> const Mat<Num_T> operator-(const Mat<Num_T> &m);
00073 
00075   template<class Num_T> const Mat<Num_T> operator*(const Mat<Num_T> &m1, const Mat<Num_T> &m2);
00077   template<class Num_T> const Vec<Num_T> operator*(const Mat<Num_T> &m, const Vec<Num_T> &v);
00079   template<class Num_T> const Mat<Num_T> operator*(const Vec<Num_T> &v, const Mat<Num_T> &m);
00081   template<class Num_T> const Mat<Num_T> operator*(const Mat<Num_T> &m, Num_T t);
00083   template<class Num_T> const Mat<Num_T> operator*(Num_T t, const Mat<Num_T> &m);
00084 
00086   template<class Num_T> const Mat<Num_T> elem_mult(const Mat<Num_T> &A, const Mat<Num_T> &B);
00088   template<class Num_T> void elem_mult_out(const Mat<Num_T> &A, const Mat<Num_T> &B, Mat<Num_T> &out);
00090   template<class Num_T> void elem_mult_out(const Mat<Num_T> &A, const Mat<Num_T> &B, const Mat<Num_T> &C, Mat<Num_T> &out);
00092   template<class Num_T> void elem_mult_out(const Mat<Num_T> &A, const Mat<Num_T> &B, const Mat<Num_T> &C, const Mat<Num_T> &D, Mat<Num_T> &out);
00094   template<class Num_T> void elem_mult_inplace(const Mat<Num_T> &A, Mat<Num_T> &B);
00096   template<class Num_T> Num_T elem_mult_sum(const Mat<Num_T> &A, const Mat<Num_T> &B);
00097 
00099   template<class Num_T> const Mat<Num_T> operator/(const Mat<Num_T> &m, Num_T t);
00100 
00102   template<class Num_T> const Mat<Num_T> elem_div(const Mat<Num_T> &A, const Mat<Num_T> &B);
00104   template<class Num_T> void elem_div_out(const Mat<Num_T> &A, const Mat<Num_T> &B, Mat<Num_T> &out);
00106   template<class Num_T> Num_T elem_div_sum(const Mat<Num_T> &A, const Mat<Num_T> &B);
00107 
00108   // -------------------------------------------------------------------------------------
00109   // Declaration of Mat
00110   // -------------------------------------------------------------------------------------
00111 
00177   template<class Num_T>
00178   class Mat {
00179   public:
00181     typedef Num_T value_type;
00182 
00184     explicit Mat(const Factory &f = DEFAULT_FACTORY);
00186     Mat(int inrow, int incol, const Factory &f = DEFAULT_FACTORY);
00188     Mat(const Mat<Num_T> &m);
00190     Mat(const Mat<Num_T> &m, const Factory &f);
00192     Mat(const Vec<Num_T> &invector, const Factory &f = DEFAULT_FACTORY);
00194     Mat(const std::string &str, const Factory &f = DEFAULT_FACTORY);
00196     Mat(const char *str, const Factory &f = DEFAULT_FACTORY);
00203     Mat(const Num_T *c_array, int rows, int cols, bool RowMajor = true, const Factory &f = DEFAULT_FACTORY);
00204 
00206     ~Mat();
00207 
00209     int cols() const { return no_cols; }
00211     int rows() const { return no_rows; }
00213     int size() const { return datasize; }
00215     void set_size(int inrow, int incol, bool copy = false);
00217     void zeros();
00219     void clear() { zeros(); }
00221     void ones();
00223     void set(const char *str);
00225     void set(const std::string &str);
00226 
00228     const Num_T &operator()(int R, int C) const;
00230     Num_T &operator()(int R, int C);
00232     const Num_T &operator()(int index) const;
00234     Num_T &operator()(int index);
00236     const Num_T &get(int R,int C) const;
00238     void set(int R,int C, const Num_T &v);
00239 
00245     const Mat<Num_T> operator()(int r1, int r2, int c1, int c2) const;
00251     const Mat<Num_T> get(int r1, int r2, int c1, int c2) const;
00252 
00254     const Vec<Num_T> get_row(int Index) const ;
00256     const Mat<Num_T> get_rows(int r1, int r2) const;
00258     const Mat<Num_T> get_rows(const Vec<int> &indexlist) const;
00260     const Vec<Num_T> get_col(int Index) const ;
00262     const Mat<Num_T> get_cols(int c1, int c2) const;
00264     const Mat<Num_T> get_cols(const Vec<int> &indexlist) const;
00266     void set_row(int Index, const Vec<Num_T> &invector);
00268     void set_col(int Index, const Vec<Num_T> &invector);
00270     void set_rows(int r, const Mat<Num_T> &m);
00272     void set_cols(int c, const Mat<Num_T> &m);
00274     void copy_row(int to, int from);
00276     void copy_col(int to, int from);
00278     void swap_rows(int r1, int r2);
00280     void swap_cols(int c1, int c2);
00281 
00283     void set_submatrix(int r1, int r2, int c1, int c2, const Mat<Num_T> &m);
00285     void set_submatrix(int r, int c, const Mat<Num_T> &m);
00287     void set_submatrix(int r1, int r2, int c1, int c2, const Num_T t);
00288 
00290     void del_row(int r);
00292     void del_rows(int r1, int r2);
00294     void del_col(int c);
00296     void del_cols(int c1, int c2);
00298     void ins_row(int r, const Vec<Num_T> &in);
00300     void ins_col(int c, const Vec<Num_T> &in);
00302     void append_row(const Vec<Num_T> &in);
00304     void append_col(const Vec<Num_T> &in);
00305 
00307     const Mat<Num_T> transpose() const;
00309     const Mat<Num_T> T() const { return this->transpose(); }
00311     const Mat<Num_T> hermitian_transpose() const;
00313     const Mat<Num_T> H() const { return this->hermitian_transpose(); }
00314 
00316     friend const Mat<Num_T> concat_horizontal <>(const Mat<Num_T> &m1, const Mat<Num_T> &m2);
00318     friend const Mat<Num_T> concat_vertical <>(const Mat<Num_T> &m1, const Mat<Num_T> &m2);
00319 
00321     Mat<Num_T>& operator=(Num_T t);
00323     Mat<Num_T>& operator=(const Mat<Num_T> &m);
00325     Mat<Num_T>& operator=(const Vec<Num_T> &v);
00327     Mat<Num_T>& operator=(const char *values);
00328 
00330     Mat<Num_T>& operator+=(const Mat<Num_T> &m);
00332     Mat<Num_T>& operator+=(Num_T t);
00334     friend const  Mat<Num_T> operator+<>(const Mat<Num_T> &m1, const Mat<Num_T> &m2);
00336     friend const Mat<Num_T> operator+<>(const Mat<Num_T> &m, Num_T t);
00338     friend const Mat<Num_T> operator+<>(Num_T t, const Mat<Num_T> &m);
00339 
00341     Mat<Num_T>& operator-=(const Mat<Num_T> &m);
00343     Mat<Num_T>& operator-=(Num_T t);
00345     friend const Mat<Num_T> operator-<>(const Mat<Num_T> &m1, const Mat<Num_T> &m2);
00347     friend const Mat<Num_T> operator-<>(const Mat<Num_T> &m, Num_T t);
00349     friend const Mat<Num_T> operator-<>(Num_T t, const Mat<Num_T> &m);
00351     friend const Mat<Num_T> operator-<>(const Mat<Num_T> &m);
00352 
00354     Mat<Num_T>& operator*=(const Mat<Num_T> &m);
00356     Mat<Num_T>& operator*=(Num_T t);
00358     friend const Mat<Num_T> operator*<>(const Mat<Num_T> &m1, const Mat<Num_T> &m2);
00360     friend const Vec<Num_T> operator*<>(const Mat<Num_T> &m, const Vec<Num_T> &v);
00372     friend const Mat<Num_T> operator*<>(const Vec<Num_T> &v,
00373                                         const Mat<Num_T> &m);
00375     friend const Mat<Num_T> operator*<>(const Mat<Num_T> &m, Num_T t);
00377     friend const Mat<Num_T> operator*<>(Num_T t, const Mat<Num_T> &m);
00378 
00380     friend const Mat<Num_T> elem_mult <>(const Mat<Num_T> &A, const Mat<Num_T> &B);
00382     friend void elem_mult_out <>(const Mat<Num_T> &A, const Mat<Num_T> &B, Mat<Num_T> &out);
00384     friend void elem_mult_out <>(const Mat<Num_T> &A, const Mat<Num_T> &B, const Mat<Num_T> &C, Mat<Num_T> &out);
00386     friend void elem_mult_out <>(const Mat<Num_T> &A, const Mat<Num_T> &B, const Mat<Num_T> &C, const Mat<Num_T> &D, Mat<Num_T> &out);
00388     friend void elem_mult_inplace <>(const Mat<Num_T> &A, Mat<Num_T> &B);
00390     friend Num_T elem_mult_sum <>(const Mat<Num_T> &A, const Mat<Num_T> &B);
00391 
00393     Mat<Num_T>& operator/=(Num_T t);
00395     friend const Mat<Num_T> operator/<>(const Mat<Num_T> &m, Num_T t);
00397     Mat<Num_T>& operator/=(const Mat<Num_T> &m);
00398 
00400     friend const Mat<Num_T> elem_div <>(const Mat<Num_T> &A, const Mat<Num_T> &B);
00402     friend void elem_div_out <>(const Mat<Num_T> &A, const Mat<Num_T> &B, Mat<Num_T> &out);
00404     friend Num_T elem_div_sum <>(const Mat<Num_T> &A, const Mat<Num_T> &B);
00405 
00407     bool operator==(const Mat<Num_T> &m) const;
00409     bool operator!=(const Mat<Num_T> &m) const;
00410 
00412     Num_T &_elem(int R,int C) { return data[R+C*no_rows]; }
00414     const Num_T &_elem(int R,int C) const { return data[R+C*no_rows]; }
00416     Num_T &_elem(int index) { return data[index]; }
00418     const Num_T &_elem(int index) const { return data[index]; }
00419 
00421     Num_T *_data() { return data; }
00423     const Num_T *_data() const { return data; }
00425     int _datasize() const { return datasize; }
00426 
00427   protected:
00429     void alloc(int rows, int cols);
00431     void free();
00432 
00435     int datasize, no_rows, no_cols;
00437 
00438     Num_T *data;
00440     const Factory &factory;
00441   };
00442 
00443   // -------------------------------------------------------------------------------------
00444   // Type definitions of mat, cmat, imat, smat, and bmat
00445   // -------------------------------------------------------------------------------------
00446 
00451   typedef Mat<double> mat;
00452 
00457   typedef Mat<std::complex<double> > cmat;
00458 
00463   typedef Mat<int> imat;
00464 
00469   typedef Mat<short int> smat;
00470 
00477   typedef Mat<bin> bmat;
00478 
00479 } //namespace itpp
00480 
00481 
00482 #include <itpp/base/vec.h>
00483 
00484 namespace itpp {
00485 
00486   // ----------------------------------------------------------------------
00487   // Declaration of input and output streams for Mat
00488   // ----------------------------------------------------------------------
00489 
00494   template <class Num_T>
00495   std::ostream &operator<<(std::ostream &os, const Mat<Num_T> &m);
00496 
00508   template <class Num_T>
00509   std::istream &operator>>(std::istream &is, Mat<Num_T> &m);
00510 
00511   // ----------------------------------------------------------------------
00512   // Implementation of templated Mat members and friends
00513   // ----------------------------------------------------------------------
00514 
00515   template<class Num_T> inline
00516   void Mat<Num_T>::alloc(int rows, int cols)
00517   {
00518     if ((rows > 0) && (cols > 0)) {
00519       datasize = rows * cols;
00520       no_rows = rows;
00521       no_cols = cols;
00522       create_elements(data, datasize, factory);
00523     }
00524     else {
00525       data = 0;
00526       datasize = 0;
00527       no_rows = 0;
00528       no_cols = 0;
00529     }
00530   }
00531 
00532   template<class Num_T> inline
00533   void Mat<Num_T>::free()
00534   {
00535     destroy_elements(data, datasize);
00536     datasize = 0;
00537     no_rows = 0;
00538     no_cols = 0;
00539   }
00540 
00541 
00542   template<class Num_T> inline
00543   Mat<Num_T>::Mat(const Factory &f) :
00544     datasize(0), no_rows(0), no_cols(0), data(0), factory(f) {}
00545 
00546   template<class Num_T> inline
00547   Mat<Num_T>::Mat(int inrow, int incol, const Factory &f) :
00548     datasize(0), no_rows(0), no_cols(0), data(0), factory(f)
00549   {
00550     it_assert_debug((inrow >= 0) && (incol >= 0), "The rows and columns must be >= 0");
00551     alloc(inrow, incol);
00552   }
00553 
00554   template<class Num_T> inline
00555   Mat<Num_T>::Mat(const Mat<Num_T> &m) :
00556     datasize(0), no_rows(0), no_cols(0), data(0), factory(m.factory)
00557   {
00558     alloc(m.no_rows, m.no_cols);
00559     copy_vector(m.datasize, m.data, data);
00560   }
00561 
00562   template<class Num_T> inline
00563   Mat<Num_T>::Mat(const Mat<Num_T> &m, const Factory &f) :
00564     datasize(0), no_rows(0), no_cols(0), data(0), factory(f)
00565   {
00566     alloc(m.no_rows, m.no_cols);
00567     copy_vector(m.datasize, m.data, data);
00568   }
00569 
00570   template<class Num_T> inline
00571   Mat<Num_T>::Mat(const Vec<Num_T> &invector, const Factory &f) :
00572     datasize(0), no_rows(0), no_cols(0), data(0), factory(f)
00573   {
00574     int size = invector.size();
00575     alloc(size, 1);
00576     copy_vector(size, invector._data(), data);
00577   }
00578 
00579   template<class Num_T> inline
00580   Mat<Num_T>::Mat(const std::string &str, const Factory &f) :
00581     datasize(0), no_rows(0), no_cols(0), data(0), factory(f)
00582   {
00583     set(str);
00584   }
00585 
00586   template<class Num_T> inline
00587   Mat<Num_T>::Mat(const char *str, const Factory &f) :
00588     datasize(0), no_rows(0), no_cols(0), data(0), factory(f)
00589   {
00590     set(str);
00591   }
00592 
00593   template<class Num_T> inline
00594   Mat<Num_T>::Mat(const Num_T *c_array, int rows, int cols, bool RowMajor, const Factory &f) :
00595     datasize(0), no_rows(0), no_cols(0), data(0), factory(f)
00596   {
00597     alloc(rows, cols);
00598 
00599     if (!RowMajor)
00600       copy_vector(datasize, c_array, data);
00601     else { // Row Major
00602       for (int i=0; i<rows; i++) {
00603         for (int j=0; j<cols; j++)
00604           data[i+j*no_rows] = c_array[i*no_cols+j];
00605       }
00606     }
00607   }
00608 
00609   template<class Num_T> inline
00610   Mat<Num_T>::~Mat()
00611   {
00612     free();
00613   }
00614 
00615 
00616   template<class Num_T>
00617   void Mat<Num_T>::set_size(int inrow, int incol, bool copy)
00618   {
00619     it_assert_debug((inrow >= 0) && (incol >= 0), "Mat<Num_T>::set_size: "
00620                     "The number of rows and columns must be >= 0");
00621     // check if we have to resize the current matrix
00622     if ((no_rows == inrow) && (no_cols == incol))
00623       return;
00624     // conditionally copy previous matrix content
00625     if (copy) {
00626       // create a temporary pointer to the allocated data
00627       Num_T* tmp = data;
00628       // store the current number of elements and number of rows
00629       int old_datasize = datasize;
00630       int old_rows = no_rows;
00631       // check the boundaries of the copied data
00632       int min_r = (no_rows < inrow) ? no_rows : inrow;
00633       int min_c = (no_cols < incol) ? no_cols : incol;
00634       // allocate new memory
00635       alloc(inrow, incol);
00636       // copy the previous data into the allocated memory
00637       for (int i = 0; i < min_c; ++i) {
00638         copy_vector(min_r, &tmp[i*old_rows], &data[i*no_rows]);
00639       }
00640       // fill-in the rest of matrix with zeros
00641       for (int i = min_r; i < inrow; ++i)
00642         for (int j = 0; j < incol; ++j)
00643           data[i+j*inrow] = Num_T(0);
00644       for (int j = min_c; j < incol; ++j)
00645         for (int i = 0; i < min_r; ++i)
00646           data[i+j*inrow] = Num_T(0);
00647       // delete old elements
00648       destroy_elements(tmp, old_datasize);
00649     }
00650     // if possible, reuse the allocated memory
00651     else if (datasize == inrow * incol) {
00652       no_rows = inrow;
00653       no_cols = incol;
00654     }
00655     // finally release old memory and allocate a new one
00656     else {
00657       free();
00658       alloc(inrow, incol);
00659     }
00660   }
00661 
00662   template<class Num_T> inline
00663   void Mat<Num_T>::zeros()
00664   {
00665     for(int i=0; i<datasize; i++)
00666       data[i] = Num_T(0);
00667   }
00668 
00669   template<class Num_T> inline
00670   void Mat<Num_T>::ones()
00671   {
00672     for(int i=0; i<datasize; i++)
00673       data[i] = Num_T(1);
00674   }
00675 
00676   template<class Num_T> inline
00677   const Num_T& Mat<Num_T>::operator()(int R,int C) const
00678   {
00679     it_assert_debug((R >= 0) && (R < no_rows) && (C >= 0) && (C < no_cols),
00680                "Mat<Num_T>::operator(): index out of range");
00681     return data[R+C*no_rows];
00682   }
00683 
00684   template<class Num_T> inline
00685   Num_T& Mat<Num_T>::operator()(int R,int C)
00686   {
00687     it_assert_debug((R >= 0) && (R < no_rows) && (C >= 0) && (C < no_cols),
00688                "Mat<Num_T>::operator(): index out of range");
00689     return data[R+C*no_rows];
00690   }
00691 
00692   template<class Num_T> inline
00693   Num_T& Mat<Num_T>::operator()(int index)
00694   {
00695     it_assert_debug((index < no_rows*no_cols) && (index >= 0),
00696                "Mat<Num_T>::operator(): index out of range");
00697     return data[index];
00698   }
00699 
00700   template<class Num_T> inline
00701   const Num_T& Mat<Num_T>::operator()(int index) const
00702   {
00703     it_assert_debug((index < no_rows*no_cols) && (index >= 0),
00704                "Mat<Num_T>::operator(): index out of range");
00705     return data[index];
00706   }
00707 
00708   template<class Num_T> inline
00709   const Num_T& Mat<Num_T>::get(int R,int C) const
00710   {
00711     it_assert_debug((R >= 0) && (R < no_rows) && (C >= 0) && (C < no_cols),
00712                "Mat<Num_T>::get(): index out of range");
00713     return data[R+C*no_rows];
00714   }
00715 
00716   template<class Num_T> inline
00717   void Mat<Num_T>::set(int R, int C, const Num_T &v)
00718   {
00719     it_assert_debug((R >= 0) && (R < no_rows) && (C >= 0) && (C < no_cols),
00720                "Mat<Num_T>::set(): index out of range");
00721     data[R+C*no_rows] = v;
00722   }
00723 
00724 
00725   template<class Num_T>
00726   void Mat<Num_T>::set(const std::string &str)
00727   {
00728     // actual row counter
00729     int rows = 0;
00730     // number of rows to allocate next time (8, 16, 32, 64, etc.)
00731     int maxrows = 8;
00732 
00733     // clean the current matrix content
00734     free();
00735 
00736     // variable to store the start of a current vector
00737     std::string::size_type beg = 0;
00738     std::string::size_type end = 0;
00739     while (end != std::string::npos) {
00740       // find next occurence of a semicolon in string str
00741       end = str.find(';', beg);
00742       // parse first row into a vector v
00743       Vec<Num_T> v(str.substr(beg, end-beg));
00744       int v_size = v.size();
00745 
00746       // this check is necessary to parse the following two strings as the
00747       // same matrix: "1 0 1; ; 1 1; " and "1 0 1; 0 0 0; 1 1 0"
00748       if ((end != std::string::npos) || (v_size > 0)) {
00749         // matrix empty -> insert v as a first row and allocate maxrows
00750         if (rows == 0) {
00751           set_size(maxrows, v_size, true);
00752           set_row(rows++, v);
00753         }
00754         else {
00755           // check if we need to resize the matrix
00756           if ((rows == maxrows) || (v_size != no_cols)) {
00757             // we need to add new rows
00758             if (rows == maxrows) {
00759               maxrows *= 2;
00760             }
00761             // check if ne need to add new colums
00762             if (v_size > no_cols) {
00763               set_size(maxrows, v_size, true);
00764             }
00765             else {
00766               set_size(maxrows, no_cols, true);
00767               // set the size of the parsed vector to the number of colums
00768               v.set_size(no_cols, true);
00769             }
00770           }
00771           // set the parsed vector as the next row
00772           set_row(rows++, v);
00773         }
00774       }
00775       // update the starting position of the next vector in the parsed
00776       // string
00777       beg = end + 1;
00778     } // if ((end != std::string::npos) || (v.size > 0))
00779 
00780     set_size(rows, no_cols, true);
00781   }
00782 
00783   template<class Num_T>
00784   void Mat<Num_T>::set(const char *str)
00785   {
00786     set(std::string(str));
00787   }
00788 
00789   template<class Num_T> inline
00790   const Mat<Num_T> Mat<Num_T>::operator()(int r1, int r2, int c1, int c2) const
00791   {
00792     if (r1 == -1) r1 = no_rows-1;
00793     if (r2 == -1) r2 = no_rows-1;
00794     if (c1 == -1) c1 = no_cols-1;
00795     if (c2 == -1) c2 = no_cols-1;
00796 
00797     it_assert_debug(r1>=0 && r2>=0 && r1<no_rows && r2<no_rows &&
00798                c1>=0 && c2>=0 && c1<no_cols && c2<no_cols, "operator()(r1,r2,c1,c2)");
00799 
00800     it_assert_debug(r2>=r1 && c2>=c1, "Mat<Num_T>::op(): r2>=r1 or c2>=c1 not fulfilled");
00801 
00802     Mat<Num_T> s(r2-r1+1, c2-c1+1);
00803 
00804     for (int i=0;i<s.no_cols;i++)
00805       copy_vector(s.no_rows, data+r1+(c1+i)*no_rows, s.data+i*s.no_rows);
00806 
00807     return s;
00808   }
00809 
00810   template<class Num_T> inline
00811   const Mat<Num_T> Mat<Num_T>::get(int r1, int r2, int c1, int c2) const
00812   {
00813     return (*this)(r1, r2, c1, c2);
00814   }
00815 
00816   template<class Num_T> inline
00817   const Vec<Num_T> Mat<Num_T>::get_row(int Index) const
00818   {
00819     it_assert_debug(Index>=0 && Index<no_rows, "Mat<Num_T>::get_row: index out of range");
00820     Vec<Num_T> a(no_cols);
00821 
00822     copy_vector(no_cols, data+Index, no_rows, a._data(), 1);
00823     return a;
00824   }
00825 
00826   template<class Num_T>
00827   const Mat<Num_T> Mat<Num_T>::get_rows(int r1, int r2) const
00828   {
00829     it_assert_debug(r1>=0 && r2<no_rows && r1<=r2, "Mat<Num_T>::get_rows: index out of range");
00830     Mat<Num_T> m(r2-r1+1, no_cols);
00831 
00832     for (int i=0; i<m.rows(); i++)
00833       copy_vector(no_cols, data+i+r1, no_rows, m.data+i, m.no_rows);
00834 
00835     return m;
00836   }
00837 
00838   template<class Num_T>
00839   const Mat<Num_T> Mat<Num_T>::get_rows(const Vec<int> &indexlist) const
00840   {
00841     Mat<Num_T> m(indexlist.size(),no_cols);
00842 
00843     for (int i=0;i<indexlist.size();i++) {
00844       it_assert_debug(indexlist(i)>=0 && indexlist(i)<no_rows, "Mat<Num_T>::get_rows: index out of range");
00845       copy_vector(no_cols, data+indexlist(i), no_rows, m.data+i, m.no_rows);
00846     }
00847 
00848     return m;
00849   }
00850 
00851   template<class Num_T> inline
00852   const Vec<Num_T> Mat<Num_T>::get_col(int Index) const
00853   {
00854     it_assert_debug(Index>=0 && Index<no_cols, "Mat<Num_T>::get_col: index out of range");
00855     Vec<Num_T> a(no_rows);
00856 
00857     copy_vector(no_rows, data+Index*no_rows, a._data());
00858 
00859     return a;
00860   }
00861 
00862   template<class Num_T> inline
00863   const Mat<Num_T> Mat<Num_T>::get_cols(int c1, int c2) const
00864   {
00865     it_assert_debug(c1>=0 && c2<no_cols && c1<=c2, "Mat<Num_T>::get_cols: index out of range");
00866     Mat<Num_T> m(no_rows, c2-c1+1);
00867 
00868     for (int i=0; i<m.cols(); i++)
00869       copy_vector(no_rows, data+(i+c1)*no_rows, m.data+i*m.no_rows);
00870 
00871     return m;
00872   }
00873 
00874   template<class Num_T> inline
00875   const Mat<Num_T> Mat<Num_T>::get_cols(const Vec<int> &indexlist) const
00876   {
00877     Mat<Num_T> m(no_rows,indexlist.size());
00878 
00879     for (int i=0; i<indexlist.size(); i++) {
00880       it_assert_debug(indexlist(i)>=0 && indexlist(i)<no_cols, "Mat<Num_T>::get_cols: index out of range");
00881       copy_vector(no_rows, data+indexlist(i)*no_rows, m.data+i*m.no_rows);
00882     }
00883 
00884     return m;
00885   }
00886 
00887   template<class Num_T> inline
00888   void Mat<Num_T>::set_row(int Index, const Vec<Num_T> &v)
00889   {
00890     it_assert_debug(Index>=0 && Index<no_rows, "Mat<Num_T>::set_row: index out of range");
00891     it_assert_debug(v.length() == no_cols, "Mat<Num_T>::set_row: lengths doesn't match");
00892 
00893     copy_vector(v.size(), v._data(), 1, data+Index, no_rows);
00894   }
00895 
00896   template<class Num_T> inline
00897   void Mat<Num_T>::set_col(int Index, const Vec<Num_T> &v)
00898   {
00899     it_assert_debug(Index>=0 && Index<no_cols, "Mat<Num_T>::set_col: index out of range");
00900     it_assert_debug(v.length() == no_rows, "Mat<Num_T>::set_col: lengths doesn't match");
00901 
00902     copy_vector(v.size(), v._data(), data+Index*no_rows);
00903   }
00904 
00905 
00906   template<class Num_T> inline
00907   void Mat<Num_T>::set_rows(int r, const Mat<Num_T> &m)
00908   {
00909     it_assert_debug((r >= 0) && (r < no_rows),
00910                     "Mat<>::set_rows(): Index out of range");
00911     it_assert_debug(no_cols == m.cols(),
00912                     "Mat<>::set_rows(): Column sizes do not match");
00913     it_assert_debug(m.rows() + r <= no_rows,
00914                     "Mat<>::set_rows(): Not enough rows");
00915 
00916     for (int i = 0; i < m.rows(); ++i) {
00917       copy_vector(no_cols, m.data+i, m.no_rows, data+i+r, no_rows);
00918     }
00919   }
00920 
00921   template<class Num_T> inline
00922   void Mat<Num_T>::set_cols(int c, const Mat<Num_T> &m)
00923   {
00924     it_assert_debug((c >= 0) && (c < no_cols),
00925                     "Mat<>::set_cols(): Index out of range");
00926     it_assert_debug(no_rows == m.rows(),
00927                     "Mat<>::set_cols(): Row sizes do not match");
00928     it_assert_debug(m.cols() + c <= no_cols,
00929                     "Mat<>::set_cols(): Not enough colums");
00930 
00931     for (int i = 0; i < m.cols(); ++i) {
00932       copy_vector(no_rows, m.data+i*no_rows, data+(i+c)*no_rows);
00933     }
00934   }
00935 
00936 
00937   template<class Num_T> inline
00938   void Mat<Num_T>::copy_row(int to, int from)
00939   {
00940     it_assert_debug(to>=0 && from>=0 && to<no_rows && from<no_rows,
00941                "Mat<Num_T>::copy_row: index out of range");
00942 
00943     if (from == to)
00944       return;
00945 
00946     copy_vector(no_cols, data+from, no_rows, data+to, no_rows);
00947   }
00948 
00949   template<class Num_T> inline
00950   void Mat<Num_T>::copy_col(int to, int from)
00951   {
00952     it_assert_debug(to>=0 && from>=0 && to<no_cols && from<no_cols,
00953                "Mat<Num_T>::copy_col: index out of range");
00954 
00955     if (from == to)
00956       return;
00957 
00958     copy_vector(no_rows, data+from*no_rows, data+to*no_rows);
00959   }
00960 
00961   template<class Num_T> inline
00962   void Mat<Num_T>::swap_rows(int r1, int r2)
00963   {
00964     it_assert_debug(r1>=0 && r2>=0 && r1<no_rows && r2<no_rows,
00965                "Mat<Num_T>::swap_rows: index out of range");
00966 
00967     if (r1 == r2)
00968       return;
00969 
00970     swap_vector(no_cols, data+r1, no_rows, data+r2, no_rows);
00971   }
00972 
00973   template<class Num_T> inline
00974   void Mat<Num_T>::swap_cols(int c1, int c2)
00975   {
00976     it_assert_debug(c1>=0 && c2>=0 && c1<no_cols && c2<no_cols,
00977                "Mat<Num_T>::swap_cols: index out of range");
00978 
00979     if (c1 == c2)
00980       return;
00981 
00982     swap_vector(no_rows, data+c1*no_rows, data+c2*no_rows);
00983   }
00984 
00985   template<class Num_T> inline
00986   void Mat<Num_T>::set_submatrix(int r1, int r2, int c1, int c2, const Mat<Num_T> &m)
00987   {
00988 
00989     if (r1 == -1) r1 = no_rows-1;
00990     if (r2 == -1) r2 = no_rows-1;
00991     if (c1 == -1) c1 = no_cols-1;
00992     if (c2 == -1) c2 = no_cols-1;
00993 
00994     it_assert_debug(r1>=0 && r2>=0 && r1<no_rows && r2<no_rows &&
00995                c1>=0 && c2>=0 && c1<no_cols && c2<no_cols, "Mat<Num_T>::set_submatrix(): index out of range");
00996 
00997     it_assert_debug(r2>=r1 && c2>=c1, "Mat<Num_T>::set_submatrix: r2<r1 or c2<c1");
00998     it_assert_debug(m.no_rows == r2-r1+1 && m.no_cols == c2-c1+1, "Mat<Num_T>::set_submatrix(): sizes don't match");
00999 
01000     for (int i=0; i<m.no_cols; i++)
01001       copy_vector(m.no_rows, m.data+i*m.no_rows, data+(c1+i)*no_rows+r1);
01002   }
01003 
01004 
01005 
01006   template<class Num_T> inline
01007   void Mat<Num_T>::set_submatrix(int r, int c, const Mat<Num_T> &m)
01008   {
01009 
01010     it_assert_debug(r>=0 && r+m.no_rows<=no_rows &&
01011                c>=0 && c+m.no_cols<=no_cols, "Mat<Num_T>::set_submatrix(): index out of range");
01012 
01013     for (int i=0; i<m.no_cols; i++)
01014       copy_vector(m.no_rows, m.data+i*m.no_rows, data+(c+i)*no_rows+r);
01015   }
01016 
01017 
01018 
01019   template<class Num_T> inline
01020   void Mat<Num_T>::set_submatrix(int r1, int r2, int c1, int c2, const Num_T t)
01021   {
01022 
01023     if (r1 == -1) r1 = no_rows-1;
01024     if (r2 == -1) r2 = no_rows-1;
01025     if (c1 == -1) c1 = no_cols-1;
01026     if (c2 == -1) c2 = no_cols-1;
01027 
01028     it_assert_debug(r1>=0 && r2>=0 && r1<no_rows && r2<no_rows &&
01029                c1>=0 && c2>=0 && c1<no_cols && c2<no_cols, "Mat<Num_T>::set_submatrix(): index out of range");
01030 
01031     it_assert_debug(r2>=r1 && c2>=c1, "Mat<Num_T>::set_submatrix: r2<r1 or c2<c1");
01032 
01033     int i, j, pos, rows = r2-r1+1;
01034 
01035     for (i=c1; i<=c2; i++) {
01036       pos = i*no_rows+r1;
01037       for (j=0; j<rows; j++) {
01038         data[pos++] = t;
01039       }
01040     }
01041   }
01042 
01043   template<class Num_T> inline
01044   void Mat<Num_T>::del_row(int r)
01045   {
01046     it_assert_debug(r>=0 && r<no_rows, "Mat<Num_T>::del_row(): index out of range");
01047     Mat<Num_T> Temp(*this);
01048     set_size(no_rows-1, no_cols, false);
01049     for (int i=0 ; i < r ; i++) {
01050       copy_vector(no_cols, &Temp.data[i], no_rows+1, &data[i], no_rows);
01051     }
01052     for (int i=r ; i < no_rows ; i++) {
01053       copy_vector(no_cols, &Temp.data[i+1], no_rows+1, &data[i], no_rows);
01054     }
01055 
01056   }
01057 
01058   template<class Num_T> inline
01059   void Mat<Num_T>::del_rows(int r1, int r2)
01060   {
01061     it_assert_debug((r1 >= 0) && (r2 < no_rows) && (r1 <= r2),
01062                     "Mat<Num_T>::del_rows(): indices out of range");
01063 
01064     Mat<Num_T> Temp(*this);
01065     int no_del_rows = r2-r1+1;
01066     set_size(no_rows-no_del_rows, no_cols, false);
01067     for (int i = 0; i < r1 ; ++i) {
01068       copy_vector(no_cols, &Temp.data[i], Temp.no_rows, &data[i], no_rows);
01069     }
01070     for (int i = r2+1; i < Temp.no_rows; ++i) {
01071       copy_vector(no_cols, &Temp.data[i], Temp.no_rows, &data[i-no_del_rows],
01072                   no_rows);
01073     }
01074   }
01075 
01076   template<class Num_T> inline
01077   void Mat<Num_T>::del_col(int c)
01078   {
01079     it_assert_debug(c>=0 && c<no_cols, "Mat<Num_T>::del_col(): index out of range");
01080     Mat<Num_T> Temp(*this);
01081 
01082     set_size(no_rows, no_cols-1, false);
01083     copy_vector(c*no_rows, Temp.data, data);
01084     copy_vector((no_cols - c)*no_rows, &Temp.data[(c+1)*no_rows], &data[c*no_rows]);
01085   }
01086 
01087   template<class Num_T> inline
01088   void Mat<Num_T>::del_cols(int c1, int c2)
01089   {
01090     it_assert_debug(c1>=0 && c2<no_cols && c1<=c2, "Mat<Num_T>::del_cols(): index out of range");
01091     Mat<Num_T> Temp(*this);
01092     int n_deleted_cols = c2-c1+1;
01093     set_size(no_rows, no_cols-n_deleted_cols, false);
01094     copy_vector(c1*no_rows, Temp.data, data);
01095     copy_vector((no_cols-c1)*no_rows, &Temp.data[(c2+1)*no_rows], &data[c1*no_rows]);
01096   }
01097 
01098   template<class Num_T> inline
01099   void Mat<Num_T>::ins_row(int r, const Vec<Num_T> &in)
01100   {
01101     it_assert_debug(r>=0 && r<=no_rows, "Mat<Num_T>::ins_row(): index out of range");
01102     it_assert_debug((in.size() == no_cols) || no_rows==0, "Mat<Num_T>::ins_row(): vector size does not match");
01103 
01104     if (no_cols==0) {
01105       no_cols = in.size();
01106     }
01107 
01108     Mat<Num_T> Temp(*this);
01109     set_size(no_rows+1, no_cols, false);
01110 
01111     for (int i=0 ; i < r ; i++) {
01112       copy_vector(no_cols, &Temp.data[i], no_rows-1, &data[i], no_rows);
01113     }
01114     copy_vector(no_cols, in._data(), 1, &data[r], no_rows);
01115     for (int i=r+1 ; i < no_rows ; i++) {
01116       copy_vector(no_cols, &Temp.data[i-1], no_rows-1, &data[i], no_rows);
01117     }
01118   }
01119 
01120   template<class Num_T> inline
01121   void Mat<Num_T>::ins_col(int c, const Vec<Num_T> &in)
01122   {
01123     it_assert_debug(c>=0 && c<=no_cols, "Mat<Num_T>::ins_col(): index out of range");
01124     it_assert_debug(in.size() == no_rows || no_cols==0, "Mat<Num_T>::ins_col(): vector size does not match");
01125 
01126     if (no_rows==0) {
01127       no_rows = in.size();
01128     }
01129 
01130     Mat<Num_T> Temp(*this);
01131     set_size(no_rows, no_cols+1, false);
01132 
01133     copy_vector(c*no_rows, Temp.data, data);
01134     copy_vector(no_rows, in._data(), &data[c*no_rows]);
01135     copy_vector((no_cols-c-1)*no_rows, &Temp.data[c*no_rows], &data[(c+1)*no_rows]);
01136   }
01137 
01138   template<class Num_T> inline
01139   void Mat<Num_T>::append_row(const Vec<Num_T> &in)
01140   {
01141     ins_row(no_rows, in);
01142   }
01143 
01144   template<class Num_T> inline
01145   void Mat<Num_T>::append_col(const Vec<Num_T> &in)
01146   {
01147     ins_col(no_cols, in);
01148   }
01149 
01150   template<class Num_T>
01151   const Mat<Num_T> Mat<Num_T>::transpose() const
01152   {
01153     Mat<Num_T> temp(no_cols, no_rows);
01154     for (int i = 0; i < no_rows; ++i) {
01155       copy_vector(no_cols, &data[i], no_rows, &temp.data[i * no_cols], 1);
01156     }
01157     return temp;
01158   }
01159 
01161   template<>
01162   const cmat Mat<std::complex<double> >::hermitian_transpose() const;
01164 
01165   template<class Num_T>
01166   const Mat<Num_T> Mat<Num_T>::hermitian_transpose() const
01167   {
01168     Mat<Num_T> temp(no_cols, no_rows);
01169     for (int i = 0; i < no_rows; ++i) {
01170       copy_vector(no_cols, &data[i], no_rows, &temp.data[i * no_cols], 1);
01171     }
01172     return temp;
01173   }
01174 
01175   template<class Num_T>
01176   const Mat<Num_T> concat_horizontal(const Mat<Num_T> &m1, const Mat<Num_T> &m2)
01177   {
01178     it_assert_debug(m1.no_rows == m2.no_rows,
01179                     "Mat<Num_T>::concat_horizontal(): wrong sizes");
01180     int no_rows = m1.no_rows;
01181     Mat<Num_T> temp(no_rows, m1.no_cols + m2.no_cols);
01182     for (int i = 0; i < m1.no_cols; ++i) {
01183       copy_vector(no_rows, &m1.data[i * no_rows], &temp.data[i * no_rows]);
01184     }
01185     for (int i = 0; i < m2.no_cols; ++i) {
01186       copy_vector(no_rows, &m2.data[i * no_rows], &temp.data[(m1.no_cols + i)
01187                                                              * no_rows]);
01188     }
01189     return temp;
01190   }
01191 
01192   template<class Num_T>
01193   const Mat<Num_T> concat_vertical(const Mat<Num_T> &m1, const Mat<Num_T> &m2)
01194   {
01195     it_assert_debug(m1.no_cols == m2.no_cols,
01196                     "Mat<Num_T>::concat_vertical; wrong sizes");
01197     int no_cols = m1.no_cols;
01198     Mat<Num_T> temp(m1.no_rows + m2.no_rows, no_cols);
01199     for (int i = 0; i < no_cols; ++i) {
01200       copy_vector(m1.no_rows, &m1.data[i * m1.no_rows],
01201                   &temp.data[i * temp.no_rows]);
01202       copy_vector(m2.no_rows, &m2.data[i * m2.no_rows],
01203                   &temp.data[i * temp.no_rows + m1.no_rows]);
01204     }
01205     return temp;
01206   }
01207 
01208   template<class Num_T> inline
01209   Mat<Num_T>& Mat<Num_T>::operator=(Num_T t)
01210   {
01211     for (int i=0; i<datasize; i++)
01212       data[i] = t;
01213     return *this;
01214   }
01215 
01216   template<class Num_T> inline
01217   Mat<Num_T>& Mat<Num_T>::operator=(const Mat<Num_T> &m)
01218   {
01219     if (this != &m) {
01220       set_size(m.no_rows,m.no_cols, false);
01221       if (m.datasize != 0)
01222         copy_vector(m.datasize, m.data, data);
01223     }
01224     return *this;
01225   }
01226 
01227   template<class Num_T> inline
01228   Mat<Num_T>& Mat<Num_T>::operator=(const Vec<Num_T> &v)
01229   {
01230     it_assert_debug((no_rows == 1 && no_cols == v.size())
01231                     || (no_cols == 1 && no_rows == v.size()),
01232                     "Mat<Num_T>::operator=(): Wrong size of the argument");
01233     set_size(v.size(), 1, false);
01234     copy_vector(v.size(), v._data(), data);
01235     return *this;
01236   }
01237 
01238   template<class Num_T> inline
01239   Mat<Num_T>& Mat<Num_T>::operator=(const char *values)
01240   {
01241     set(values);
01242     return *this;
01243   }
01244 
01245   //-------------------- Templated friend functions --------------------------
01246 
01247   template<class Num_T> inline
01248   Mat<Num_T>& Mat<Num_T>::operator+=(const Mat<Num_T> &m)
01249   {
01250     if (datasize == 0)
01251       operator=(m);
01252     else {
01253       int i, j, m_pos=0, pos=0;
01254       it_assert_debug(m.no_rows==no_rows && m.no_cols==no_cols,"Mat<Num_T>::operator+=: wrong sizes");
01255       for (i=0; i<no_cols; i++) {
01256         for (j=0; j<no_rows; j++)
01257           data[pos+j] += m.data[m_pos+j];
01258         pos += no_rows;
01259         m_pos += m.no_rows;
01260       }
01261     }
01262     return *this;
01263   }
01264 
01265   template<class Num_T> inline
01266   Mat<Num_T>& Mat<Num_T>::operator+=(Num_T t)
01267   {
01268     for (int i=0; i<datasize; i++)
01269       data[i] += t;
01270     return *this;
01271   }
01272 
01273   template<class Num_T> inline
01274   const Mat<Num_T> operator+(const Mat<Num_T> &m1, const Mat<Num_T> &m2)
01275   {
01276     Mat<Num_T> r(m1.no_rows, m1.no_cols);
01277     int i, j, m1_pos=0, m2_pos=0, r_pos=0;
01278 
01279     it_assert_debug(m1.no_rows==m2.no_rows && m1.no_cols == m2.no_cols, "Mat<Num_T>::operator+: wrong sizes");
01280 
01281     for (i=0; i<r.no_cols; i++) {
01282       for (j=0; j<r.no_rows; j++)
01283         r.data[r_pos+j] = m1.data[m1_pos+j] + m2.data[m2_pos+j];
01284       // next column
01285       m1_pos += m1.no_rows;
01286       m2_pos += m2.no_rows;
01287       r_pos += r.no_rows;
01288     }
01289 
01290     return r;
01291   }
01292 
01293 
01294   template<class Num_T> inline
01295   const Mat<Num_T> operator+(const Mat<Num_T> &m, Num_T t)
01296   {
01297     Mat<Num_T> r(m.no_rows, m.no_cols);
01298 
01299     for (int i=0; i<r.datasize; i++)
01300       r.data[i] = m.data[i] + t;
01301 
01302     return r;
01303   }
01304 
01305   template<class Num_T> inline
01306   const Mat<Num_T> operator+(Num_T t, const Mat<Num_T> &m)
01307   {
01308     Mat<Num_T> r(m.no_rows, m.no_cols);
01309 
01310     for (int i=0; i<r.datasize; i++)
01311       r.data[i] = t + m.data[i];
01312 
01313     return r;
01314   }
01315 
01316   template<class Num_T> inline
01317   Mat<Num_T>& Mat<Num_T>::operator-=(const Mat<Num_T> &m)
01318   {
01319     int i,j, m_pos=0, pos=0;
01320 
01321     if (datasize == 0) {
01322       set_size(m.no_rows, m.no_cols, false);
01323       for (i=0; i<no_cols; i++) {
01324         for (j=0; j<no_rows; j++)
01325           data[pos+j] = -m.data[m_pos+j];
01326         // next column
01327         m_pos += m.no_rows;
01328         pos += no_rows;
01329       }
01330     } else {
01331       it_assert_debug(m.no_rows==no_rows && m.no_cols==no_cols,"Mat<Num_T>::operator-=: wrong sizes");
01332       for (i=0; i<no_cols; i++) {
01333         for (j=0; j<no_rows; j++)
01334           data[pos+j] -= m.data[m_pos+j];
01335         // next column
01336         m_pos += m.no_rows;
01337         pos += no_rows;
01338       }
01339     }
01340     return *this;
01341   }
01342 
01343   template<class Num_T> inline
01344   const Mat<Num_T> operator-(const Mat<Num_T> &m1, const Mat<Num_T> &m2)
01345   {
01346     Mat<Num_T> r(m1.no_rows, m1.no_cols);
01347     int i, j, m1_pos=0, m2_pos=0, r_pos=0;
01348 
01349     it_assert_debug(m1.no_rows==m2.no_rows && m1.no_cols == m2.no_cols, "Mat<Num_T>::operator-: wrong sizes");
01350 
01351     for (i=0; i<r.no_cols; i++) {
01352       for (j=0; j<r.no_rows; j++)
01353         r.data[r_pos+j] = m1.data[m1_pos+j] - m2.data[m2_pos+j];
01354       // next column
01355       m1_pos += m1.no_rows;
01356       m2_pos += m2.no_rows;
01357       r_pos += r.no_rows;
01358     }
01359 
01360     return r;
01361   }
01362 
01363   template<class Num_T> inline
01364   Mat<Num_T>& Mat<Num_T>::operator-=(Num_T t)
01365   {
01366     for (int i=0; i<datasize; i++)
01367       data[i] -= t;
01368     return *this;
01369   }
01370 
01371   template<class Num_T> inline
01372   const Mat<Num_T> operator-(const Mat<Num_T> &m, Num_T t)
01373   {
01374     Mat<Num_T> r(m.no_rows, m.no_cols);
01375     int i, j, m_pos=0, r_pos=0;
01376 
01377     for (i=0; i<r.no_cols; i++) {
01378       for (j=0; j<r.no_rows; j++)
01379         r.data[r_pos+j] = m.data[m_pos+j] - t;
01380       // next column
01381       m_pos += m.no_rows;
01382       r_pos += r.no_rows;
01383     }
01384 
01385     return r;
01386   }
01387 
01388   template<class Num_T> inline
01389   const Mat<Num_T> operator-(Num_T t, const Mat<Num_T> &m)
01390   {
01391     Mat<Num_T> r(m.no_rows, m.no_cols);
01392     int i, j, m_pos=0, r_pos=0;
01393 
01394     for (i=0; i<r.no_cols; i++) {
01395       for (j=0; j<r.no_rows; j++)
01396         r.data[r_pos+j] = t - m.data[m_pos+j];
01397       // next column
01398       m_pos += m.no_rows;
01399       r_pos += r.no_rows;
01400     }
01401 
01402     return r;
01403   }
01404 
01405   template<class Num_T> inline
01406   const Mat<Num_T> operator-(const Mat<Num_T> &m)
01407   {
01408     Mat<Num_T> r(m.no_rows, m.no_cols);
01409     int i, j, m_pos=0, r_pos=0;
01410 
01411     for (i=0; i<r.no_cols; i++) {
01412       for (j=0; j<r.no_rows; j++)
01413         r.data[r_pos+j] = -m.data[m_pos+j];
01414       // next column
01415       m_pos += m.no_rows;
01416       r_pos += r.no_rows;
01417     }
01418 
01419     return r;
01420   }
01421 
01422 #if defined(HAVE_BLAS)
01423   template<> mat& Mat<double>::operator*=(const mat &m);
01424   template<> cmat& Mat<std::complex<double> >::operator*=(const cmat &m);
01425 #endif
01426 
01427   template<class Num_T> inline
01428   Mat<Num_T>& Mat<Num_T>::operator*=(const Mat<Num_T> &m)
01429   {
01430     it_assert_debug(no_cols == m.no_rows,"Mat<Num_T>::operator*=: wrong sizes");
01431     Mat<Num_T> r(no_rows, m.no_cols);
01432 
01433     Num_T tmp;
01434 
01435     int i,j,k, r_pos=0, pos=0, m_pos=0;
01436 
01437     for (i=0; i<r.no_cols; i++) {
01438       for (j=0; j<r.no_rows; j++) {
01439         tmp = Num_T(0);
01440         pos = 0;
01441         for (k=0; k<no_cols; k++) {
01442           tmp += data[pos+j] * m.data[m_pos+k];
01443           pos += no_rows;
01444         }
01445         r.data[r_pos+j] = tmp;
01446       }
01447       r_pos += r.no_rows;
01448       m_pos += m.no_rows;
01449     }
01450     operator=(r); // time consuming
01451     return *this;
01452   }
01453 
01454   template<class Num_T> inline
01455   Mat<Num_T>& Mat<Num_T>::operator*=(Num_T t)
01456   {
01457     scal_vector(datasize, t, data);
01458     return *this;
01459   }
01460 
01461 #if defined(HAVE_BLAS)
01462   template<> const mat operator*(const mat &m1, const mat &m2);
01463   template<> const cmat operator*(const cmat &m1, const cmat &m2);
01464 #endif
01465 
01466 
01467   template<class Num_T>
01468   const Mat<Num_T> operator*(const Mat<Num_T> &m1, const Mat<Num_T> &m2)
01469   {
01470     it_assert_debug(m1.no_cols == m2.no_rows,"Mat<Num_T>::operator*: wrong sizes");
01471     Mat<Num_T> r(m1.no_rows, m2.no_cols);
01472 
01473     Num_T tmp;
01474     int i, j, k;
01475     Num_T *tr=r.data, *t1, *t2=m2.data;
01476 
01477     for (i=0; i<r.no_cols; i++) {
01478       for (j=0; j<r.no_rows; j++) {
01479         tmp = Num_T(0); t1 = m1.data+j;
01480         for (k=m1.no_cols; k>0; k--) {
01481           tmp += *(t1) * *(t2++);
01482           t1 += m1.no_rows;
01483         }
01484         *(tr++) = tmp; t2 -= m2.no_rows;
01485       }
01486       t2 += m2.no_rows;
01487     }
01488 
01489     return r;
01490   }
01491 
01492 #if defined(HAVE_BLAS)
01493   template<> const vec operator*(const mat &m, const vec &v);
01494   template<> const cvec operator*(const cmat &m, const cvec &v);
01495 #endif
01496 
01497   template<class Num_T>
01498   const Vec<Num_T> operator*(const Mat<Num_T> &m, const Vec<Num_T> &v)
01499   {
01500     it_assert_debug(m.no_cols == v.size(),"Mat<Num_T>::operator*: wrong sizes");
01501     Vec<Num_T> r(m.no_rows);
01502     int i, k, m_pos;
01503 
01504     for (i=0; i<m.no_rows; i++) {
01505       r(i) = Num_T(0);
01506       m_pos = 0;
01507       for (k=0; k<m.no_cols; k++) {
01508         r(i) += m.data[m_pos+i] * v(k);
01509         m_pos += m.no_rows;
01510       }
01511     }
01512 
01513     return r;
01514   }
01515 
01516   template<class Num_T> inline
01517   const Mat<Num_T> operator*(const Vec<Num_T> &v, const Mat<Num_T> &m)
01518   {
01519     it_assert((m.no_rows == 1),"Mat<Num_T>::operator*(): wrong sizes");
01520     it_warning("Mat<Num_T>::operator*(v, m): This operator is deprecated. "
01521                "Please use outer_product(v, m.get_row(0)) instead.");
01522     return outer_product(v, m.get_row(0));
01523   }
01524 
01525   template<class Num_T> inline
01526   const Mat<Num_T> operator*(const Mat<Num_T> &m, Num_T t)
01527   {
01528     Mat<Num_T> r(m.no_rows, m.no_cols);
01529 
01530     for (int i=0; i<r.datasize; i++)
01531       r.data[i] = m.data[i] * t;
01532 
01533     return r;
01534   }
01535 
01536   template<class Num_T> inline
01537   const Mat<Num_T> operator*(Num_T t, const Mat<Num_T> &m)
01538   {
01539     return operator*(m, t);
01540   }
01541 
01542   template<class Num_T> inline
01543   const Mat<Num_T> elem_mult(const Mat<Num_T> &A, const Mat<Num_T> &B)
01544   {
01545     Mat<Num_T> out;
01546     elem_mult_out(A,B,out);
01547     return out;
01548   }
01549 
01550   template<class Num_T> inline
01551   void elem_mult_out(const Mat<Num_T> &A, const Mat<Num_T> &B, Mat<Num_T> &out)
01552   {
01553     it_assert_debug( (A.no_rows==B.no_rows) && (A.no_cols==B.no_cols), "Mat<Num_T>::elem_mult_out: wrong sizes");
01554 
01555     if( (out.no_rows != A.no_rows) || (out.no_cols != A.no_cols) )
01556       out.set_size(A.no_rows, A.no_cols);
01557 
01558     for(int i=0; i<out.datasize; i++)
01559       out.data[i] = A.data[i] * B.data[i];
01560   }
01561 
01562   template<class Num_T> inline
01563   void elem_mult_out(const Mat<Num_T> &A, const Mat<Num_T> &B, const Mat<Num_T> &C, Mat<Num_T> &out)
01564   {
01565     it_assert_debug( (A.no_rows==B.no_rows==C.no_rows) \
01566                 && (A.no_cols==B.no_cols==C.no_cols), \
01567                 "Mat<Num_T>::elem_mult_out: wrong sizes" );
01568 
01569     if( (out.no_rows != A.no_rows) || (out.no_cols != A.no_cols) )
01570       out.set_size(A.no_rows, A.no_cols);
01571 
01572     for(int i=0; i<out.datasize; i++)
01573       out.data[i] = A.data[i] * B.data[i] * C.data[i];
01574   }
01575 
01576   template<class Num_T> inline
01577   void elem_mult_out(const Mat<Num_T> &A, const Mat<Num_T> &B, const Mat<Num_T> &C, const Mat<Num_T> &D, Mat<Num_T> &out)
01578   {
01579     it_assert_debug( (A.no_rows==B.no_rows==C.no_rows==D.no_rows) \
01580                 && (A.no_cols==B.no_cols==C.no_cols==D.no_cols), \
01581                 "Mat<Num_T>::elem_mult_out: wrong sizes" );
01582     if( (out.no_rows != A.no_rows) || (out.no_cols != A.no_cols) )
01583       out.set_size(A.no_rows, A.no_cols);
01584 
01585     for(int i=0; i<out.datasize; i++)
01586       out.data[i] = A.data[i] * B.data[i] * C.data[i] * D.data[i];
01587   }
01588 
01589   template<class Num_T> inline
01590   void elem_mult_inplace(const Mat<Num_T> &A, Mat<Num_T> &B)
01591   {
01592     it_assert_debug( (A.no_rows==B.no_rows) && (A.no_cols==B.no_cols), \
01593                 "Mat<Num_T>::elem_mult_inplace: wrong sizes" );
01594 
01595     for(int i=0; i<B.datasize; i++)
01596       B.data[i] *= A.data[i];
01597   }
01598 
01599   template<class Num_T> inline
01600   Num_T elem_mult_sum(const Mat<Num_T> &A, const Mat<Num_T> &B)
01601   {
01602     it_assert_debug( (A.no_rows==B.no_rows) && (A.no_cols==B.no_cols), "Mat<Num_T>::elem_mult_sum: wrong sizes" );
01603 
01604     Num_T acc = 0;
01605 
01606     for(int i=0; i<A.datasize; i++)
01607       acc += A.data[i] * B.data[i];
01608 
01609     return acc;
01610   }
01611 
01612   template<class Num_T> inline
01613   Mat<Num_T>& Mat<Num_T>::operator/=(Num_T t)
01614   {
01615     for (int i=0; i<datasize; i++)
01616       data[i] /= t;
01617     return *this;
01618   }
01619 
01620   template<class Num_T> inline
01621   const Mat<Num_T> operator/(const Mat<Num_T> &m, Num_T t)
01622   {
01623     Mat<Num_T> r(m.no_rows, m.no_cols);
01624 
01625     for (int i=0; i<r.datasize; i++)
01626       r.data[i] = m.data[i] / t;
01627 
01628     return r;
01629   }
01630 
01631   template<class Num_T> inline
01632   Mat<Num_T>& Mat<Num_T>::operator/=(const Mat<Num_T> &m)
01633   {
01634     it_assert_debug(m.no_rows==no_rows && m.no_cols==no_cols, "Mat<Num_T>::operator/=: wrong sizes");
01635 
01636     for (int i=0; i<datasize; i++)
01637       data[i] /= m.data[i];
01638     return *this;
01639   }
01640 
01641   template<class Num_T> inline
01642   const Mat<Num_T> elem_div(const Mat<Num_T> &A, const Mat<Num_T> &B)
01643   {
01644     Mat<Num_T> out;
01645     elem_div_out(A,B,out);
01646     return out;
01647   }
01648 
01649   template<class Num_T> inline
01650   void elem_div_out(const Mat<Num_T> &A, const Mat<Num_T> &B, Mat<Num_T> &out)
01651   {
01652     it_assert_debug( (A.no_rows==B.no_rows) && (A.no_cols==B.no_cols), "Mat<Num_T>::elem_div_out: wrong sizes");
01653 
01654     if( (out.no_rows != A.no_rows) || (out.no_cols != A.no_cols) )
01655       out.set_size(A.no_rows, A.no_cols);
01656 
01657     for(int i=0; i<out.datasize; i++)
01658       out.data[i] = A.data[i] / B.data[i];
01659   }
01660 
01661   template<class Num_T> inline
01662   Num_T elem_div_sum(const Mat<Num_T> &A, const Mat<Num_T> &B)
01663   {
01664     it_assert_debug( (A.no_rows==B.no_rows) && (A.no_cols==B.no_cols), "Mat<Num_T>::elem_div_sum: wrong sizes" );
01665 
01666     Num_T acc = 0;
01667 
01668     for(int i=0; i<A.datasize; i++)
01669       acc += A.data[i] / B.data[i];
01670 
01671     return acc;
01672   }
01673 
01674   template<class Num_T>
01675   bool Mat<Num_T>::operator==(const Mat<Num_T> &m) const
01676   {
01677     if (no_rows!=m.no_rows || no_cols != m.no_cols) return false;
01678     for (int i=0;i<datasize;i++) {
01679       if (data[i]!=m.data[i]) return false;
01680     }
01681     return true;
01682   }
01683 
01684   template<class Num_T>
01685   bool Mat<Num_T>::operator!=(const Mat<Num_T> &m) const
01686   {
01687     if (no_rows != m.no_rows || no_cols != m.no_cols) return true;
01688     for (int i=0;i<datasize;i++) {
01689       if (data[i]!=m.data[i]) return true;
01690     }
01691     return false;
01692   }
01693 
01694   template <class Num_T>
01695   std::ostream &operator<<(std::ostream &os, const Mat<Num_T> &m)
01696   {
01697     int i;
01698 
01699     switch (m.rows()) {
01700     case 0 :
01701       os << "[]";
01702       break;
01703     case 1 :
01704       os << '[' << m.get_row(0) << ']';
01705       break;
01706     default:
01707       os << '[' << m.get_row(0) << std::endl;
01708       for (i=1; i<m.rows()-1; i++)
01709               os << ' ' << m.get_row(i) << std::endl;
01710       os << ' ' << m.get_row(m.rows()-1) << ']';
01711     }
01712 
01713     return os;
01714   }
01715 
01716   template <class Num_T>
01717   std::istream &operator>>(std::istream &is, Mat<Num_T> &m)
01718   {
01719     std::ostringstream buffer;
01720     bool started = false;
01721     bool finished = false;
01722     bool brackets = false;
01723     bool within_double_brackets = false;
01724     char c;
01725 
01726     while (!finished) {
01727       if (is.eof()) {
01728         finished = true;
01729       } else {
01730         c = is.get();
01731 
01732         if (is.eof() || (c == '\n')) {
01733           if (brackets) {
01734             // Right bracket missing
01735             is.setstate(std::ios_base::failbit);
01736             finished = true;
01737           } else if (!((c == '\n') && !started)) {
01738             finished = true;
01739           }
01740         } else if ((c == ' ') || (c == '\t')) {
01741           if (started) {
01742             buffer << ' ';
01743           }
01744         } else if (c == '[') {
01745           if ((started && !brackets) || within_double_brackets) {
01746             // Unexpected left bracket
01747             is.setstate(std::ios_base::failbit);
01748             finished = true;
01749           } else if (!started) {
01750             started = true;
01751             brackets = true;
01752           } else {
01753             within_double_brackets = true;
01754           }
01755         } else if (c == ']') {
01756           if (!started || !brackets) {
01757             // Unexpected right bracket
01758             is.setstate(std::ios_base::failbit);
01759             finished = true;
01760           } else if (within_double_brackets) {
01761             within_double_brackets = false;
01762             buffer << ';';
01763           } else {
01764             finished = true;
01765           }
01766           while (!is.eof() && (((c = is.peek()) == ' ') || (c == '\t'))) {
01767             is.get();
01768           }
01769           if (!is.eof() && (c == '\n')) {
01770             is.get();
01771           }
01772         } else {
01773           started = true;
01774           buffer << c;
01775         }
01776       }
01777     }
01778 
01779     if (!started) {
01780       m.set_size(0, false);
01781     } else {
01782       m.set(buffer.str());
01783     }
01784 
01785     return is;
01786   }
01787 
01789 
01790   // ---------------------------------------------------------------------
01791   // Instantiations
01792   // ---------------------------------------------------------------------
01793 
01794 #ifdef HAVE_EXTERN_TEMPLATE
01795 
01796   // class instantiations
01797 
01798   extern template class Mat<double>;
01799   extern template class Mat<std::complex<double> >;
01800   extern template class Mat<int>;
01801   extern template class Mat<short int>;
01802   extern template class Mat<bin>;
01803 
01804   // addition operators
01805 
01806   extern template const mat operator+(const mat &m1, const mat &m2);
01807   extern template const cmat operator+(const cmat &m1, const cmat &m2);
01808   extern template const imat operator+(const imat &m1, const imat &m2);
01809   extern template const smat operator+(const smat &m1, const smat &m2);
01810   extern template const bmat operator+(const bmat &m1, const bmat &m2);
01811 
01812   extern template const mat operator+(const mat &m, double t);
01813   extern template const cmat operator+(const cmat &m, std::complex<double> t);
01814   extern template const imat operator+(const imat &m, int t);
01815   extern template const smat operator+(const smat &m, short t);
01816   extern template const bmat operator+(const bmat &m, bin t);
01817 
01818   extern template const mat operator+(double t, const mat &m);
01819   extern template const cmat operator+(std::complex<double> t, const cmat &m);
01820   extern template const imat operator+(int t, const imat &m);
01821   extern template const smat operator+(short t, const smat &m);
01822   extern template const bmat operator+(bin t, const bmat &m);
01823 
01824   // subraction operators
01825 
01826   extern template const mat operator-(const mat &m1, const mat &m2);
01827   extern template const cmat operator-(const cmat &m1, const cmat &m2);
01828   extern template const imat operator-(const imat &m1, const imat &m2);
01829   extern template const smat operator-(const smat &m1, const smat &m2);
01830   extern template const bmat operator-(const bmat &m1, const bmat &m2);
01831 
01832   extern template const mat operator-(const mat &m, double t);
01833   extern template const cmat operator-(const cmat &m, std::complex<double> t);
01834   extern template const imat operator-(const imat &m, int t);
01835   extern template const smat operator-(const smat &m, short t);
01836   extern template const bmat operator-(const bmat &m, bin t);
01837 
01838   extern template const mat operator-(double t, const mat &m);
01839   extern template const cmat operator-(std::complex<double> t, const cmat &m);
01840   extern template const imat operator-(int t, const imat &m);
01841   extern template const smat operator-(short t, const smat &m);
01842   extern template const bmat operator-(bin t, const bmat &m);
01843 
01844   // unary minus
01845 
01846   extern template const mat operator-(const mat &m);
01847   extern template const cmat operator-(const cmat &m);
01848   extern template const imat operator-(const imat &m);
01849   extern template const smat operator-(const smat &m);
01850   extern template const bmat operator-(const bmat &m);
01851 
01852   // multiplication operators
01853 
01854 #if !defined(HAVE_BLAS)
01855   extern template const mat operator*(const mat &m1, const mat &m2);
01856   extern template const cmat operator*(const cmat &m1, const cmat &m2);
01857 #endif
01858   extern template const imat operator*(const imat &m1, const imat &m2);
01859   extern template const smat operator*(const smat &m1, const smat &m2);
01860   extern template const bmat operator*(const bmat &m1, const bmat &m2);
01861 
01862 #if !defined(HAVE_BLAS)
01863   extern template const vec operator*(const mat &m, const vec &v);
01864   extern template const cvec operator*(const cmat &m, const cvec &v);
01865 #endif
01866   extern template const ivec operator*(const imat &m, const ivec &v);
01867   extern template const svec operator*(const smat &m, const svec &v);
01868   extern template const bvec operator*(const bmat &m, const bvec &v);
01869 
01870   extern template const mat operator*(const vec &v, const mat &m);
01871   extern template const cmat operator*(const cvec &v, const cmat &m);
01872   extern template const imat operator*(const ivec &v, const imat &m);
01873   extern template const smat operator*(const svec &v, const smat &m);
01874   extern template const bmat operator*(const bvec &v, const bmat &m);
01875 
01876   extern template const mat operator*(const mat &m, double t);
01877   extern template const cmat operator*(const cmat &m, std::complex<double> t);
01878   extern template const imat operator*(const imat &m, int t);
01879   extern template const smat operator*(const smat &m, short t);
01880   extern template const bmat operator*(const bmat &m, bin t);
01881 
01882   extern template const mat operator*(double t, const mat &m);
01883   extern template const cmat operator*(std::complex<double> t, const cmat &m);
01884   extern template const imat operator*(int t, const imat &m);
01885   extern template const smat operator*(short t, const smat &m);
01886   extern template const bmat operator*(bin t, const bmat &m);
01887 
01888   // elementwise multiplication
01889 
01890   extern template const mat elem_mult(const mat &A, const mat &B);
01891   extern template const cmat elem_mult(const cmat &A, const cmat &B);
01892   extern template const imat elem_mult(const imat &A, const imat &B);
01893   extern template const smat elem_mult(const smat &A, const smat &B);
01894   extern template const bmat elem_mult(const bmat &A, const bmat &B);
01895 
01896   extern template void elem_mult_out(const mat &A, const mat &B, mat &out);
01897   extern template void elem_mult_out(const cmat &A, const cmat &B, cmat &out);
01898   extern template void elem_mult_out(const imat &A, const imat &B, imat &out);
01899   extern template void elem_mult_out(const smat &A, const smat &B, smat &out);
01900   extern template void elem_mult_out(const bmat &A, const bmat &B, bmat &out);
01901 
01902   extern template void elem_mult_out(const mat &A, const mat &B, const mat &C,
01903                                      mat &out);
01904   extern template void elem_mult_out(const cmat &A, const cmat &B,
01905                                      const cmat &C, cmat &out);
01906   extern template void elem_mult_out(const imat &A, const imat &B,
01907                                      const imat &C, imat &out);
01908   extern template void elem_mult_out(const smat &A, const smat &B,
01909                                      const smat &C, smat &out);
01910   extern template void elem_mult_out(const bmat &A, const bmat &B,
01911                                      const bmat &C, bmat &out);
01912 
01913   extern template void elem_mult_out(const mat &A, const mat &B, const mat &C,
01914                                      const mat &D, mat &out);
01915   extern template void elem_mult_out(const cmat &A, const cmat &B,
01916                                      const cmat &C, const cmat &D, cmat &out);
01917   extern template void elem_mult_out(const imat &A, const imat &B,
01918                                      const imat &C, const imat &D, imat &out);
01919   extern template void elem_mult_out(const smat &A, const smat &B,
01920                                      const smat &C, const smat &D, smat &out);
01921   extern template void elem_mult_out(const bmat &A, const bmat &B,
01922                                      const bmat &C, const bmat &D, bmat &out);
01923 
01924   extern template void elem_mult_inplace(const mat &A, mat &B);
01925   extern template void elem_mult_inplace(const cmat &A, cmat &B);
01926   extern template void elem_mult_inplace(const imat &A, imat &B);
01927   extern template void elem_mult_inplace(const smat &A, smat &B);
01928   extern template void elem_mult_inplace(const bmat &A, bmat &B);
01929 
01930   extern template double elem_mult_sum(const mat &A, const mat &B);
01931   extern template std::complex<double> elem_mult_sum(const cmat &A,
01932                                                      const cmat &B);
01933   extern template int elem_mult_sum(const imat &A, const imat &B);
01934   extern template short elem_mult_sum(const smat &A, const smat &B);
01935   extern template bin elem_mult_sum(const bmat &A, const bmat &B);
01936 
01937   // division operator
01938 
01939   extern template const mat operator/(const mat &m, double t);
01940   extern template const cmat operator/(const cmat &m, std::complex<double> t);
01941   extern template const imat operator/(const imat &m, int t);
01942   extern template const smat operator/(const smat &m, short t);
01943   extern template const bmat operator/(const bmat &m, bin t);
01944 
01945   // elementwise division
01946 
01947   extern template const mat elem_div(const mat &A, const mat &B);
01948   extern template const cmat elem_div(const cmat &A, const cmat &B);
01949   extern template const imat elem_div(const imat &A, const imat &B);
01950   extern template const smat elem_div(const smat &A, const smat &B);
01951   extern template const bmat elem_div(const bmat &A, const bmat &B);
01952 
01953   extern template void elem_div_out(const mat &A, const mat &B, mat &out);
01954   extern template void elem_div_out(const cmat &A, const cmat &B, cmat &out);
01955   extern template void elem_div_out(const imat &A, const imat &B, imat &out);
01956   extern template void elem_div_out(const smat &A, const smat &B, smat &out);
01957   extern template void elem_div_out(const bmat &A, const bmat &B, bmat &out);
01958 
01959   extern template double elem_div_sum(const mat &A, const mat &B);
01960   extern template std::complex<double> elem_div_sum(const cmat &A,
01961                                                     const cmat &B);
01962   extern template int elem_div_sum(const imat &A, const imat &B);
01963   extern template short elem_div_sum(const smat &A, const smat &B);
01964   extern template bin elem_div_sum(const bmat &A, const bmat &B);
01965 
01966   // concatenation
01967 
01968   extern template const mat concat_horizontal(const mat &m1, const mat &m2);
01969   extern template const cmat concat_horizontal(const cmat &m1, const cmat &m2);
01970   extern template const imat concat_horizontal(const imat &m1, const imat &m2);
01971   extern template const smat concat_horizontal(const smat &m1, const smat &m2);
01972   extern template const bmat concat_horizontal(const bmat &m1, const bmat &m2);
01973 
01974   extern template const mat concat_vertical(const mat &m1, const mat &m2);
01975   extern template const cmat concat_vertical(const cmat &m1, const cmat &m2);
01976   extern template const imat concat_vertical(const imat &m1, const imat &m2);
01977   extern template const smat concat_vertical(const smat &m1, const smat &m2);
01978   extern template const bmat concat_vertical(const bmat &m1, const bmat &m2);
01979 
01980   // I/O streams
01981 
01982   extern template std::ostream &operator<<(std::ostream &os, const mat  &m);
01983   extern template std::ostream &operator<<(std::ostream &os, const cmat &m);
01984   extern template std::ostream &operator<<(std::ostream &os, const imat  &m);
01985   extern template std::ostream &operator<<(std::ostream &os, const smat  &m);
01986   extern template std::ostream &operator<<(std::ostream &os, const bmat  &m);
01987 
01988   extern template std::istream &operator>>(std::istream &is, mat  &m);
01989   extern template std::istream &operator>>(std::istream &is, cmat &m);
01990   extern template std::istream &operator>>(std::istream &is, imat  &m);
01991   extern template std::istream &operator>>(std::istream &is, smat  &m);
01992   extern template std::istream &operator>>(std::istream &is, bmat  &m);
01993 
01994 #endif // HAVE_EXTERN_TEMPLATE
01995 
01997 
01998 } // namespace itpp
01999 
02000 #endif // #ifndef MAT_H
SourceForge Logo

Generated on Thu Apr 24 13:38:58 2008 for IT++ by Doxygen 1.5.5