dune-common  2.7.0
fmatrix.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=4 sw=2 sts=2:
3 #ifndef DUNE_FMATRIX_HH
4 #define DUNE_FMATRIX_HH
5 
6 #include <cmath>
7 #include <cstddef>
8 #include <iostream>
9 #include <algorithm>
10 #include <initializer_list>
11 
14 #include <dune/common/fvector.hh>
16 #include <dune/common/precision.hh>
19 
20 namespace Dune
21 {
22 
34  template< class K, int ROWS, int COLS = ROWS > class FieldMatrix;
35 
36  template< class K, int ROWS, int COLS >
37  struct DenseMatVecTraits< FieldMatrix<K,ROWS,COLS> >
38  {
40 
41  // each row is implemented by a field vector
43 
45  typedef const row_type &const_row_reference;
46 
47  typedef std::array<row_type,ROWS> container_type;
48  typedef K value_type;
49  typedef typename container_type::size_type size_type;
50  };
51 
52  template< class K, int ROWS, int COLS >
53  struct FieldTraits< FieldMatrix<K,ROWS,COLS> >
54  {
57  };
58 
67  template<class K, int ROWS, int COLS>
68  class FieldMatrix : public DenseMatrix< FieldMatrix<K,ROWS,COLS> >
69  {
70  std::array< FieldVector<K,COLS>, ROWS > _data;
72  public:
73 
75  enum {
77  rows = ROWS,
79  cols = COLS
80  };
81 
82  typedef typename Base::size_type size_type;
83  typedef typename Base::row_type row_type;
84 
87 
88  //===== constructors
91  constexpr FieldMatrix() = default;
92 
95  FieldMatrix(std::initializer_list<Dune::FieldVector<K, cols> > const &l) {
96  assert(l.size() == rows); // Actually, this is not needed any more!
97  std::copy_n(l.begin(), std::min(static_cast<std::size_t>(ROWS),
98  l.size()),
99  _data.begin());
100  }
101 
102  template <class T,
103  typename = std::enable_if_t<HasDenseMatrixAssigner<FieldMatrix, T>::value>>
104  FieldMatrix(T const& rhs)
105  {
106  *this = rhs;
107  }
108 
109  using Base::operator=;
110 
112  FieldMatrix& operator=(const FieldMatrix&) = default;
113 
115  template<typename T>
117  {
118  _data = x._data;
119  return *this;
120  }
121 
123  template <typename T, int rows, int cols>
125 
127  template <class OtherScalar>
128  friend auto operator+ ( const FieldMatrix& matrixA,
129  const FieldMatrix<OtherScalar,ROWS,COLS>& matrixB)
130  {
132 
133  for (size_type i = 0; i < ROWS; ++i)
134  for (size_type j = 0; j < COLS; ++j)
135  result[i][j] = matrixA[i][j] + matrixB[i][j];
136 
137  return result;
138  }
139 
141  template <class OtherScalar>
142  friend auto operator- ( const FieldMatrix& matrixA,
143  const FieldMatrix<OtherScalar,ROWS,COLS>& matrixB)
144  {
146 
147  for (size_type i = 0; i < ROWS; ++i)
148  for (size_type j = 0; j < COLS; ++j)
149  result[i][j] = matrixA[i][j] - matrixB[i][j];
150 
151  return result;
152  }
153 
155  template <class Scalar,
156  std::enable_if_t<IsNumber<Scalar>::value, int> = 0>
157  friend auto operator* ( const FieldMatrix& matrix, Scalar scalar)
158  {
160 
161  for (size_type i = 0; i < ROWS; ++i)
162  for (size_type j = 0; j < COLS; ++j)
163  result[i][j] = matrix[i][j] * scalar;
164 
165  return result;
166  }
167 
169  template <class Scalar,
170  std::enable_if_t<IsNumber<Scalar>::value, int> = 0>
171  friend auto operator* ( Scalar scalar, const FieldMatrix& matrix)
172  {
174 
175  for (size_type i = 0; i < ROWS; ++i)
176  for (size_type j = 0; j < COLS; ++j)
177  result[i][j] = scalar * matrix[i][j];
178 
179  return result;
180  }
181 
183  template <class Scalar,
184  std::enable_if_t<IsNumber<Scalar>::value, int> = 0>
185  friend auto operator/ ( const FieldMatrix& matrix, Scalar scalar)
186  {
188 
189  for (size_type i = 0; i < ROWS; ++i)
190  for (size_type j = 0; j < COLS; ++j)
191  result[i][j] = matrix[i][j] / scalar;
192 
193  return result;
194  }
195 
198  template <class OtherScalar, int otherCols>
199  friend auto operator* ( const FieldMatrix& matrixA,
201  {
203 
204  for (size_type i = 0; i < matrixA.mat_rows(); ++i)
205  for (size_type j = 0; j < matrixB.mat_cols(); ++j)
206  {
207  result[i][j] = 0;
208  for (size_type k = 0; k < matrixA.mat_cols(); ++k)
209  result[i][j] += matrixA[i][k] * matrixB[k][j];
210  }
211 
212  return result;
213  }
214 
216  template<int l>
218  {
220 
221  for (size_type i=0; i<l; i++) {
222  for (size_type j=0; j<cols; j++) {
223  C[i][j] = 0;
224  for (size_type k=0; k<rows; k++)
225  C[i][j] += M[i][k]*(*this)[k][j];
226  }
227  }
228  return C;
229  }
230 
231  using Base::rightmultiply;
232 
234  template <int r, int c>
236  {
237  static_assert(r == c, "Cannot rightmultiply with non-square matrix");
238  static_assert(r == cols, "Size mismatch");
239  FieldMatrix<K,rows,cols> C(*this);
240 
241  for (size_type i=0; i<rows; i++)
242  for (size_type j=0; j<cols; j++) {
243  (*this)[i][j] = 0;
244  for (size_type k=0; k<cols; k++)
245  (*this)[i][j] += C[i][k]*M[k][j];
246  }
247  return *this;
248  }
249 
251  template<int l>
253  {
255 
256  for (size_type i=0; i<rows; i++) {
257  for (size_type j=0; j<l; j++) {
258  C[i][j] = 0;
259  for (size_type k=0; k<cols; k++)
260  C[i][j] += (*this)[i][k]*M[k][j];
261  }
262  }
263  return C;
264  }
265 
266  // make this thing a matrix
267  constexpr size_type mat_rows() const { return ROWS; }
268  constexpr size_type mat_cols() const { return COLS; }
269 
271  {
272  DUNE_ASSERT_BOUNDS(i < ROWS);
273  return _data[i];
274  }
275 
277  {
278  DUNE_ASSERT_BOUNDS(i < ROWS);
279  return _data[i];
280  }
281  };
282 
283 #ifndef DOXYGEN // hide specialization
284 
286  template<class K>
287  class FieldMatrix<K,1,1> : public DenseMatrix< FieldMatrix<K,1,1> >
288  {
289  FieldVector<K,1> _data;
290  typedef DenseMatrix< FieldMatrix<K,1,1> > Base;
291  public:
292  // standard constructor and everything is sufficient ...
293 
294  //===== type definitions and constants
295 
297  typedef typename Base::size_type size_type;
298 
300  enum {
303  blocklevel = 1
304  };
305 
306  typedef typename Base::row_type row_type;
307 
308  typedef typename Base::row_reference row_reference;
310 
312  enum {
315  rows = 1,
318  cols = 1
319  };
320 
321  //===== constructors
324  FieldMatrix() = default;
325 
328  FieldMatrix(std::initializer_list<Dune::FieldVector<K, 1>> const &l)
329  {
330  std::copy_n(l.begin(), std::min(static_cast< std::size_t >( 1 ), l.size()), &_data);
331  }
332 
333  template <class T,
334  typename = std::enable_if_t<HasDenseMatrixAssigner<FieldMatrix, T>::value>>
335  FieldMatrix(T const& rhs)
336  {
337  *this = rhs;
338  }
339 
340  using Base::operator=;
341 
343  template <class OtherScalar>
344  friend auto operator+ ( const FieldMatrix& matrixA,
345  const FieldMatrix<OtherScalar,1,1>& matrixB)
346  {
347  return FieldMatrix<typename PromotionTraits<K,OtherScalar>::PromotedType,1,1>{matrixA[0][0] + matrixB[0][0]};
348  }
349 
351  template <class Scalar,
352  std::enable_if_t<IsNumber<Scalar>::value, int> = 0>
353  friend auto operator+ ( const FieldMatrix& matrix,
354  const Scalar& scalar)
355  {
356  return FieldMatrix<typename PromotionTraits<K,Scalar>::PromotedType,1,1>{matrix[0][0] + scalar};
357  }
358 
360  template <class Scalar,
361  std::enable_if_t<IsNumber<Scalar>::value, int> = 0>
362  friend auto operator+ ( const Scalar& scalar,
363  const FieldMatrix& matrix)
364  {
365  return FieldMatrix<typename PromotionTraits<Scalar,K>::PromotedType,1,1>{scalar + matrix[0][0]};
366  }
367 
369  template <class OtherScalar>
370  friend auto operator- ( const FieldMatrix& matrixA,
371  const FieldMatrix<OtherScalar,1,1>& matrixB)
372  {
373  return FieldMatrix<typename PromotionTraits<K,OtherScalar>::PromotedType,1,1>{matrixA[0][0] - matrixB[0][0]};
374  }
375 
377  template <class Scalar,
378  std::enable_if_t<IsNumber<Scalar>::value, int> = 0>
379  friend auto operator- ( const FieldMatrix& matrix,
380  const Scalar& scalar)
381  {
382  return FieldMatrix<typename PromotionTraits<K,Scalar>::PromotedType,1,1>{matrix[0][0] - scalar};
383  }
384 
386  template <class Scalar,
387  std::enable_if_t<IsNumber<Scalar>::value, int> = 0>
388  friend auto operator- ( const Scalar& scalar,
389  const FieldMatrix& matrix)
390  {
391  return FieldMatrix<typename PromotionTraits<Scalar,K>::PromotedType,1,1>{scalar - matrix[0][0]};
392  }
393 
395  template <class Scalar,
396  std::enable_if_t<IsNumber<Scalar>::value, int> = 0>
397  friend auto operator* ( const FieldMatrix& matrix, Scalar scalar)
398  {
399  return FieldMatrix<typename PromotionTraits<K,Scalar>::PromotedType,1,1> {matrix[0][0] * scalar};
400  }
401 
403  template <class Scalar,
404  std::enable_if_t<IsNumber<Scalar>::value, int> = 0>
405  friend auto operator* ( Scalar scalar, const FieldMatrix& matrix)
406  {
407  return FieldMatrix<typename PromotionTraits<K,Scalar>::PromotedType,1,1> {scalar * matrix[0][0]};
408  }
409 
411  template <class Scalar,
412  std::enable_if_t<IsNumber<Scalar>::value, int> = 0>
413  friend auto operator/ ( const FieldMatrix& matrix, Scalar scalar)
414  {
415  return FieldMatrix<typename PromotionTraits<K,Scalar>::PromotedType,1,1> {matrix[0][0] / scalar};
416  }
417 
418  //===== solve
419 
422  template <class OtherScalar, int otherCols>
423  friend auto operator* ( const FieldMatrix& matrixA,
424  const FieldMatrix<OtherScalar, 1, otherCols>& matrixB)
425  {
426  FieldMatrix<typename PromotionTraits<K,OtherScalar>::PromotedType,1,otherCols> result;
427 
428  for (size_type j = 0; j < matrixB.mat_cols(); ++j)
429  result[0][j] = matrixA[0][0] * matrixB[0][j];
430 
431  return result;
432  }
433 
435  template<int l>
436  FieldMatrix<K,l,1> leftmultiplyany (const FieldMatrix<K,l,1>& M) const
437  {
438  FieldMatrix<K,l,1> C;
439  for (size_type j=0; j<l; j++)
440  C[j][0] = M[j][0]*(*this)[0][0];
441  return C;
442  }
443 
446  {
447  _data[0] *= M[0][0];
448  return *this;
449  }
450 
452  template<int l>
453  FieldMatrix<K,1,l> rightmultiplyany (const FieldMatrix<K,1,l>& M) const
454  {
455  FieldMatrix<K,1,l> C;
456 
457  for (size_type j=0; j<l; j++)
458  C[0][j] = M[0][j]*_data[0];
459  return C;
460  }
461 
462  // make this thing a matrix
463  constexpr size_type mat_rows() const { return 1; }
464  constexpr size_type mat_cols() const { return 1; }
465 
467  {
469  DUNE_ASSERT_BOUNDS(i == 0);
470  return _data;
471  }
472 
474  {
476  DUNE_ASSERT_BOUNDS(i == 0);
477  return _data;
478  }
479 
481  FieldMatrix& operator+= (const K& k)
482  {
483  _data[0] += k;
484  return (*this);
485  }
486 
488  FieldMatrix& operator-= (const K& k)
489  {
490  _data[0] -= k;
491  return (*this);
492  }
493 
495  FieldMatrix& operator*= (const K& k)
496  {
497  _data[0] *= k;
498  return (*this);
499  }
500 
502  FieldMatrix& operator/= (const K& k)
503  {
504  _data[0] /= k;
505  return (*this);
506  }
507 
508  //===== conversion operator
509 
510  operator const K& () const { return _data[0]; }
511 
512  };
513 
515  template<typename K>
516  std::ostream& operator<< (std::ostream& s, const FieldMatrix<K,1,1>& a)
517  {
518  s << a[0][0];
519  return s;
520  }
521 
522 #endif // DOXYGEN
523 
524  namespace FMatrixHelp {
525 
527  template <typename K>
528  static inline K invertMatrix (const FieldMatrix<K,1,1> &matrix, FieldMatrix<K,1,1> &inverse)
529  {
530  using real_type = typename FieldTraits<K>::real_type;
531  inverse[0][0] = real_type(1.0)/matrix[0][0];
532  return matrix[0][0];
533  }
534 
536  template <typename K>
537  static inline K invertMatrix_retTransposed (const FieldMatrix<K,1,1> &matrix, FieldMatrix<K,1,1> &inverse)
538  {
539  return invertMatrix(matrix,inverse);
540  }
541 
542 
544  template <typename K>
545  static inline K invertMatrix (const FieldMatrix<K,2,2> &matrix, FieldMatrix<K,2,2> &inverse)
546  {
547  using real_type = typename FieldTraits<K>::real_type;
548  // code generated by maple
549  K det = (matrix[0][0]*matrix[1][1] - matrix[0][1]*matrix[1][0]);
550  K det_1 = real_type(1.0)/det;
551  inverse[0][0] = matrix[1][1] * det_1;
552  inverse[0][1] = - matrix[0][1] * det_1;
553  inverse[1][0] = - matrix[1][0] * det_1;
554  inverse[1][1] = matrix[0][0] * det_1;
555  return det;
556  }
557 
560  template <typename K>
561  static inline K invertMatrix_retTransposed (const FieldMatrix<K,2,2> &matrix, FieldMatrix<K,2,2> &inverse)
562  {
563  using real_type = typename FieldTraits<K>::real_type;
564  // code generated by maple
565  K det = (matrix[0][0]*matrix[1][1] - matrix[0][1]*matrix[1][0]);
566  K det_1 = real_type(1.0)/det;
567  inverse[0][0] = matrix[1][1] * det_1;
568  inverse[1][0] = - matrix[0][1] * det_1;
569  inverse[0][1] = - matrix[1][0] * det_1;
570  inverse[1][1] = matrix[0][0] * det_1;
571  return det;
572  }
573 
575  template <typename K>
576  static inline K invertMatrix (const FieldMatrix<K,3,3> &matrix, FieldMatrix<K,3,3> &inverse)
577  {
578  using real_type = typename FieldTraits<K>::real_type;
579  // code generated by maple
580  K t4 = matrix[0][0] * matrix[1][1];
581  K t6 = matrix[0][0] * matrix[1][2];
582  K t8 = matrix[0][1] * matrix[1][0];
583  K t10 = matrix[0][2] * matrix[1][0];
584  K t12 = matrix[0][1] * matrix[2][0];
585  K t14 = matrix[0][2] * matrix[2][0];
586 
587  K det = (t4*matrix[2][2]-t6*matrix[2][1]-t8*matrix[2][2]+
588  t10*matrix[2][1]+t12*matrix[1][2]-t14*matrix[1][1]);
589  K t17 = real_type(1.0)/det;
590 
591  inverse[0][0] = (matrix[1][1] * matrix[2][2] - matrix[1][2] * matrix[2][1])*t17;
592  inverse[0][1] = -(matrix[0][1] * matrix[2][2] - matrix[0][2] * matrix[2][1])*t17;
593  inverse[0][2] = (matrix[0][1] * matrix[1][2] - matrix[0][2] * matrix[1][1])*t17;
594  inverse[1][0] = -(matrix[1][0] * matrix[2][2] - matrix[1][2] * matrix[2][0])*t17;
595  inverse[1][1] = (matrix[0][0] * matrix[2][2] - t14) * t17;
596  inverse[1][2] = -(t6-t10) * t17;
597  inverse[2][0] = (matrix[1][0] * matrix[2][1] - matrix[1][1] * matrix[2][0]) * t17;
598  inverse[2][1] = -(matrix[0][0] * matrix[2][1] - t12) * t17;
599  inverse[2][2] = (t4-t8) * t17;
600 
601  return det;
602  }
603 
605  template <typename K>
606  static inline K invertMatrix_retTransposed (const FieldMatrix<K,3,3> &matrix, FieldMatrix<K,3,3> &inverse)
607  {
608  using real_type = typename FieldTraits<K>::real_type;
609  // code generated by maple
610  K t4 = matrix[0][0] * matrix[1][1];
611  K t6 = matrix[0][0] * matrix[1][2];
612  K t8 = matrix[0][1] * matrix[1][0];
613  K t10 = matrix[0][2] * matrix[1][0];
614  K t12 = matrix[0][1] * matrix[2][0];
615  K t14 = matrix[0][2] * matrix[2][0];
616 
617  K det = (t4*matrix[2][2]-t6*matrix[2][1]-t8*matrix[2][2]+
618  t10*matrix[2][1]+t12*matrix[1][2]-t14*matrix[1][1]);
619  K t17 = real_type(1.0)/det;
620 
621  inverse[0][0] = (matrix[1][1] * matrix[2][2] - matrix[1][2] * matrix[2][1])*t17;
622  inverse[1][0] = -(matrix[0][1] * matrix[2][2] - matrix[0][2] * matrix[2][1])*t17;
623  inverse[2][0] = (matrix[0][1] * matrix[1][2] - matrix[0][2] * matrix[1][1])*t17;
624  inverse[0][1] = -(matrix[1][0] * matrix[2][2] - matrix[1][2] * matrix[2][0])*t17;
625  inverse[1][1] = (matrix[0][0] * matrix[2][2] - t14) * t17;
626  inverse[2][1] = -(t6-t10) * t17;
627  inverse[0][2] = (matrix[1][0] * matrix[2][1] - matrix[1][1] * matrix[2][0]) * t17;
628  inverse[1][2] = -(matrix[0][0] * matrix[2][1] - t12) * t17;
629  inverse[2][2] = (t4-t8) * t17;
630 
631  return det;
632  }
633 
635  template< class K, int m, int n, int p >
636  static inline void multMatrix ( const FieldMatrix< K, m, n > &A,
637  const FieldMatrix< K, n, p > &B,
639  {
640  typedef typename FieldMatrix< K, m, p > :: size_type size_type;
641 
642  for( size_type i = 0; i < m; ++i )
643  {
644  for( size_type j = 0; j < p; ++j )
645  {
646  ret[ i ][ j ] = K( 0 );
647  for( size_type k = 0; k < n; ++k )
648  ret[ i ][ j ] += A[ i ][ k ] * B[ k ][ j ];
649  }
650  }
651  }
652 
654  template <typename K, int rows, int cols>
656  {
657  typedef typename FieldMatrix<K,rows,cols>::size_type size_type;
658 
659  for(size_type i=0; i<cols; i++)
660  for(size_type j=0; j<cols; j++)
661  {
662  ret[i][j]=0.0;
663  for(size_type k=0; k<rows; k++)
664  ret[i][j]+=matrix[k][i]*matrix[k][j];
665  }
666  }
667 
669 
671  template <typename K, int rows, int cols>
672  static inline void multAssignTransposed( const FieldMatrix<K,rows,cols> &matrix, const FieldVector<K,rows> & x, FieldVector<K,cols> & ret)
673  {
674  typedef typename FieldMatrix<K,rows,cols>::size_type size_type;
675 
676  for(size_type i=0; i<cols; ++i)
677  {
678  ret[i] = 0.0;
679  for(size_type j=0; j<rows; ++j)
680  ret[i] += matrix[j][i]*x[j];
681  }
682  }
683 
685  template <typename K, int rows, int cols>
687  {
689  multAssign(matrix,x,ret);
690  return ret;
691  }
692 
694  template <typename K, int rows, int cols>
696  {
698  multAssignTransposed( matrix, x, ret );
699  return ret;
700  }
701 
702  } // end namespace FMatrixHelp
703 
706 } // end namespace
707 
708 #include "fmatrixev.hh"
709 #endif
Dune::FieldMatrix::mat_cols
constexpr size_type mat_cols() const
Definition: fmatrix.hh:268
Dune::FMatrixHelp::multTransposed
static FieldVector< K, cols > multTransposed(const FieldMatrix< K, rows, cols > &matrix, const FieldVector< K, rows > &x)
calculates ret = matrix^T * x
Definition: fmatrix.hh:695
exceptions.hh
A few common exception classes.
Dune::DenseMatVecTraits< FieldMatrix< K, ROWS, COLS > >::value_type
K value_type
Definition: fmatrix.hh:48
Dune::FieldTraits::real_type
T real_type
export the type representing the real type of the field
Definition: ftraits.hh:28
Dune::FieldMatrix::FieldMatrix
FieldMatrix(std::initializer_list< Dune::FieldVector< K, cols > > const &l)
Constructor initializing the matrix from a list of vector.
Definition: fmatrix.hh:95
Dune::FieldMatrix::leftmultiplyany
FieldMatrix< K, l, cols > leftmultiplyany(const FieldMatrix< K, l, rows > &M) const
Multiplies M from the left to this matrix, this matrix is not modified.
Definition: fmatrix.hh:217
Dune::DenseMatVecTraits< FieldMatrix< K, ROWS, COLS > >::size_type
container_type::size_type size_type
Definition: fmatrix.hh:49
Dune::DenseMatVecTraits< FieldMatrix< K, ROWS, COLS > >::const_row_reference
const typedef row_type & const_row_reference
Definition: fmatrix.hh:45
promotiontraits.hh
Compute type of the result of an arithmetic operation involving two different number types.
Dune::AlignedNumberImpl::min
auto min(const AlignedNumber< T, align > &a, const AlignedNumber< T, align > &b)
Definition: debugalign.hh:434
Dune::FieldTraits< FieldMatrix< K, ROWS, COLS > >::real_type
FieldTraits< K >::real_type real_type
Definition: fmatrix.hh:56
Dune::FieldMatrix::mat_access
row_reference mat_access(size_type i)
Definition: fmatrix.hh:270
Dune::FieldMatrix::size_type
Base::size_type size_type
Definition: fmatrix.hh:82
boundschecking.hh
Macro for wrapping boundary checks.
DUNE_UNUSED_PARAMETER
#define DUNE_UNUSED_PARAMETER(parm)
A macro to mark intentionally unused function parameters with.
Definition: unused.hh:25
Dune::FieldMatrix::mat_access
const_row_reference mat_access(size_type i) const
Definition: fmatrix.hh:276
fmatrixev.hh
Eigenvalue computations for the FieldMatrix class.
Dune::FieldTraits< FieldMatrix< K, ROWS, COLS > >::field_type
FieldTraits< K >::field_type field_type
Definition: fmatrix.hh:55
Dune::FMatrixHelp::multTransposedMatrix
static void multTransposedMatrix(const FieldMatrix< K, rows, cols > &matrix, FieldMatrix< K, cols, cols > &ret)
calculates ret= A_t*A
Definition: fmatrix.hh:655
Dune::DenseMatrix::rightmultiply
MAT & rightmultiply(const DenseMatrix< M2 > &M)
Multiplies M from the right to this matrix.
Definition: densematrix.hh:672
Dune::FieldMatrix::cols
@ cols
The number of columns.
Definition: fmatrix.hh:79
Dune::FieldMatrix
A dense n x m matrix.
Definition: densematrix.hh:41
Dune::DenseMatrix::const_row_reference
Traits::const_row_reference const_row_reference
The type used to represent a reference to a constant row (usually const row_type &)
Definition: densematrix.hh:199
Dune::FieldTraits::field_type
T field_type
export the type representing the field
Definition: ftraits.hh:26
fvector.hh
Implements a vector constructed from a given type representing a field and a compile-time given size.
Dune::DenseMatrix
A dense n x m matrix.
Definition: densematrix.hh:28
Dune::FieldMatrix::rightmultiply
FieldMatrix & rightmultiply(const FieldMatrix< K, r, c > &M)
Multiplies M from the right to this matrix.
Definition: fmatrix.hh:235
typetraits.hh
Traits for type conversions and type information.
Dune::FieldMatrix::FieldMatrix
FieldMatrix(T const &rhs)
Definition: fmatrix.hh:104
Dune::FieldMatrix::operator+
friend auto operator+(const FieldMatrix &matrixA, const FieldMatrix< OtherScalar, ROWS, COLS > &matrixB)
vector space addition – two-argument version
Definition: fmatrix.hh:128
Dune::DenseMatrixHelp::multAssign
static void multAssign(const DenseMatrix< MAT > &matrix, const DenseVector< V1 > &x, DenseVector< V2 > &ret)
calculates ret = matrix * x
Definition: densematrix.hh:1198
Dune::FieldMatrix::rows
@ rows
The number of rows.
Definition: fmatrix.hh:77
Dune::FieldMatrix::operator=
FieldMatrix & operator=(const FieldMatrix &)=default
copy assignment operator
Dune::operator<<
std::ostream & operator<<(std::ostream &s, const bigunsignedint< k > &x)
Definition: bigunsignedint.hh:273
Dune::DenseMatrix::row_type
Traits::row_type row_type
The type used to represent a row (must fulfill the Dune::DenseVector interface)
Definition: densematrix.hh:193
Dune::FMatrixHelp::multMatrix
static void multMatrix(const FieldMatrix< K, m, n > &A, const FieldMatrix< K, n, p > &B, FieldMatrix< K, m, p > &ret)
calculates ret = A * B
Definition: fmatrix.hh:636
Dune::FieldMatrix::rightmultiplyany
FieldMatrix< K, rows, l > rightmultiplyany(const FieldMatrix< K, cols, l > &M) const
Multiplies M from the right to this matrix, this matrix is not modified.
Definition: fmatrix.hh:252
Dune::DenseMatrix::row_reference
Traits::row_reference row_reference
The type used to represent a reference to a row (usually row_type &)
Definition: densematrix.hh:196
Dune::FieldTraits
Definition: ftraits.hh:23
Dune::FMatrixHelp::mult
static FieldVector< K, rows > mult(const FieldMatrix< K, rows, cols > &matrix, const FieldVector< K, cols > &x)
calculates ret = matrix * x
Definition: fmatrix.hh:686
Dune::FMatrixHelp::invertMatrix_retTransposed
static K invertMatrix_retTransposed(const FieldMatrix< K, 1, 1 > &matrix, FieldMatrix< K, 1, 1 > &inverse)
invert scalar without changing the original matrix
Definition: fmatrix.hh:537
Dune::DenseMatrix::size_type
Traits::size_type size_type
The type used for the index access and size operation.
Definition: densematrix.hh:190
Dune::FieldMatrix::row_type
Base::row_type row_type
Definition: fmatrix.hh:83
Dune::DenseMatVecTraits< FieldMatrix< K, ROWS, COLS > >::row_reference
row_type & row_reference
Definition: fmatrix.hh:44
Dune::FieldMatrix::operator/
friend auto operator/(const FieldMatrix &matrix, Scalar scalar)
vector space division by scalar
Definition: fmatrix.hh:185
Dune::FieldMatrix::operator*
friend auto operator*(const FieldMatrix &matrix, Scalar scalar)
vector space multiplication with scalar
Definition: fmatrix.hh:157
Dune::DenseMatVecTraits
Definition: matvectraits.hh:29
Dune::FieldVector< K, COLS >
precision.hh
Various precision settings for calculations with FieldMatrix and FieldVector.
Dune::FieldMatrix::FieldMatrix
constexpr FieldMatrix()=default
Default constructor.
Dune::DenseMatVecTraits< FieldMatrix< K, ROWS, COLS > >::derived_type
FieldMatrix< K, ROWS, COLS > derived_type
Definition: fmatrix.hh:39
Dune::FieldMatrix::row_reference
Base::row_reference row_reference
Definition: fmatrix.hh:85
Dune::DenseMatVecTraits< FieldMatrix< K, ROWS, COLS > >::container_type
std::array< row_type, ROWS > container_type
Definition: fmatrix.hh:47
Dune::FMatrixHelp::multAssignTransposed
static void multAssignTransposed(const FieldMatrix< K, rows, cols > &matrix, const FieldVector< K, rows > &x, FieldVector< K, cols > &ret)
calculates ret = matrix^T * x
Definition: fmatrix.hh:672
DUNE_ASSERT_BOUNDS
#define DUNE_ASSERT_BOUNDS(cond)
If DUNE_CHECK_BOUNDS is defined: check if condition cond holds; otherwise, do nothing.
Definition: boundschecking.hh:28
Dune::FieldMatrix::mat_rows
constexpr size_type mat_rows() const
Definition: fmatrix.hh:267
Dune::FMatrixHelp::invertMatrix
static K invertMatrix(const FieldMatrix< K, 1, 1 > &matrix, FieldMatrix< K, 1, 1 > &inverse)
invert scalar without changing the original matrix
Definition: fmatrix.hh:528
Dune::Simd::Scalar
typename Overloads::ScalarType< std::decay_t< V > >::type Scalar
Element type of some SIMD type.
Definition: simd/interface.hh:233
densematrix.hh
Implements a matrix constructed from a given type representing a field and a compile-time given numbe...
Dune::FieldMatrix::operator-
friend auto operator-(const FieldMatrix &matrixA, const FieldMatrix< OtherScalar, ROWS, COLS > &matrixB)
vector space subtraction – two-argument version
Definition: fmatrix.hh:142
Dune::DenseMatVecTraits< FieldMatrix< K, ROWS, COLS > >::row_type
FieldVector< K, COLS > row_type
Definition: fmatrix.hh:42
Dune::FieldMatrix::operator=
FieldMatrix & operator=(const FieldMatrix< T, ROWS, COLS > &x)
copy assignment from FieldMatrix over a different field
Definition: fmatrix.hh:116
Dune
Dune namespace.
Definition: alignedallocator.hh:13
Dune::FieldMatrix::const_row_reference
Base::const_row_reference const_row_reference
Definition: fmatrix.hh:86