src/math/RectMatrix.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 RectMatrix.hpp
 * @author Teng Lin
 * @date 10/11/2004
 * @version 1.0
 */

#ifndef MATH_RECTMATRIX_HPP
#define MATH_RECTMATRIX_HPP
#include <math.h>
#include <cmath>
#include "Vector.hpp"

namespace oopse {

  /**
   * @class RectMatrix RectMatrix.hpp "math/RectMatrix.hpp"
   * @brief rectangular matrix class
   */
  template<typename Real, unsigned int Row, unsigned int Col> 
  class RectMatrix {
  public:
    typedef Real ElemType;
    typedef Real* ElemPoinerType;
            
    /** default constructor */
    RectMatrix() {
      for (unsigned int i = 0; i < Row; i++)
        for (unsigned int j = 0; j < Col; j++)
          this->data_[i][j] = 0.0;
    }

    /** Constructs and initializes every element of this matrix to a scalar */ 
    RectMatrix(Real s) {
      for (unsigned int i = 0; i < Row; i++)
        for (unsigned int j = 0; j < Col; j++)
          this->data_[i][j] = s;
    }

    RectMatrix(Real* array) {
      for (unsigned int i = 0; i < Row; i++)
        for (unsigned int j = 0; j < Col; j++)
          this->data_[i][j] = array[i * Row + j];
    }

    /** copy constructor */
    RectMatrix(const RectMatrix<Real, Row, Col>& m) {
      *this = m;
    }
            
    /** destructor*/
    ~RectMatrix() {}

    /** copy assignment operator */
    RectMatrix<Real, Row, Col>& operator =(const RectMatrix<Real, Row, Col>& m) {
      if (this == &m)
        return *this;
                
      for (unsigned int i = 0; i < Row; i++)
        for (unsigned int j = 0; j < Col; 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) {
      //assert( i < Row && j < Col);
      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 < Row; i++) {
        for (unsigned int j = 0; j < Col; j++) {
          array[i * Row + j] = this->data_[i][j];
        }
      }
    }


    /** Returns the pointer of internal array */
    Real* getArrayPointer() {
      return &this->data_[0][0];
    }

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

      for (unsigned int i = 0; i < Col; 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 Vector<Real, Row>& v) {

      for (unsigned int i = 0; i < Col; 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
     */                
    Vector<Real, Col> getColumn(unsigned int col) {
      Vector<Real, Col> v;

      for (unsigned int j = 0; j < Row; 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 Vector<Real, Col>& v){

      for (unsigned int j = 0; j < Row; 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 < Row && j < Row);

      for (unsigned int k = 0; k < Col; 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 < Col && j < Col);
                    
      for (unsigned int k = 0; k < Row; 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 RectMatrix<Real, Row, Col>& m) {
      for (unsigned int i = 0; i < Row; i++)
        for (unsigned int j = 0; j < Col; 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 RectMatrix<Real, Row, Col>& m) {
      return !(*this == m);
    }

    /** Negates the value of this matrix in place. */           
    inline void negate() {
      for (unsigned int i = 0; i < Row; i++)
        for (unsigned int j = 0; j < Col; 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 RectMatrix<Real, Row, Col>& m) {
      for (unsigned int i = 0; i < Row; i++)
        for (unsigned int j = 0; j < Col; 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 RectMatrix<Real, Row, Col>& m ) {
      for (unsigned int i = 0; i < Row; i++)
        for (unsigned int j = 0; j < Col; 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 RectMatrix<Real, Row, Col>& m1, const RectMatrix<Real, Row, Col>& m2 ) {
      for (unsigned int i = 0; i < Row; i++)
        for (unsigned int j = 0; j < Col; 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 RectMatrix<Real, Row, Col>& m ) {
      for (unsigned int i = 0; i < Row; i++)
        for (unsigned int j = 0; j < Col; 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 RectMatrix<Real, Row, Col>& m1, const RectMatrix<Real, Row, Col>& m2){
      for (unsigned int i = 0; i < Row; i++)
        for (unsigned int j = 0; j < Col; 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 < Row; i++)
        for (unsigned int j = 0; j < Col; 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 RectMatrix<Real, Row, Col>& m ) {
      for (unsigned int i = 0; i < Row; i++)
        for (unsigned int j = 0; j < Col; 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 < Row; i++)
        for (unsigned int j = 0; j < Col; 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 RectMatrix<Real, Row, Col>& m ) {
      for (unsigned int i = 0; i < Row; i++)
        for (unsigned int j = 0; j < Col; j++)  
          this->data_[i][j] = m.data_[i][j] / s;
    }

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

    /**
     *  Divides every element of this matrix by a scalar.
     * @param s the scalar value
     */
    RectMatrix<Real, Row, Col>& 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
     */
    RectMatrix<Real, Row, Col>& operator += (const RectMatrix<Real, Row, Col>& m) {
      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
     */
    RectMatrix<Real, Row, Col>& operator -= (const RectMatrix<Real, Row, Col>& m){
      sub(m);
      return *this;
    }

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

      return result;
    }

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

        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 <= getNRow());
        assert(beginCol + m.getNCol() - 1 <= getNCol());

        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];
    }
    
    unsigned int getNRow() const {return Row;}
    unsigned int getNCol() const {return Col;}        

  protected:
    Real data_[Row][Col];
  };

  /** Negate the value of every element of this matrix. */
  template<typename Real, unsigned int Row, unsigned int Col> 
  inline RectMatrix<Real, Row, Col> operator -(const RectMatrix<Real, Row, Col>& m) {
    RectMatrix<Real, Row, Col> 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, unsigned int Row, unsigned int Col> 
  inline RectMatrix<Real, Row, Col> operator + (const RectMatrix<Real, Row, Col>& m1,const RectMatrix<Real, Row, Col>& m2) {
    RectMatrix<Real, Row, Col> result;

    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, unsigned int Row, unsigned int Col> 
  inline RectMatrix<Real, Row, Col> operator - (const RectMatrix<Real, Row, Col>& m1, const RectMatrix<Real, Row, Col>& m2) {
    RectMatrix<Real, Row, Col> result;

    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, unsigned int Row, unsigned int Col> 
  inline RectMatrix<Real, Row, Col> operator *(const RectMatrix<Real, Row, Col>& m, Real s) {
    RectMatrix<Real, Row, Col> result;

    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, unsigned int Row, unsigned int Col> 
  inline RectMatrix<Real, Row, Col> operator *(Real s, const RectMatrix<Real, Row, Col>& m) {
    RectMatrix<Real, Row, Col> result;

    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, unsigned int Row, unsigned int Col, unsigned int SameDim> 
  inline RectMatrix<Real, Row, Col> operator *(const RectMatrix<Real, Row, SameDim>& m1, const RectMatrix<Real, SameDim, Col>& m2) {
    RectMatrix<Real, Row, Col> result;

    for (unsigned int i = 0; i < Row; i++)
      for (unsigned int j = 0; j < Col; 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, unsigned int Row, unsigned int Col>
  inline Vector<Real, Row> operator *(const RectMatrix<Real, Row, Col>& m, const Vector<Real, Col>& v) {
    Vector<Real, Row> result;

    for (unsigned int i = 0; i < Row ; i++)
      for (unsigned int j = 0; j < Col ; 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, unsigned int Row, unsigned int Col> 
  inline RectMatrix<Real, Row, Col> operator /(const RectMatrix<Real, Row, Col>& m, Real s) {
    RectMatrix<Real, Row, Col> result;

    result.div(s, m);

    return result;
  }    

  /**
   * Write to an output stream
   */
  template<typename Real,  unsigned int Row, unsigned int Col>
  std::ostream &operator<< ( std::ostream& o, const RectMatrix<Real, Row, Col>& m) {
    for (unsigned int i = 0; i < Row ; i++) {
      o << "(";
      for (unsigned int j = 0; j < Col ; j++) {
        o << m(i, j);
        if (j != Col -1)
          o << "\t";
      }
      o << ")" << std::endl;
    }
    return o;        
  }    
}
#endif //MATH_RECTMATRIX_HPP
Revision:	2596
Committed:	Wed Feb 22 20:35:16 2006 UTC (18 years, 4 months ago) by tim
File size:	17311 byte(s)
Log Message:	Adding Hydrodynamics Module
#	User	Rev	Content
1	gezelter	2204	/*
2	gezelter	1930	* Copyright (c) 2005 The University of Notre Dame. All Rights Reserved.
3	tim	1564	*
4	gezelter	1930	* 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	tim	1564	*/
41	gezelter	1930
42	tim	1564	/**
43			* @file RectMatrix.hpp
44			* @author Teng Lin
45			* @date 10/11/2004
46			* @version 1.0
47			*/
48
49			#ifndef MATH_RECTMATRIX_HPP
50			#define MATH_RECTMATRIX_HPP
51	tim	1937	#include <math.h>
52	tim	1567	#include <cmath>
53	tim	1564	#include "Vector.hpp"
54
55			namespace oopse {
56
57	gezelter	2204	/**
58			* @class RectMatrix RectMatrix.hpp "math/RectMatrix.hpp"
59			* @brief rectangular matrix class
60			*/
61			template<typename Real, unsigned int Row, unsigned int Col>
62			class RectMatrix {
63			public:
64			typedef Real ElemType;
65			typedef Real* ElemPoinerType;
66	tim	1630
67	gezelter	2204	/** default constructor */
68			RectMatrix() {
69			for (unsigned int i = 0; i < Row; i++)
70			for (unsigned int j = 0; j < Col; j++)
71			this->data_[i][j] = 0.0;
72			}
73	tim	1564
74	gezelter	2204	/** Constructs and initializes every element of this matrix to a scalar */
75			RectMatrix(Real s) {
76			for (unsigned int i = 0; i < Row; i++)
77			for (unsigned int j = 0; j < Col; j++)
78			this->data_[i][j] = s;
79			}
80	tim	1564
81	gezelter	2204	RectMatrix(Real* array) {
82			for (unsigned int i = 0; i < Row; i++)
83			for (unsigned int j = 0; j < Col; j++)
84			this->data_[i][j] = array[i * Row + j];
85			}
86	tim	1644
87	gezelter	2204	/** copy constructor */
88			RectMatrix(const RectMatrix<Real, Row, Col>& m) {
89			*this = m;
90			}
91	tim	1630
92	gezelter	2204	/** destructor*/
93			~RectMatrix() {}
94	tim	1564
95	gezelter	2204	/** copy assignment operator */
96			RectMatrix<Real, Row, Col>& operator =(const RectMatrix<Real, Row, Col>& m) {
97			if (this == &m)
98			return *this;
99	tim	1630
100	gezelter	2204	for (unsigned int i = 0; i < Row; i++)
101			for (unsigned int j = 0; j < Col; j++)
102			this->data_[i][j] = m.data_[i][j];
103			return *this;
104			}
105	tim	1564
106	gezelter	2204	/**
107			* Return the reference of a single element of this matrix.
108			* @return the reference of a single element of this matrix
109			* @param i row index
110			* @param j Column index
111			*/
112			Real& operator()(unsigned int i, unsigned int j) {
113			//assert( i < Row && j < Col);
114			return this->data_[i][j];
115			}
116	tim	1564
117	gezelter	2204	/**
118			* Return the value of a single element of this matrix.
119			* @return the value of a single element of this matrix
120			* @param i row index
121			* @param j Column index
122			*/
123			Real operator()(unsigned int i, unsigned int j) const {
124	tim	1630
125	gezelter	2204	return this->data_[i][j];
126			}
127	tim	1564
128	gezelter	2204	/**
129			* Copy the internal data to an array
130			* @param array the pointer of destination array
131			*/
132			void getArray(Real* array) {
133			for (unsigned int i = 0; i < Row; i++) {
134			for (unsigned int j = 0; j < Col; j++) {
135			array[i * Row + j] = this->data_[i][j];
136			}
137			}
138			}
139	tim	1644
140
141	gezelter	2204	/** Returns the pointer of internal array */
142			Real* getArrayPointer() {
143			return &this->data_[0][0];
144			}
145	tim	1564
146	gezelter	2204	/**
147			* Returns a row of this matrix as a vector.
148			* @return a row of this matrix as a vector
149			* @param row the row index
150			*/
151			Vector<Real, Row> getRow(unsigned int row) {
152			Vector<Real, Row> v;
153	tim	1564
154	tim	2596	for (unsigned int i = 0; i < Col; i++)
155	gezelter	2204	v[i] = this->data_[row][i];
156	tim	1564
157	gezelter	2204	return v;
158			}
159	tim	1564
160	gezelter	2204	/**
161			* Sets a row of this matrix
162			* @param row the row index
163			* @param v the vector to be set
164			*/
165			void setRow(unsigned int row, const Vector<Real, Row>& v) {
166	tim	1564
167	tim	2596	for (unsigned int i = 0; i < Col; i++)
168	gezelter	2204	this->data_[row][i] = v[i];
169			}
170	tim	1564
171	gezelter	2204	/**
172			* Returns a column of this matrix as a vector.
173			* @return a column of this matrix as a vector
174			* @param col the column index
175			*/
176			Vector<Real, Col> getColumn(unsigned int col) {
177			Vector<Real, Col> v;
178	tim	1564
179	tim	2596	for (unsigned int j = 0; j < Row; j++)
180	gezelter	2204	v[j] = this->data_[j][col];
181	tim	1564
182	gezelter	2204	return v;
183			}
184	tim	1564
185	gezelter	2204	/**
186			* Sets a column of this matrix
187			* @param col the column index
188			* @param v the vector to be set
189			*/
190			void setColumn(unsigned int col, const Vector<Real, Col>& v){
191	tim	1564
192	tim	2596	for (unsigned int j = 0; j < Row; j++)
193	gezelter	2204	this->data_[j][col] = v[j];
194			}
195	tim	1594
196	gezelter	2204	/**
197			* swap two rows of this matrix
198			* @param i the first row
199			* @param j the second row
200			*/
201			void swapRow(unsigned int i, unsigned int j){
202			assert(i < Row && j < Row);
203	tim	1594
204	gezelter	2204	for (unsigned int k = 0; k < Col; k++)
205			std::swap(this->data_[i][k], this->data_[j][k]);
206			}
207	tim	1594
208	gezelter	2204	/**
209			* swap two Columns of this matrix
210			* @param i the first Column
211			* @param j the second Column
212			*/
213			void swapColumn(unsigned int i, unsigned int j){
214			assert(i < Col && j < Col);
215	tim	1630
216	gezelter	2204	for (unsigned int k = 0; k < Row; k++)
217			std::swap(this->data_[k][i], this->data_[k][j]);
218			}
219	tim	1564
220	gezelter	2204	/**
221			* Tests if this matrix is identical to matrix m
222			* @return true if this matrix is equal to the matrix m, return false otherwise
223			* @m matrix to be compared
224			*
225			* @todo replace operator == by template function equal
226			*/
227			bool operator ==(const RectMatrix<Real, Row, Col>& m) {
228			for (unsigned int i = 0; i < Row; i++)
229			for (unsigned int j = 0; j < Col; j++)
230			if (!equal(this->data_[i][j], m.data_[i][j]))
231			return false;
232	tim	1564
233	gezelter	2204	return true;
234			}
235	tim	1564
236	gezelter	2204	/**
237			* Tests if this matrix is not equal to matrix m
238			* @return true if this matrix is not equal to the matrix m, return false otherwise
239			* @m matrix to be compared
240			*/
241			bool operator !=(const RectMatrix<Real, Row, Col>& m) {
242			return !(*this == m);
243			}
244	tim	1564
245	gezelter	2204	/** Negates the value of this matrix in place. */
246			inline void negate() {
247			for (unsigned int i = 0; i < Row; i++)
248			for (unsigned int j = 0; j < Col; j++)
249			this->data_[i][j] = -this->data_[i][j];
250			}
251	tim	1630
252	gezelter	2204	/**
253			* Sets the value of this matrix to the negation of matrix m.
254			* @param m the source matrix
255			*/
256			inline void negate(const RectMatrix<Real, Row, Col>& m) {
257			for (unsigned int i = 0; i < Row; i++)
258			for (unsigned int j = 0; j < Col; j++)
259			this->data_[i][j] = -m.data_[i][j];
260			}
261	tim	1630
262	gezelter	2204	/**
263			* Sets the value of this matrix to the sum of itself and m (*this += m).
264			* @param m the other matrix
265			*/
266			inline void add( const RectMatrix<Real, Row, Col>& m ) {
267			for (unsigned int i = 0; i < Row; i++)
268			for (unsigned int j = 0; j < Col; j++)
269			this->data_[i][j] += m.data_[i][j];
270			}
271	tim	1630
272	gezelter	2204	/**
273			* Sets the value of this matrix to the sum of m1 and m2 (*this = m1 + m2).
274			* @param m1 the first matrix
275			* @param m2 the second matrix
276			*/
277			inline void add( const RectMatrix<Real, Row, Col>& m1, const RectMatrix<Real, Row, Col>& m2 ) {
278			for (unsigned int i = 0; i < Row; i++)
279			for (unsigned int j = 0; j < Col; j++)
280			this->data_[i][j] = m1.data_[i][j] + m2.data_[i][j];
281			}
282	tim	1630
283	gezelter	2204	/**
284			* Sets the value of this matrix to the difference of itself and m (*this -= m).
285			* @param m the other matrix
286			*/
287			inline void sub( const RectMatrix<Real, Row, Col>& m ) {
288			for (unsigned int i = 0; i < Row; i++)
289			for (unsigned int j = 0; j < Col; j++)
290			this->data_[i][j] -= m.data_[i][j];
291			}
292	tim	1630
293	gezelter	2204	/**
294			* Sets the value of this matrix to the difference of matrix m1 and m2 (*this = m1 - m2).
295			* @param m1 the first matrix
296			* @param m2 the second matrix
297			*/
298			inline void sub( const RectMatrix<Real, Row, Col>& m1, const RectMatrix<Real, Row, Col>& m2){
299			for (unsigned int i = 0; i < Row; i++)
300			for (unsigned int j = 0; j < Col; j++)
301			this->data_[i][j] = m1.data_[i][j] - m2.data_[i][j];
302			}
303	tim	1564
304	gezelter	2204	/**
305			* Sets the value of this matrix to the scalar multiplication of itself (this = s).
306			* @param s the scalar value
307			*/
308			inline void mul( Real s ) {
309			for (unsigned int i = 0; i < Row; i++)
310			for (unsigned int j = 0; j < Col; j++)
311			this->data_[i][j] *= s;
312			}
313	tim	1564
314	gezelter	2204	/**
315			* Sets the value of this matrix to the scalar multiplication of matrix m (this = s m).
316			* @param s the scalar value
317			* @param m the matrix
318			*/
319			inline void mul( Real s, const RectMatrix<Real, Row, Col>& m ) {
320			for (unsigned int i = 0; i < Row; i++)
321			for (unsigned int j = 0; j < Col; j++)
322			this->data_[i][j] = s * m.data_[i][j];
323			}
324	tim	1564
325	gezelter	2204	/**
326			* Sets the value of this matrix to the scalar division of itself (*this /= s ).
327			* @param s the scalar value
328			*/
329			inline void div( Real s) {
330			for (unsigned int i = 0; i < Row; i++)
331			for (unsigned int j = 0; j < Col; j++)
332			this->data_[i][j] /= s;
333			}
334	tim	1564
335	gezelter	2204	/**
336			* Sets the value of this matrix to the scalar division of matrix m (*this = m /s).
337			* @param s the scalar value
338			* @param m the matrix
339			*/
340			inline void div( Real s, const RectMatrix<Real, Row, Col>& m ) {
341			for (unsigned int i = 0; i < Row; i++)
342			for (unsigned int j = 0; j < Col; j++)
343			this->data_[i][j] = m.data_[i][j] / s;
344			}
345	tim	1564
346	gezelter	2204	/**
347			* Multiples a scalar into every element of this matrix.
348			* @param s the scalar value
349			*/
350			RectMatrix<Real, Row, Col>& operator *=(const Real s) {
351			this->mul(s);
352			return *this;
353			}
354	tim	1564
355	gezelter	2204	/**
356			* Divides every element of this matrix by a scalar.
357			* @param s the scalar value
358			*/
359			RectMatrix<Real, Row, Col>& operator /=(const Real s) {
360			this->div(s);
361			return *this;
362			}
363	tim	1564
364	gezelter	2204	/**
365			* Sets the value of this matrix to the sum of the other matrix and itself (*this += m).
366			* @param m the other matrix
367			*/
368			RectMatrix<Real, Row, Col>& operator += (const RectMatrix<Real, Row, Col>& m) {
369			add(m);
370			return *this;
371			}
372	tim	1564
373	gezelter	2204	/**
374			* Sets the value of this matrix to the differerence of itself and the other matrix (*this -= m)
375			* @param m the other matrix
376			*/
377			RectMatrix<Real, Row, Col>& operator -= (const RectMatrix<Real, Row, Col>& m){
378			sub(m);
379			return *this;
380			}
381	tim	1564
382	gezelter	2204	/** Return the transpose of this matrix */
383			RectMatrix<Real, Col, Row> transpose() const{
384			RectMatrix<Real, Col, Row> result;
385	tim	1630
386	gezelter	2204	for (unsigned int i = 0; i < Row; i++)
387			for (unsigned int j = 0; j < Col; j++)
388			result(j, i) = this->data_[i][j];
389	tim	1630
390	gezelter	2204	return result;
391			}
392	tim	2596
393			template<class MatrixType>
394			void setSubMatrix(unsigned int beginRow, unsigned int beginCol, const MatrixType& m) {
395			assert(beginRow + m.getNRow() -1 <= getNRow());
396			assert(beginCol + m.getNCol() -1 <= getNCol());
397
398			for (unsigned int i = 0; i < m.getNRow(); ++i)
399			for (unsigned int j = 0; j < m.getNCol(); ++j)
400			this->data_[beginRow+i][beginCol+j] = m(i, j);
401			}
402
403			template<class MatrixType>
404			void getSubMatrix(unsigned int beginRow, unsigned int beginCol, MatrixType& m) {
405			assert(beginRow + m.getNRow() -1 <= getNRow());
406			assert(beginCol + m.getNCol() - 1 <= getNCol());
407
408			for (unsigned int i = 0; i < m.getNRow(); ++i)
409			for (unsigned int j = 0; j < m.getNCol(); ++j)
410			m(i, j) = this->data_[beginRow+i][beginCol+j];
411			}
412
413			unsigned int getNRow() const {return Row;}
414			unsigned int getNCol() const {return Col;}
415
416	gezelter	2204	protected:
417			Real data_[Row][Col];
418			};
419	tim	1564
420	gezelter	2204	/** Negate the value of every element of this matrix. */
421			template<typename Real, unsigned int Row, unsigned int Col>
422			inline RectMatrix<Real, Row, Col> operator -(const RectMatrix<Real, Row, Col>& m) {
423			RectMatrix<Real, Row, Col> result(m);
424	tim	1564
425	gezelter	2204	result.negate();
426	tim	1564
427	gezelter	2204	return result;
428			}
429	tim	1564
430	gezelter	2204	/**
431			* Return the sum of two matrixes (m1 + m2).
432			* @return the sum of two matrixes
433			* @param m1 the first matrix
434			* @param m2 the second matrix
435			*/
436			template<typename Real, unsigned int Row, unsigned int Col>
437			inline RectMatrix<Real, Row, Col> operator + (const RectMatrix<Real, Row, Col>& m1,const RectMatrix<Real, Row, Col>& m2) {
438			RectMatrix<Real, Row, Col> result;
439	tim	1564
440	gezelter	2204	result.add(m1, m2);
441	tim	1564
442	gezelter	2204	return result;
443			}
444	tim	1564
445	gezelter	2204	/**
446			* Return the difference of two matrixes (m1 - m2).
447			* @return the sum of two matrixes
448			* @param m1 the first matrix
449			* @param m2 the second matrix
450			*/
451			template<typename Real, unsigned int Row, unsigned int Col>
452			inline RectMatrix<Real, Row, Col> operator - (const RectMatrix<Real, Row, Col>& m1, const RectMatrix<Real, Row, Col>& m2) {
453			RectMatrix<Real, Row, Col> result;
454	tim	1564
455	gezelter	2204	result.sub(m1, m2);
456	tim	1564
457	gezelter	2204	return result;
458			}
459	tim	1564
460	gezelter	2204	/**
461			* Return the multiplication of scalra and matrix (m * s).
462			* @return the multiplication of a scalra and a matrix
463			* @param m the matrix
464			* @param s the scalar
465			*/
466			template<typename Real, unsigned int Row, unsigned int Col>
467			inline RectMatrix<Real, Row, Col> operator *(const RectMatrix<Real, Row, Col>& m, Real s) {
468			RectMatrix<Real, Row, Col> result;
469	tim	1564
470	gezelter	2204	result.mul(s, m);
471	tim	1564
472	gezelter	2204	return result;
473			}
474	tim	1564
475	gezelter	2204	/**
476			* Return the multiplication of a scalra and a matrix (s * m).
477			* @return the multiplication of a scalra and a matrix
478			* @param s the scalar
479			* @param m the matrix
480			*/
481			template<typename Real, unsigned int Row, unsigned int Col>
482			inline RectMatrix<Real, Row, Col> operator *(Real s, const RectMatrix<Real, Row, Col>& m) {
483			RectMatrix<Real, Row, Col> result;
484	tim	1564
485	gezelter	2204	result.mul(s, m);
486	tim	1564
487	gezelter	2204	return result;
488			}
489	tim	1564
490	gezelter	2204	/**
491			* Return the multiplication of two matrixes (m1 * m2).
492			* @return the multiplication of two matrixes
493			* @param m1 the first matrix
494			* @param m2 the second matrix
495			*/
496			template<typename Real, unsigned int Row, unsigned int Col, unsigned int SameDim>
497			inline RectMatrix<Real, Row, Col> operator *(const RectMatrix<Real, Row, SameDim>& m1, const RectMatrix<Real, SameDim, Col>& m2) {
498			RectMatrix<Real, Row, Col> result;
499	tim	1564
500	gezelter	2204	for (unsigned int i = 0; i < Row; i++)
501			for (unsigned int j = 0; j < Col; j++)
502			for (unsigned int k = 0; k < SameDim; k++)
503			result(i, j) += m1(i, k) * m2(k, j);
504	tim	1564
505	gezelter	2204	return result;
506			}
507	tim	1564
508	gezelter	2204	/**
509			* Return the multiplication of a matrix and a vector (m * v).
510			* @return the multiplication of a matrix and a vector
511			* @param m the matrix
512			* @param v the vector
513			*/
514			template<typename Real, unsigned int Row, unsigned int Col>
515			inline Vector<Real, Row> operator *(const RectMatrix<Real, Row, Col>& m, const Vector<Real, Col>& v) {
516			Vector<Real, Row> result;
517	tim	1564
518	gezelter	2204	for (unsigned int i = 0; i < Row ; i++)
519			for (unsigned int j = 0; j < Col ; j++)
520			result[i] += m(i, j) * v[j];
521	tim	1564
522	gezelter	2204	return result;
523			}
524	tim	1564
525	gezelter	2204	/**
526			* Return the scalar division of matrix (m / s).
527			* @return the scalar division of matrix
528			* @param m the matrix
529			* @param s the scalar
530			*/
531			template<typename Real, unsigned int Row, unsigned int Col>
532			inline RectMatrix<Real, Row, Col> operator /(const RectMatrix<Real, Row, Col>& m, Real s) {
533			RectMatrix<Real, Row, Col> result;
534	tim	1564
535	gezelter	2204	result.div(s, m);
536	tim	1564
537	gezelter	2204	return result;
538			}
539	tim	1586
540	gezelter	2204	/**
541			* Write to an output stream
542			*/
543			template<typename Real, unsigned int Row, unsigned int Col>
544			std::ostream &operator<< ( std::ostream& o, const RectMatrix<Real, Row, Col>& m) {
545			for (unsigned int i = 0; i < Row ; i++) {
546			o << "(";
547			for (unsigned int j = 0; j < Col ; j++) {
548			o << m(i, j);
549			if (j != Col -1)
550			o << "\t";
551			}
552			o << ")" << std::endl;
553			}
554			return o;
555			}
556	tim	1564	}
557			#endif //MATH_RECTMATRIX_HPP