ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/trunk/src/math/RectMatrix.hpp
Revision: 891
Committed: Wed Feb 22 20:35:16 2006 UTC (19 years, 8 months ago) by tim
File size: 17311 byte(s)
Log Message: 
Adding Hydrodynamics Module

File Contents

# User Rev Content
1 gezelter 507 /*
2 gezelter 246 * Copyright (c) 2005 The University of Notre Dame. All Rights Reserved.
3 tim 71 *
4 gezelter 246 * 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 71 */
41 gezelter 246
42 tim 71 /**
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 253 #include <math.h>
52 tim 74 #include <cmath>
53 tim 71 #include "Vector.hpp"
54    
55     namespace oopse {
56    
57 gezelter 507 /**
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 137
67 gezelter 507 /** 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 71
74 gezelter 507 /** 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 71
81 gezelter 507 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 151
87 gezelter 507 /** copy constructor */
88     RectMatrix(const RectMatrix<Real, Row, Col>& m) {
89     *this = m;
90     }
91 tim 137
92 gezelter 507 /** destructor*/
93     ~RectMatrix() {}
94 tim 71
95 gezelter 507 /** copy assignment operator */
96     RectMatrix<Real, Row, Col>& operator =(const RectMatrix<Real, Row, Col>& m) {
97     if (this == &m)
98     return *this;
99 tim 137
100 gezelter 507 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 71
106 gezelter 507 /**
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 71
117 gezelter 507 /**
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 137
125 gezelter 507 return this->data_[i][j];
126     }
127 tim 71
128 gezelter 507 /**
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 151
140    
141 gezelter 507 /** Returns the pointer of internal array */
142     Real* getArrayPointer() {
143     return &this->data_[0][0];
144     }
145 tim 71
146 gezelter 507 /**
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 71
154 tim 891 for (unsigned int i = 0; i < Col; i++)
155 gezelter 507 v[i] = this->data_[row][i];
156 tim 71
157 gezelter 507 return v;
158     }
159 tim 71
160 gezelter 507 /**
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 71
167 tim 891 for (unsigned int i = 0; i < Col; i++)
168 gezelter 507 this->data_[row][i] = v[i];
169     }
170 tim 71
171 gezelter 507 /**
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 71
179 tim 891 for (unsigned int j = 0; j < Row; j++)
180 gezelter 507 v[j] = this->data_[j][col];
181 tim 71
182 gezelter 507 return v;
183     }
184 tim 71
185 gezelter 507 /**
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 71
192 tim 891 for (unsigned int j = 0; j < Row; j++)
193 gezelter 507 this->data_[j][col] = v[j];
194     }
195 tim 101
196 gezelter 507 /**
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 101
204 gezelter 507 for (unsigned int k = 0; k < Col; k++)
205     std::swap(this->data_[i][k], this->data_[j][k]);
206     }
207 tim 101
208 gezelter 507 /**
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 137
216 gezelter 507 for (unsigned int k = 0; k < Row; k++)
217     std::swap(this->data_[k][i], this->data_[k][j]);
218     }
219 tim 71
220 gezelter 507 /**
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 71
233 gezelter 507 return true;
234     }
235 tim 71
236 gezelter 507 /**
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 71
245 gezelter 507 /** 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 137
252 gezelter 507 /**
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 137
262 gezelter 507 /**
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 137
272 gezelter 507 /**
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 137
283 gezelter 507 /**
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 137
293 gezelter 507 /**
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 71
304 gezelter 507 /**
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 71
314 gezelter 507 /**
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 71
325 gezelter 507 /**
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 71
335 gezelter 507 /**
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 71
346 gezelter 507 /**
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 71
355 gezelter 507 /**
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 71
364 gezelter 507 /**
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 71
373 gezelter 507 /**
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 71
382 gezelter 507 /** Return the transpose of this matrix */
383     RectMatrix<Real, Col, Row> transpose() const{
384     RectMatrix<Real, Col, Row> result;
385 tim 137
386 gezelter 507 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 137
390 gezelter 507 return result;
391     }
392 tim 891
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 507 protected:
417     Real data_[Row][Col];
418     };
419 tim 71
420 gezelter 507 /** 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 71
425 gezelter 507 result.negate();
426 tim 71
427 gezelter 507 return result;
428     }
429 tim 71
430 gezelter 507 /**
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 71
440 gezelter 507 result.add(m1, m2);
441 tim 71
442 gezelter 507 return result;
443     }
444 tim 71
445 gezelter 507 /**
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 71
455 gezelter 507 result.sub(m1, m2);
456 tim 71
457 gezelter 507 return result;
458     }
459 tim 71
460 gezelter 507 /**
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 71
470 gezelter 507 result.mul(s, m);
471 tim 71
472 gezelter 507 return result;
473     }
474 tim 71
475 gezelter 507 /**
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 71
485 gezelter 507 result.mul(s, m);
486 tim 71
487 gezelter 507 return result;
488     }
489 tim 71
490 gezelter 507 /**
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 71
500 gezelter 507 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 71
505 gezelter 507 return result;
506     }
507 tim 71
508 gezelter 507 /**
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 71
518 gezelter 507 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 71
522 gezelter 507 return result;
523     }
524 tim 71
525 gezelter 507 /**
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 71
535 gezelter 507 result.div(s, m);
536 tim 71
537 gezelter 507 return result;
538     }
539 tim 93
540 gezelter 507 /**
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 71 }
557     #endif //MATH_RECTMATRIX_HPP

Properties

Name Value
svn:executable *