Loading...
Searching...
No Matches
fmatrix.hh
Go to the documentation of this file.
// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:
// SPDX-FileCopyrightInfo: Copyright © DUNE Project contributors, see file LICENSE.md in module root
// SPDX-License-Identifier: LicenseRef-GPL-2.0-only-with-DUNE-exception
#ifndef DUNE_FMATRIX_HH
#define DUNE_FMATRIX_HH
#include <cmath>
#include <cstddef>
#include <iostream>
#include <algorithm>
#include <initializer_list>
#include <dune/common/boundschecking.hh>
#include <dune/common/exceptions.hh>
#include <dune/common/fvector.hh>
#include <dune/common/densematrix.hh>
#include <dune/common/precision.hh>
#include <dune/common/promotiontraits.hh>
#include <dune/common/typetraits.hh>
#include <dune/common/matrixconcepts.hh>
namespace Dune
{
namespace Impl
{
template<class M>
class ColumnVectorView
{
public:
using value_type = typename M::value_type;
using size_type = typename M::size_type;
constexpr ColumnVectorView(M& matrix, size_type col) :
matrix_(matrix),
col_(col)
{}
constexpr size_type N () const {
return matrix_.N();
}
template<class M_ = M,
std::enable_if_t<std::is_same_v<M_,M> and not std::is_const_v<M_>, int> = 0>
constexpr value_type& operator[] (size_type row) {
return matrix_[row][col_];
}
constexpr const value_type& operator[] (size_type row) const {
return matrix_[row][col_];
}
protected:
M& matrix_;
const size_type col_;
};
}
template<typename M>
struct FieldTraits< Impl::ColumnVectorView<M> >
{
using field_type = typename FieldTraits<M>::field_type;
using real_type = typename FieldTraits<M>::real_type;
};
template< class K, int ROWS, int COLS = ROWS > class FieldMatrix;
template< class K, int ROWS, int COLS >
struct DenseMatVecTraits< FieldMatrix<K,ROWS,COLS> >
{
typedef FieldMatrix<K,ROWS,COLS> derived_type;
// each row is implemented by a field vector
typedef FieldVector<K,COLS> row_type;
typedef row_type &row_reference;
typedef const row_type &const_row_reference;
typedef std::array<row_type,ROWS> container_type;
typedef K value_type;
typedef typename container_type::size_type size_type;
};
template< class K, int ROWS, int COLS >
struct FieldTraits< FieldMatrix<K,ROWS,COLS> >
{
typedef typename FieldTraits<K>::field_type field_type;
typedef typename FieldTraits<K>::real_type real_type;
};
template<class K, int ROWS, int COLS>
class FieldMatrix : public DenseMatrix< FieldMatrix<K,ROWS,COLS> >
{
template<class,int,int> friend class FieldMatrix;
std::array< FieldVector<K,COLS>, ROWS > _data;
typedef DenseMatrix< FieldMatrix<K,ROWS,COLS> > Base;
public:
constexpr static int rows = ROWS;
constexpr static int cols = COLS;
typedef typename Base::size_type size_type;
typedef typename Base::row_type row_type;
typedef typename Base::row_reference row_reference;
typedef typename Base::const_row_reference const_row_reference;
//===== constructors
constexpr FieldMatrix() = default;
constexpr FieldMatrix(std::initializer_list<Dune::FieldVector<K, cols> > const &l) {
assert(l.size() == rows); // Actually, this is not needed any more!
for(std::size_t i=0; i<std::min<std::size_t>(ROWS, l.size()); ++i)
_data[i] = std::data(l)[i];
}
FieldMatrix(const FieldMatrix&) = default;
template <class T,
typename = std::enable_if_t<HasDenseMatrixAssigner<FieldMatrix, T>::value>>
constexpr FieldMatrix(T const& rhs)
: _data{}
{
*this = rhs;
}
using Base::operator=;
constexpr FieldMatrix& operator=(const FieldMatrix&) = default;
template<typename T>
constexpr FieldMatrix& operator=(const FieldMatrix<T, ROWS, COLS>& x)
{
// The copy must be done element-by-element since a std::array does not have
// a converting assignment operator.
for (std::size_t i = 0; i < _data.size(); ++i)
_data[i] = x._data[i];
return *this;
}
template <typename T, int rows, int cols>
FieldMatrix& operator=(FieldMatrix<T,rows,cols> const&) = delete;
constexpr FieldMatrix<K, COLS, ROWS> transposed() const
{
Dune::FieldMatrix<K, COLS, ROWS> AT;
for( int i = 0; i < ROWS; ++i )
for( int j = 0; j < COLS; ++j )
AT[j][i] = (*this)[i][j];
return AT;
}
template <class OtherScalar>
friend constexpr auto operator+ ( const FieldMatrix& matrixA,
const FieldMatrix<OtherScalar,ROWS,COLS>& matrixB)
{
FieldMatrix<typename PromotionTraits<K,OtherScalar>::PromotedType,ROWS,COLS> result;
for (size_type i = 0; i < ROWS; ++i)
for (size_type j = 0; j < COLS; ++j)
result[i][j] = matrixA[i][j] + matrixB[i][j];
return result;
}
template <class OtherScalar>
friend constexpr auto operator- ( const FieldMatrix& matrixA,
const FieldMatrix<OtherScalar,ROWS,COLS>& matrixB)
{
FieldMatrix<typename PromotionTraits<K,OtherScalar>::PromotedType,ROWS,COLS> result;
for (size_type i = 0; i < ROWS; ++i)
for (size_type j = 0; j < COLS; ++j)
result[i][j] = matrixA[i][j] - matrixB[i][j];
return result;
}
template <class Scalar,
std::enable_if_t<IsNumber<Scalar>::value, int> = 0>
friend constexpr auto operator* ( const FieldMatrix& matrix, Scalar scalar)
{
FieldMatrix<typename PromotionTraits<K,Scalar>::PromotedType,ROWS,COLS> result;
for (size_type i = 0; i < ROWS; ++i)
for (size_type j = 0; j < COLS; ++j)
result[i][j] = matrix[i][j] * scalar;
return result;
}
template <class Scalar,
std::enable_if_t<IsNumber<Scalar>::value, int> = 0>
friend constexpr auto operator* ( Scalar scalar, const FieldMatrix& matrix)
{
FieldMatrix<typename PromotionTraits<K,Scalar>::PromotedType,ROWS,COLS> result;
for (size_type i = 0; i < ROWS; ++i)
for (size_type j = 0; j < COLS; ++j)
result[i][j] = scalar * matrix[i][j];
return result;
}
template <class Scalar,
std::enable_if_t<IsNumber<Scalar>::value, int> = 0>
friend constexpr auto operator/ ( const FieldMatrix& matrix, Scalar scalar)
{
FieldMatrix<typename PromotionTraits<K,Scalar>::PromotedType,ROWS,COLS> result;
for (size_type i = 0; i < ROWS; ++i)
for (size_type j = 0; j < COLS; ++j)
result[i][j] = matrix[i][j] / scalar;
return result;
}
template <class OtherScalar, int otherCols>
friend constexpr auto operator* ( const FieldMatrix& matrixA,
const FieldMatrix<OtherScalar, COLS, otherCols>& matrixB)
{
FieldMatrix<typename PromotionTraits<K,OtherScalar>::PromotedType,ROWS,otherCols> result;
for (size_type i = 0; i < matrixA.mat_rows(); ++i)
for (size_type j = 0; j < matrixB.mat_cols(); ++j)
{
result[i][j] = 0;
for (size_type k = 0; k < matrixA.mat_cols(); ++k)
result[i][j] += matrixA[i][k] * matrixB[k][j];
}
return result;
}
template <class OtherMatrix, std::enable_if_t<
Impl::IsStaticSizeMatrix_v<OtherMatrix>
and not Impl::IsFieldMatrix_v<OtherMatrix>
, int> = 0>
friend constexpr auto operator* ( const FieldMatrix& matrixA,
const OtherMatrix& matrixB)
{
using Field = typename PromotionTraits<K, typename OtherMatrix::field_type>::PromotedType;
Dune::FieldMatrix<Field, rows ,OtherMatrix::cols> result;
for (std::size_t j=0; j<rows; ++j)
matrixB.mtv(matrixA[j], result[j]);
return result;
}
template <class OtherMatrix, std::enable_if_t<
Impl::IsStaticSizeMatrix_v<OtherMatrix>
and not Impl::IsFieldMatrix_v<OtherMatrix>
, int> = 0>
friend constexpr auto operator* ( const OtherMatrix& matrixA,
const FieldMatrix& matrixB)
{
using Field = typename PromotionTraits<K, typename OtherMatrix::field_type>::PromotedType;
Dune::FieldMatrix<Field, OtherMatrix::rows, cols> result;
for (std::size_t j=0; j<cols; ++j)
{
auto B_j = Impl::ColumnVectorView(matrixB, j);
auto result_j = Impl::ColumnVectorView(result, j);
matrixA.mv(B_j, result_j);
}
return result;
}
template<int l>
constexpr FieldMatrix<K,l,cols> leftmultiplyany (const FieldMatrix<K,l,rows>& M) const
{
FieldMatrix<K,l,cols> C;
for (size_type i=0; i<l; i++) {
for (size_type j=0; j<cols; j++) {
C[i][j] = 0;
for (size_type k=0; k<rows; k++)
C[i][j] += M[i][k]*(*this)[k][j];
}
}
return C;
}
using Base::rightmultiply;
template <int r, int c>
constexpr FieldMatrix& rightmultiply (const FieldMatrix<K,r,c>& M)
{
static_assert(r == c, "Cannot rightmultiply with non-square matrix");
static_assert(r == cols, "Size mismatch");
FieldMatrix<K,rows,cols> C(*this);
for (size_type i=0; i<rows; i++)
for (size_type j=0; j<cols; j++) {
(*this)[i][j] = 0;
for (size_type k=0; k<cols; k++)
(*this)[i][j] += C[i][k]*M[k][j];
}
return *this;
}
template<int l>
constexpr FieldMatrix<K,rows,l> rightmultiplyany (const FieldMatrix<K,cols,l>& M) const
{
FieldMatrix<K,rows,l> C;
for (size_type i=0; i<rows; i++) {
for (size_type j=0; j<l; j++) {
C[i][j] = 0;
for (size_type k=0; k<cols; k++)
C[i][j] += (*this)[i][k]*M[k][j];
}
}
return C;
}
// make this thing a matrix
static constexpr size_type mat_rows() { return ROWS; }
static constexpr size_type mat_cols() { return COLS; }
constexpr row_reference mat_access ( size_type i )
{
DUNE_ASSERT_BOUNDS(i < ROWS);
return _data[i];
}
constexpr const_row_reference mat_access ( size_type i ) const
{
DUNE_ASSERT_BOUNDS(i < ROWS);
return _data[i];
}
};
#ifndef DOXYGEN // hide specialization
template<class K>
class FieldMatrix<K,1,1> : public DenseMatrix< FieldMatrix<K,1,1> >
{
FieldVector<K,1> _data;
typedef DenseMatrix< FieldMatrix<K,1,1> > Base;
public:
// standard constructor and everything is sufficient ...
//===== type definitions and constants
typedef typename Base::size_type size_type;
constexpr static int blocklevel = 1;
typedef typename Base::row_type row_type;
typedef typename Base::row_reference row_reference;
typedef typename Base::const_row_reference const_row_reference;
constexpr static int rows = 1;
constexpr static int cols = 1;
//===== constructors
constexpr FieldMatrix() = default;
FieldMatrix(std::initializer_list<Dune::FieldVector<K, 1>> const &l)
{
std::copy_n(l.begin(), std::min<std::size_t>(1, l.size()), &_data);
}
template <class T,
typename = std::enable_if_t<HasDenseMatrixAssigner<FieldMatrix, T>::value>>
constexpr FieldMatrix(T const& rhs)
{
*this = rhs;
}
using Base::operator=;
constexpr FieldMatrix<K, 1, 1> transposed() const
{
return *this;
}
template <class OtherScalar>
friend constexpr auto operator+ ( const FieldMatrix& matrixA,
const FieldMatrix<OtherScalar,1,1>& matrixB)
{
return FieldMatrix<typename PromotionTraits<K,OtherScalar>::PromotedType,1,1>{matrixA[0][0] + matrixB[0][0]};
}
template <class Scalar,
std::enable_if_t<IsNumber<Scalar>::value, int> = 0>
friend constexpr auto operator+ ( const FieldMatrix& matrix,
const Scalar& scalar)
{
return FieldMatrix<typename PromotionTraits<K,Scalar>::PromotedType,1,1>{matrix[0][0] + scalar};
}
template <class Scalar,
std::enable_if_t<IsNumber<Scalar>::value, int> = 0>
friend constexpr auto operator+ ( const Scalar& scalar,
const FieldMatrix& matrix)
{
return FieldMatrix<typename PromotionTraits<Scalar,K>::PromotedType,1,1>{scalar + matrix[0][0]};
}
template <class OtherScalar>
friend constexpr auto operator- ( const FieldMatrix& matrixA,
const FieldMatrix<OtherScalar,1,1>& matrixB)
{
return FieldMatrix<typename PromotionTraits<K,OtherScalar>::PromotedType,1,1>{matrixA[0][0] - matrixB[0][0]};
}
template <class Scalar,
std::enable_if_t<IsNumber<Scalar>::value, int> = 0>
friend constexpr auto operator- ( const FieldMatrix& matrix,
const Scalar& scalar)
{
return FieldMatrix<typename PromotionTraits<K,Scalar>::PromotedType,1,1>{matrix[0][0] - scalar};
}
template <class Scalar,
std::enable_if_t<IsNumber<Scalar>::value, int> = 0>
friend constexpr auto operator- ( const Scalar& scalar,
const FieldMatrix& matrix)
{
return FieldMatrix<typename PromotionTraits<Scalar,K>::PromotedType,1,1>{scalar - matrix[0][0]};
}
template <class Scalar,
std::enable_if_t<IsNumber<Scalar>::value, int> = 0>
friend constexpr auto operator* ( const FieldMatrix& matrix, Scalar scalar)
{
return FieldMatrix<typename PromotionTraits<K,Scalar>::PromotedType,1,1> {matrix[0][0] * scalar};
}
template <class Scalar,
std::enable_if_t<IsNumber<Scalar>::value, int> = 0>
friend constexpr auto operator* ( Scalar scalar, const FieldMatrix& matrix)
{
return FieldMatrix<typename PromotionTraits<K,Scalar>::PromotedType,1,1> {scalar * matrix[0][0]};
}
template <class Scalar,
std::enable_if_t<IsNumber<Scalar>::value, int> = 0>
friend constexpr auto operator/ ( const FieldMatrix& matrix, Scalar scalar)
{
return FieldMatrix<typename PromotionTraits<K,Scalar>::PromotedType,1,1> {matrix[0][0] / scalar};
}
//===== solve
template <class OtherScalar, int otherCols>
friend constexpr auto operator* ( const FieldMatrix& matrixA,
const FieldMatrix<OtherScalar, 1, otherCols>& matrixB)
{
FieldMatrix<typename PromotionTraits<K,OtherScalar>::PromotedType,1,otherCols> result;
for (size_type j = 0; j < matrixB.mat_cols(); ++j)
result[0][j] = matrixA[0][0] * matrixB[0][j];
return result;
}
template <class OtherMatrix, std::enable_if_t<
Impl::IsStaticSizeMatrix_v<OtherMatrix>
and not Impl::IsFieldMatrix_v<OtherMatrix>
and (OtherMatrix::rows==1)
, int> = 0>
friend constexpr auto operator* ( const FieldMatrix& matrixA,
const OtherMatrix& matrixB)
{
using Field = typename PromotionTraits<K, typename OtherMatrix::field_type>::PromotedType;
Dune::FieldMatrix<Field, rows ,OtherMatrix::cols> result;
for (std::size_t j=0; j<rows; ++j)
matrixB.mtv(matrixA[j], result[j]);
return result;
}
template <class OtherMatrix, std::enable_if_t<
Impl::IsStaticSizeMatrix_v<OtherMatrix>
and not Impl::IsFieldMatrix_v<OtherMatrix>
and (OtherMatrix::cols==1)
, int> = 0>
friend constexpr auto operator* ( const OtherMatrix& matrixA,
const FieldMatrix& matrixB)
{
using Field = typename PromotionTraits<K, typename OtherMatrix::field_type>::PromotedType;
Dune::FieldMatrix<Field, OtherMatrix::rows, cols> result;
for (std::size_t j=0; j<cols; ++j)
{
auto B_j = Impl::ColumnVectorView(matrixB, j);
auto result_j = Impl::ColumnVectorView(result, j);
matrixA.mv(B_j, result_j);
}
return result;
}
template<int l>
constexpr FieldMatrix<K,l,1> leftmultiplyany (const FieldMatrix<K,l,1>& M) const
{
FieldMatrix<K,l,1> C;
for (size_type j=0; j<l; j++)
C[j][0] = M[j][0]*(*this)[0][0];
return C;
}
constexpr FieldMatrix& rightmultiply (const FieldMatrix& M)
{
_data[0] *= M[0][0];
return *this;
}
template<int l>
constexpr FieldMatrix<K,1,l> rightmultiplyany (const FieldMatrix<K,1,l>& M) const
{
FieldMatrix<K,1,l> C;
for (size_type j=0; j<l; j++)
C[0][j] = M[0][j]*_data[0];
return C;
}
// make this thing a matrix
static constexpr size_type mat_rows() { return 1; }
static constexpr size_type mat_cols() { return 1; }
constexpr row_reference mat_access ([[maybe_unused]] size_type i)
{
DUNE_ASSERT_BOUNDS(i == 0);
return _data;
}
constexpr const_row_reference mat_access ([[maybe_unused]] size_type i) const
{
DUNE_ASSERT_BOUNDS(i == 0);
return _data;
}
constexpr FieldMatrix& operator+= (const K& k)
{
_data[0] += k;
return (*this);
}
constexpr FieldMatrix& operator-= (const K& k)
{
_data[0] -= k;
return (*this);
}
constexpr FieldMatrix& operator*= (const K& k)
{
_data[0] *= k;
return (*this);
}
constexpr FieldMatrix& operator/= (const K& k)
{
_data[0] /= k;
return (*this);
}
//===== conversion operator
constexpr operator const K& () const { return _data[0]; }
};
template<typename K>
std::ostream& operator<< (std::ostream& s, const FieldMatrix<K,1,1>& a)
{
s << a[0][0];
return s;
}
#endif // DOXYGEN
namespace FMatrixHelp {
template <typename K>
static constexpr K invertMatrix (const FieldMatrix<K,1,1> &matrix, FieldMatrix<K,1,1> &inverse)
{
using real_type = typename FieldTraits<K>::real_type;
inverse[0][0] = real_type(1.0)/matrix[0][0];
return matrix[0][0];
}
template <typename K>
static constexpr K invertMatrix_retTransposed (const FieldMatrix<K,1,1> &matrix, FieldMatrix<K,1,1> &inverse)
{
return invertMatrix(matrix,inverse);
}
template <typename K>
static constexpr K invertMatrix (const FieldMatrix<K,2,2> &matrix, FieldMatrix<K,2,2> &inverse)
{
using real_type = typename FieldTraits<K>::real_type;
// code generated by maple
K det = (matrix[0][0]*matrix[1][1] - matrix[0][1]*matrix[1][0]);
K det_1 = real_type(1.0)/det;
inverse[0][0] = matrix[1][1] * det_1;
inverse[0][1] = - matrix[0][1] * det_1;
inverse[1][0] = - matrix[1][0] * det_1;
inverse[1][1] = matrix[0][0] * det_1;
return det;
}
template <typename K>
static constexpr K invertMatrix_retTransposed (const FieldMatrix<K,2,2> &matrix, FieldMatrix<K,2,2> &inverse)
{
using real_type = typename FieldTraits<K>::real_type;
// code generated by maple
K det = (matrix[0][0]*matrix[1][1] - matrix[0][1]*matrix[1][0]);
K det_1 = real_type(1.0)/det;
inverse[0][0] = matrix[1][1] * det_1;
inverse[1][0] = - matrix[0][1] * det_1;
inverse[0][1] = - matrix[1][0] * det_1;
inverse[1][1] = matrix[0][0] * det_1;
return det;
}
template <typename K>
static constexpr K invertMatrix (const FieldMatrix<K,3,3> &matrix, FieldMatrix<K,3,3> &inverse)
{
using real_type = typename FieldTraits<K>::real_type;
// code generated by maple
K t4 = matrix[0][0] * matrix[1][1];
K t6 = matrix[0][0] * matrix[1][2];
K t8 = matrix[0][1] * matrix[1][0];
K t10 = matrix[0][2] * matrix[1][0];
K t12 = matrix[0][1] * matrix[2][0];
K t14 = matrix[0][2] * matrix[2][0];
K det = (t4*matrix[2][2]-t6*matrix[2][1]-t8*matrix[2][2]+
t10*matrix[2][1]+t12*matrix[1][2]-t14*matrix[1][1]);
K t17 = real_type(1.0)/det;
inverse[0][0] = (matrix[1][1] * matrix[2][2] - matrix[1][2] * matrix[2][1])*t17;
inverse[0][1] = -(matrix[0][1] * matrix[2][2] - matrix[0][2] * matrix[2][1])*t17;
inverse[0][2] = (matrix[0][1] * matrix[1][2] - matrix[0][2] * matrix[1][1])*t17;
inverse[1][0] = -(matrix[1][0] * matrix[2][2] - matrix[1][2] * matrix[2][0])*t17;
inverse[1][1] = (matrix[0][0] * matrix[2][2] - t14) * t17;
inverse[1][2] = -(t6-t10) * t17;
inverse[2][0] = (matrix[1][0] * matrix[2][1] - matrix[1][1] * matrix[2][0]) * t17;
inverse[2][1] = -(matrix[0][0] * matrix[2][1] - t12) * t17;
inverse[2][2] = (t4-t8) * t17;
return det;
}
template <typename K>
static constexpr K invertMatrix_retTransposed (const FieldMatrix<K,3,3> &matrix, FieldMatrix<K,3,3> &inverse)
{
using real_type = typename FieldTraits<K>::real_type;
// code generated by maple
K t4 = matrix[0][0] * matrix[1][1];
K t6 = matrix[0][0] * matrix[1][2];
K t8 = matrix[0][1] * matrix[1][0];
K t10 = matrix[0][2] * matrix[1][0];
K t12 = matrix[0][1] * matrix[2][0];
K t14 = matrix[0][2] * matrix[2][0];
K det = (t4*matrix[2][2]-t6*matrix[2][1]-t8*matrix[2][2]+
t10*matrix[2][1]+t12*matrix[1][2]-t14*matrix[1][1]);
K t17 = real_type(1.0)/det;
inverse[0][0] = (matrix[1][1] * matrix[2][2] - matrix[1][2] * matrix[2][1])*t17;
inverse[1][0] = -(matrix[0][1] * matrix[2][2] - matrix[0][2] * matrix[2][1])*t17;
inverse[2][0] = (matrix[0][1] * matrix[1][2] - matrix[0][2] * matrix[1][1])*t17;
inverse[0][1] = -(matrix[1][0] * matrix[2][2] - matrix[1][2] * matrix[2][0])*t17;
inverse[1][1] = (matrix[0][0] * matrix[2][2] - t14) * t17;
inverse[2][1] = -(t6-t10) * t17;
inverse[0][2] = (matrix[1][0] * matrix[2][1] - matrix[1][1] * matrix[2][0]) * t17;
inverse[1][2] = -(matrix[0][0] * matrix[2][1] - t12) * t17;
inverse[2][2] = (t4-t8) * t17;
return det;
}
template< class K, int m, int n, int p >
static constexpr void multMatrix ( const FieldMatrix< K, m, n > &A,
const FieldMatrix< K, n, p > &B,
FieldMatrix< K, m, p > &ret )
{
typedef typename FieldMatrix< K, m, p > :: size_type size_type;
for( size_type i = 0; i < m; ++i )
{
for( size_type j = 0; j < p; ++j )
{
ret[ i ][ j ] = K( 0 );
for( size_type k = 0; k < n; ++k )
ret[ i ][ j ] += A[ i ][ k ] * B[ k ][ j ];
}
}
}
template <typename K, int rows, int cols>
static constexpr void multTransposedMatrix(const FieldMatrix<K,rows,cols> &matrix, FieldMatrix<K,cols,cols>& ret)
{
typedef typename FieldMatrix<K,rows,cols>::size_type size_type;
for(size_type i=0; i<cols; i++)
for(size_type j=0; j<cols; j++)
{
ret[i][j]=0.0;
for(size_type k=0; k<rows; k++)
ret[i][j]+=matrix[k][i]*matrix[k][j];
}
}
using Dune::DenseMatrixHelp::multAssign;
template <typename K, int rows, int cols>
static constexpr void multAssignTransposed( const FieldMatrix<K,rows,cols> &matrix, const FieldVector<K,rows> & x, FieldVector<K,cols> & ret)
{
typedef typename FieldMatrix<K,rows,cols>::size_type size_type;
for(size_type i=0; i<cols; ++i)
{
ret[i] = 0.0;
for(size_type j=0; j<rows; ++j)
ret[i] += matrix[j][i]*x[j];
}
}
template <typename K, int rows, int cols>
static constexpr FieldVector<K,rows> mult(const FieldMatrix<K,rows,cols> &matrix, const FieldVector<K,cols> & x)
{
FieldVector<K,rows> ret;
multAssign(matrix,x,ret);
return ret;
}
template <typename K, int rows, int cols>
static constexpr FieldVector<K,cols> multTransposed(const FieldMatrix<K,rows,cols> &matrix, const FieldVector<K,rows> & x)
{
FieldVector<K,cols> ret;
multAssignTransposed( matrix, x, ret );
return ret;
}
} // end namespace FMatrixHelp
} // end namespace
#include "fmatrixev.hh"
#endif
typetraits.hh
Traits for type conversions and type information.
precision.hh
Various precision settings for calculations with FieldMatrix and FieldVector.
exceptions.hh
A few common exception classes.
densematrix.hh
Implements a matrix constructed from a given type representing a field and a compile-time given numbe...
boundschecking.hh
Macro for wrapping boundary checks.
fmatrixev.hh
Eigenvalue computations for the FieldMatrix class.
matrixconcepts.hh
fvector.hh
Implements a vector constructed from a given type representing a field and a compile-time given size.
promotiontraits.hh
Compute type of the result of an arithmetic operation involving two different number types.
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:30
Dune::operator<<
std::ostream & operator<<(std::ostream &s, const bigunsignedint< k > &x)
Definition bigunsignedint.hh:301
Dune::operator*
constexpr decltype(auto) operator*() const
Dereferencing operator.
Definition iteratorfacades.hh:1119
Dune::Simd::Scalar
typename Overloads::ScalarType< std::decay_t< V > >::type Scalar
Element type of some SIMD type.
Definition simd/interface.hh:235
Dune
Dune namespace
Definition alignedallocator.hh:13
Dune::DenseMatrixHelp::multAssign
static void multAssign(const DenseMatrix< MAT > &matrix, const DenseVector< V1 > &x, DenseVector< V2 > &ret)
calculates ret = matrix * x
Definition densematrix.hh:1174
Dune::FMatrixHelp
Definition fmatrix.hh:675
Dune::FMatrixHelp::multAssign
static void multAssign(const DenseMatrix< MAT > &matrix, const DenseVector< V1 > &x, DenseVector< V2 > &ret)
calculates ret = matrix * x
Definition densematrix.hh:1174
Dune::FMatrixHelp::multAssignTransposed
static constexpr 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:823
Dune::FMatrixHelp::multTransposed
static constexpr FieldVector< K, cols > multTransposed(const FieldMatrix< K, rows, cols > &matrix, const FieldVector< K, rows > &x)
calculates ret = matrix^T * x
Definition fmatrix.hh:846
Dune::FMatrixHelp::invertMatrix_retTransposed
static constexpr K invertMatrix_retTransposed(const FieldMatrix< K, 1, 1 > &matrix, FieldMatrix< K, 1, 1 > &inverse)
invert scalar without changing the original matrix
Definition fmatrix.hh:688
Dune::FMatrixHelp::multTransposedMatrix
static constexpr void multTransposedMatrix(const FieldMatrix< K, rows, cols > &matrix, FieldMatrix< K, cols, cols > &ret)
calculates ret= A_t*A
Definition fmatrix.hh:806
Dune::FMatrixHelp::multMatrix
static constexpr 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:787
Dune::FMatrixHelp::invertMatrix
static constexpr K invertMatrix(const FieldMatrix< K, 1, 1 > &matrix, FieldMatrix< K, 1, 1 > &inverse)
invert scalar without changing the original matrix
Definition fmatrix.hh:679
Dune::FMatrixHelp::mult
static constexpr FieldVector< K, rows > mult(const FieldMatrix< K, rows, cols > &matrix, const FieldVector< K, cols > &x)
calculates ret = matrix * x
Definition fmatrix.hh:837
Dune::DenseMatrix< FieldMatrix< K, ROWS, COLS > >::operator+=
constexpr derived_type & operator+=(const DenseMatrix< Other > &x)
Definition densematrix.hh:294
Dune::DenseMatrix< FieldMatrix< K, ROWS, COLS > >::operator*=
constexpr derived_type & operator*=(const field_type &k)
Definition densematrix.hh:326
Dune::DenseMatrix< FieldMatrix< K, ROWS, COLS > >::operator-=
constexpr derived_type & operator-=(const DenseMatrix< Other > &x)
Definition densematrix.hh:317
Dune::FieldMatrix< K, dim, dim >::M
constexpr size_type M() const
Definition densematrix.hh:708
Dune::DenseMatrix< FieldMatrix< K, ROWS, COLS > >::rightmultiply
FieldMatrix< K, ROWS, COLS > & rightmultiply(const DenseMatrix< M2 > &M)
Definition densematrix.hh:650
Dune::DenseMatrix< FieldMatrix< K, ROWS, COLS > >::operator/=
constexpr derived_type & operator/=(const field_type &k)
Definition densematrix.hh:334
Dune::DenseMatrix< FieldMatrix< K, ROWS, COLS > >::operator-
constexpr derived_type operator-() const
Definition densematrix.hh:303
Dune::DenseMatrix::mtv
constexpr void mtv(const X &x, Y &y) const
y = A^T x
Definition densematrix.hh:392
Dune::DenseMatrix< FieldMatrix< K, ROWS, COLS > >::blocklevel
static constexpr int blocklevel
Definition densematrix.hh:183
Dune::DenseMatrix< FieldMatrix< K, ROWS, COLS > >::row_type
Traits::row_type row_type
Definition densematrix.hh:174
Dune::DenseMatrix< FieldMatrix< K, ROWS, COLS > >::size_type
Traits::size_type size_type
Definition densematrix.hh:171
Dune::DenseMatrix< FieldMatrix< K, ROWS, COLS > >::const_row_reference
Traits::const_row_reference const_row_reference
Definition densematrix.hh:180
Dune::DenseMatrix< FieldMatrix< K, ROWS, COLS > >::row_reference
Traits::row_reference row_reference
Definition densematrix.hh:177
Dune::DenseMatrix< FieldMatrix< K, ROWS, COLS > >::DenseMatrix
friend class DenseMatrix
Definition densematrix.hh:153
Dune::FieldMatrix
A dense n x m matrix.
Definition fmatrix.hh:117
Dune::FieldMatrix::mat_cols
static constexpr size_type mat_cols()
Definition fmatrix.hh:375
Dune::FieldMatrix::FieldMatrix
constexpr FieldMatrix()=default
Default constructor.
Dune::FieldMatrix::rightmultiplyany
constexpr 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:359
Dune::FieldMatrix::const_row_reference
Base::const_row_reference const_row_reference
Definition fmatrix.hh:132
Dune::FieldMatrix::FieldMatrix
friend class FieldMatrix
Definition fmatrix.hh:118
Dune::FieldMatrix::rightmultiply
constexpr FieldMatrix & rightmultiply(const FieldMatrix< K, r, c > &M)
Multiplies M from the right to this matrix.
Definition fmatrix.hh:342
Dune::FieldMatrix::row_type
Base::row_type row_type
Definition fmatrix.hh:129
Dune::FieldMatrix::size_type
Base::size_type size_type
Definition fmatrix.hh:128
Dune::FieldMatrix::mat_access
constexpr const_row_reference mat_access(size_type i) const
Definition fmatrix.hh:383
Dune::FieldMatrix::operator=
constexpr FieldMatrix & operator=(const FieldMatrix< T, ROWS, COLS > &x)
copy assignment from FieldMatrix over a different field
Definition fmatrix.hh:166
Dune::FieldMatrix::FieldMatrix
constexpr FieldMatrix(T const &rhs)
copy constructor from assignable type T
Definition fmatrix.hh:153
Dune::FieldMatrix::operator=
FieldMatrix & operator=(FieldMatrix< T, rows, cols > const &)=delete
no copy assignment from FieldMatrix of different size
Dune::FieldMatrix::row_reference
Base::row_reference row_reference
Definition fmatrix.hh:131
Dune::FieldMatrix::FieldMatrix
constexpr FieldMatrix(std::initializer_list< Dune::FieldVector< K, cols > > const &l)
Constructor initializing the matrix from a list of vector.
Definition fmatrix.hh:141
Dune::FieldMatrix< K, dim, dim >::rows
static constexpr int rows
Definition fmatrix.hh:124
Dune::FieldMatrix::mat_rows
static constexpr size_type mat_rows()
Definition fmatrix.hh:374
Dune::FieldMatrix::operator=
constexpr FieldMatrix & operator=(const FieldMatrix &)=default
copy assignment operator
Dune::FieldMatrix< K, dim, dim >::cols
static constexpr int cols
Definition fmatrix.hh:126
Dune::FieldMatrix::transposed
constexpr FieldMatrix< K, COLS, ROWS > transposed() const
Return transposed of the matrix as FieldMatrix.
Definition fmatrix.hh:180
Dune::FieldMatrix::operator/
friend constexpr auto operator/(const FieldMatrix &matrix, Scalar scalar)
vector space division by scalar
Definition fmatrix.hh:248
Dune::FieldMatrix::operator+
friend constexpr auto operator+(const FieldMatrix &matrixA, const FieldMatrix< OtherScalar, ROWS, COLS > &matrixB)
vector space addition – two-argument version
Definition fmatrix.hh:191
Dune::FieldMatrix::leftmultiplyany
constexpr 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:324
Dune::FieldMatrix::mat_access
constexpr row_reference mat_access(size_type i)
Definition fmatrix.hh:377
Dune::FieldMatrix::FieldMatrix
FieldMatrix(const FieldMatrix &)=default
copy constructor
Dune::FieldVector
vector space out of a tensor product of fields.
Definition fvector.hh:97
Dune::DenseMatVecTraits< FieldMatrix< K, ROWS, COLS > >::container_type
std::array< row_type, ROWS > container_type
Definition fmatrix.hh:95
Dune::DenseMatVecTraits< FieldMatrix< K, ROWS, COLS > >::derived_type
FieldMatrix< K, ROWS, COLS > derived_type
Definition fmatrix.hh:87
Dune::DenseMatVecTraits< FieldMatrix< K, ROWS, COLS > >::value_type
K value_type
Definition fmatrix.hh:96
Dune::DenseMatVecTraits< FieldMatrix< K, ROWS, COLS > >::const_row_reference
const row_type & const_row_reference
Definition fmatrix.hh:93
Dune::DenseMatVecTraits< FieldMatrix< K, ROWS, COLS > >::row_type
FieldVector< K, COLS > row_type
Definition fmatrix.hh:90
Dune::DenseMatVecTraits< FieldMatrix< K, ROWS, COLS > >::size_type
container_type::size_type size_type
Definition fmatrix.hh:97
Dune::DenseMatVecTraits< FieldMatrix< K, ROWS, COLS > >::row_reference
row_type & row_reference
Definition fmatrix.hh:92
Dune::FieldTraits< FieldMatrix< K, ROWS, COLS > >::field_type
FieldTraits< K >::field_type field_type
Definition fmatrix.hh:103
Dune::FieldTraits< FieldMatrix< K, ROWS, COLS > >::real_type
FieldTraits< K >::real_type real_type
Definition fmatrix.hh:104
Dune::FieldTraits
Definition ftraits.hh:26
Dune::FieldTraits::field_type
T field_type
export the type representing the field
Definition ftraits.hh:28
Dune::FieldTraits::real_type
T real_type
export the type representing the real type of the field
Definition ftraits.hh:30
Dune::DenseMatVecTraits
Definition matvectraits.hh:31
Dune::PromotionTraits::PromotedType
decltype(std::declval< T1 >()+std::declval< T2 >()) PromotedType
Definition promotiontraits.hh:28