| 1 | /* | 
|---|
| 2 | * RealSpaceMatrix.hpp | 
|---|
| 3 | * | 
|---|
| 4 | *  Created on: Jun 25, 2010 | 
|---|
| 5 | *      Author: crueger | 
|---|
| 6 | */ | 
|---|
| 7 |  | 
|---|
| 8 | #ifndef REALSPACEMATRIX_HPP_ | 
|---|
| 9 | #define REALSPACEMATRIX_HPP_ | 
|---|
| 10 |  | 
|---|
| 11 | #include <iosfwd> | 
|---|
| 12 | #include "Helpers/defs.hpp" | 
|---|
| 13 |  | 
|---|
| 14 | /** | 
|---|
| 15 | * Simple class to store matrices with a few extra functions | 
|---|
| 16 | */ | 
|---|
| 17 |  | 
|---|
| 18 | class Vector; | 
|---|
| 19 | class MatrixContent; | 
|---|
| 20 |  | 
|---|
| 21 | /** 3x3 Matrix class. | 
|---|
| 22 | * This class has some specific features for 3D space such as rotations or | 
|---|
| 23 | * hard-coded determinants. | 
|---|
| 24 | */ | 
|---|
| 25 | class RealSpaceMatrix | 
|---|
| 26 | { | 
|---|
| 27 | friend Vector operator*(const RealSpaceMatrix&,const Vector&); | 
|---|
| 28 | public: | 
|---|
| 29 | RealSpaceMatrix(); | 
|---|
| 30 |  | 
|---|
| 31 | /** | 
|---|
| 32 | * construct a matrix from a 3x3 double array that contains all Elements. | 
|---|
| 33 | * | 
|---|
| 34 | * The elements are laid out in the following way: | 
|---|
| 35 | * array -> matrix | 
|---|
| 36 | * 0 -> (0,0) | 
|---|
| 37 | * 1 -> (1,0) | 
|---|
| 38 | * 2 -> (2,0) | 
|---|
| 39 | * 3 -> (0,1) | 
|---|
| 40 | * 4 -> (1,1) | 
|---|
| 41 | * 5 -> (2,1) | 
|---|
| 42 | * 6 -> (0,2) | 
|---|
| 43 | * 7 -> (1,2) | 
|---|
| 44 | * 8 -> (2,2) | 
|---|
| 45 | * | 
|---|
| 46 | */ | 
|---|
| 47 | RealSpaceMatrix(const double*); | 
|---|
| 48 | RealSpaceMatrix(const RealSpaceMatrix&); | 
|---|
| 49 | RealSpaceMatrix(const MatrixContent&); | 
|---|
| 50 | virtual ~RealSpaceMatrix(); | 
|---|
| 51 |  | 
|---|
| 52 | /** | 
|---|
| 53 | * Set the matrix to a unit matrix. | 
|---|
| 54 | */ | 
|---|
| 55 | void setIdentity(); | 
|---|
| 56 |  | 
|---|
| 57 | /** | 
|---|
| 58 | * Set all matrix entries to zero. | 
|---|
| 59 | */ | 
|---|
| 60 | void setZero(); | 
|---|
| 61 |  | 
|---|
| 62 | /** | 
|---|
| 63 | * Sets all matrix corresponding to a rotation matrix, | 
|---|
| 64 | * first around the x-, then the y-, finally the z-axis | 
|---|
| 65 | * with given angles. | 
|---|
| 66 | */ | 
|---|
| 67 | void setRotation(const double x, const double y, const double z); | 
|---|
| 68 |  | 
|---|
| 69 | /** | 
|---|
| 70 | * Access the matrix at index (i,j) | 
|---|
| 71 | */ | 
|---|
| 72 | double &at(size_t i, size_t j); | 
|---|
| 73 | /** | 
|---|
| 74 | * Access the matrix at index (i,j) | 
|---|
| 75 | */ | 
|---|
| 76 | const double at(size_t i, size_t j) const; | 
|---|
| 77 |  | 
|---|
| 78 | /** | 
|---|
| 79 | * Set the matrix at index (i,j). | 
|---|
| 80 | * | 
|---|
| 81 | * Slightly faster than at(i,j)=x | 
|---|
| 82 | */ | 
|---|
| 83 | void set(size_t i, size_t j, const double value); | 
|---|
| 84 |  | 
|---|
| 85 | /** | 
|---|
| 86 | * get the ith row of the matrix as a vector | 
|---|
| 87 | */ | 
|---|
| 88 | Vector &row(size_t); | 
|---|
| 89 | const Vector &row(size_t) const; | 
|---|
| 90 |  | 
|---|
| 91 | /** | 
|---|
| 92 | * get the ith column of the matrix as a vector | 
|---|
| 93 | */ | 
|---|
| 94 | Vector &column(size_t); | 
|---|
| 95 | const Vector &column(size_t) const; | 
|---|
| 96 |  | 
|---|
| 97 | /** | 
|---|
| 98 | * get the diagonal of the matrix as a vector | 
|---|
| 99 | */ | 
|---|
| 100 | Vector &diagonal(); | 
|---|
| 101 | const Vector &diagonal() const; | 
|---|
| 102 |  | 
|---|
| 103 | /** | 
|---|
| 104 | * Calculate the determinant of the matrix | 
|---|
| 105 | */ | 
|---|
| 106 | double determinant() const; | 
|---|
| 107 |  | 
|---|
| 108 | /** | 
|---|
| 109 | * Calculate the inverse of the matrix. | 
|---|
| 110 | * | 
|---|
| 111 | * Rather costly, so use precomputation as often as possible. | 
|---|
| 112 | */ | 
|---|
| 113 | RealSpaceMatrix invert() const; | 
|---|
| 114 |  | 
|---|
| 115 | /** | 
|---|
| 116 | * Diagonalizes a matrix and sets its rows to the resulting eigenvalues. | 
|---|
| 117 | * The eigenvalues are returned as a vector. | 
|---|
| 118 | * | 
|---|
| 119 | * Rather costly, so use precomputation as often as possible. | 
|---|
| 120 | */ | 
|---|
| 121 | Vector transformToEigenbasis(); | 
|---|
| 122 |  | 
|---|
| 123 | /** | 
|---|
| 124 | * Calculate the transpose of the matrix. | 
|---|
| 125 | */ | 
|---|
| 126 | RealSpaceMatrix transpose() const; | 
|---|
| 127 | RealSpaceMatrix &transpose(); | 
|---|
| 128 |  | 
|---|
| 129 | // operators | 
|---|
| 130 | RealSpaceMatrix &operator=(const RealSpaceMatrix&); | 
|---|
| 131 |  | 
|---|
| 132 | const RealSpaceMatrix &operator+=(const RealSpaceMatrix&); | 
|---|
| 133 | const RealSpaceMatrix &operator-=(const RealSpaceMatrix&); | 
|---|
| 134 | const RealSpaceMatrix &operator*=(const RealSpaceMatrix&); | 
|---|
| 135 |  | 
|---|
| 136 | const RealSpaceMatrix &operator*=(const double); | 
|---|
| 137 |  | 
|---|
| 138 | const RealSpaceMatrix operator+(const RealSpaceMatrix&) const; | 
|---|
| 139 | const RealSpaceMatrix operator-(const RealSpaceMatrix&) const; | 
|---|
| 140 | const RealSpaceMatrix operator*(const RealSpaceMatrix&) const; | 
|---|
| 141 |  | 
|---|
| 142 | bool operator==(const RealSpaceMatrix&) const; | 
|---|
| 143 |  | 
|---|
| 144 | private: | 
|---|
| 145 | RealSpaceMatrix(MatrixContent*); | 
|---|
| 146 | void createViews(); | 
|---|
| 147 | MatrixContent *content; | 
|---|
| 148 | // we keep around some Vector views of the matrix, to return references | 
|---|
| 149 | Vector* rows_ptr[NDIM]; | 
|---|
| 150 | Vector* columns_ptr[NDIM]; | 
|---|
| 151 | Vector* diagonal_ptr; | 
|---|
| 152 | }; | 
|---|
| 153 |  | 
|---|
| 154 | const RealSpaceMatrix operator*(const double,const RealSpaceMatrix&); | 
|---|
| 155 | const RealSpaceMatrix operator*(const RealSpaceMatrix&,const double); | 
|---|
| 156 |  | 
|---|
| 157 | /** | 
|---|
| 158 | * Takes a symmetric matrix that stores the lower diagonal and produces a | 
|---|
| 159 | * full matrix. | 
|---|
| 160 | * | 
|---|
| 161 | * The array is laid out as follows: | 
|---|
| 162 | * | 
|---|
| 163 | * array -> matrix | 
|---|
| 164 | * 0 -> (0,0) | 
|---|
| 165 | * 1 -> (1,0);(0,1) | 
|---|
| 166 | * 2 -> (1,1) | 
|---|
| 167 | * 3 -> (2,0);(0,2) | 
|---|
| 168 | * 4 -> (2,1);(1,2) | 
|---|
| 169 | * 5 -> (2,2) | 
|---|
| 170 | */ | 
|---|
| 171 | RealSpaceMatrix ReturnFullMatrixforSymmetric(const double * const cell_size); | 
|---|
| 172 |  | 
|---|
| 173 | std::ostream &operator<<(std::ostream&,const RealSpaceMatrix&); | 
|---|
| 174 | Vector operator*(const RealSpaceMatrix&,const Vector&); | 
|---|
| 175 | Vector& operator*=(Vector&,const RealSpaceMatrix&); | 
|---|
| 176 |  | 
|---|
| 177 | #endif /* REALSPACEMATRIX_HPP_ */ | 
|---|