Reference documentation for deal.II version 9.5.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
Public Types | Public Member Functions | Static Public Attributes | Protected Attributes | Private Attributes | List of all members
TensorProductMatrixSymmetricSum< dim, Number, n_rows_1d > Class Template Reference

#include <deal.II/lac/tensor_product_matrix.h>

Public Types

using value_type = Number
 

Public Member Functions

 TensorProductMatrixSymmetricSum ()=default
 
template<typename T >
 TensorProductMatrixSymmetricSum (const T &mass_matrix, const T &derivative_matrix)
 
template<typename T >
void reinit (const T &mass_matrix, const T &derivative_matrix)
 
unsigned int m () const
 
unsigned int n () const
 
void vmult (const ArrayView< Number > &dst, const ArrayView< const Number > &src) const
 
void vmult (const ArrayView< Number > &dst, const ArrayView< const Number > &src, AlignedVector< Number > &tmp) const
 
void apply_inverse (const ArrayView< Number > &dst, const ArrayView< const Number > &src) const
 
std::size_t memory_consumption () const
 

Static Public Attributes

static constexpr int n_rows_1d_static = n_rows_1d
 

Protected Attributes

std::array< Table< 2, Number >, dim > mass_matrix
 
std::array< Table< 2, Number >, dim > derivative_matrix
 
std::array< AlignedVector< Number >, dim > eigenvalues
 
std::array< Table< 2, Number >, dim > eigenvectors
 

Private Attributes

AlignedVector< Number > tmp_array
 
Threads::Mutex mutex
 

Detailed Description

template<int dim, typename Number, int n_rows_1d = -1>
class TensorProductMatrixSymmetricSum< dim, Number, n_rows_1d >

This is a special matrix class defined as the tensor product (or Kronecker product) of 1d matrices of the type

\begin{align*} L &= A_1 \otimes M_0 + M_1 \otimes A_0 \end{align*}

in 2d and

\begin{align*} L &= A_2 \otimes M_1 \otimes M_0 + M_2 \otimes A_1 \otimes M_0 + M_2 \otimes M_1 \otimes A_0 \end{align*}

in 3d. The typical application setting is a discretization of the Laplacian \(L\) on a Cartesian (axis-aligned) geometry, where it can be exactly represented by the Kronecker or tensor product of a 1d mass matrix \(M\) and a 1d Laplace matrix \(A\) in each tensor direction (due to symmetry \(M\) and \(A\) are the same in each dimension). The dimension of the resulting class is the product of the one-dimensional matrices.

This class implements two basic operations, namely the usual multiplication by a vector and the inverse. For both operations, fast tensorial techniques can be applied that implement the operator evaluation in \(\text{size}(M)^{d+1}\) arithmetic operations, considerably less than \(\text{size}(M)^{2d}\) for the naive forward transformation and \(\text{size}(M)^{3d}\) for setting up the inverse of \(L\).

Interestingly, the exact inverse of the matrix \(L\) can be found through tensor products due to 1964's work by Lynch et al. [126],

\begin{align*} L^{-1} &= S_1 \otimes S_0 (\Lambda_1 \otimes I + I \otimes \Lambda_0)^{-1} S_1^\mathrm T \otimes S_0^\mathrm T, \end{align*}

where \(S_d\) is the matrix of eigenvectors to the generalized eigenvalue problem in the given tensor direction \(d\):

\begin{align*} A_d s &= \lambda M_d s, d = 0, \quad \ldots,\mathrm{dim}, \end{align*}

and \(\Lambda_d\) is the diagonal matrix representing the generalized eigenvalues \(\lambda\). Note that the vectors \(s\) are such that they simultaneously diagonalize \(A_d\) and \(M_d\), i.e. \(S_d^{\mathrm T} A_d S_d = \Lambda_d\) and \(S_d^{\mathrm T} M_d S_d = I\). This method of matrix inversion is called fast diagonalization method.

This class requires LAPACK support.

Note
This class allows for two modes of usage. The first is a use case with run time constants for the matrix dimensions that is achieved by setting the optional template parameter n_rows_1d to -1. The second mode of usage that is faster allows to set the template parameter as a compile time constant, giving significantly faster code in particular for small sizes of the matrix.
This class can work with scalar types (float, double) and VectorizedArray types.
Template Parameters
dimDimension of the problem. Currently, 1D, 2d, and 3d codes are implemented.
NumberArithmetic type of the underlying array elements. Note that the underlying LAPACK implementation supports only float and double numbers, so only these two types are currently supported by the generic class. Nevertheless, a template specialization for the vectorized types VectorizedArray<float> and VectorizedArray<double> exists. This is necessary to perform LAPACK calculations for each vectorization lane, i.e. for the supported float and double numbers.
n_rows_1dCompile-time number of rows of 1d matrices (only valid if the number of rows and columns coincide for each dimension). By default at -1, which means that the number of rows is determined at run-time by means of the matrices passed to the reinit() function.

Definition at line 113 of file tensor_product_matrix.h.

Member Typedef Documentation

◆ value_type

template<int dim, typename Number , int n_rows_1d = -1>
using TensorProductMatrixSymmetricSum< dim, Number, n_rows_1d >::value_type = Number

Type of matrix entries. This alias is analogous to value_type in the standard library containers.

Definition at line 120 of file tensor_product_matrix.h.

Constructor & Destructor Documentation

◆ TensorProductMatrixSymmetricSum() [1/2]

template<int dim, typename Number , int n_rows_1d = -1>
TensorProductMatrixSymmetricSum< dim, Number, n_rows_1d >::TensorProductMatrixSymmetricSum ( )
default

Default constructor.

◆ TensorProductMatrixSymmetricSum() [2/2]

template<int dim, typename Number , int n_rows_1d = -1>
template<typename T >
TensorProductMatrixSymmetricSum< dim, Number, n_rows_1d >::TensorProductMatrixSymmetricSum ( const T &  mass_matrix,
const T &  derivative_matrix 
)

Constructor that is equivalent to the empty constructor and immediately calling reinit(mass_matrix, derivative_matrix).

Member Function Documentation

◆ reinit()

template<int dim, typename Number , int n_rows_1d = -1>
template<typename T >
void TensorProductMatrixSymmetricSum< dim, Number, n_rows_1d >::reinit ( const T &  mass_matrix,
const T &  derivative_matrix 
)

Initializes the tensor product matrix by copying the arrays of 1d mass matrices mass_matrix and 1d derivative matrices derivative_matrix into its base class counterparts, respectively, and by assembling the regarding generalized eigenvalues and eigenvectors in eigenvalues and eigenvectors, respectively. Note that the current implementation requires each \(M_{d}\) to be symmetric and positive definite and every \(A_{d}\) to be symmetric and invertible but not necessarily positive definite. Columns and rows filled with zero are ignored.

Warning
This class accepts the following types: "std::array<Table<2, Number>, dim>", "std::array<FullMatrix<Number>, dim>", and "Table<2, Number>". In the latter case, we consider the same 1d mass matrix mass_matrix and the same 1d derivative matrix derivative_matrix for each tensor direction.

◆ m()

template<int dim, typename Number , int n_rows_1d = -1>
unsigned int TensorProductMatrixSymmetricSum< dim, Number, n_rows_1d >::m ( ) const

Return the number of rows of the tensor product matrix resulting from the Kronecker product of 1d matrices, which is described in the main documentation of TensorProductMatrixSymmetricSum.

◆ n()

template<int dim, typename Number , int n_rows_1d = -1>
unsigned int TensorProductMatrixSymmetricSum< dim, Number, n_rows_1d >::n ( ) const

Return the number of columns of the tensor product matrix resulting from the Kronecker product of 1d matrices, which is described in the main documentation of TensorProductMatrixSymmetricSum.

◆ vmult() [1/2]

template<int dim, typename Number , int n_rows_1d = -1>
void TensorProductMatrixSymmetricSum< dim, Number, n_rows_1d >::vmult ( const ArrayView< Number > &  dst,
const ArrayView< const Number > &  src 
) const

Implements a matrix-vector product with the underlying matrix as described in the main documentation of TensorProductMatrixSymmetricSum. This function is operating on ArrayView to allow checks of array bounds with respect to dst and src.

Warning
This function works on an internal temporal array, leading to increased memory consumption if many instances of this class are created, e.g., a different object on every cell with different underlying coefficients each. Furthermore, only one thread runs this function at a time (ensured internally with a mutex). If these two limitations are an issue, please consider the other version of this function.

◆ vmult() [2/2]

template<int dim, typename Number , int n_rows_1d = -1>
void TensorProductMatrixSymmetricSum< dim, Number, n_rows_1d >::vmult ( const ArrayView< Number > &  dst,
const ArrayView< const Number > &  src,
AlignedVector< Number > &  tmp 
) const

Same as above but letting the user provide a user-owned temporary array, resolving the two issues described above. This array is resized internally to the needed size.

◆ apply_inverse()

template<int dim, typename Number , int n_rows_1d = -1>
void TensorProductMatrixSymmetricSum< dim, Number, n_rows_1d >::apply_inverse ( const ArrayView< Number > &  dst,
const ArrayView< const Number > &  src 
) const

Implements a matrix-vector product with the underlying matrix as described in the main documentation of TensorProductMatrixSymmetricSum. This function is operating on ArrayView to allow checks of array bounds with respect to dst and src.

◆ memory_consumption()

template<int dim, typename Number , int n_rows_1d = -1>
std::size_t TensorProductMatrixSymmetricSum< dim, Number, n_rows_1d >::memory_consumption ( ) const

Return the memory consumption of the allocated memory in this class.

Member Data Documentation

◆ n_rows_1d_static

template<int dim, typename Number , int n_rows_1d = -1>
constexpr int TensorProductMatrixSymmetricSum< dim, Number, n_rows_1d >::n_rows_1d_static = n_rows_1d
staticconstexpr

The static number of rows of the 1d matrices. For more details, see the description of the template parameter n_rows_1d.

Definition at line 126 of file tensor_product_matrix.h.

◆ mass_matrix

template<int dim, typename Number , int n_rows_1d = -1>
std::array<Table<2, Number>, dim> TensorProductMatrixSymmetricSum< dim, Number, n_rows_1d >::mass_matrix
protected

An array containing a mass matrix for each tensor direction.

Definition at line 224 of file tensor_product_matrix.h.

◆ derivative_matrix

template<int dim, typename Number , int n_rows_1d = -1>
std::array<Table<2, Number>, dim> TensorProductMatrixSymmetricSum< dim, Number, n_rows_1d >::derivative_matrix
protected

An array containing a derivative matrix for each tensor direction.

Definition at line 229 of file tensor_product_matrix.h.

◆ eigenvalues

template<int dim, typename Number , int n_rows_1d = -1>
std::array<AlignedVector<Number>, dim> TensorProductMatrixSymmetricSum< dim, Number, n_rows_1d >::eigenvalues
protected

An array storing the generalized eigenvalues for each tensor direction.

Definition at line 235 of file tensor_product_matrix.h.

◆ eigenvectors

template<int dim, typename Number , int n_rows_1d = -1>
std::array<Table<2, Number>, dim> TensorProductMatrixSymmetricSum< dim, Number, n_rows_1d >::eigenvectors
protected

An array storing the generalized eigenvectors for each tensor direction.

Definition at line 241 of file tensor_product_matrix.h.

◆ tmp_array

template<int dim, typename Number , int n_rows_1d = -1>
AlignedVector<Number> TensorProductMatrixSymmetricSum< dim, Number, n_rows_1d >::tmp_array
mutableprivate

An array for temporary data.

Definition at line 247 of file tensor_product_matrix.h.

◆ mutex

template<int dim, typename Number , int n_rows_1d = -1>
Threads::Mutex TensorProductMatrixSymmetricSum< dim, Number, n_rows_1d >::mutex
mutableprivate

A mutex that guards access to the array tmp_array.

Definition at line 252 of file tensor_product_matrix.h.


The documentation for this class was generated from the following file: