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

Generated on Sun Dec 9 17:31:00 2007 for IT++ by Doxygen 1.5.4