src/math/DynamicRectMatrix.hpp

/*
 * Copyright (c) 2005 The University of Notre Dame. All Rights Reserved.
 *
 * The University of Notre Dame grants you ("Licensee") a
 * non-exclusive, royalty free, license to use, modify and
 * redistribute this software in source and binary code form, provided
 * that the following conditions are met:
 *
 * 1. Acknowledgement of the program authors must be made in any
 *    publication of scientific results based in part on use of the
 *    program.  An acceptable form of acknowledgement is citation of
 *    the article in which the program was described (Matthew
 *    A. Meineke, Charles F. Vardeman II, Teng Lin, Christopher
 *    J. Fennell and J. Daniel Gezelter, "OOPSE: An Object-Oriented
 *    Parallel Simulation Engine for Molecular Dynamics,"
 *    J. Comput. Chem. 26, pp. 252-271 (2005))
 *
 * 2. Redistributions of source code must retain the above copyright
 *    notice, this list of conditions and the following disclaimer.
 *
 * 3. Redistributions in binary form must reproduce the above copyright
 *    notice, this list of conditions and the following disclaimer in the
 *    documentation and/or other materials provided with the
 *    distribution.
 *
 * This software is provided "AS IS," without a warranty of any
 * kind. All express or implied conditions, representations and
 * warranties, including any implied warranty of merchantability,
 * fitness for a particular purpose or non-infringement, are hereby
 * excluded.  The University of Notre Dame and its licensors shall not
 * be liable for any damages suffered by licensee as a result of
 * using, modifying or distributing the software or its
 * derivatives. In no event will the University of Notre Dame or its
 * licensors be liable for any lost revenue, profit or data, or for
 * direct, indirect, special, consequential, incidental or punitive
 * damages, however caused and regardless of the theory of liability,
 * arising out of the use of or inability to use software, even if the
 * University of Notre Dame has been advised of the possibility of
 * such damages.
 */
 
/**
 * @file DynamicRectMatrix.hpp
 * @author Teng Lin
 * @date 10/11/2004
 * @version 1.0
 */

#ifndef MATH_DYNAMICRECTMATRIX_HPP
#define MATH_DYNAMICRECTMATRIX_HPP
#include <math.h>
#include <cmath>
#include "math/DynamicVector.hpp"

namespace oopse {

  /**
   * @class DynamicRectMatrix DynamicRectMatrix.hpp "math/DynamicRectMatrix.hpp"
   * @brief rectangular matrix class
   */
  template<typename Real> 
  class DynamicRectMatrix {
  public:
    typedef Real ElemType;
    typedef Real* ElemPoinerType;
    typedef DynamicRectMatrix<Real> SelfType;
            
    /** default constructor */
    DynamicRectMatrix(){}
    DynamicRectMatrix(int nrow, int ncol) {
      allocate(nrow, ncol);

      for (unsigned int i = 0; i < nrow_; i++)
        for (unsigned int j = 0; j < ncol_; j++)
          this->data_[i][j] = 0.0;
    }

    /** Constructs and initializes every element of this matrix to a scalar */ 
    DynamicRectMatrix(int nrow, int ncol, Real s) {
      allocate(nrow, ncol);
      for (unsigned int i = 0; i < nrow_; i++)
        for (unsigned int j = 0; j < ncol_; j++)
          this->data_[i][j] = s;
    }

    DynamicRectMatrix(int nrow, int ncol, Real* array) {
      allocate(nrow, ncol);
      for (unsigned int i = 0; i < nrow_; i++)
        for (unsigned int j = 0; j < ncol_; j++)
          this->data_[i][j] = array[i * nrow_ + j];
    }

    /** copy constructor */
    DynamicRectMatrix(const SelfType& m) {
      allocate(m.getNRow(), m.getNCol());
      
      for (unsigned int i = 0; i < nrow_; i++)
        for (unsigned int j = 0; j < ncol_; j++)
            this->data_[i][j] = m.data_[i][j];
    }
            
    /** destructor*/
    ~DynamicRectMatrix() { deallocate();}

    /** copy assignment operator */
    DynamicRectMatrix<Real> operator =(const DynamicRectMatrix<Real> m) {
      if (this == &m)
          return *this;
      if (nrow_ != m.getNRow() || ncol_ != m.getNCol()) {
        deallocate();
        allocate(m.getNRow(), m.getNCol());
      }
      
      for (unsigned int i = 0; i < nrow_; i++)
        for (unsigned int j = 0; j < ncol_; j++)
          this->data_[i][j] = m.data_[i][j];
      return *this;
    }
            
    /**
     * Return the reference of a single element of this matrix.
     * @return the reference of a single element of this matrix 
     * @param i row index
     * @param j Column index
     */
    Real& operator()(unsigned int i, unsigned int j) {
      return this->data_[i][j];
    }

    /**
     * Return the value of a single element of this matrix.
     * @return the value of a single element of this matrix 
     * @param i row index
     * @param j Column index
     */        
    Real operator()(unsigned int i, unsigned int j) const  {
                
      return this->data_[i][j];  
    }

    /** 
     * Copy the internal data to an array
     * @param array the pointer of destination array
     */
    void getArray(Real* array) {
      for (unsigned int i = 0; i < nrow_; i++) {
        for (unsigned int j = 0; j < ncol_; j++) {
          array[i * nrow_ + j] = this->data_[i][j];
        }
      }
    }

    /**
     * Returns a row of  this matrix as a vector.
     * @return a row of  this matrix as a vector 
     * @param row the row index
     */                
    DynamicVector<Real> getRow(unsigned int row) {
      DynamicVector<Real> v;

      for (unsigned int i = 0; i < ncol_; i++)
        v[i] = this->data_[row][i];

      return v;
    }

    /**
     * Sets a row of  this matrix
     * @param row the row index
     * @param v the vector to be set
     */                
    void setRow(unsigned int row, const DynamicVector<Real>& v) {
      assert(v.size() == nrow_);
      for (unsigned int i = 0; i < ncol_; i++)
        this->data_[row][i] = v[i];
    }

    /**
     * Returns a column of  this matrix as a vector.
     * @return a column of  this matrix as a vector 
     * @param col the column index
     */                
    DynamicVector<Real> getColumn(unsigned int col) {
      DynamicVector<Real> v(ncol_);

      for (unsigned int j = 0; j < nrow_; j++)
        v[j] = this->data_[j][col];

      return v;
    }

    /**
     * Sets a column of  this matrix
     * @param col the column index
     * @param v the vector to be set
     */                
    void setColumn(unsigned int col, const DynamicVector<Real>& v){

      for (unsigned int j = 0; j < nrow_; j++)
        this->data_[j][col] = v[j];
    }         

    /**
     * swap two rows of this matrix
     * @param i the first row
     * @param j the second row
     */
    void swapRow(unsigned int i, unsigned int j){
      assert(i < nrow_ && j < nrow_);

      for (unsigned int k = 0; k < ncol_; k++)
        std::swap(this->data_[i][k], this->data_[j][k]);
    }

    /**
     * swap two Columns of this matrix
     * @param i the first Column
     * @param j the second Column
     */
    void swapColumn(unsigned int i, unsigned int j){
      assert(i < ncol_ && j < ncol_);
                    
      for (unsigned int k = 0; k < nrow_; k++)
        std::swap(this->data_[k][i], this->data_[k][j]);
    }

    /**
     * Tests if this matrix is identical to matrix m
     * @return true if this matrix is equal to the matrix m, return false otherwise
     * @m matrix to be compared
     *
     * @todo replace operator == by template function equal
     */
    bool operator ==(const DynamicRectMatrix<Real> m) {
      assert(nrow_ == m.getNRow() && ncol_ == m.getNCol());
      for (unsigned int i = 0; i < nrow_; i++)
        for (unsigned int j = 0; j < ncol_; j++)
          if (!equal(this->data_[i][j], m.data_[i][j]))
            return false;

      return true;
    }

    /**
     * Tests if this matrix is not equal to matrix m
     * @return true if this matrix is not equal to the matrix m, return false otherwise
     * @m matrix to be compared
     */
    bool operator !=(const DynamicRectMatrix<Real> m) {
      return !(*this == m);
    }

    /** Negates the value of this matrix in place. */           
    inline void negate() {
      for (unsigned int i = 0; i < nrow_; i++)
        for (unsigned int j = 0; j < ncol_; j++)
          this->data_[i][j] = -this->data_[i][j];
    }
            
    /**
     * Sets the value of this matrix to the negation of matrix m.
     * @param m the source matrix
     */
    inline void negate(const DynamicRectMatrix<Real> m) {
      for (unsigned int i = 0; i < nrow_; i++)
        for (unsigned int j = 0; j < ncol_; j++)
          this->data_[i][j] = -m.data_[i][j];        
    }
            
    /**
     * Sets the value of this matrix to the sum of itself and m (*this += m).
     * @param m the other matrix
     */
    inline void add( const DynamicRectMatrix<Real> m ) {
      assert(nrow_ == m.getNRow() && ncol_ == m.getNCol());
      for (unsigned int i = 0; i < nrow_; i++)
        for (unsigned int j = 0; j < ncol_; j++)        
          this->data_[i][j] += m.data_[i][j];
    }
            
    /**
     * Sets the value of this matrix to the sum of m1 and m2 (*this = m1 + m2).
     * @param m1 the first matrix
     * @param m2 the second matrix
     */
    inline void add( const DynamicRectMatrix<Real> m1, const DynamicRectMatrix<Real> m2 ) {
      assert(m1.getNRow() == m2.getNRow() && m1.getNCol() == m2.getNCol());
      for (unsigned int i = 0; i < nrow_; i++)
        for (unsigned int j = 0; j < ncol_; j++)        
          this->data_[i][j] = m1.data_[i][j] + m2.data_[i][j];
    }
            
    /**
     * Sets the value of this matrix to the difference  of itself and m (*this -= m).
     * @param m the other matrix
     */
    inline void sub( const DynamicRectMatrix<Real> m ) {
      assert(nrow_ == m.getNRow() && ncol_ == m.getNCol());
      for (unsigned int i = 0; i < nrow_; i++)
        for (unsigned int j = 0; j < ncol_; j++)        
          this->data_[i][j] -= m.data_[i][j];
    }
            
    /**
     * Sets the value of this matrix to the difference of matrix m1 and m2 (*this = m1 - m2).
     * @param m1 the first matrix
     * @param m2 the second matrix
     */
    inline void sub( const DynamicRectMatrix<Real> m1, const DynamicRectMatrix<Real> m2){
      assert(m1.getNRow() == m2.getNRow() && m1.getNCol() == m2.getNCol());
      for (unsigned int i = 0; i < nrow_; i++)
        for (unsigned int j = 0; j < ncol_; j++)        
          this->data_[i][j] = m1.data_[i][j] - m2.data_[i][j];
    }

    /**
     * Sets the value of this matrix to the scalar multiplication of itself (*this *= s).
     * @param s the scalar value
     */
    inline void mul( Real s ) {
      for (unsigned int i = 0; i < nrow_; i++)
        for (unsigned int j = 0; j < ncol_; j++)  
          this->data_[i][j] *= s;
    }

    /**
     * Sets the value of this matrix to the scalar multiplication of matrix m  (*this = s * m).
     * @param s the scalar value
     * @param m the matrix
     */
    inline void mul( Real s, const DynamicRectMatrix<Real> m ) {
      assert(nrow_ == m.getNRow() && ncol_ == m.getNCol());    
      for (unsigned int i = 0; i < nrow_; i++)
        for (unsigned int j = 0; j < ncol_; j++)  
          this->data_[i][j] = s * m.data_[i][j];
    }

    /**
     * Sets the value of this matrix to the scalar division of itself  (*this /= s ).
     * @param s the scalar value
     */             
    inline void div( Real s) {
      for (unsigned int i = 0; i < nrow_; i++)
        for (unsigned int j = 0; j < ncol_; j++)  
          this->data_[i][j] /= s;
    }

    /**
     * Sets the value of this matrix to the scalar division of matrix m  (*this = m /s).
     * @param s the scalar value
     * @param m the matrix
     */
    inline void div( Real s, const DynamicRectMatrix<Real> m ) {
      assert(nrow_ == m.getNRow() && ncol_ == m.getNCol());
      for (unsigned int i = 0; i < nrow_; i++)
        for (unsigned int j = 0; j < ncol_; j++)  
          this->data_[i][j] = m.data_[i][j] / s;
    }

    /**
     *  Multiples a scalar into every element of this matrix.
     * @param s the scalar value
     */
    DynamicRectMatrix<Real> operator *=(const Real s) {
      this->mul(s);
      return *this;
    }

    /**
     *  Divides every element of this matrix by a scalar.
     * @param s the scalar value
     */
    DynamicRectMatrix<Real> operator /=(const Real s) {
      this->div(s);
      return *this;
    }

    /**
     * Sets the value of this matrix to the sum of the other matrix and itself (*this += m).
     * @param m the other matrix
     */
    DynamicRectMatrix<Real> operator += (const DynamicRectMatrix<Real> m) {
      assert(nrow_ == m.getNRow() && ncol_ == m.getNCol());
      add(m);
      return *this;
    }

    /**
     * Sets the value of this matrix to the differerence of itself and the other matrix (*this -= m) 
     * @param m the other matrix
     */
    DynamicRectMatrix<Real> operator -= (const DynamicRectMatrix<Real> m){
      assert(nrow_ == m.getNRow() && ncol_ == m.getNCol());
      sub(m);
      return *this;
    }

    /** Return the transpose of this matrix */
    DynamicRectMatrix<Real> transpose() const{
      DynamicRectMatrix<Real> result(ncol_,nrow_);
                
      for (unsigned int i = 0; i < nrow_; i++)
        for (unsigned int j = 0; j < ncol_; j++)              
          result(j, i) = this->data_[i][j];

      return result;
    }

    unsigned int getNRow() const {return nrow_;}
    unsigned int getNCol() const {return ncol_;}

    template<class MatrixType>
    void setSubMatrix(unsigned int beginRow, unsigned int beginCol, const MatrixType& m) {
        assert(beginRow + m.getNRow() -1 <= nrow_);
        assert(beginCol + m.getNCol() -1 <= ncol_);

        for (unsigned int i = 0; i < m.getNRow(); ++i)
            for (unsigned int j = 0; j < m.getNCol(); ++j)
                this->data_[beginRow+i][beginCol+j] = m(i, j);
    }

    template<class MatrixType>
    void getSubMatrix(unsigned int beginRow, unsigned int beginCol, MatrixType& m) {
        assert(beginRow + m.getNRow() -1 <= nrow_);
        assert(beginCol + m.getNCol() - 1 <= ncol_);

        for (unsigned int i = 0; i < m.getNRow(); ++i)
            for (unsigned int j = 0; j < m.getNCol(); ++j)
                m(i, j) = this->data_[beginRow+i][beginCol+j];
    }
    
  protected:
    Real** data_;
    unsigned int nrow_;
    unsigned int ncol_;
  private:
    void allocate(int nrow, int ncol) {
        nrow_ = nrow;
        ncol_ = ncol;
        data_ = new Real*[nrow_];
        for (int i = 0; i < nrow_; ++i)
            data_[i] = new Real[ncol_];
    }

    void deallocate() {
      for (int i = 0; i < nrow_; ++i)
        delete data_[i];
      delete []data_;

      nrow_ = 0;
      ncol_ = 0;
      data_ = NULL;
    }
    
  };

  /** Negate the value of every element of this matrix. */
  template<typename Real> 
  inline DynamicRectMatrix<Real> operator -(const DynamicRectMatrix<Real> m) {
    DynamicRectMatrix<Real> result(m);

    result.negate();

    return result;
  }
    
  /**
   * Return the sum of two matrixes  (m1 + m2). 
   * @return the sum of two matrixes
   * @param m1 the first matrix
   * @param m2 the second matrix
   */ 
  template<typename Real> 
  inline DynamicRectMatrix<Real> operator + (const DynamicRectMatrix<Real> m1,const DynamicRectMatrix<Real> m2) {
    
    DynamicRectMatrix<Real> result(m1.getNRow(), m1.getNCol());

    result.add(m1, m2);

    return result;
  }
    
  /**
   * Return the difference of two matrixes  (m1 - m2). 
   * @return the sum of two matrixes
   * @param m1 the first matrix
   * @param m2 the second matrix
   */
  template<typename Real> 
  inline DynamicRectMatrix<Real> operator - (const DynamicRectMatrix<Real> m1, const DynamicRectMatrix<Real> m2) {
    DynamicRectMatrix<Real> result(m1.getNRow(), m1.getNCol());

    result.sub(m1, m2);

    return result;
  }

  /**
   * Return the multiplication of scalra and  matrix  (m * s). 
   * @return the multiplication of a scalra and  a matrix 
   * @param m the matrix
   * @param s the scalar
   */
  template<typename Real> 
  inline DynamicRectMatrix<Real> operator *(const DynamicRectMatrix<Real> m, Real s) {
    DynamicRectMatrix<Real> result(m.getNRow(), m.getNCol());

    result.mul(s, m);

    return result;
  }

  /**
   * Return the multiplication of a scalra and  a matrix  (s * m). 
   * @return the multiplication of a scalra and  a matrix 
   * @param s the scalar
   * @param m the matrix
   */
  template<typename Real> 
  inline DynamicRectMatrix<Real> operator *(Real s, const DynamicRectMatrix<Real> m) {
    DynamicRectMatrix<Real> result(m.getNRow(), m.getNCol());

    result.mul(s, m);

    return result;
  }
    
  /**
   * Return the multiplication of two matrixes  (m1 * m2). 
   * @return the multiplication of two matrixes
   * @param m1 the first matrix
   * @param m2 the second matrix
   */
  template<typename Real> 
  inline DynamicRectMatrix<Real> operator *(const DynamicRectMatrix<Real>& m1, const DynamicRectMatrix<Real>& m2) {
    assert(m1.getNCol() == m2.getNRow());
    unsigned int sameDim = m1.getNCol();
    int nrow = m1.getNRow();
    int ncol = m2.getNCol();
    DynamicRectMatrix<Real> result(nrow, ncol );
    for (unsigned int i = 0; i < nrow; i++)
      for (unsigned int j = 0; j < ncol; j++)
        for (unsigned int k = 0; k < sameDim; k++)
          result(i, j)  += m1(i, k) * m2(k, j);                

    return result;
  }
    
  /**
   * Return the multiplication of  a matrix and a vector  (m * v). 
   * @return the multiplication of a matrix and a vector
   * @param m the matrix
   * @param v the vector
   */
  template<typename Real>
  inline DynamicVector<Real> operator *(const DynamicRectMatrix<Real> m, const DynamicVector<Real>& v) {
    int nrow = m.getNRow();
    int ncol = m.getNCol();
    assert(ncol = v.size());
    DynamicVector<Real> result(nrow);
    
    for (unsigned int i = 0; i < nrow ; i++)
      for (unsigned int j = 0; j < ncol ; j++)            
        result[i] += m(i, j) * v[j];
            
    return result;                                                                 
  }

  /**
   * Return the scalar division of matrix   (m / s). 
   * @return the scalar division of matrix  
   * @param m the matrix
   * @param s the scalar
   */
  template<typename Real> 
  inline DynamicRectMatrix<Real> operator /(const DynamicRectMatrix<Real> m, Real s) {
    DynamicRectMatrix<Real> result(m.getNRow(), m.getNCol());

    result.div(s, m);

    return result;
  }    

  /**
   * Write to an output stream
   */
  template<typename Real>
  std::ostream &operator<< ( std::ostream& o, const DynamicRectMatrix<Real> m) {
    for (unsigned int i = 0; i < m.getNRow() ; i++) {
      o << "(";
      for (unsigned int j = 0; j < m.getNCol() ; j++) {
        o << m(i, j);
        if (j != m.getNCol() -1)
          o << "\t";
      }
      o << ")" << std::endl;
    }
    return o;        
  }    
}
#endif //MATH_RECTMATRIX_HPP

Revision:	1349
Committed:	Wed May 20 19:35:05 2009 UTC (15 years, 11 months ago) by cli2
Original Path:	*trunk/src/math/DynamicRectMatrix.hpp*
File size:	18755 byte(s)
Log Message:	Modifications to support RMSD calculations, removing classes we aren't using.
#	User	Rev	Content
1	tim	891	/*
2			* Copyright (c) 2005 The University of Notre Dame. All Rights Reserved.
3			*
4			* The University of Notre Dame grants you ("Licensee") a
5			* non-exclusive, royalty free, license to use, modify and
6			* redistribute this software in source and binary code form, provided
7			* that the following conditions are met:
8			*
9			* 1. Acknowledgement of the program authors must be made in any
10			* publication of scientific results based in part on use of the
11			* program. An acceptable form of acknowledgement is citation of
12			* the article in which the program was described (Matthew
13			* A. Meineke, Charles F. Vardeman II, Teng Lin, Christopher
14			* J. Fennell and J. Daniel Gezelter, "OOPSE: An Object-Oriented
15			* Parallel Simulation Engine for Molecular Dynamics,"
16			* J. Comput. Chem. 26, pp. 252-271 (2005))
17			*
18			* 2. Redistributions of source code must retain the above copyright
19			* notice, this list of conditions and the following disclaimer.
20			*
21			* 3. Redistributions in binary form must reproduce the above copyright
22			* notice, this list of conditions and the following disclaimer in the
23			* documentation and/or other materials provided with the
24			* distribution.
25			*
26			* This software is provided "AS IS," without a warranty of any
27			* kind. All express or implied conditions, representations and
28			* warranties, including any implied warranty of merchantability,
29			* fitness for a particular purpose or non-infringement, are hereby
30			* excluded. The University of Notre Dame and its licensors shall not
31			* be liable for any damages suffered by licensee as a result of
32			* using, modifying or distributing the software or its
33			* derivatives. In no event will the University of Notre Dame or its
34			* licensors be liable for any lost revenue, profit or data, or for
35			* direct, indirect, special, consequential, incidental or punitive
36			* damages, however caused and regardless of the theory of liability,
37			* arising out of the use of or inability to use software, even if the
38			* University of Notre Dame has been advised of the possibility of
39			* such damages.
40			*/
41
42			/**
43			* @file DynamicRectMatrix.hpp
44			* @author Teng Lin
45			* @date 10/11/2004
46			* @version 1.0
47			*/
48
49			#ifndef MATH_DYNAMICRECTMATRIX_HPP
50			#define MATH_DYNAMICRECTMATRIX_HPP
51			#include <math.h>
52			#include <cmath>
53			#include "math/DynamicVector.hpp"
54
55			namespace oopse {
56
57			/**
58			* @class DynamicRectMatrix DynamicRectMatrix.hpp "math/DynamicRectMatrix.hpp"
59			* @brief rectangular matrix class
60			*/
61			template<typename Real>
62			class DynamicRectMatrix {
63			public:
64			typedef Real ElemType;
65			typedef Real* ElemPoinerType;
66			typedef DynamicRectMatrix<Real> SelfType;
67
68			/** default constructor */
69	cli2	1349	DynamicRectMatrix(){}
70	tim	891	DynamicRectMatrix(int nrow, int ncol) {
71			allocate(nrow, ncol);
72
73			for (unsigned int i = 0; i < nrow_; i++)
74			for (unsigned int j = 0; j < ncol_; j++)
75			this->data_[i][j] = 0.0;
76			}
77
78			/** Constructs and initializes every element of this matrix to a scalar */
79			DynamicRectMatrix(int nrow, int ncol, Real s) {
80			allocate(nrow, ncol);
81			for (unsigned int i = 0; i < nrow_; i++)
82			for (unsigned int j = 0; j < ncol_; j++)
83			this->data_[i][j] = s;
84			}
85
86			DynamicRectMatrix(int nrow, int ncol, Real* array) {
87			allocate(nrow, ncol);
88			for (unsigned int i = 0; i < nrow_; i++)
89			for (unsigned int j = 0; j < ncol_; j++)
90			this->data_[i][j] = array[i * nrow_ + j];
91			}
92
93			/** copy constructor */
94			DynamicRectMatrix(const SelfType& m) {
95			allocate(m.getNRow(), m.getNCol());
96
97			for (unsigned int i = 0; i < nrow_; i++)
98			for (unsigned int j = 0; j < ncol_; j++)
99			this->data_[i][j] = m.data_[i][j];
100			}
101
102			/** destructor*/
103			~DynamicRectMatrix() { deallocate();}
104
105			/** copy assignment operator */
106			DynamicRectMatrix<Real> operator =(const DynamicRectMatrix<Real> m) {
107			if (this == &m)
108			return *this;
109			if (nrow_ != m.getNRow() \|\| ncol_ != m.getNCol()) {
110			deallocate();
111			allocate(m.getNRow(), m.getNCol());
112			}
113
114			for (unsigned int i = 0; i < nrow_; i++)
115			for (unsigned int j = 0; j < ncol_; j++)
116			this->data_[i][j] = m.data_[i][j];
117			return *this;
118			}
119
120			/**
121			* Return the reference of a single element of this matrix.
122			* @return the reference of a single element of this matrix
123			* @param i row index
124			* @param j Column index
125			*/
126			Real& operator()(unsigned int i, unsigned int j) {
127			return this->data_[i][j];
128			}
129
130			/**
131			* Return the value of a single element of this matrix.
132			* @return the value of a single element of this matrix
133			* @param i row index
134			* @param j Column index
135			*/
136			Real operator()(unsigned int i, unsigned int j) const {
137
138			return this->data_[i][j];
139			}
140
141			/**
142			* Copy the internal data to an array
143			* @param array the pointer of destination array
144			*/
145			void getArray(Real* array) {
146			for (unsigned int i = 0; i < nrow_; i++) {
147			for (unsigned int j = 0; j < ncol_; j++) {
148			array[i * nrow_ + j] = this->data_[i][j];
149			}
150			}
151			}
152
153			/**
154			* Returns a row of this matrix as a vector.
155			* @return a row of this matrix as a vector
156			* @param row the row index
157			*/
158			DynamicVector<Real> getRow(unsigned int row) {
159			DynamicVector<Real> v;
160
161			for (unsigned int i = 0; i < ncol_; i++)
162			v[i] = this->data_[row][i];
163
164			return v;
165			}
166
167			/**
168			* Sets a row of this matrix
169			* @param row the row index
170			* @param v the vector to be set
171			*/
172			void setRow(unsigned int row, const DynamicVector<Real>& v) {
173			assert(v.size() == nrow_);
174			for (unsigned int i = 0; i < ncol_; i++)
175			this->data_[row][i] = v[i];
176			}
177
178			/**
179			* Returns a column of this matrix as a vector.
180			* @return a column of this matrix as a vector
181			* @param col the column index
182			*/
183			DynamicVector<Real> getColumn(unsigned int col) {
184			DynamicVector<Real> v(ncol_);
185
186			for (unsigned int j = 0; j < nrow_; j++)
187			v[j] = this->data_[j][col];
188
189			return v;
190			}
191
192			/**
193			* Sets a column of this matrix
194			* @param col the column index
195			* @param v the vector to be set
196			*/
197			void setColumn(unsigned int col, const DynamicVector<Real>& v){
198
199			for (unsigned int j = 0; j < nrow_; j++)
200			this->data_[j][col] = v[j];
201			}
202
203			/**
204			* swap two rows of this matrix
205			* @param i the first row
206			* @param j the second row
207			*/
208			void swapRow(unsigned int i, unsigned int j){
209			assert(i < nrow_ && j < nrow_);
210
211			for (unsigned int k = 0; k < ncol_; k++)
212			std::swap(this->data_[i][k], this->data_[j][k]);
213			}
214
215			/**
216			* swap two Columns of this matrix
217			* @param i the first Column
218			* @param j the second Column
219			*/
220			void swapColumn(unsigned int i, unsigned int j){
221			assert(i < ncol_ && j < ncol_);
222
223			for (unsigned int k = 0; k < nrow_; k++)
224			std::swap(this->data_[k][i], this->data_[k][j]);
225			}
226
227			/**
228			* Tests if this matrix is identical to matrix m
229			* @return true if this matrix is equal to the matrix m, return false otherwise
230			* @m matrix to be compared
231			*
232			* @todo replace operator == by template function equal
233			*/
234			bool operator ==(const DynamicRectMatrix<Real> m) {
235			assert(nrow_ == m.getNRow() && ncol_ == m.getNCol());
236			for (unsigned int i = 0; i < nrow_; i++)
237			for (unsigned int j = 0; j < ncol_; j++)
238			if (!equal(this->data_[i][j], m.data_[i][j]))
239			return false;
240
241			return true;
242			}
243
244			/**
245			* Tests if this matrix is not equal to matrix m
246			* @return true if this matrix is not equal to the matrix m, return false otherwise
247			* @m matrix to be compared
248			*/
249			bool operator !=(const DynamicRectMatrix<Real> m) {
250			return !(*this == m);
251			}
252
253			/** Negates the value of this matrix in place. */
254			inline void negate() {
255			for (unsigned int i = 0; i < nrow_; i++)
256			for (unsigned int j = 0; j < ncol_; j++)
257			this->data_[i][j] = -this->data_[i][j];
258			}
259
260			/**
261			* Sets the value of this matrix to the negation of matrix m.
262			* @param m the source matrix
263			*/
264			inline void negate(const DynamicRectMatrix<Real> m) {
265			for (unsigned int i = 0; i < nrow_; i++)
266			for (unsigned int j = 0; j < ncol_; j++)
267			this->data_[i][j] = -m.data_[i][j];
268			}
269
270			/**
271			* Sets the value of this matrix to the sum of itself and m (*this += m).
272			* @param m the other matrix
273			*/
274			inline void add( const DynamicRectMatrix<Real> m ) {
275			assert(nrow_ == m.getNRow() && ncol_ == m.getNCol());
276			for (unsigned int i = 0; i < nrow_; i++)
277			for (unsigned int j = 0; j < ncol_; j++)
278			this->data_[i][j] += m.data_[i][j];
279			}
280
281			/**
282			* Sets the value of this matrix to the sum of m1 and m2 (*this = m1 + m2).
283			* @param m1 the first matrix
284			* @param m2 the second matrix
285			*/
286			inline void add( const DynamicRectMatrix<Real> m1, const DynamicRectMatrix<Real> m2 ) {
287			assert(m1.getNRow() == m2.getNRow() && m1.getNCol() == m2.getNCol());
288			for (unsigned int i = 0; i < nrow_; i++)
289			for (unsigned int j = 0; j < ncol_; j++)
290			this->data_[i][j] = m1.data_[i][j] + m2.data_[i][j];
291			}
292
293			/**
294			* Sets the value of this matrix to the difference of itself and m (*this -= m).
295			* @param m the other matrix
296			*/
297			inline void sub( const DynamicRectMatrix<Real> m ) {
298			assert(nrow_ == m.getNRow() && ncol_ == m.getNCol());
299			for (unsigned int i = 0; i < nrow_; i++)
300			for (unsigned int j = 0; j < ncol_; j++)
301			this->data_[i][j] -= m.data_[i][j];
302			}
303
304			/**
305			* Sets the value of this matrix to the difference of matrix m1 and m2 (*this = m1 - m2).
306			* @param m1 the first matrix
307			* @param m2 the second matrix
308			*/
309			inline void sub( const DynamicRectMatrix<Real> m1, const DynamicRectMatrix<Real> m2){
310			assert(m1.getNRow() == m2.getNRow() && m1.getNCol() == m2.getNCol());
311			for (unsigned int i = 0; i < nrow_; i++)
312			for (unsigned int j = 0; j < ncol_; j++)
313			this->data_[i][j] = m1.data_[i][j] - m2.data_[i][j];
314			}
315
316			/**
317			* Sets the value of this matrix to the scalar multiplication of itself (this = s).
318			* @param s the scalar value
319			*/
320			inline void mul( Real s ) {
321			for (unsigned int i = 0; i < nrow_; i++)
322			for (unsigned int j = 0; j < ncol_; j++)
323			this->data_[i][j] *= s;
324			}
325
326			/**
327			* Sets the value of this matrix to the scalar multiplication of matrix m (this = s m).
328			* @param s the scalar value
329			* @param m the matrix
330			*/
331			inline void mul( Real s, const DynamicRectMatrix<Real> m ) {
332			assert(nrow_ == m.getNRow() && ncol_ == m.getNCol());
333			for (unsigned int i = 0; i < nrow_; i++)
334			for (unsigned int j = 0; j < ncol_; j++)
335			this->data_[i][j] = s * m.data_[i][j];
336			}
337
338			/**
339			* Sets the value of this matrix to the scalar division of itself (*this /= s ).
340			* @param s the scalar value
341			*/
342			inline void div( Real s) {
343			for (unsigned int i = 0; i < nrow_; i++)
344			for (unsigned int j = 0; j < ncol_; j++)
345			this->data_[i][j] /= s;
346			}
347
348			/**
349			* Sets the value of this matrix to the scalar division of matrix m (*this = m /s).
350			* @param s the scalar value
351			* @param m the matrix
352			*/
353			inline void div( Real s, const DynamicRectMatrix<Real> m ) {
354			assert(nrow_ == m.getNRow() && ncol_ == m.getNCol());
355			for (unsigned int i = 0; i < nrow_; i++)
356			for (unsigned int j = 0; j < ncol_; j++)
357			this->data_[i][j] = m.data_[i][j] / s;
358			}
359
360			/**
361			* Multiples a scalar into every element of this matrix.
362			* @param s the scalar value
363			*/
364			DynamicRectMatrix<Real> operator *=(const Real s) {
365			this->mul(s);
366			return *this;
367			}
368
369			/**
370			* Divides every element of this matrix by a scalar.
371			* @param s the scalar value
372			*/
373			DynamicRectMatrix<Real> operator /=(const Real s) {
374			this->div(s);
375			return *this;
376			}
377
378			/**
379			* Sets the value of this matrix to the sum of the other matrix and itself (*this += m).
380			* @param m the other matrix
381			*/
382			DynamicRectMatrix<Real> operator += (const DynamicRectMatrix<Real> m) {
383			assert(nrow_ == m.getNRow() && ncol_ == m.getNCol());
384			add(m);
385			return *this;
386			}
387
388			/**
389			* Sets the value of this matrix to the differerence of itself and the other matrix (*this -= m)
390			* @param m the other matrix
391			*/
392			DynamicRectMatrix<Real> operator -= (const DynamicRectMatrix<Real> m){
393			assert(nrow_ == m.getNRow() && ncol_ == m.getNCol());
394			sub(m);
395			return *this;
396			}
397
398			/** Return the transpose of this matrix */
399			DynamicRectMatrix<Real> transpose() const{
400			DynamicRectMatrix<Real> result(ncol_,nrow_);
401
402			for (unsigned int i = 0; i < nrow_; i++)
403			for (unsigned int j = 0; j < ncol_; j++)
404			result(j, i) = this->data_[i][j];
405
406			return result;
407			}
408
409			unsigned int getNRow() const {return nrow_;}
410			unsigned int getNCol() const {return ncol_;}
411
412			template<class MatrixType>
413			void setSubMatrix(unsigned int beginRow, unsigned int beginCol, const MatrixType& m) {
414			assert(beginRow + m.getNRow() -1 <= nrow_);
415			assert(beginCol + m.getNCol() -1 <= ncol_);
416
417			for (unsigned int i = 0; i < m.getNRow(); ++i)
418			for (unsigned int j = 0; j < m.getNCol(); ++j)
419			this->data_[beginRow+i][beginCol+j] = m(i, j);
420			}
421
422			template<class MatrixType>
423			void getSubMatrix(unsigned int beginRow, unsigned int beginCol, MatrixType& m) {
424			assert(beginRow + m.getNRow() -1 <= nrow_);
425			assert(beginCol + m.getNCol() - 1 <= ncol_);
426
427			for (unsigned int i = 0; i < m.getNRow(); ++i)
428			for (unsigned int j = 0; j < m.getNCol(); ++j)
429			m(i, j) = this->data_[beginRow+i][beginCol+j];
430			}
431
432			protected:
433			Real** data_;
434			unsigned int nrow_;
435			unsigned int ncol_;
436			private:
437			void allocate(int nrow, int ncol) {
438			nrow_ = nrow;
439			ncol_ = ncol;
440			data_ = new Real*[nrow_];
441			for (int i = 0; i < nrow_; ++i)
442			data_[i] = new Real[ncol_];
443			}
444
445			void deallocate() {
446			for (int i = 0; i < nrow_; ++i)
447			delete data_[i];
448			delete []data_;
449
450			nrow_ = 0;
451			ncol_ = 0;
452			data_ = NULL;
453			}
454
455			};
456
457			/** Negate the value of every element of this matrix. */
458			template<typename Real>
459			inline DynamicRectMatrix<Real> operator -(const DynamicRectMatrix<Real> m) {
460			DynamicRectMatrix<Real> result(m);
461
462			result.negate();
463
464			return result;
465			}
466
467			/**
468			* Return the sum of two matrixes (m1 + m2).
469			* @return the sum of two matrixes
470			* @param m1 the first matrix
471			* @param m2 the second matrix
472			*/
473			template<typename Real>
474			inline DynamicRectMatrix<Real> operator + (const DynamicRectMatrix<Real> m1,const DynamicRectMatrix<Real> m2) {
475
476			DynamicRectMatrix<Real> result(m1.getNRow(), m1.getNCol());
477
478			result.add(m1, m2);
479
480			return result;
481			}
482
483			/**
484			* Return the difference of two matrixes (m1 - m2).
485			* @return the sum of two matrixes
486			* @param m1 the first matrix
487			* @param m2 the second matrix
488			*/
489			template<typename Real>
490			inline DynamicRectMatrix<Real> operator - (const DynamicRectMatrix<Real> m1, const DynamicRectMatrix<Real> m2) {
491			DynamicRectMatrix<Real> result(m1.getNRow(), m1.getNCol());
492
493			result.sub(m1, m2);
494
495			return result;
496			}
497
498			/**
499			* Return the multiplication of scalra and matrix (m * s).
500			* @return the multiplication of a scalra and a matrix
501			* @param m the matrix
502			* @param s the scalar
503			*/
504			template<typename Real>
505			inline DynamicRectMatrix<Real> operator *(const DynamicRectMatrix<Real> m, Real s) {
506			DynamicRectMatrix<Real> result(m.getNRow(), m.getNCol());
507
508			result.mul(s, m);
509
510			return result;
511			}
512
513			/**
514			* Return the multiplication of a scalra and a matrix (s * m).
515			* @return the multiplication of a scalra and a matrix
516			* @param s the scalar
517			* @param m the matrix
518			*/
519			template<typename Real>
520			inline DynamicRectMatrix<Real> operator *(Real s, const DynamicRectMatrix<Real> m) {
521			DynamicRectMatrix<Real> result(m.getNRow(), m.getNCol());
522
523			result.mul(s, m);
524
525			return result;
526			}
527
528			/**
529			* Return the multiplication of two matrixes (m1 * m2).
530			* @return the multiplication of two matrixes
531			* @param m1 the first matrix
532			* @param m2 the second matrix
533			*/
534			template<typename Real>
535			inline DynamicRectMatrix<Real> operator *(const DynamicRectMatrix<Real>& m1, const DynamicRectMatrix<Real>& m2) {
536			assert(m1.getNCol() == m2.getNRow());
537			unsigned int sameDim = m1.getNCol();
538			int nrow = m1.getNRow();
539			int ncol = m2.getNCol();
540			DynamicRectMatrix<Real> result(nrow, ncol );
541			for (unsigned int i = 0; i < nrow; i++)
542			for (unsigned int j = 0; j < ncol; j++)
543			for (unsigned int k = 0; k < sameDim; k++)
544			result(i, j) += m1(i, k) * m2(k, j);
545
546			return result;
547			}
548
549			/**
550			* Return the multiplication of a matrix and a vector (m * v).
551			* @return the multiplication of a matrix and a vector
552			* @param m the matrix
553			* @param v the vector
554			*/
555			template<typename Real>
556			inline DynamicVector<Real> operator *(const DynamicRectMatrix<Real> m, const DynamicVector<Real>& v) {
557			int nrow = m.getNRow();
558			int ncol = m.getNCol();
559			assert(ncol = v.size());
560			DynamicVector<Real> result(nrow);
561
562			for (unsigned int i = 0; i < nrow ; i++)
563			for (unsigned int j = 0; j < ncol ; j++)
564			result[i] += m(i, j) * v[j];
565
566			return result;
567			}
568
569			/**
570			* Return the scalar division of matrix (m / s).
571			* @return the scalar division of matrix
572			* @param m the matrix
573			* @param s the scalar
574			*/
575			template<typename Real>
576			inline DynamicRectMatrix<Real> operator /(const DynamicRectMatrix<Real> m, Real s) {
577			DynamicRectMatrix<Real> result(m.getNRow(), m.getNCol());
578
579			result.div(s, m);
580
581			return result;
582			}
583
584			/**
585			* Write to an output stream
586			*/
587			template<typename Real>
588			std::ostream &operator<< ( std::ostream& o, const DynamicRectMatrix<Real> m) {
589			for (unsigned int i = 0; i < m.getNRow() ; i++) {
590			o << "(";
591			for (unsigned int j = 0; j < m.getNCol() ; j++) {
592			o << m(i, j);
593			if (j != m.getNCol() -1)
594			o << "\t";
595			}
596			o << ")" << std::endl;
597			}
598			return o;
599			}
600			}
601			#endif //MATH_RECTMATRIX_HPP
602