waves
Basic FE playground
tMatrix.h
Go to the documentation of this file.
1 #ifndef KATOPTRON_TMATRIX_H
2 #define KATOPTRON_TMATRIX_H
3 
4 #include "Stokhos_Tpetra_MP_Vector.hpp"
5 
6 template <class T, int rows>
7 class tVector;
8 
9 template <class T, int rows, int cols>
10 class tMatrix
11 {
12 public:
14  {
15  for (auto i = 0; i < rows; ++i)
16  for (auto j = 0; j < cols; ++j)
17  data[i * cols + j] = (T)0;
18  }
19  T &operator()(int row, int col);
21  void operator+=(tMatrix rhs);
22  T det(); //only for 3x3
23  tMatrix inv(); //only for 3x3
26  template <class U>
28  void print();
29  void clean()
30  {
31  for (auto i = 0; i < rows; ++i)
32  for (auto j = 0; j < cols; ++j)
33  data[i * cols + j] = (T)0;
34  }
35 
36 private:
37  T data[rows * cols];
38  int numRows = rows;
39  int numCols = cols;
40  int getRows() { return numRows; }
41  int getCols() { return numCols; }
42 };
43 
44 template <class T, int rows, int cols>
45 template <class U>
47 {
49  for (auto i = 0; i < rows; ++i)
50  for (auto j = 0; j < cols; ++j)
51  results(i) += data[i * cols + j] * rhs(j);
52  return results;
53 }
54 
55 template <class T, int rows, int cols>
57 {
58  tMatrix<T, rows, cols> results;
59  for (auto i = 0; i < rows; ++i)
60  for (auto j = 0; j < cols; ++j)
61  results(i, j) = data[i * cols + j] + rhs(i, j);
62  return results;
63 }
64 
65 template <class T, int rows, int cols>
67 {
68  for (auto i = 0; i < rows; ++i)
69  for (auto j = 0; j < cols; ++j)
70  data[i * cols + j] += rhs(i, j);
71 }
72 
73 template <class T, int rows, int cols>
75 {
76  tMatrix<T, rows, cols> results;
77  for (auto i = 0; i < rows; ++i)
78  for (auto j = 0; j < rows; ++j)
79  for (auto k = 0; k < rows; ++k)
80  results(i, j) += data[i * cols + k] * rhs(k, j);
81  return results;
82 }
83 
84 template <class T, int rows, int cols>
86 {
87  return data[row * cols + col];
88 }
89 
90 template <class T, int rows, int cols>
92 {
93  return data[0 * cols + 0] * (data[1 * cols + 1] * data[2 * cols + 2] - data[1 * cols + 2] * data[2 * cols + 1]) - data[0 * cols + 1] * (data[1 * cols + 0] * data[2 * cols + 2] - data[1 * cols + 2] * data[2 * cols + 0]) + data[0 * cols + 2] * (data[1 * cols + 0] * data[2 * cols + 1] - data[1 * cols + 1] * data[2 * cols + 0]);
94 }
95 
96 template <class T, int rows, int cols>
98 {
99  tMatrix<T, rows, cols> results;
100  T det = this->det();
101  results(0, 0) = (1. / det) * (data[1 * cols + 1] * data[2 * cols + 2] - data[1 * cols + 2] * data[2 * cols + 1]);
102  results(0, 1) = (1. / det) * (data[0 * cols + 2] * data[2 * cols + 1] - data[2 * cols + 2] * data[0 * cols + 1]);
103  results(0, 2) = (1. / det) * (data[0 * cols + 1] * data[1 * cols + 2] - data[1 * cols + 1] * data[0 * cols + 2]);
104 
105  results(1, 0) = (1. / det) * (data[1 * cols + 2] * data[2 * cols + 0] - data[2 * cols + 2] * data[1 * cols + 0]);
106  results(1, 1) = (1. / det) * (data[0 * cols + 0] * data[2 * cols + 2] - data[2 * cols + 0] * data[0 * cols + 2]);
107  results(1, 2) = (1. / det) * (data[0 * cols + 2] * data[1 * cols + 0] - data[1 * cols + 2] * data[0 * cols + 0]);
108 
109  results(2, 0) = (1. / det) * (data[1 * cols + 0] * data[2 * cols + 1] - data[2 * cols + 0] * data[1 * cols + 1]);
110  results(2, 1) = (1. / det) * (data[0 * cols + 1] * data[2 * cols + 0] - data[2 * cols + 1] * data[0 * cols + 0]);
111  results(2, 2) = (1. / det) * (data[0 * cols + 0] * data[1 * cols + 1] - data[1 * cols + 0] * data[0 * cols + 1]);
112  return results;
113 }
114 
115 template <class T, int rows, int cols>
117 {
118  tVector<T, rows> results;
119  for (auto i = 0; i < rows; ++i)
120  for (auto j = 0; j < cols; ++j)
121  results(i) += data[i * cols + j] * rhs(j);
122  return results;
123 }
124 
125 template <class T, int rows, int cols>
127 {
128  for (auto i = 0; i < rows; ++i)
129  {
130  for (auto j = 0; j < cols; ++j)
131  std::cout << data[i * cols + j] << " ";
132  std::cout << std::endl;
133  }
134 }
135 
136 template <class T, int rows>
137 class tVector
138 {
139 public:
141  {
142  for (auto i = 0; i < rows; ++i)
143  data[i] = (T)0;
144  }
145  T &operator()(int row);
147  void operator+=(tVector rhs);
148  T dotproduct(tVector rhs);
149  template <class U>
150  typename Sacado::Promote<T, U>::type dotproduct(tVector<U, rows> rhs);
152  void print();
153  void clean()
154  {
155  for (auto i = 0; i < rows; ++i)
156  data[i] = (T)0;
157  }
158 
159 private:
160  T data[rows];
161  int numRows = rows;
162  int getRows() { return numRows; }
163 };
164 
165 template <class T, int rows>
167 {
168  return data[row];
169 }
170 
171 template <class T, int rows>
173 {
174  tVector<T, rows> results;
175  for (auto i = 0; i < rows; ++i)
176  results(i) = data[i] + rhs(i);
177  return results;
178 }
179 
180 template <class T, int rows>
182 {
183  for (auto i = 0; i < rows; ++i)
184  data[i] += rhs(i);
185 }
186 
187 template <class T, int rows>
189 {
190  T result = (T)0;
191  for (auto i = 0; i < rows; ++i)
192  result += data[i] * rhs(i);
193  return result;
194 }
195 
196 template <class T, int rows>
197 template <class U>
198 typename Sacado::Promote<T, U>::type tVector<T, rows>::dotproduct(tVector<U, rows> rhs)
199 {
200  typename Sacado::Promote<T, U>::type result = 0;
201  for (auto i = 0; i < rows; ++i)
202  result += data[i] * rhs(i);
203  return result;
204 }
205 
206 template <class T, int rows>
208 {
209  tMatrix<T, rows, rows> results;
210  for (auto i = 0; i < rows; ++i)
211  for (auto j = 0; j < rows; ++j)
212  results(i, j) = data[i] * rhs(j);
213  return results;
214 }
215 
216 template <class T, int rows>
218 {
219  for (auto i = 0; i < rows; ++i)
220  std::cout << data[i] << std::endl;
221 }
222 
223 #endif //KATOPTRON_TMATRIX_H
tMatrix::tMatrix
tMatrix()
Definition: tMatrix.h:13
tVector::ddotproduct
tMatrix< T, rows, rows > ddotproduct(tVector rhs)
Definition: tMatrix.h:207
tMatrix::numCols
int numCols
Definition: tMatrix.h:39
tVector::numRows
int numRows
Definition: tMatrix.h:161
tMatrix::operator+
tMatrix operator+(tMatrix rhs)
Definition: tMatrix.h:56
tMatrix::det
T det()
Definition: tMatrix.h:91
tVector::operator+=
void operator+=(tVector rhs)
Definition: tMatrix.h:181
tVector::data
T data[rows]
Definition: tMatrix.h:160
tMatrix::print
void print()
Definition: tMatrix.h:126
tMatrix::clean
void clean()
Definition: tMatrix.h:29
tVector::getRows
int getRows()
Definition: tMatrix.h:162
tVector::dotproduct
T dotproduct(tVector rhs)
Definition: tMatrix.h:188
tMatrix::operator()
T & operator()(int row, int col)
Definition: tMatrix.h:85
tVector::operator+
tVector operator+(tVector rhs)
Definition: tMatrix.h:172
tVector::clean
void clean()
Definition: tMatrix.h:153
tMatrix::data
T data[rows *cols]
Definition: tMatrix.h:37
tVector
Definition: tMatrix.h:7
tVector::print
void print()
Definition: tMatrix.h:217
tMatrix
Definition: tMatrix.h:10
tMatrix::operator+=
void operator+=(tMatrix rhs)
Definition: tMatrix.h:66
tMatrix::numRows
int numRows
Definition: tMatrix.h:38
tMatrix::inv
tMatrix inv()
Definition: tMatrix.h:97
tMatrix::getRows
int getRows()
Definition: tMatrix.h:40
tMatrix::getCols
int getCols()
Definition: tMatrix.h:41
tVector::operator()
T & operator()(int row)
Definition: tMatrix.h:166
tMatrix::operator*
tMatrix operator*(tMatrix rhs)
Definition: tMatrix.h:74
tVector::tVector
tVector()
Definition: tMatrix.h:140