Dune Core Modules (2.8.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
19
20namespace 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 {
39 typedef FieldMatrix<K,ROWS,COLS> derived_type;
40
41 // each row is implemented by a field vector
42 typedef FieldVector<K,COLS> row_type;
43
44 typedef row_type &row_reference;
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 {
55 typedef typename FieldTraits<K>::field_type field_type;
56 typedef typename FieldTraits<K>::real_type real_type;
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
85 typedef typename Base::row_reference row_reference;
86 typedef typename Base::const_row_reference const_row_reference;
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
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,
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,
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
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");
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 static constexpr size_type mat_rows() { return ROWS; }
268 static constexpr size_type mat_cols() { return COLS; }
269
270 row_reference mat_access ( size_type i )
271 {
272 DUNE_ASSERT_BOUNDS(i < ROWS);
273 return _data[i];
274 }
275
276 const_row_reference mat_access ( size_type i ) const
277 {
278 DUNE_ASSERT_BOUNDS(i < ROWS);
279 return _data[i];
280 }
281 };
282
283#ifndef DOXYGEN // hide specialization
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;
309 typedef typename Base::const_row_reference const_row_reference;
310
312 enum {
315 rows = 1,
318 cols = 1
319 };
320
321 //===== constructors
324 constexpr 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 static constexpr size_type mat_rows() { return 1; }
464 static constexpr size_type mat_cols() { return 1; }
465
466 row_reference mat_access ([[maybe_unused]] size_type i)
467 {
468 DUNE_ASSERT_BOUNDS(i == 0);
469 return _data;
470 }
471
472 const_row_reference mat_access ([[maybe_unused]] size_type i) const
473 {
474 DUNE_ASSERT_BOUNDS(i == 0);
475 return _data;
476 }
477
479 FieldMatrix& operator+= (const K& k)
480 {
481 _data[0] += k;
482 return (*this);
483 }
484
486 FieldMatrix& operator-= (const K& k)
487 {
488 _data[0] -= k;
489 return (*this);
490 }
491
493 FieldMatrix& operator*= (const K& k)
494 {
495 _data[0] *= k;
496 return (*this);
497 }
498
500 FieldMatrix& operator/= (const K& k)
501 {
502 _data[0] /= k;
503 return (*this);
504 }
505
506 //===== conversion operator
507
508 operator const K& () const { return _data[0]; }
509
510 };
511
513 template<typename K>
514 std::ostream& operator<< (std::ostream& s, const FieldMatrix<K,1,1>& a)
515 {
516 s << a[0][0];
517 return s;
518 }
519
520#endif // DOXYGEN
521
522 namespace FMatrixHelp {
523
525 template <typename K>
526 static inline K invertMatrix (const FieldMatrix<K,1,1> &matrix, FieldMatrix<K,1,1> &inverse)
527 {
528 using real_type = typename FieldTraits<K>::real_type;
529 inverse[0][0] = real_type(1.0)/matrix[0][0];
530 return matrix[0][0];
531 }
532
534 template <typename K>
535 static inline K invertMatrix_retTransposed (const FieldMatrix<K,1,1> &matrix, FieldMatrix<K,1,1> &inverse)
536 {
537 return invertMatrix(matrix,inverse);
538 }
539
540
542 template <typename K>
543 static inline K invertMatrix (const FieldMatrix<K,2,2> &matrix, FieldMatrix<K,2,2> &inverse)
544 {
545 using real_type = typename FieldTraits<K>::real_type;
546 // code generated by maple
547 K det = (matrix[0][0]*matrix[1][1] - matrix[0][1]*matrix[1][0]);
548 K det_1 = real_type(1.0)/det;
549 inverse[0][0] = matrix[1][1] * det_1;
550 inverse[0][1] = - matrix[0][1] * det_1;
551 inverse[1][0] = - matrix[1][0] * det_1;
552 inverse[1][1] = matrix[0][0] * det_1;
553 return det;
554 }
555
558 template <typename K>
559 static inline K invertMatrix_retTransposed (const FieldMatrix<K,2,2> &matrix, FieldMatrix<K,2,2> &inverse)
560 {
561 using real_type = typename FieldTraits<K>::real_type;
562 // code generated by maple
563 K det = (matrix[0][0]*matrix[1][1] - matrix[0][1]*matrix[1][0]);
564 K det_1 = real_type(1.0)/det;
565 inverse[0][0] = matrix[1][1] * det_1;
566 inverse[1][0] = - matrix[0][1] * det_1;
567 inverse[0][1] = - matrix[1][0] * det_1;
568 inverse[1][1] = matrix[0][0] * det_1;
569 return det;
570 }
571
573 template <typename K>
574 static inline K invertMatrix (const FieldMatrix<K,3,3> &matrix, FieldMatrix<K,3,3> &inverse)
575 {
576 using real_type = typename FieldTraits<K>::real_type;
577 // code generated by maple
578 K t4 = matrix[0][0] * matrix[1][1];
579 K t6 = matrix[0][0] * matrix[1][2];
580 K t8 = matrix[0][1] * matrix[1][0];
581 K t10 = matrix[0][2] * matrix[1][0];
582 K t12 = matrix[0][1] * matrix[2][0];
583 K t14 = matrix[0][2] * matrix[2][0];
584
585 K det = (t4*matrix[2][2]-t6*matrix[2][1]-t8*matrix[2][2]+
586 t10*matrix[2][1]+t12*matrix[1][2]-t14*matrix[1][1]);
587 K t17 = real_type(1.0)/det;
588
589 inverse[0][0] = (matrix[1][1] * matrix[2][2] - matrix[1][2] * matrix[2][1])*t17;
590 inverse[0][1] = -(matrix[0][1] * matrix[2][2] - matrix[0][2] * matrix[2][1])*t17;
591 inverse[0][2] = (matrix[0][1] * matrix[1][2] - matrix[0][2] * matrix[1][1])*t17;
592 inverse[1][0] = -(matrix[1][0] * matrix[2][2] - matrix[1][2] * matrix[2][0])*t17;
593 inverse[1][1] = (matrix[0][0] * matrix[2][2] - t14) * t17;
594 inverse[1][2] = -(t6-t10) * t17;
595 inverse[2][0] = (matrix[1][0] * matrix[2][1] - matrix[1][1] * matrix[2][0]) * t17;
596 inverse[2][1] = -(matrix[0][0] * matrix[2][1] - t12) * t17;
597 inverse[2][2] = (t4-t8) * t17;
598
599 return det;
600 }
601
603 template <typename K>
604 static inline K invertMatrix_retTransposed (const FieldMatrix<K,3,3> &matrix, FieldMatrix<K,3,3> &inverse)
605 {
606 using real_type = typename FieldTraits<K>::real_type;
607 // code generated by maple
608 K t4 = matrix[0][0] * matrix[1][1];
609 K t6 = matrix[0][0] * matrix[1][2];
610 K t8 = matrix[0][1] * matrix[1][0];
611 K t10 = matrix[0][2] * matrix[1][0];
612 K t12 = matrix[0][1] * matrix[2][0];
613 K t14 = matrix[0][2] * matrix[2][0];
614
615 K det = (t4*matrix[2][2]-t6*matrix[2][1]-t8*matrix[2][2]+
616 t10*matrix[2][1]+t12*matrix[1][2]-t14*matrix[1][1]);
617 K t17 = real_type(1.0)/det;
618
619 inverse[0][0] = (matrix[1][1] * matrix[2][2] - matrix[1][2] * matrix[2][1])*t17;
620 inverse[1][0] = -(matrix[0][1] * matrix[2][2] - matrix[0][2] * matrix[2][1])*t17;
621 inverse[2][0] = (matrix[0][1] * matrix[1][2] - matrix[0][2] * matrix[1][1])*t17;
622 inverse[0][1] = -(matrix[1][0] * matrix[2][2] - matrix[1][2] * matrix[2][0])*t17;
623 inverse[1][1] = (matrix[0][0] * matrix[2][2] - t14) * t17;
624 inverse[2][1] = -(t6-t10) * t17;
625 inverse[0][2] = (matrix[1][0] * matrix[2][1] - matrix[1][1] * matrix[2][0]) * t17;
626 inverse[1][2] = -(matrix[0][0] * matrix[2][1] - t12) * t17;
627 inverse[2][2] = (t4-t8) * t17;
628
629 return det;
630 }
631
633 template< class K, int m, int n, int p >
634 static inline void multMatrix ( const FieldMatrix< K, m, n > &A,
635 const FieldMatrix< K, n, p > &B,
637 {
638 typedef typename FieldMatrix< K, m, p > :: size_type size_type;
639
640 for( size_type i = 0; i < m; ++i )
641 {
642 for( size_type j = 0; j < p; ++j )
643 {
644 ret[ i ][ j ] = K( 0 );
645 for( size_type k = 0; k < n; ++k )
646 ret[ i ][ j ] += A[ i ][ k ] * B[ k ][ j ];
647 }
648 }
649 }
650
652 template <typename K, int rows, int cols>
654 {
655 typedef typename FieldMatrix<K,rows,cols>::size_type size_type;
656
657 for(size_type i=0; i<cols; i++)
658 for(size_type j=0; j<cols; j++)
659 {
660 ret[i][j]=0.0;
661 for(size_type k=0; k<rows; k++)
662 ret[i][j]+=matrix[k][i]*matrix[k][j];
663 }
664 }
665
666 using Dune::DenseMatrixHelp::multAssign;
667
669 template <typename K, int rows, int cols>
671 {
672 typedef typename FieldMatrix<K,rows,cols>::size_type size_type;
673
674 for(size_type i=0; i<cols; ++i)
675 {
676 ret[i] = 0.0;
677 for(size_type j=0; j<rows; ++j)
678 ret[i] += matrix[j][i]*x[j];
679 }
680 }
681
683 template <typename K, int rows, int cols>
685 {
687 multAssign(matrix,x,ret);
688 return ret;
689 }
690
692 template <typename K, int rows, int cols>
694 {
696 multAssignTransposed( matrix, x, ret );
697 return ret;
698 }
699
700 } // end namespace FMatrixHelp
701
704} // end namespace
705
706#include "fmatrixev.hh"
707#endif
Macro for wrapping boundary checks.
A dense n x m matrix.
Definition: densematrix.hh:165
derived_type operator-() const
Matrix negation.
Definition: densematrix.hh:326
constexpr size_type M() const
number of columns
Definition: densematrix.hh:731
FieldMatrix< K, ROWS, COLS > & rightmultiply(const DenseMatrix< M2 > &M)
Multiplies M from the right to this matrix.
Definition: densematrix.hh:673
derived_type & operator/=(const field_type &k)
vector space division by scalar
Definition: densematrix.hh:357
derived_type & operator*=(const field_type &k)
vector space multiplication with scalar
Definition: densematrix.hh:349
derived_type & operator-=(const DenseMatrix< Other > &x)
vector space subtraction
Definition: densematrix.hh:340
@ blocklevel
The number of block levels we contain. This is 1.
Definition: densematrix.hh:205
Traits::row_type row_type
The type used to represent a row (must fulfill the Dune::DenseVector interface)
Definition: densematrix.hh:194
Traits::size_type size_type
The type used for the index access and size operation.
Definition: densematrix.hh:191
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:200
Traits::row_reference row_reference
The type used to represent a reference to a row (usually row_type &)
Definition: densematrix.hh:197
derived_type & operator+=(const DenseMatrix< Other > &x)
vector space addition
Definition: densematrix.hh:317
A dense n x m matrix.
Definition: fmatrix.hh:69
constexpr FieldMatrix()=default
Default constructor.
FieldMatrix & operator=(const FieldMatrix< T, ROWS, COLS > &x)
copy assignment from FieldMatrix over a different field
Definition: fmatrix.hh:116
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
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
FieldMatrix & rightmultiply(const FieldMatrix< K, r, c > &M)
Multiplies M from the right to this matrix.
Definition: fmatrix.hh:235
@ rows
The number of rows.
Definition: fmatrix.hh:77
@ cols
The number of columns.
Definition: fmatrix.hh:79
friend auto operator*(const FieldMatrix &matrix, Scalar scalar)
vector space multiplication with scalar
Definition: fmatrix.hh:157
FieldMatrix & operator=(FieldMatrix< T, rows, cols > const &)=delete
no copy assignment from FieldMatrix of different size
FieldMatrix(std::initializer_list< Dune::FieldVector< K, cols > > const &l)
Constructor initializing the matrix from a list of vector.
Definition: fmatrix.hh:95
friend auto operator/(const FieldMatrix &matrix, Scalar scalar)
vector space division by scalar
Definition: fmatrix.hh:185
friend auto operator+(const FieldMatrix &matrixA, const FieldMatrix< OtherScalar, ROWS, COLS > &matrixB)
vector space addition – two-argument version
Definition: fmatrix.hh:128
FieldMatrix & operator=(const FieldMatrix &)=default
copy assignment operator
vector space out of a tensor product of fields.
Definition: fvector.hh:95
Implements a matrix constructed from a given type representing a field and a compile-time given numbe...
A few common exception classes.
Traits for type conversions and type information.
static FieldVector< K, cols > multTransposed(const FieldMatrix< K, rows, cols > &matrix, const FieldVector< K, rows > &x)
calculates ret = matrix^T * x
Definition: fmatrix.hh:693
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:535
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:634
static K invertMatrix(const FieldMatrix< K, 1, 1 > &matrix, FieldMatrix< K, 1, 1 > &inverse)
invert scalar without changing the original matrix
Definition: fmatrix.hh:526
static FieldVector< K, rows > mult(const FieldMatrix< K, rows, cols > &matrix, const FieldVector< K, cols > &x)
calculates ret = matrix * x
Definition: fmatrix.hh:684
static void multTransposedMatrix(const FieldMatrix< K, rows, cols > &matrix, FieldMatrix< K, cols, cols > &ret)
calculates ret= A_t*A
Definition: fmatrix.hh:653
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:670
Eigenvalue computations for the FieldMatrix class.
Implements a vector constructed from a given type representing a field and a compile-time given size.
#define DUNE_ASSERT_BOUNDS(cond)
If DUNE_CHECK_BOUNDS is defined: check if condition cond holds; otherwise, do nothing.
Definition: boundschecking.hh:28
auto min(ADLTag< 0 >, const V &v1, const V &v2)
implements binary Simd::min()
Definition: defaults.hh:87
typename Overloads::ScalarType< std::decay_t< V > >::type Scalar
Element type of some SIMD type.
Definition: interface.hh:233
Dune namespace.
Definition: alignedallocator.hh:11
Various precision settings for calculations with FieldMatrix and FieldVector.
Compute type of the result of an arithmetic operation involving two different number types.
Creative Commons License   |  Legal Statements / Impressum  |  Hosted by TU Dresden  |  generated with Hugo v0.111.3 (Nov 21, 23:30, 2024)