RDKit
Open-source cheminformatics and machine learning.
Matrix.h
Go to the documentation of this file.
1 //
2 // Copyright (C) 2004-2006 Rational Discovery LLC
3 //
4 // @@ All Rights Reserved @@
5 // This file is part of the RDKit.
6 // The contents are covered by the terms of the BSD license
7 // which is included in the file license.txt, found at the root
8 // of the RDKit source tree.
9 //
10 #ifndef __RD_MATRIX_H__
11 #define __RD_MATRIX_H__
12 
13 #include <RDGeneral/Invariant.h>
14 #include "Vector.h"
15 #include <iostream>
16 #include <iomanip>
17 #include <cstring>
18 #include <boost/smart_ptr.hpp>
19 
20 //#ifndef INVARIANT_SILENT_METHOD
21 //#define INVARIANT_SILENT_METHOD
22 //#endif
23 
24 namespace RDNumeric {
25 
26 //! A matrix class for general, non-square matrices
27 template <class TYPE>
28 class Matrix {
29  public:
30  typedef boost::shared_array<TYPE> DATA_SPTR;
31 
32  //! Initialize with a size.
33  Matrix(unsigned int nRows, unsigned int nCols)
34  : d_nRows(nRows), d_nCols(nCols), d_dataSize(nRows * nCols) {
35  TYPE *data = new TYPE[d_dataSize];
36  memset(static_cast<void *>(data), 0, d_dataSize * sizeof(TYPE));
37  d_data.reset(data);
38  };
39 
40  //! Initialize with a size and default value.
41  Matrix(unsigned int nRows, unsigned int nCols, TYPE val)
42  : d_nRows(nRows), d_nCols(nCols), d_dataSize(nRows * nCols) {
43  TYPE *data = new TYPE[d_dataSize];
44  unsigned int i;
45  for (i = 0; i < d_dataSize; i++) {
46  data[i] = val;
47  }
48  d_data.reset(data);
49  }
50 
51  //! Initialize from a pointer.
52  /*!
53  <b>NOTE:</b> this does not take ownership of the data,
54  if you delete the data externally, this Matrix will be sad.
55  */
56  Matrix(unsigned int nRows, unsigned int nCols, DATA_SPTR data)
57  : d_nRows(nRows), d_nCols(nCols), d_dataSize(nRows * nCols) {
58  d_data = data;
59  }
60 
61  //! copy constructor
62  /*! We make a copy of the other vector's data.
63  */
64  Matrix(const Matrix<TYPE> &other)
65  : d_nRows(other.numRows()),
66  d_nCols(other.numCols()),
68  TYPE *data = new TYPE[d_dataSize];
69  const TYPE *otherData = other.getData();
70  memcpy(static_cast<void *>(data), static_cast<const void *>(otherData),
71  d_dataSize * sizeof(TYPE));
72  d_data.reset(data);
73  }
74 
75  virtual ~Matrix() {}
76 
77  //! returns the number of rows
78  inline unsigned int numRows() const { return d_nRows; }
79 
80  //! returns the number of columns
81  inline unsigned int numCols() const { return d_nCols; }
82 
83  inline unsigned int getDataSize() const { return d_dataSize; }
84 
85  //! returns a particular element of the matrix
86  inline virtual TYPE getVal(unsigned int i, unsigned int j) const {
87  PRECONDITION(i < d_nRows, "bad index");
88  PRECONDITION(j < d_nCols, "bad index");
89  unsigned int id = i * d_nCols + j;
90  return d_data[id];
91  }
92 
93  //! sets a particular element of the matrix
94  inline virtual void setVal(unsigned int i, unsigned int j, TYPE val) {
95  PRECONDITION(i < d_nRows, "bad index");
96  PRECONDITION(j < d_nCols, "bad index");
97  unsigned int id = i * d_nCols + j;
98 
99  d_data[id] = val;
100  }
101 
102  //! returns a copy of a row of the matrix
103  inline virtual void getRow(unsigned int i, Vector<TYPE> &row) const {
104  PRECONDITION(i < d_nRows, "bad index");
105  PRECONDITION(d_nCols == row.size(), "");
106  unsigned int id = i * d_nCols;
107  TYPE *rData = row.getData();
108  TYPE *data = d_data.get();
109  memcpy(static_cast<void *>(rData), static_cast<void *>(&data[id]),
110  d_nCols * sizeof(TYPE));
111  }
112 
113  //! returns a copy of a column of the matrix
114  inline virtual void getCol(unsigned int i, Vector<TYPE> &col) const {
115  PRECONDITION(i < d_nCols, "bad index");
116  PRECONDITION(d_nRows == col.size(), "");
117  unsigned int j, id;
118  TYPE *rData = col.getData();
119  TYPE *data = d_data.get();
120  for (j = 0; j < d_nRows; j++) {
121  id = j * d_nCols + i;
122  rData[j] = data[id];
123  }
124  }
125 
126  //! returns a pointer to our data array
127  inline TYPE *getData() { return d_data.get(); }
128 
129  //! returns a const pointer to our data array
130  inline const TYPE *getData() const { return d_data.get(); }
131 
132  //! Copy operator.
133  /*! We make a copy of the other Matrix's data.
134  */
135 
137  PRECONDITION(d_nRows == other.numRows(),
138  "Num rows mismatch in matrix copying");
139  PRECONDITION(d_nCols == other.numCols(),
140  "Num cols mismatch in matrix copying");
141  const TYPE *otherData = other.getData();
142  TYPE *data = d_data.get();
143  memcpy(static_cast<void *>(data), static_cast<const void *>(otherData),
144  d_dataSize * sizeof(TYPE));
145  return *this;
146  }
147 
148  //! Matrix addition.
149  /*! Perform a element by element addition of other Matrix to this Matrix
150  */
151  virtual Matrix<TYPE> &operator+=(const Matrix<TYPE> &other) {
152  PRECONDITION(d_nRows == other.numRows(),
153  "Num rows mismatch in matrix addition");
154  PRECONDITION(d_nCols == other.numCols(),
155  "Num cols mismatch in matrix addition");
156  const TYPE *oData = other.getData();
157  unsigned int i;
158  TYPE *data = d_data.get();
159  for (i = 0; i < d_dataSize; i++) {
160  data[i] += oData[i];
161  }
162  return *this;
163  }
164 
165  //! Matrix subtraction.
166  /*! Perform a element by element subtraction of other Matrix from this Matrix
167  */
168  virtual Matrix<TYPE> &operator-=(const Matrix<TYPE> &other) {
169  PRECONDITION(d_nRows == other.numRows(),
170  "Num rows mismatch in matrix addition");
171  PRECONDITION(d_nCols == other.numCols(),
172  "Num cols mismatch in matrix addition");
173  const TYPE *oData = other.getData();
174  unsigned int i;
175  TYPE *data = d_data.get();
176  for (i = 0; i < d_dataSize; i++) {
177  data[i] -= oData[i];
178  }
179  return *this;
180  }
181 
182  //! Multiplication by a scalar
183  virtual Matrix<TYPE> &operator*=(TYPE scale) {
184  unsigned int i;
185  TYPE *data = d_data.get();
186  for (i = 0; i < d_dataSize; i++) {
187  data[i] *= scale;
188  }
189  return *this;
190  }
191 
192  //! division by a scalar
193  virtual Matrix<TYPE> &operator/=(TYPE scale) {
194  unsigned int i;
195  TYPE *data = d_data.get();
196  for (i = 0; i < d_dataSize; i++) {
197  data[i] /= scale;
198  }
199  return *this;
200  }
201 
202  //! copies the transpose of this Matrix into another, returns the result
203  /*!
204  \param transpose the Matrix to store the results
205 
206  \return the transpose of this matrix.
207  This is just a reference to the argument.
208 
209  */
211  unsigned int tRows = transpose.numRows();
212  unsigned int tCols = transpose.numCols();
213  PRECONDITION(d_nCols == tRows, "Size mismatch during transposing");
214  PRECONDITION(d_nRows == tCols, "Size mismatch during transposing");
215  unsigned int i, j;
216  unsigned int idA, idAt, idT;
217  TYPE *tData = transpose.getData();
218  TYPE *data = d_data.get();
219  for (i = 0; i < d_nRows; i++) {
220  idA = i * d_nCols;
221  for (j = 0; j < d_nCols; j++) {
222  idAt = idA + j;
223  idT = j * tCols + i;
224  tData[idT] = data[idAt];
225  }
226  }
227  return transpose;
228  }
229 
230  protected:
231  Matrix() : d_nRows(0), d_nCols(0), d_dataSize(0), d_data(){};
232  unsigned int d_nRows;
233  unsigned int d_nCols;
234  unsigned int d_dataSize;
235  DATA_SPTR d_data;
236 
237  private:
238  Matrix<TYPE> &operator=(const Matrix<TYPE> &other);
239 };
240 
241 //! Matrix multiplication
242 /*!
243  Multiply a Matrix A with a second Matrix B
244  so the result is C = A*B
245 
246  \param A the the first Matrix used in the multiplication
247  \param B the Matrix by which to multiply
248  \param C Matrix to use for the results
249 
250  \return the results of multiplying A by B.
251  This is just a reference to C.
252 */
253 template <class TYPE>
255  Matrix<TYPE> &C) {
256  unsigned int aRows = A.numRows();
257  unsigned int aCols = A.numCols();
258  unsigned int cRows = C.numRows();
259  unsigned int cCols = C.numCols();
260  unsigned int bRows = B.numRows();
261  unsigned int bCols = B.numCols();
262  CHECK_INVARIANT(aCols == bRows, "Size mismatch during multiplication");
263  CHECK_INVARIANT(aRows == cRows, "Size mismatch during multiplication");
264  CHECK_INVARIANT(bCols == cCols, "Size mismatch during multiplication");
265 
266  // we have the sizes check do do the multiplication
267  TYPE *cData = C.getData();
268  const TYPE *bData = B.getData();
269  const TYPE *aData = A.getData();
270  unsigned int i, j, k;
271  unsigned int idA, idAt, idB, idC, idCt;
272  for (i = 0; i < aRows; i++) {
273  idC = i * cCols;
274  idA = i * aCols;
275  for (j = 0; j < cCols; j++) {
276  idCt = idC + j;
277  cData[idCt] = (TYPE)0.0;
278  for (k = 0; k < aCols; k++) {
279  idAt = idA + k;
280  idB = k * bCols + j;
281  cData[idCt] += (aData[idAt] * bData[idB]);
282  }
283  }
284  }
285  return C;
286 };
287 
288 //! Matrix-Vector multiplication
289 /*!
290  Multiply a Matrix A with a Vector x
291  so the result is y = A*x
292 
293  \param A the matrix used in the multiplication
294  \param x the Vector by which to multiply
295  \param y Vector to use for the results
296 
297  \return the results of multiplying x by this
298  This is just a reference to y.
299 */
300 template <class TYPE>
302  Vector<TYPE> &y) {
303  unsigned int aRows = A.numRows();
304  unsigned int aCols = A.numCols();
305  unsigned int xSiz = x.size();
306  unsigned int ySiz = y.size();
307  CHECK_INVARIANT(aCols == xSiz, "Size mismatch during multiplication");
308  CHECK_INVARIANT(aRows == ySiz, "Size mismatch during multiplication");
309  unsigned int i, j;
310  unsigned int idA, idAt;
311  const TYPE *xData = x.getData();
312  const TYPE *aData = A.getData();
313  TYPE *yData = y.getData();
314  for (i = 0; i < aRows; i++) {
315  idA = i * aCols;
316  yData[i] = (TYPE)(0.0);
317  for (j = 0; j < aCols; j++) {
318  idAt = idA + j;
319  yData[i] += (aData[idAt] * xData[j]);
320  }
321  }
322  return y;
323 };
324 
326 };
327 
328 //! ostream operator for Matrix's
329 template <class TYPE>
330 std::ostream &operator<<(std::ostream &target,
331  const RDNumeric::Matrix<TYPE> &mat) {
332  unsigned int nr = mat.numRows();
333  unsigned int nc = mat.numCols();
334  target << "Rows: " << mat.numRows() << " Columns: " << mat.numCols() << "\n";
335 
336  unsigned int i, j;
337  for (i = 0; i < nr; i++) {
338  for (j = 0; j < nc; j++) {
339  target << std::setw(7) << std::setprecision(3) << mat.getVal(i, j);
340  }
341  target << "\n";
342  }
343  return target;
344 }
345 
346 #endif
Matrix(unsigned int nRows, unsigned int nCols)
Initialize with a size.
Definition: Matrix.h:33
DATA_SPTR d_data
Definition: Matrix.h:235
virtual Matrix< TYPE > & operator-=(const Matrix< TYPE > &other)
Matrix subtraction.
Definition: Matrix.h:168
unsigned int numRows() const
returns the number of rows
Definition: Matrix.h:78
std::ostream & operator<<(std::ostream &target, const RDNumeric::Matrix< TYPE > &mat)
ostream operator for Matrix&#39;s
Definition: Matrix.h:330
#define CHECK_INVARIANT(expr, mess)
Definition: Invariant.h:99
virtual Matrix< TYPE > & operator*=(TYPE scale)
Multiplication by a scalar.
Definition: Matrix.h:183
virtual void setVal(unsigned int i, unsigned int j, TYPE val)
sets a particular element of the matrix
Definition: Matrix.h:94
Matrix(unsigned int nRows, unsigned int nCols, TYPE val)
Initialize with a size and default value.
Definition: Matrix.h:41
Matrix< TYPE > & assign(const Matrix< TYPE > &other)
Copy operator.
Definition: Matrix.h:136
virtual Matrix< TYPE > & transpose(Matrix< TYPE > &transpose) const
copies the transpose of this Matrix into another, returns the result
Definition: Matrix.h:210
unsigned int size() const
return the size (dimension) of the vector
Definition: Vector.h:77
virtual ~Matrix()
Definition: Matrix.h:75
unsigned int d_nCols
Definition: Matrix.h:233
Matrix(unsigned int nRows, unsigned int nCols, DATA_SPTR data)
Initialize from a pointer.
Definition: Matrix.h:56
A matrix class for general, non-square matrices.
Definition: Matrix.h:28
boost::shared_array< TYPE > DATA_SPTR
Definition: Matrix.h:30
const TYPE * getData() const
returns a const pointer to our data array
Definition: Matrix.h:130
unsigned int numCols() const
returns the number of columns
Definition: Matrix.h:81
Matrix< double > DoubleMatrix
Definition: Matrix.h:323
Matrix< TYPE > & multiply(const Matrix< TYPE > &A, const Matrix< TYPE > &B, Matrix< TYPE > &C)
Matrix multiplication.
Definition: Matrix.h:254
unsigned int d_nRows
Definition: Matrix.h:231
virtual TYPE getVal(unsigned int i, unsigned int j) const
returns a particular element of the matrix
Definition: Matrix.h:86
virtual Matrix< TYPE > & operator+=(const Matrix< TYPE > &other)
Matrix addition.
Definition: Matrix.h:151
virtual Matrix< TYPE > & operator/=(TYPE scale)
division by a scalar
Definition: Matrix.h:193
#define PRECONDITION(expr, mess)
Definition: Invariant.h:107
virtual void getCol(unsigned int i, Vector< TYPE > &col) const
returns a copy of a column of the matrix
Definition: Matrix.h:114
virtual void getRow(unsigned int i, Vector< TYPE > &row) const
returns a copy of a row of the matrix
Definition: Matrix.h:103
TYPE * getData()
returns a pointer to our data array
Definition: Vector.h:102
A class to represent vectors of numbers.
Definition: Vector.h:28
unsigned int d_dataSize
Definition: Matrix.h:234
TYPE * getData()
returns a pointer to our data array
Definition: Matrix.h:127
unsigned int getDataSize() const
Definition: Matrix.h:83
Matrix(const Matrix< TYPE > &other)
copy constructor
Definition: Matrix.h:64