KASKADE 7 development version
|
Classes and functions for basic vector and matrix arithmetics. More...
Classes | |
class | Kaskade::LinearProductSpace< Scalar_, Seq > |
class | Kaskade::DynamicMatrix< K > |
A LAPACK-compatible dense matrix class with shape specified at runtime. More... | |
class | Kaskade::LocalMatrix< Entry, diagonal, SortedRowIdx, SortedColIdx > |
Providing a matrix or array interface to LAPACK-ordered entries. More... | |
class | Kaskade::LocalMatrices< Entry, diagonal, SortedRowIdx, SortedColIdx, IT, IndexType > |
A structure for holding a sequence of several local matrices to be filled sequentially and to be scattered jointly. More... | |
class | Kaskade::LocalVector< Entry, SortedIdx > |
Providing a vector interface for contiguously stored entries. More... | |
class | Kaskade::LocalVectors< Entry, SortedIdx, Vector > |
A structure for holding a sequence of several local vectors to be filled sequentially and to be scattered jointly. More... | |
struct | Kaskade::MatrixTraits< Matrix > |
Defines domain and range types for matrix classes. More... | |
class | Kaskade::NumaDenseMatrix< Entry > |
A dense matrix class tailored towards NUMA architectures. More... | |
class | Kaskade::NumaVector< Entry > |
A vector class tailored towards NUMA architectures. This vector distributes its entries in blocks of approximately equal size to the different nodes of NUMA machines in order to exploit the larger accumulated memory bandwidth for fast vector operations. This is most beneficial for large vectors exceeding the cache of the CPUs. Modern processors (2015-01-01) have up to 6MB 3rd level cache, so this is intended for vectors of more than 50,000 doubles, say. More... | |
class | Kaskade::SLAPMatrix< Num, offset > |
Very simple dense matrix class for interfacing with LAPACK routines and the optimization tool uncmin. More... | |
class | Kaskade::DualPairing< X, Xstar > |
Abstract base class for dual pairing of \( X \) and its dual space \( X^* \). More... | |
class | Kaskade::DefaultDualPairing< X, Xstar > |
Default implementation of dual pairing (relies on scalar product operator * being defined) More... | |
class | Kaskade::ZeroDualPairing< X, Xstar > |
Zero dual pairing implementation (mainly useful whenever a dual pairing is needed for formal language reasons, but never used) More... | |
class | Kaskade::SymmetricLinearOperator< X, Xstar > |
Interface for symmetric linear operators. More... | |
class | Kaskade::SymmetricMatrixOperator< X, Matrix > |
A symmetric operator represented by a matrix. More... | |
class | Kaskade::SymmetricLinearOperatorWrapper< X, Xstar > |
Wrapper class to present (hopefully) symmetric linear operators in the SymmetricLinearOperator interface. More... | |
class | Kaskade::SymmetricPreconditioner< X, Xstar > |
Interface for symmetric preconditioners. More... | |
class | Kaskade::IdentityPreconditioner< X > |
The trivial (identity) preconditioner. More... | |
class | Kaskade::SymmetricPreconditionerWrapper< X, Xstar > |
Wrapper class presenting a symmetric preconditioner interface for any preconditioner. More... | |
struct | Kaskade::StaticIndexRange< first, last > |
A compile-time index range. More... | |
class | Kaskade::Tensor< Entry, Size0, Sizes > |
A class for representing tensors of arbitrary static rank and extents. More... | |
class | Kaskade::ThreadedMatrixDetail::CRSChunkPatternInfo< Index > |
A base class representing basic meta information about sparsity patterns of NUMA matrix chunks. More... | |
class | Kaskade::ThreadedMatrixDetail::CRSChunkPatternCreator< Index > |
A class supporting two-stage construction of sparsity patterns of NUMA matrix chunks. More... | |
class | Kaskade::ThreadedMatrixDetail::CRSChunkPattern< Index > |
This class maintains the sparsity structure of a couple of matrix rows (a NUMA matrix chunk). More... | |
class | Kaskade::ThreadedMatrixDetail::CRSChunk< Entry, Index > |
This class stores a couple of compressed row storage rows in memory allocated locally on a NUMA node. More... | |
class | Kaskade::NumaCRSPatternCreator< Index > |
A NUMA-aware creator for matrix sparsity patterns. More... | |
class | Kaskade::NumaCRSPattern< Index > |
A NUMA-aware compressed row storage sparsity pattern. More... | |
class | Kaskade::NumaBCRSMatrix< Entry, Index > |
A NUMA-aware compressed row storage matrix adhering mostly to the Dune ISTL interface (to complete...) More... | |
class | Kaskade::MatrixAsTriplet< Scalar, SparseIndexInt > |
class | Kaskade::Determinant< dim, Scalar, byValue > |
A class for computing determinants and their derivatives. More... | |
Functions | |
template<class CoordType , int dim, std::enable_if_t< std::is_floating_point_v< CoordType >, int > = 0> | |
Dune::FieldVector< CoordType, dim > | Kaskade::projectToUnitSimplex (Dune::FieldVector< CoordType, dim > x) |
Projects a vector onto the unit simplex. More... | |
template<class T , int n, class S , class enable = typename std::enable_if_t<std::is_arithmetic_v<S>>> | |
Dune::FieldVector< T, n > | Dune::operator* (Dune::FieldVector< T, n > x, S s) |
Scalar-vector multiplication \( (s,x) \mapsto sx \). More... | |
template<class T , int n, class S , class enable = typename std::enable_if_t<std::is_arithmetic_v<S>>> | |
Dune::FieldVector< T, n > | Dune::operator* (S s, Dune::FieldVector< T, n > x) |
Scalar-vector multiplication \( (x,s) \mapsto sx \). More... | |
template<class T , int n> | |
Dune::FieldVector< T, n > | Dune::operator+ (Dune::FieldVector< T, n > x, Dune::FieldVector< T, n > const &y) |
Vector addition \( (x,y) \mapsto x+y \). More... | |
template<class T , int n> | |
Dune::FieldVector< T, n > | Dune::operator+ (Dune::FieldVector< T, n > x, T const &y) |
Vector addition \( (x,y) \mapsto x+y \). More... | |
template<class T , int n> | |
Dune::FieldVector< T, n > | Dune::operator- (Dune::FieldVector< T, n > x) |
Vector negation \( x \mapsto - x \). More... | |
template<class T , int n> | |
Dune::FieldVector< T, n > | Dune::max (Dune::FieldVector< T, n > x, Dune::FieldVector< T, n > const &y) |
Componentwise maximum. More... | |
template<class T , int n> | |
Dune::FieldVector< T, n > | Dune::min (Dune::FieldVector< T, n > x, Dune::FieldVector< T, n > const &y) |
Componentwise minimum. More... | |
template<class T , int n, int m> | |
Dune::FieldVector< T, n *m > | Dune::asVector (Dune::FieldMatrix< T, n, m > const &x) |
Converts a matrix of size n x m to a vector of size n*m by concatenating all the columns. More... | |
template<int n, int m, class T > | |
Dune::FieldMatrix< T, n, m > | Dune::asMatrix (Dune::FieldVector< T, n *m > const &x) |
Converts a vector of size nm to a matrix of size n x m by filling columns successively. More... | |
template<class T , int n, int m> | |
Dune::FieldMatrix< T, n, m > | Dune::outerProduct (Dune::FieldVector< T, n > const &x, Dune::FieldVector< T, m > const &y) |
outer vector product \( (x,y) \mapsto x y^T \). More... | |
template<class T , int n> | |
Dune::FieldVector< T, n > | Dune::normalize (Dune::FieldVector< T, n > x) |
returns the normalized vector More... | |
template<class T , int n, int m> | |
Dune::FieldVector< T, n+m > | Dune::vertcat (Dune::FieldVector< T, n > const &x, Dune::FieldVector< T, m > const &y) |
Concatenation of vectors. More... | |
template<class T > | |
T | Dune::vectorProduct (Dune::FieldVector< T, 2 > const &x, Dune::FieldVector< T, 2 > const &y) |
vector product \( (x,y) \mapsto x \times y \in \mathbb{R} \). More... | |
template<class T > | |
Dune::FieldMatrix< T, 1, 2 > | Dune::vectorProductMatrix (Dune::FieldVector< T, 2 > const &x) |
the matrix \( A(x) \in \mathbb{R}^{1\times 2}\) satisfying \( A(x) b = x\times b \) for all \( b \) More... | |
template<class T > | |
Dune::FieldVector< T, 3 > | Dune::vectorProduct (Dune::FieldVector< T, 3 > const &x, Dune::FieldVector< T, 3 > const &y) |
vector product \( (x,y) \mapsto x \times y \). More... | |
template<class T > | |
Dune::FieldMatrix< T, 3, 3 > | Dune::vectorProductMatrix (Dune::FieldVector< T, 3 > const &x) |
the skew-symmetric matrix \( A(x)\in \mathbb{R}^{3\times 3} \) satisfying \( A(x) b = x\times b \) for all \( b \) More... | |
template<class T , class S , int n, int m> | |
Dune::FieldVector< S, n > | Dune::operator* (Dune::FieldMatrix< T, n, m > const &A, Dune::FieldVector< S, m > const &x) |
Matrix-vector multiplication \( (A,x) \mapsto Ax \). More... | |
template<class T , int n, int m> | |
Dune::FieldMatrix< T, n, m > | Dune::operator- (Dune::FieldMatrix< T, n, m > const &A) |
Matrix negation \( A \mapsto - A \). More... | |
template<class Scalar , int n, class ... Rest> | |
auto | Kaskade::horzcat (Dune::FieldVector< Scalar, n > const &a, Rest ... rest) |
Concatenates an arbitrary number of column vectors and matrices with coinciding number of rows to a larger matrix. More... | |
template<class Scalar , int n> | |
Dune::FieldVector< Scalar, n *n > | Kaskade::vectorize (Dune::FieldMatrix< Scalar, n, n > const &A) |
matrix entries as vector, concatenating columns More... | |
template<int n, class Scalar > | |
Dune::FieldMatrix< Scalar, n, n > | Kaskade::unvectorize (Dune::FieldVector< Scalar, n *n > const &v) |
forming quadratic matrix columns of vector segments More... | |
template<class T , int n, int m> | |
T | Kaskade::contraction (Dune::FieldMatrix< T, n, m > const &A, Dune::FieldMatrix< T, n, m > const &B) |
Matrix contraction (Frobenius product) \( (A,B) \mapsto A:B = \sum_{i,j} A_{ij} B_{ij} \). More... | |
template<class T , int n> | |
T | Dune::trace (Dune::FieldMatrix< T, n, n > const &A) |
Matrix trace \( A \mapsto \mathrm{tr}\,A = \sum_{i} A_{ii} \). More... | |
template<class T , int n> | |
T | determinant (Dune::FieldMatrix< T, n, n > const &A) |
Matrix determinant \( A \mapsto \mathrm{det}\,A \). More... | |
template<class T , int n> | |
Dune::FieldVector< T, n > | diag (Dune::FieldMatrix< T, n, n > const &A) |
Matrix diagonal as a vector. More... | |
template<class T , int n> | |
Dune::FieldMatrix< T, n, n > | unitMatrix () |
Returns the identity matrix of size n \( \mapsto I \). More... | |
template<class T > | |
Dune::FieldMatrix< T, 3, 3 > | skewMatrix (Dune::FieldVector< T, 3 > v) |
Returns the skew symmetric matrix related to the cross product with (3d) vector v. More... | |
template<class T , int n> | |
Dune::FieldVector< T, n > | unitVector (int i) |
Returns the unit vector \( e_i \) with \( (e_i)_k = \delta_{ik} \). More... | |
template<class MatrixZ , class Matrix > | |
void | MatMult (MatrixZ &z, Matrix const &x, Matrix const &y) |
Computes Z = X*Y. X and Y need to be compatible, i.e. X.M()==Y.N() has to hold. No aliasing may occur, i.e. &Z != &X and &Z != &Y must hold. Z is reshaped as needed. More... | |
template<class B , class A , class OutIter > | |
OutIter | Kaskade::vectorToSequence (Dune::BlockVector< B, A > const &v, OutIter iter) |
copies the entries of a vector sequentially to the output iterator. More... | |
template<class B , class A , class InIter > | |
InIter | Kaskade::vectorFromSequence (Dune::BlockVector< B, A > &v, InIter in) |
copies the elements obtained from the input iterator sequentially to the entries of the vector. More... | |
void | Kaskade::DynamicMatrixDetail::gemm (bool transposeA, bool transposeB, int m, int n, int k, double alpha, double const *A, int lda, double const *B, int ldb, double beta, double *C, int ldc) |
computes \( C \leftarrow \alpha A B + \beta C \) More... | |
template<class Scalar , int d> | |
Dune::FieldVector< Scalar, d > | Kaskade::eigenvalues (Dune::FieldMatrix< Scalar, d, d > A) |
Computes the eigenvalues of a symmetric matrix \( A \). More... | |
template<class Scalar , int d> | |
Dune::FieldVector< Scalar, d > | Kaskade::eig (Dune::FieldMatrix< Scalar, d, d > &A, bool computeEigenvectors=true) |
Computes the eigenvalues and eigenvectors of a symmetric matrix \( A \). More... | |
template<class Entry , int n> | |
Dune::FieldMatrix< Entry, n, n > | Kaskade::diag (Dune::FieldVector< Entry, n > const &d) |
Creates a fixed-size diagonal matrix from the given diagonal vector. More... | |
template<int n, class Entry , int m, class A > | |
Dune::BlockVector< Dune::FieldVector< Entry, n >, A > | Kaskade::reshapeBlocks (Dune::BlockVector< Dune::FieldVector< Entry, m >, A > const &x) |
Reshapes a nested block vector to differen entry sizes. More... | |
void | Kaskade::MatrixMultiply (SLAPMatrix< double > A, std::vector< double > &in, std::vector< double > &out) |
Matrix multiplication. More... | |
template<class Scalar , int n, class Index = size_t> | |
NumaBCRSMatrix< Dune::FieldMatrix< Scalar, n, n >, Index > | Kaskade::sparseUnitMatrix (Index N) |
creates a unit matrix in NUMA block compressed row storage format More... | |
template<class Scalar , int n, int m, class Index = size_t> | |
NumaBCRSMatrix< Dune::FieldMatrix< Scalar, n, m >, Index > | Kaskade::sparseZeroMatrix (Index N, Index M) |
creates a zero matrix in NUMA block compressed row storage format More... | |
template<class Scalar > | |
void | Kaskade::deleteLowerTriangle (std::vector< int > &ridx, std::vector< int > &cidx, std::vector< Scalar > &data) |
removes subdiagonal entries from the matrix entries stored in elementary triplet format More... | |
void | umfpack_triplet_to_col (int rows, int cols, std::vector< long > const &ridx, std::vector< long > const &cidx, std::vector< double > const &values, std::vector< long > &Ap, std::vector< long > &Ai, std::vector< double > &Az) |
Triplet to compressed column storage conversion . This is just a frontent to DirectType::UMFPACK utility routines. More... | |
void | umfpack_triplet_to_col (int rows, int cols, std::vector< int > const &ridx, std::vector< int > const &cidx, std::vector< double > const &values, std::vector< int > &Ap, std::vector< int > &Ai, std::vector< double > &Az) |
Triplet to compressed column storage conversion . This is just a frontent to DirectType::UMFPACK utility routines. More... | |
void | umfpack_col_to_triplet (std::vector< long > const &Ap, std::vector< long > &Ti) |
Compressed column storage to triplet storage conversion. This is just a frontent to DirectType::UMFPACK utility routines. More... | |
void | umfpack_col_to_triplet (std::vector< int > const &Ap, std::vector< int > &Ti) |
Compressed column storage to triplet storage conversion. This is just a frontent to DirectType::UMFPACK utility routines. More... | |
template<int m, class Scalar , class Sequence > | |
auto const & | component (LinearProductSpace< Scalar, Sequence > const &x) |
Provides access to the m-th component of a product space. More... | |
template<class X0 , class ... Xs> | |
auto | cartesianProduct (X0 const &x0, Xs const &... xs) |
Creates a cartesian product of the given vectors, i.e. a LinearProductSpace object. More... | |
template<class Scalar , class Seq , class OutIter > | |
OutIter | vectorToSequence (LinearProductSpace< Scalar, Seq > const &v, OutIter i) |
writes the coefficients of a vector to a flat scalar sequence <Scalar,Seq> More... | |
template<class Scalar , class Seq , class InIter > | |
InIter | vectorFromSequence (LinearProductSpace< Scalar, Seq > &v, InIter i) |
reads the coefficients of a vector from a flat scalar sequence <Scalar,Seq> More... | |
template<class IndexP , class EntryP , class IndexA , class EntryA > | |
auto | conjugation (NumaBCRSMatrix< EntryP, IndexP > const &P, NumaBCRSMatrix< EntryA, IndexA > const &A, bool onlyLowerTriangle=false, bool createDiagonal=false) |
Creates the conjugation product \( P^T A P\). More... | |
template<class Entry > | |
DynamicMatrix< typename EntryTraits< Entry >::transpose_type > | transpose (DynamicMatrix< Entry > const &A) |
Computes the transpose of a matrix. More... | |
template<class Entry > | |
std::tuple< DynamicMatrix< Dune::FieldMatrix< typename EntryTraits< Entry >::field_type, 1, 1 > >, std::vector< typename EntryTraits< Entry >::real_type >, DynamicMatrix< Dune::FieldMatrix< typename EntryTraits< Entry >::field_type, 1, 1 > > > | svd (DynamicMatrix< Entry > const &A) |
Computes the singular value decomposition \( A = U \Sigma V^T \). More... | |
template<class Entry > | |
EntryTraits< Entry >::real_type | two_norm (DynamicMatrix< Entry > const &A) |
Computes the spectral norm of the given matrix. More... | |
template<class Entry > | |
EntryTraits< Entry >::real_type | condition (DynamicMatrix< Entry > const &A) |
Computes the condition number of the given square matrix. More... | |
template<class MEntry , class Vector > | |
Vector | gesv (DynamicMatrix< MEntry > const &A, Vector const &b) |
Computes the solution of \( Ax = b \) for general invertible matrices. More... | |
template<class MEntry , class Vector > | |
Vector | posv (DynamicMatrix< MEntry > const &A, Vector const &b) |
Computes the solution of \( Ax = b \) if \( A \) is symmetric positive definite. More... | |
template<class Scalar > | |
void | invert (DynamicMatrix< Scalar > &A) |
Computes the inverse \( A^{-1} \). More... | |
template<class Entry > | |
DynamicMatrix< Entry > & | invertSpd (DynamicMatrix< Entry > &A, bool symmetrize=true) |
In place computation of the inverse \( A^{-1} \). More... | |
template<class Entry > | |
DynamicMatrix< typename EntryTraits< Entry >::transpose_type > | pseudoinverse (DynamicMatrix< Entry > A) |
Computes the right pseudoinverse \( A^+ \) of a wide matrix such that \( A A^+ = I \). More... | |
template<class Entry > | |
std::vector< typename EntryTraits< Entry >::real_type > | eigenvalues (DynamicMatrix< Entry > A) |
Computes the eigenvalues of a symmetric matrix \( A \). More... | |
template<class Entry > | |
std::vector< typename EntryTraits< Entry >::real_type > | eig (DynamicMatrix< Entry > &A, bool computeEigenvectors=true) |
Computes the eigenvalues and eigenvectors of a symmetric matrix \( A \). More... | |
template<class Scalar , int d> | |
bool | makePositiveSemiDefinite (Dune::FieldMatrix< Scalar, d, d > &A) |
Computes a positive semi-definite approximation to \( A \) by cutting off negative eigenvalues. More... | |
template<class Entry > | |
bool | makePositiveSemiDefinite (DynamicMatrix< Entry > &A, typename ScalarTraits< typename EntryTraits< Entry >::field_type >::Real threshold=0) |
Computes a positive semi-definite approximation to \( A \) by cutting off negative eigenvalues. More... | |
template<class Entry > | |
std::tuple< std::vector< typename EntryTraits< Entry >::field_type >, DynamicMatrix< typename EntryTraits< Entry >::field_type >, std::vector< int > > | qr (DynamicMatrix< Entry > const &A) |
Computes a rank-revealing QR decomposition of \( A \) using column pivoting. More... | |
template<class MEntry , class VEntry > | |
auto | operator* (DynamicMatrix< MEntry > const &A, Dune::DynamicVector< VEntry > const &x) |
matrix-vector product \( Ax \) More... | |
template<class MEntry , class VEntry > | |
auto | operator* (DynamicMatrix< MEntry > const &A, Dune::BlockVector< VEntry > const &x) |
matrix-vector product \( Ax \) More... | |
template<class MEntry > | |
auto | operator* (typename Dune::FieldTraits< DynamicMatrix< MEntry > >::field_type x, DynamicMatrix< MEntry > A) |
multiplication with scalar More... | |
template<class MEntry > | |
auto | operator- (DynamicMatrix< MEntry > const &A) |
Multiplication with -1. More... | |
template<class MEntry > | |
auto | operator- (DynamicMatrix< MEntry > A, DynamicMatrix< MEntry > const &B) |
matrix subtraction \( A-B \) More... | |
template<class MEntry > | |
auto | operator+ (DynamicMatrix< MEntry > A, DynamicMatrix< MEntry > const &B) |
matrix addition \( A+B \) More... | |
template<class EntryA , class EntryB > | |
void | matrixMatrixProduct (DynamicMatrix< EntryA > const &A, DynamicMatrix< EntryB > const &B, DynamicMatrix< typename DynamicMatrixDetail::ProductTraits< EntryA, EntryB >::type > &C) |
matrix-matrix product \( A B \) More... | |
template<class EntryA , class EntryB > | |
DynamicMatrix< typename DynamicMatrixDetail::ProductTraits< EntryA, EntryB >::type > | operator* (DynamicMatrix< EntryA > const &A, DynamicMatrix< EntryB > const &B) |
matrix-matrix product \( A B \) More... | |
template<class MEntry > | |
DynamicMatrix< MEntry > | submatrix (DynamicMatrix< MEntry > const &A, int rowstart, int rowend, int colstart, int colend) |
extracts contiguous submatrices, copying the data More... | |
template<class EntryA , class EntryB > | |
auto | kron (DynamicMatrix< EntryA > const &A, DynamicMatrix< EntryB > const &B) |
Computes the Kronecker product of two matrices. More... | |
template<class EntryA , class EntryB > | |
auto | horzcat (DynamicMatrix< EntryA > const &A, DynamicMatrix< EntryB > const &B) |
Concatenates two matrices horizontally. More... | |
template<class EntryX , class EntryY > | |
auto | vertcat (Dune::DynamicVector< EntryX > const &x, Dune::DynamicVector< EntryY > const &y) |
Concatenes two (column) vectors vertically. More... | |
template<class Entry > | |
DynamicMatrix< Entry > | vertcat (DynamicMatrix< Entry > const &A, DynamicMatrix< Entry > const &B) |
Concatenates two matrices vertically. More... | |
template<class Entry > | |
auto | dynamicUnitMatrix (Entry const &d, int n) |
Returns the diagonal matrix with the given entry replicated on the diagonal. More... | |
template<class Target , class Source , class RowIndices , class ColIndices > | |
Target | submatrix (Source const &A, RowIndices const &ri, ColIndices const &ci) |
extracts a sparse or dense submatrix of \( A \) More... | |
template<class EntryA , class EntryB , class IndexA , class IndexB > | |
auto | operator* (NumaBCRSMatrix< EntryA, IndexA > const &A, NumaBCRSMatrix< EntryB, IndexB > const &B) |
Computes the matrix-matrix product \( (A,B) \mapsto AB \). More... | |
template<class Entry , int Size0, int ... Sizes> | |
Tensor< Entry, Size0, Sizes... > | operator* (typename EntryTraits< Entry >::field_type a, Tensor< Entry, Size0, Sizes... > A) |
Scalar times tensor multiplication. More... | |
template<class Entry , int n> | |
Entry | trace (Tensor< Entry, n, n > const &A) |
Trace of a square matrix-shaped tensor. More... | |
template<class Entry , class Index > | |
std::ostream & | operator<< (std::ostream &out, NumaBCRSMatrix< Entry, Index > const &A) |
Writes a NumaBCRSMatrix to an output stream. More... | |
template<class SparseMatrix > | |
DynamicMatrix< typename SparseMatrix::block_type > | full (SparseMatrix const &A) |
Converts a sparse NumaBCRSMatrix to a dense matrix. More... | |
template<class SparseMatrix , class RowRange , class ColRange > | |
DynamicMatrix< typename SparseMatrix::block_type > | full (SparseMatrix const &A, RowRange const &rows, ColRange const &cols) |
Converts a subrange of a sparse matrix to a dense matrix. More... | |
template<class Entry , class Index , class Index2 > | |
NumaBCRSMatrix< Entry, Index > | operator+ (NumaBCRSMatrix< Entry, Index > const &A, NumaBCRSMatrix< Entry, Index2 > const &B) |
Matrix addition \( (A,B) \mapsto A+B \). The sparsity patterns of both matrices can be different. The size of the matrices have to be the same. Both must not have symmetric storage. More... | |
template<class Entry , class Index > | |
auto | transpose (NumaBCRSMatrix< Entry, Index > const &A) |
Creates the transposed sparse matrix \( A^T \). More... | |
template<class Target , class Source , class RowIndices , class ColIndices > | |
Target | eraseRowsNCols (Source const &A, RowIndices const &ri, ColIndices const &ci) |
"deletes" rows/columns by extracting a copy of the matrix keeping the non-deleted rows/columns. More... | |
template<class Target , class Source , class RowIndices > | |
Target | eraseRows (Source const &A, RowIndices const &ri) |
"deletes" rows by extracting a copy of the matrix keeping only the non-deleted rows. More... | |
template<class Target , class Source , class ColIndices > | |
Target | eraseCols (Source const &A, ColIndices const &ci) |
"deletes" columns by extracting a copy of the matrix keeping only the non-deleted columns. More... | |
template<class Source > | |
std::vector< size_t > | nonZeroColumns (Source const &A) |
returns all indices of nonzero columns. More... | |
template<int row2, int col2, int row1, int col1, class Scalar , class Index > | |
NumaBCRSMatrix< Dune::FieldMatrix< Scalar, row2, col2 >, Index > | reshapeBlocks (NumaBCRSMatrix< Dune::FieldMatrix< Scalar, row1, col1 >, Index > const &A) |
reshapes NumaBRCSMAtrix entry block structure More... | |
template<int blockrows, int blockcols, class Scalar , class Index > | |
NumaBCRSMatrix< Dune::FieldMatrix< Scalar, blockrows, blockcols >, Index > | horzcat (NumaBCRSMatrix< Dune::FieldMatrix< Scalar, blockrows, blockcols >, Index > const &A, NumaBCRSMatrix< Dune::FieldMatrix< Scalar, blockrows, blockcols >, Index > const &B) |
concatenate two matrices horizontally More... | |
template<int blockrows, int blockcols, class Scalar , class Index > | |
NumaBCRSMatrix< Dune::FieldMatrix< Scalar, blockrows, blockcols >, Index > | vertcat (NumaBCRSMatrix< Dune::FieldMatrix< Scalar, blockrows, blockcols >, Index > const &A, NumaBCRSMatrix< Dune::FieldMatrix< Scalar, blockrows, blockcols >, Index > const &B) |
concatenate two matrices vertically More... | |
template<int blockrows, int blockcols, class Index = size_t> | |
NumaBCRSMatrix< Dune::FieldMatrix< double, blockrows, blockcols >, Index > | diagcat (NumaBCRSMatrix< Dune::FieldMatrix< double, blockrows, blockcols >, Index > const &A, NumaBCRSMatrix< Dune::FieldMatrix< double, blockrows, blockcols >, Index > const &B) |
concatenate two matrices diagonally, resulting matrix is zero on the off block-diagonal More... | |
Variables | |
constexpr StaticIndexRange< 0, std::numeric_limits< int >::max()> | _ |
Convenience static range denoting the full available range, analogous to Matlab's :. More... | |
template<class T > | |
constexpr int | Kaskade::rank = ScalarDetail::Rank<T>::value |
Reports the rank of vector, matrix, and tensor types of static size. More... | |
Classes and functions for basic vector and matrix arithmetics.
Dune::FieldMatrix< T, n, m > Dune::asMatrix | ( | Dune::FieldVector< T, n *m > const & | x | ) |
Converts a vector of size nm to a matrix of size n x m by filling columns successively.
Definition at line 149 of file fixdune.hh.
Dune::FieldVector< T, n *m > Dune::asVector | ( | Dune::FieldMatrix< T, n, m > const & | x | ) |
Converts a matrix of size n x m to a vector of size n*m by concatenating all the columns.
Definition at line 135 of file fixdune.hh.
Referenced by Kaskade::NedelecSimplexShapeFunctionSet< ctype, dimension, T >::interpolate().
|
related |
Creates a cartesian product of the given vectors, i.e. a LinearProductSpace object.
Definition at line 478 of file linearspace.hh.
|
related |
Provides access to the m-th component of a product space.
This simplifies the access to individual components of a product space, e.g., a variable set of finite element functions. Instead of
one can write
m | the index of the component (nonnegative) |
Scalar | a scalar type |
Sequence | a boost::fusion sequence type |
Definition at line 457 of file linearspace.hh.
|
related |
Computes the condition number of the given square matrix.
For invertible matrices, the condition number is given as \( \kappa(A) = \|A\| \, \|A^{-1}\| \). For singular matrices, \( +\infty \) is returned.
Currently, the condition is computed via the signular value decomposition, and takes \( O(n^3) \) time.
Definition at line 844 of file dynamicMatrix.hh.
|
related |
Creates the conjugation product \( P^T A P\).
Note that for typical multigrid Galerkin projections the memory access patterns of the conjugation are quite bad for performance. Consider assembling the projected matrix directly on the coarser discretization.
onlyLowerTriangle. | If true, we will only touch the lower triangle of A and only create the lower triangle of the resulting matrix P^T A P |
createDiagonal. | If true, the diagonal of the result will be structurally nonzero, even if not resulting from the sparsity patterns of P and A. |
Definition at line 123 of file conjugation.hh.
Referenced by MlStack< Grd >::MlStack().
T Kaskade::contraction | ( | Dune::FieldMatrix< T, n, m > const & | A, |
Dune::FieldMatrix< T, n, m > const & | B | ||
) |
Matrix contraction (Frobenius product) \( (A,B) \mapsto A:B = \sum_{i,j} A_{ij} B_{ij} \).
Definition at line 494 of file fixdune.hh.
Referenced by Kaskade::Elastomechanics::MaterialLaws::StVenantKirchhoff< dim, Scalar >::d0(), Kaskade::Elastomechanics::MaterialLaws::StVenantKirchhoff< dim, Scalar >::d1(), and Kaskade::Elastomechanics::MaterialLaws::StVenantKirchhoff< dim, Scalar >::d2().
void Kaskade::deleteLowerTriangle | ( | std::vector< int > & | ridx, |
std::vector< int > & | cidx, | ||
std::vector< Scalar > & | data | ||
) |
removes subdiagonal entries from the matrix entries stored in elementary triplet format
Definition at line 747 of file triplet.hh.
T determinant | ( | Dune::FieldMatrix< T, n, n > const & | A | ) |
Matrix determinant \( A \mapsto \mathrm{det}\,A \).
Definition at line 531 of file fixdune.hh.
Dune::FieldVector< T, n > diag | ( | Dune::FieldMatrix< T, n, n > const & | A | ) |
Dune::FieldMatrix< Entry, n, n > Kaskade::diag | ( | Dune::FieldVector< Entry, n > const & | d | ) |
Creates a fixed-size diagonal matrix from the given diagonal vector.
Definition at line 351 of file dynamicMatrixOps.hh.
|
related |
concatenate two matrices diagonally, resulting matrix is zero on the off block-diagonal
Requires the matrices to have same number of entry types of Dune::FieldMatrix<double,blockrows,blockcols>
A | NumaBRCSMatrix, first submatrix (top left) |
B | NumaBRCSMatrix, second matrix (bottom right) |
Definition at line 3228 of file threadedMatrix.hh.
|
related |
Returns the diagonal matrix with the given entry replicated on the diagonal.
In order to obtain a unit matrix with m by m unit matrix entries on the diagonal, use
Definition at line 334 of file dynamicMatrixOps.hh.
Dune::FieldVector< Scalar, d > Kaskade::eig | ( | Dune::FieldMatrix< Scalar, d, d > & | A, |
bool | computeEigenvectors = true |
||
) |
Computes the eigenvalues and eigenvectors of a symmetric matrix \( A \).
Computes a vector of eigenvalues \( d \) and, optionally, an orthonormal matrix \( V \) of column eigenvectors, such that \( A = V D V^T \) with \( D = \mathrm{diag}(d) \).
A | on entry, a symmetric square matrix, on exit, contains the eigenvectors in columns if computeEigenvectors is true |
computeEigenvectors | if true, compute eigenvectors in addition to eigenvalues |
|
related |
Computes the eigenvalues and eigenvectors of a symmetric matrix \( A \).
A | on entry, a symmetric square matrix, on exit, contains the eigenvectors in columns if computeEigenvectors is true |
computeEigenvectors | if true, compute eigenvectors in addition to eigenvalues |
Dune::FieldVector< Scalar, d > Kaskade::eigenvalues | ( | Dune::FieldMatrix< Scalar, d, d > | A | ) |
Computes the eigenvalues of a symmetric matrix \( A \).
|
related |
Computes the eigenvalues of a symmetric matrix \( A \).
|
related |
"deletes" columns by extracting a copy of the matrix keeping only the non-deleted columns.
Target | the target matrix format , usually Dune::BCRSMatrix. |
Source | the source matrix, usually Dune::BCRSMatrix or NumaBCRSMatrix. |
ColIndices | an STL sequence type |
Definition at line 3030 of file threadedMatrix.hh.
|
related |
"deletes" rows by extracting a copy of the matrix keeping only the non-deleted rows.
Target | the target matrix format , usually Dune::BCRSMatrix. |
Source | the source matrix, usually Dune::BCRSMatrix or NumaBCRSMatrix. |
RowIndices | an STL sequence type |
Definition at line 3010 of file threadedMatrix.hh.
|
related |
"deletes" rows/columns by extracting a copy of the matrix keeping the non-deleted rows/columns.
Target | the target matrix format , usually Dune::BCRSMatrix. |
Source | the source matrix, usually Dune::BCRSMatrix or NumaBCRSMatrix. |
RowIndices | an STL sequence type |
ColIndices | an STL sequence type |
Definition at line 2979 of file threadedMatrix.hh.
|
related |
Converts a sparse NumaBCRSMatrix to a dense matrix.
Be aware that, depending on the number of nonzeros in the sparse matrix, the memory requirement of the dense matrix may be much higher. Use this only for matrices up to 1000x1000, and rather less.
Definition at line 2814 of file threadedMatrix.hh.
|
related |
Converts a subrange of a sparse matrix to a dense matrix.
This is equivalent to Matlab's full(A(rows,cols)) (except that we index 0-based).
SparseMatrix | a sparse matrix type adhering to the Dune::BCRSMatrix interface (e.g., NumaBCRSMatrix) |
rows | an STL range of row indices |
cols | an STL range of column indices |
The row and column index ranges need not be sorted, but cols must not contain duplicates.
Be aware that, depending on the number of nonzeros in the sparse matrix, the memory requirement of the dense matrix may be much higher. Use this only for matrices up to 1000x1000, and rather less.
Definition at line 2841 of file threadedMatrix.hh.
void Kaskade::DynamicMatrixDetail::gemm | ( | bool | transposeA, |
bool | transposeB, | ||
int | m, | ||
int | n, | ||
int | k, | ||
double | alpha, | ||
double const * | A, | ||
int | lda, | ||
double const * | B, | ||
int | ldb, | ||
double | beta, | ||
double * | C, | ||
int | ldc | ||
) |
computes \( C \leftarrow \alpha A B + \beta C \)
\( A\in\R^{m\times k}, B\in\R^{k\times n}, C\in\R^{m\times n} \)
Referenced by Kaskade::DynamicMatrix< K >::matrixMatrixProduct().
|
related |
Computes the solution of \( Ax = b \) for general invertible matrices.
Vector | a vector type, usually Dune::DynamicVector or Dune::BlockVector |
The field types of MEntry and VEntry must be the same.
auto Kaskade::horzcat | ( | Dune::FieldVector< Scalar, n > const & | a, |
Rest ... | rest | ||
) |
Concatenates an arbitrary number of column vectors and matrices with coinciding number of rows to a larger matrix.
This works as Matlab's [x,y,z] construction.
Definition at line 341 of file fixdune.hh.
|
related |
Concatenates two matrices horizontally.
A and B must have the same number of scalar rows.
Definition at line 285 of file dynamicMatrixOps.hh.
Referenced by Kaskade::contactLinesearch().
|
related |
concatenate two matrices horizontally
Requires the matrices to have same number of rows and entry types of Dune::FieldMatrix<Scalar,blockrows,blockcols>
A | NumaBRCSMatrix, first submatrix (left) |
B | NumaBRCSMatrix, second matrix (right) |
Definition at line 3135 of file threadedMatrix.hh.
|
related |
Computes the inverse \( A^{-1} \).
A | the input matrix to be inverted, is overwritten with its inverse |
|
related |
In place computation of the inverse \( A^{-1} \).
Entry | the type of matrix entries, either a scalar or a FieldMatrix<scalar,n,n> |
A | the quadratic input matrix to be inverted (has to be spd), is overwritten with its inverse |
symmetrize | if true, the complete inverse is computed, otherwise only the lower triangular part |
Can throw NonpositiveMatrixException or SingularMatrixException.
|
related |
Computes the Kronecker product of two matrices.
kron(A,B) with \( A\in \mathbf{R}^{n_a\times m_a} \) and \( B\in \mathbf{R}^{n_b\times m_b} \) computes a block matrix \( C \in\mathbf{R}^{n_a n_b \times m_a m_b}\) with \( n_b \times m_b \) blocks such that the \( (i,j)\)-block is \( A_{ij} B \). This is compatible with Matlab's kron function.
The entry types of A and B must be multiplicable.
Definition at line 206 of file dynamicMatrixOps.hh.
|
related |
Computes a positive semi-definite approximation to \( A \) by cutting off negative eigenvalues.
\( A \) has to be symmetric.
|
related |
Computes a positive semi-definite approximation to \( A \) by cutting off negative eigenvalues.
\( A \) has to be symmetric.
void MatMult | ( | MatrixZ & | z, |
Matrix const & | x, | ||
Matrix const & | y | ||
) |
Computes Z = X*Y. X and Y need to be compatible, i.e. X.M()==Y.N() has to hold. No aliasing may occur, i.e. &Z != &X and &Z != &Y must hold. Z is reshaped as needed.
Definition at line 614 of file fixdune.hh.
Referenced by Kaskade::TransferData< Space, CoarseningPolicy >::transferMatrix().
|
related |
matrix-matrix product \( A B \)
Definition at line 123 of file dynamicMatrixOps.hh.
void Kaskade::MatrixMultiply | ( | SLAPMatrix< double > | A, |
std::vector< double > & | in, | ||
std::vector< double > & | out | ||
) |
Matrix multiplication.
out = A*in out is resized properly
Referenced by Kaskade::CUBThetaModelFunction::MatMultiply().
Dune::FieldVector< T, n > Dune::max | ( | Dune::FieldVector< T, n > | x, |
Dune::FieldVector< T, n > const & | y | ||
) |
Componentwise maximum.
Definition at line 110 of file fixdune.hh.
Referenced by boundingBox(), Kaskade::boundingBox(), clapack_dgels(), Kaskade::coarsening(), Bridge::FixedSystemLinearization< Functional >::cols(), Kaskade::contactLinesearch(), PointwiseCorrection< PolyEquation >::correct(), corrTangent(), BisectionPolicy< Cell, BoundingBox >::createChildren(), Kaskade::PenaltyFunction< Base >::d0(), Kaskade::PenaltyFunction< Base >::d1(), Kaskade::PenaltyFunction< Base >::d2(), Kaskade::PenaltyFunction< Base >::d3(), dgels(), dgelsy(), dgesvd(), Kaskade::ContinuousLagrangeMapperImplementation< ScalarType, GV, ShapeFunctionFilter >::dofOnEntity(), dstev(), Kaskade::edgeAveragingFull(), estimateGeometricTruncationError(), Kaskade::ThreadedMatrixDetail::CRSChunkPattern< Index >::exists(), Kaskade::BoundaryInterpolationDetail::getChildVertexDisplacements(), Kaskade::getContactConstraints(), Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::getRhs(), Kaskade::HierarchicExtensionSimplexShapeFunctionSet< ctype, dimension, Scalar >::HierarchicExtensionSimplexShapeFunctionSet(), Kaskade::HierarchicSimplexShapeFunctionSet< ctype, dimension, Scalar >::HierarchicSimplexShapeFunctionSet(), InteriorPointTangentialPredictor< IPF, DomainVector >::initializeCorrector(), Kaskade::SumFunctional< Fu1, Fu2 >::integrationOrder(), Kaskade::makeAuxiliarySpaceMultigridStack(), Dune::max(), Kaskade::CellData< Grd, T >::maxElement(), Kaskade::maxNorm(), Kaskade::maxOrder(), InteriorPointTangentialPredictor< IPF, DomainVector >::muOnSuccess(), Kaskade::StrakosTichyEnergyErrorTerminationCriterion< R >::operator bool(), Kaskade::StrakosTichyPTerminationCriterion< R >::operator bool(), Kaskade::HierarchicalBasisErrorEstimator2< Functional, VariableSetDescription, ExtensionVariableSetDescription, ExtensionSpace, NormFunctional, LinearSolverLA, LinearSolverHA, LinearSolverLM, LinearSolverHM, lumpM, RefinementStrategy >::operator()(), Kaskade::HierarchicalBasisErrorEstimator< Functional, ExtensionVariableSetDescription, ExtensionSpace, NormFunctional, AdjustRHS >::operator()(), Kaskade::ScalarProductView< Function >::order(), Kaskade::BoundaryDisplacementByInterpolation< GridView >::order(), Kaskade::FunctionViews::Difference< Function1, Function2 >::order(), Kaskade::FunctionViews::Gradient< Function, asVector >::order(), Kaskade::FunctionViews::GradientAbsSquare< Function >::order(), Kaskade::Elastomechanics::orientationPreservingStepsize(), Kaskade::ThreadedMatrixDetail::CRSChunkPattern< Index >::position(), RangeEncoder< UInt >::push(), Kaskade::qpLinesearch(), Kaskade::InverseLocalMatrixStorageBase::regularize(), Kaskade::StrakosTichyEnergyErrorTerminationCriterion< R >::relativeError(), Kaskade::StrakosTichyPTerminationCriterion< R >::relativeError(), Bridge::FixedSystemLinearization< Functional >::rows(), Kaskade::VectorialConverterBase< GridView >::solve(), Kaskade::TransferData< Space, CoarseningPolicy >::TransferMatrix::TransferMatrix(), Kaskade::uniformEmbeddedErrorEstimation(), Kaskade::UniformPartitionedMapper< Implementation, Tagger, SFData >::update(), Kaskade::UniformScalarMapper< Implementation, SFData >::update(), InteriorPointTangentialPredictor< IPF, DomainVector >::updateModelOfHomotopy(), Kaskade::CG_Detail::Regularization< X, Xstar, Derived >::updateRegularization(), and Kaskade::writeVTK().
Dune::FieldVector< T, n > Dune::min | ( | Dune::FieldVector< T, n > | x, |
Dune::FieldVector< T, n > const & | y | ||
) |
Componentwise minimum.
Definition at line 122 of file fixdune.hh.
Referenced by Kaskade::NumaCRSPatternCreator< Index >::addDiagonal(), GeometricObject::ImplementationDetail::DistanceImpl< Metric, Scalar, dim, Ball< Scalar, dim >, BoundingBox< Scalar, dim > >::apply(), Kaskade::VariationalFunctionalAssembler< F, SparseIndex, BoundaryDetector, QuadRule >::assemble(), boundingBox(), Kaskade::boundingBox(), checkInside(), clapack_dgels(), Bridge::FixedSystemLinearization< Functional >::cols(), Kaskade::PolynomialModel::computeCoefficients(), Kaskade::OmegaModel::computeCoefficients(), Kaskade::EtaModel::computeCoefficients(), Kaskade::contactLinesearch(), PointwiseCorrection< PolyEquation >::correct(), corrTangent(), dgels(), dgelsy(), dgesvd(), Kaskade::Bridge::Vector< Implementation >::doapplyAsDualTo(), Kaskade::getContactConstraints(), Bridge::IPLinearization< BarrierFunctional, VectorImpl, ImageImpl, QuadRule >::getRHSBlocks(), Kaskade::gridIterate(), Dune::min(), InteriorPointTangentialPredictor< IPF, DomainVector >::muOnSuccess(), Kaskade::HierarchicalBasisErrorEstimator2< Functional, VariableSetDescription, ExtensionVariableSetDescription, ExtensionSpace, NormFunctional, LinearSolverLA, LinearSolverHA, LinearSolverLM, LinearSolverHM, lumpM, RefinementStrategy >::operator()(), Kaskade::HierarchicalBasisErrorEstimator< Functional, ExtensionVariableSetDescription, ExtensionSpace, NormFunctional, AdjustRHS >::operator()(), Kaskade::Tensor< Entry, Size0 >::operator()(), Kaskade::Tensor< Entry, Size0, Sizes >::operator()(), Kaskade::parallelFor(), Kaskade::parallelForNodes(), Kaskade::PatchDomainDecompositionPreconditioner< Space, m, StorageTag, SparseMatrixIndex >::PatchDomainDecompositionPreconditioner(), Kaskade::qpLinesearch(), Bridge::FixedSystemLinearization< Functional >::rows(), Kaskade::PPCGAsNormalSolver< Assembler_, PrecondAssembler, Domain_, Range_, VariableSet, components >::setRelativeAccuracy(), Kaskade::TransferData< Space, CoarseningPolicy >::TransferData(), Kaskade::CG_Detail::Regularization< X, Xstar, Derived >::updateRegularization(), Kaskade::UniformSampler< Function >::writeVTK(), and Kaskade::writeVTK().
|
related |
returns all indices of nonzero columns.
Source | the source matrix, usually Dune::BCRSMatrix or NumaBCRSMatrix. |
Definition at line 3049 of file threadedMatrix.hh.
Dune::FieldVector< T, n > Dune::normalize | ( | Dune::FieldVector< T, n > | x | ) |
returns the normalized vector
Definition at line 178 of file fixdune.hh.
Dune::FieldVector< S, n > Dune::operator* | ( | Dune::FieldMatrix< T, n, m > const & | A, |
Dune::FieldVector< S, m > const & | x | ||
) |
Matrix-vector multiplication \( (A,x) \mapsto Ax \).
Definition at line 264 of file fixdune.hh.
Dune::FieldVector< T, n > Dune::operator* | ( | Dune::FieldVector< T, n > | x, |
S | s | ||
) |
Scalar-vector multiplication \( (s,x) \mapsto sx \).
Definition at line 51 of file fixdune.hh.
|
related |
matrix-matrix product \( A B \)
matrixMatrixProduct(A,B,C)
for performance reasons. Remember, though, that the product takes O(n^3) time, whereas the matrix allocation and copy takes O(n^2) time. For sufficiently large matrices, the extra allocation/copy overhead will be negligible. Definition at line 155 of file dynamicMatrixOps.hh.
|
related |
matrix-vector product \( Ax \)
A.mv(x,y)
, which offers higher performance. Definition at line 51 of file dynamicMatrixOps.hh.
|
related |
matrix-vector product \( Ax \)
A.mv(x,y)
, which offers higher performance. Definition at line 34 of file dynamicMatrixOps.hh.
|
related |
Computes the matrix-matrix product \( (A,B) \mapsto AB \).
EntryA | the type of entries of matrix A. Usually Dune::FieldMatrix<Real,na,ma>. |
EntryB | the type of entries of matrix B. Usually Dune::FieldMatrix<Real,nb,mb>. |
As matrix entries \( A_{ij} B_{jk}\) have to be multiplied, the entry size must be compatible, i.e. ma==nb. As an extension, if entries of A or B are scalar (i.e. na==ma==1), we treat them implicitly as multiples of the unit matrix of compatible size.
Definition at line 37 of file matrixProduct.hh.
Dune::FieldVector< T, n > Dune::operator* | ( | S | s, |
Dune::FieldVector< T, n > | x | ||
) |
Scalar-vector multiplication \( (x,s) \mapsto sx \).
Definition at line 63 of file fixdune.hh.
|
related |
multiplication with scalar
Definition at line 65 of file dynamicMatrixOps.hh.
|
related |
Dune::FieldVector< T, n > Dune::operator+ | ( | Dune::FieldVector< T, n > | x, |
Dune::FieldVector< T, n > const & | y | ||
) |
Vector addition \( (x,y) \mapsto x+y \).
Definition at line 74 of file fixdune.hh.
Dune::FieldVector< T, n > Dune::operator+ | ( | Dune::FieldVector< T, n > | x, |
T const & | y | ||
) |
Vector addition \( (x,y) \mapsto x+y \).
Definition at line 85 of file fixdune.hh.
|
related |
matrix addition \( A+B \)
Definition at line 109 of file dynamicMatrixOps.hh.
|
related |
Matrix addition \( (A,B) \mapsto A+B \). The sparsity patterns of both matrices can be different. The size of the matrices have to be the same. Both must not have symmetric storage.
Definition at line 2885 of file threadedMatrix.hh.
Dune::FieldMatrix< T, n, m > Dune::operator- | ( | Dune::FieldMatrix< T, n, m > const & | A | ) |
Matrix negation \( A \mapsto - A \).
Definition at line 280 of file fixdune.hh.
Dune::FieldVector< T, n > Dune::operator- | ( | Dune::FieldVector< T, n > | x | ) |
Vector negation \( x \mapsto - x \).
Definition at line 98 of file fixdune.hh.
|
related |
matrix subtraction \( A-B \)
Definition at line 93 of file dynamicMatrixOps.hh.
|
related |
Multiplication with -1.
Definition at line 82 of file dynamicMatrixOps.hh.
|
related |
Writes a NumaBCRSMatrix to an output stream.
The matrix is written as text in triplet format, one entry per row.
Definition at line 2795 of file threadedMatrix.hh.
Dune::FieldMatrix< T, n, m > Dune::outerProduct | ( | Dune::FieldVector< T, n > const & | x, |
Dune::FieldVector< T, m > const & | y | ||
) |
outer vector product \( (x,y) \mapsto x y^T \).
Definition at line 164 of file fixdune.hh.
Referenced by Kaskade::BoundaryInterpolationDisplacement< GridView >::derivative(), Kaskade::NedelecSimplexShapeFunction< ctype, dimension, T >::evaluateDerivative(), and Kaskade::Elastomechanics::MaterialLaws::OrthotropicNonLinearMaterial< dim, Scalar >::OrthotropicNonLinearMaterial().
|
related |
Computes the solution of \( Ax = b \) if \( A \) is symmetric positive definite.
Vector | a vector type, usually Dune::DynamicVector or Dune::BlockVector |
A | must be symmetric positive definite |
The field types of MEntry and VEntry must be the same.
Dune::FieldVector< CoordType, dim > Kaskade::projectToUnitSimplex | ( | Dune::FieldVector< CoordType, dim > | x | ) |
Projects a vector onto the unit simplex.
Definition at line 63 of file barycentric.hh.
|
related |
Computes the right pseudoinverse \( A^+ \) of a wide matrix such that \( A A^+ = I \).
Entry | the type of matrix entries of \( A \). This needs to be scalar, i.e. either an arithmetic type or Dune::FieldMatrix<*,1,1>. |
A | the matrix \( A \in\mathbb{R}^{n\times m}\) |
|
related |
Computes a rank-revealing QR decomposition of \( A \) using column pivoting.
This computes a unitary matrix \( Q \) and an upper triangular matrix \( R \) such that \( A P = QR \) with a permutation matrix \( P \). The permutation matrix is given as an index permutation vector \( p \), such that \( P_{k,j} = \delta_{k,p_j} \), i.e. the \( j \)-th column of \( AP \) is the \( p_j\)-th column of \( A \).
Entry | either a scalar arithmetic type or a square field matrix |
Referenced by Kaskade::getContactArea(), and Kaskade::getContactConstraints().
Dune::BlockVector< Dune::FieldVector< Entry, n >, A > Kaskade::reshapeBlocks | ( | Dune::BlockVector< Dune::FieldVector< Entry, m >, A > const & | x | ) |
Reshapes a nested block vector to differen entry sizes.
Use this to convert a vector [[1,2],[3,4]] to [1,2,3,4] or [[1],[2],[3],[4]].
n | the static block size of the nested result vector |
Entry | the type of entries in the nested vector |
m | the static length of the nested source vector |
x | a nested dynamic block vector |
n and x must be compatible, i.e. there must be some integer S such that n*S = m*N, where N is the size of the given block vector.
Definition at line 45 of file reshape.hh.
|
related |
reshapes NumaBRCSMAtrix entry block structure
Converts a NumaBCRSMAtrix with entry-type Dune::FieldMatrix<Scalar,row1,col1> to entry type Dune::FieldMatrix<Scalar,row2,col2>, i.e. the size of the contained blocks are reshaped. Row size row1 must be divisible by row2 and column size col1 must be divisible by col2.
row1 | number of rows in input block |
col1 | number of columns in input block |
row2 | desired number of block-rows |
col2 | desired number of column-rows |
A | NumaBRCSMatrix whose entry-type is to be reshaped |
Definition at line 3087 of file threadedMatrix.hh.
Dune::FieldMatrix< T, 3, 3 > skewMatrix | ( | Dune::FieldVector< T, 3 > | v | ) |
Returns the skew symmetric matrix related to the cross product with (3d) vector v.
Definition at line 568 of file fixdune.hh.
NumaBCRSMatrix< Dune::FieldMatrix< Scalar, n, n >, Index > Kaskade::sparseUnitMatrix | ( | Index | N | ) |
creates a unit matrix in NUMA block compressed row storage format
The entries are fixed-dimensional unit matrices.
Scalar | the field type of entries |
n | dimension of the entries |
N | the size of the matrix to create (N blocks of size n each) |
Definition at line 2918 of file threadedMatrix.hh.
NumaBCRSMatrix< Dune::FieldMatrix< Scalar, n, m >, Index > Kaskade::sparseZeroMatrix | ( | Index | N, |
Index | M | ||
) |
creates a zero matrix in NUMA block compressed row storage format
This is an empty sparsity pattern - all entries are structurally zero.
Scalar | the field type of entries |
n | dimension of the entries |
N,M | the size of the matrix to create (NxM blocks of size n each) |
Definition at line 2942 of file threadedMatrix.hh.
|
related |
extracts contiguous submatrices, copying the data
The submatrix is defined by half-open row and column ranges.
Definition at line 178 of file dynamicMatrixOps.hh.
|
related |
extracts a sparse or dense submatrix of \( A \)
This works as Matlab's A(ri,ci)
submatrix extraction.
Target | the target matrix format, usually Dune::BCRSMatrix, NumaBCRSMatrix, or DynamicMatrix. |
Source | the source matrix, usually Dune::BCRSMatrix or NumaBCRSMatrix. |
RowIndices | an STL sequence type |
ColIndices | an STL sequence type |
For target Dune::BCRSMatrix, the complexity is \( n_r n_c \log n_s \), where \( n_r \) is the size of ri, \( n_c \) the size of ci, and \( n_s \) the number of nonzeros per row in A.
Definition at line 173 of file matrixOps.hh.
|
related |
Computes the singular value decomposition \( A = U \Sigma V^T \).
If \( A \) is a \( n \times m \) matrix, then the SVD gives orthonormal \( U, V \) of size \( n\times n \) and \( m\times m\), respectively, and a \( n\times m\) matrix \( \Sigma \) with only nonnegative diagonal entries.
Definition at line 807 of file dynamicMatrix.hh.
T Dune::trace | ( | Dune::FieldMatrix< T, n, n > const & | A | ) |
Matrix trace \( A \mapsto \mathrm{tr}\,A = \sum_{i} A_{ii} \).
Note: Instead of trace(A*B) use the equivalent, but faster contraction(A,B).
Definition at line 515 of file fixdune.hh.
Referenced by Kaskade::Elastomechanics::MaterialLaws::OrthotropicNonLinearMaterial< dim, Scalar >::d0(), Kaskade::Elastomechanics::MaterialLaws::StVenantKirchhoff< dim, Scalar >::d1(), Kaskade::Elastomechanics::MaterialLaws::OrthotropicNonLinearMaterial< dim, Scalar >::d1(), Kaskade::Elastomechanics::LameNavier< dim_, Scalar_ >::d1(), Kaskade::Elastomechanics::MaterialLaws::StVenantKirchhoff< dim, Scalar >::d2(), Kaskade::Elastomechanics::MaterialLaws::OrthotropicNonLinearMaterial< dim, Scalar >::d2(), Kaskade::Elastomechanics::deviatoricPart(), and Kaskade::Elastomechanics::MaterialLaws::StVenantKirchhoff< dim, Scalar >::setLinearizationPoint().
|
related |
|
related |
Computes the transpose of a matrix.
This is a deep transpose, i.e. if the entries are non-scalar, they are transposed, too.
Definition at line 765 of file dynamicMatrix.hh.
Referenced by Kaskade::DynamicMatrix< K >::usmtv().
|
related |
Creates the transposed sparse matrix \( A^T \).
Definition at line 2956 of file threadedMatrix.hh.
|
related |
Computes the spectral norm of the given matrix.
This computes \( \max_{\|v\|_2=1} \|Av\|_2 \).
Definition at line 824 of file dynamicMatrix.hh.
void umfpack_col_to_triplet | ( | std::vector< int > const & | Ap, |
std::vector< int > & | Ti | ||
) |
Compressed column storage to triplet storage conversion. This is just a frontent to DirectType::UMFPACK utility routines.
void umfpack_col_to_triplet | ( | std::vector< long > const & | Ap, |
std::vector< long > & | Ti | ||
) |
Compressed column storage to triplet storage conversion. This is just a frontent to DirectType::UMFPACK utility routines.
Referenced by Kaskade::MatrixAsTriplet< Scalar, SparseIndexInt >::operator+=(), and Kaskade::MatrixAsTriplet< Scalar, SparseIndexInt >::operator-=().
void umfpack_triplet_to_col | ( | int | rows, |
int | cols, | ||
std::vector< int > const & | ridx, | ||
std::vector< int > const & | cidx, | ||
std::vector< double > const & | values, | ||
std::vector< int > & | Ap, | ||
std::vector< int > & | Ai, | ||
std::vector< double > & | Az | ||
) |
Triplet to compressed column storage conversion . This is just a frontent to DirectType::UMFPACK utility routines.
void umfpack_triplet_to_col | ( | int | rows, |
int | cols, | ||
std::vector< long > const & | ridx, | ||
std::vector< long > const & | cidx, | ||
std::vector< double > const & | values, | ||
std::vector< long > & | Ap, | ||
std::vector< long > & | Ai, | ||
std::vector< double > & | Az | ||
) |
Triplet to compressed column storage conversion . This is just a frontent to DirectType::UMFPACK utility routines.
rows | number of rows | |
cols | number of columns | |
[in] | ridx | row indices of matrix entries |
[in] | cidx | column indices of matrix entries |
[in] | values | matrix entries |
[out] | Ap | column start array |
[out] | Ai | row indices |
[out] | Az | matrix entries |
The output arrays
Referenced by Kaskade::BoomerAMG< Op >::BoomerAMG(), Kaskade::Euclid< Op >::Euclid(), Kaskade::MatrixAsTriplet< Scalar, SparseIndexInt >::operator+=(), and Kaskade::MatrixAsTriplet< Scalar, SparseIndexInt >::operator-=().
Dune::FieldMatrix< T, n, n > unitMatrix | ( | ) |
Returns the identity matrix of size n \( \mapsto I \).
Definition at line 555 of file fixdune.hh.
Dune::FieldVector< T, n > unitVector | ( | int | i | ) |
Returns the unit vector \( e_i \) with \( (e_i)_k = \delta_{ik} \).
Definition at line 585 of file fixdune.hh.
Dune::FieldMatrix< Scalar, n, n > Kaskade::unvectorize | ( | Dune::FieldVector< Scalar, n *n > const & | v | ) |
forming quadratic matrix columns of vector segments
Definition at line 475 of file fixdune.hh.
InIter Kaskade::vectorFromSequence | ( | Dune::BlockVector< B, A > & | v, |
InIter | in | ||
) |
copies the elements obtained from the input iterator sequentially to the entries of the vector.
The value type of the input iterator must be convertible to the scalar entries of the vector. The input iterator is advanced by the number of scalar entries in the vector.
Definition at line 127 of file crsutil.hh.
|
related |
reads the coefficients of a vector from a flat scalar sequence <Scalar,Seq>
Definition at line 505 of file linearspace.hh.
Referenced by Kaskade::LinearProductSpace< Scalar_, Seq >::read().
Dune::FieldVector< Scalar, n *n > Kaskade::vectorize | ( | Dune::FieldMatrix< Scalar, n, n > const & | A | ) |
matrix entries as vector, concatenating columns
Definition at line 459 of file fixdune.hh.
T Dune::vectorProduct | ( | Dune::FieldVector< T, 2 > const & | x, |
Dune::FieldVector< T, 2 > const & | y | ||
) |
vector product \( (x,y) \mapsto x \times y \in \mathbb{R} \).
Definition at line 209 of file fixdune.hh.
Referenced by GeomTools::crossProduct().
Dune::FieldVector< T, 3 > Dune::vectorProduct | ( | Dune::FieldVector< T, 3 > const & | x, |
Dune::FieldVector< T, 3 > const & | y | ||
) |
vector product \( (x,y) \mapsto x \times y \).
Definition at line 231 of file fixdune.hh.
Dune::FieldMatrix< T, 1, 2 > Dune::vectorProductMatrix | ( | Dune::FieldVector< T, 2 > const & | x | ) |
the matrix \( A(x) \in \mathbb{R}^{1\times 2}\) satisfying \( A(x) b = x\times b \) for all \( b \)
Definition at line 219 of file fixdune.hh.
Dune::FieldMatrix< T, 3, 3 > Dune::vectorProductMatrix | ( | Dune::FieldVector< T, 3 > const & | x | ) |
the skew-symmetric matrix \( A(x)\in \mathbb{R}^{3\times 3} \) satisfying \( A(x) b = x\times b \) for all \( b \)
Definition at line 244 of file fixdune.hh.
OutIter Kaskade::vectorToSequence | ( | Dune::BlockVector< B, A > const & | v, |
OutIter | iter | ||
) |
copies the entries of a vector sequentially to the output iterator.
The Vector type can be Dune::BlockVector<B,A>, Dune::FieldVector<K,size>, or a finite element function. The type of the scalar entries of the vector must be convertible to the value type of the output iterator.
Definition at line 76 of file crsutil.hh.
|
related |
writes the coefficients of a vector to a flat scalar sequence <Scalar,Seq>
Definition at line 493 of file linearspace.hh.
Referenced by Kaskade::LinearProductSpace< Scalar_, Seq >::write().
|
related |
Concatenes two (column) vectors vertically.
Definition at line 296 of file dynamicMatrixOps.hh.
Dune::FieldVector< T, n+m > Dune::vertcat | ( | Dune::FieldVector< T, n > const & | x, |
Dune::FieldVector< T, m > const & | y | ||
) |
Concatenation of vectors.
Definition at line 191 of file fixdune.hh.
|
related |
Concatenates two matrices vertically.
A and B must have the same number of columns.
Definition at line 309 of file dynamicMatrixOps.hh.
|
related |
concatenate two matrices vertically
Requires the matrices to have same number of rows and entry types of Dune::FieldMatrix<Scalar,blockrows,blockcols>
A | NumaBRCSMatrix, first submatrix (top) |
B | NumaBRCSMatrix, second matrix (bottom) |
Definition at line 3182 of file threadedMatrix.hh.
|
constexpr |
Convenience static range denoting the full available range, analogous to Matlab's :.
Use this to take complete dimensions when slicing tensors:
Definition at line 74 of file tensor.hh.
Referenced by Kaskade::writeMarc().