[ VIGRA Homepage | Class Index | Function Index | File Index | Main Page ]
![]() |
Matrix functions | ![]() |
---|
Functions | |
template<class T, class C1, class C2, class C3> bool | symmetricEigensystem (MultiArrayView< 2, T, C1 > const &a, MultiArrayView< 2, T, C2 > &ew, MultiArrayView< 2, T, C3 > &ev) |
template<class T, class C1, class C2, class C3> bool | nonsymmetricEigensystem (MultiArrayView< 2, T, C1 > const &a, MultiArrayView< 2, std::complex< T >, C2 > &ew, MultiArrayView< 2, T, C3 > &ev) |
template<class POLYNOMIAL, class VECTOR> bool | polynomialRootsEigenvalueMethod (POLYNOMIAL const &poly, VECTOR &roots, bool polishRoots) |
template<class POLYNOMIAL, class VECTOR> bool | polynomialRealRootsEigenvalueMethod (POLYNOMIAL const &p, VECTOR &roots, bool polishRoots) |
template<class T, class C1, class C2> void | inverse (const MultiArrayView< 2, T, C1 > &v, MultiArrayView< 2, T, C2 > &r) |
template<class T, class C> TemporaryMatrix< T > | inverse (const MultiArrayView< 2, T, C > &v) |
template<class T, class C1, class C2, class C3> void | qrDecomposition (MultiArrayView< 2, T, C1 > const &a, MultiArrayView< 2, T, C2 > &q, MultiArrayView< 2, T, C3 > &r) |
template<class T, class C1, class C2, class C3> void | reverseElimination (const MultiArrayView< 2, T, C1 > &r, const MultiArrayView< 2, T, C2 > &b, MultiArrayView< 2, T, C3 > &x) |
template<class T, class C1, class C2, class C3> bool | linearSolve (const MultiArrayView< 2, T, C1 > &a, const MultiArrayView< 2, T, C2 > &b, MultiArrayView< 2, T, C3 > &res) |
template<class T, class C> unsigned int | rowCount (const MultiArrayView< 2, T, C > &x) |
template<class T, class C> unsigned int | columnCount (const MultiArrayView< 2, T, C > &x) |
template<class T, class C> MultiArrayView< 2, T, C > | rowVector (MultiArrayView< 2, T, C > const &m, int d) |
template<class T, class C> MultiArrayView< 2, T, C > | columnVector (MultiArrayView< 2, T, C > const &m, int d) |
template<class T, class C> bool | isSymmetric (const MultiArrayView< 2, T, C > &v) |
template<class T, class ALLOC> Matrix< T, ALLLOC >::SquaredNormType | squaredNorm (const Matrix< T, ALLLOC > &a) |
template<class T, class ALLOC> Matrix< T, ALLLOC >::NormType | norm (const Matrix< T, ALLLOC > &a) |
template<class T, class C> void | identityMatrix (MultiArrayView< 2, T, C > &r) |
template<class T> TemporaryMatrix< T > | identityMatrix (unsigned int size) |
template<class T, class C1, class C2> void | diagonalMatrix (MultiArrayView< 2, T, C1 > const &v, MultiArrayView< 2, T, C2 > &r) |
template<class T, class C> TemporaryMatrix< T > | diagonalMatrix (MultiArrayView< 2, T, C > const &v) |
template<class T, class C1, class C2> void | transpose (const MultiArrayView< 2, T, C1 > &v, MultiArrayView< 2, T, C2 > &r) |
template<class T, class C> TemporaryMatrix< T > | transpose (MultiArrayView< 2, T, C > const &v) |
template<class T, class C1, class C2, class C3> void | add (const MultiArrayView< 2, T, C1 > &a, const MultiArrayView< 2, T, C2 > &b, MultiArrayView< 2, T, C3 > &r) |
template<class T, class C1, class C2> TemporaryMatrix< T > | operator+ (const MultiArrayView< 2, T, C1 > &a, const MultiArrayView< 2, T, C2 > &b) |
template<class T, class C1, class C2, class C3> void | sub (const MultiArrayView< 2, T, C1 > &a, const MultiArrayView< 2, T, C2 > &b, MultiArrayView< 2, T, C3 > &r) |
template<class T, class C1, class C2> TemporaryMatrix< T > | operator- (const MultiArrayView< 2, T, C1 > &a, const MultiArrayView< 2, T, C2 > &b) |
template<class T, class C> TemporaryMatrix< T > | operator- (const MultiArrayView< 2, T, C > &a) |
template<class T, class C1, class C2> T | dot (const MultiArrayView< 2, T, C1 > &x, const MultiArrayView< 2, T, C2 > &y) |
template<class T, class C1, class C2> T | dot (const MultiArrayView< 1, T, C1 > &x, const MultiArrayView< 1, T, C2 > &y) |
template<class T, class C1, class C2, class C3> void | cross (const MultiArrayView< 1, T, C1 > &x, const MultiArrayView< 1, T, C2 > &y, MultiArrayView< 1, T, C3 > &r) |
template<class T, class C1, class C2, class C3> void | cross (const MultiArrayView< 2, T, C1 > &x, const MultiArrayView< 2, T, C2 > &y, MultiArrayView< 2, T, C3 > &r) |
template<class T, class C1, class C2> TemporaryMatrix< T > | cross (const MultiArrayView< 2, T, C1 > &x, const MultiArrayView< 2, T, C2 > &y) |
template<class T, class C1, class C2, class C3> void | outer (const MultiArrayView< 2, T, C1 > &x, const MultiArrayView< 2, T, C2 > &y, MultiArrayView< 2, T, C3 > &r) |
template<class T, class C1, class C2> TemporaryMatrix< T > | outer (const MultiArrayView< 2, T, C1 > &x, const MultiArrayView< 2, T, C2 > &y) |
template<class T, class C1, class C2> void | smul (const MultiArrayView< 2, T, C1 > &a, T b, MultiArrayView< 2, T, C2 > &r) |
template<class T, class C2, class C3> void | smul (T a, const MultiArrayView< 2, T, C2 > &b, MultiArrayView< 2, T, C3 > &r) |
template<class T, class C1, class C2, class C3> void | mmul (const MultiArrayView< 2, T, C1 > &a, const MultiArrayView< 2, T, C2 > &b, MultiArrayView< 2, T, C3 > &r) |
template<class T, class C1, class C2> TemporaryMatrix< T > | mmul (const MultiArrayView< 2, T, C1 > &a, const MultiArrayView< 2, T, C2 > &b) |
template<class T, class C1, class C2, class C3> void | pmul (const MultiArrayView< 2, T, C1 > &a, const MultiArrayView< 2, T, C2 > &b, MultiArrayView< 2, T, C3 > &r) |
template<class T, class C1, class C2> TemporaryMatrix< T > | pmul (const MultiArrayView< 2, T, C1 > &a, const MultiArrayView< 2, T, C2 > &b) |
template<class T, class C> TemporaryMatrix< T > | operator * (const MultiArrayView< 2, T, C > &a, T b) |
template<class T, class C> TemporaryMatrix< T > | operator * (T a, const MultiArrayView< 2, T, C > &b) |
template<class T, class A, int N, class DATA, class DERIVED> TinyVector< T, N > | operator * (const Matrix< T, A > &a, const TinyVectorBase< T, N, DATA, DERIVED > &b) |
template<class T, int N, class DATA, class DERIVED, class A> TinyVector< T, N > | operator * (const TinyVectorBase< T, N, DATA, DERIVED > &a, const Matrix< T, A > &b) |
template<class T, class C1, class C2> TemporaryMatrix< T > | operator * (const MultiArrayView< 2, T, C1 > &a, const MultiArrayView< 2, T, C2 > &b) |
template<class T, class C1, class C2> void | sdiv (const MultiArrayView< 2, T, C1 > &a, T b, MultiArrayView< 2, T, C2 > &r) |
template<class T, class C1, class C2, class C3> void | pdiv (const MultiArrayView< 2, T, C1 > &a, const MultiArrayView< 2, T, C2 > &b, MultiArrayView< 2, T, C3 > &r) |
template<class T, class C1, class C2> TemporaryMatrix< T > | pdiv (const MultiArrayView< 2, T, C1 > &a, const MultiArrayView< 2, T, C2 > &b) |
template<class T, class C> TemporaryMatrix< T > | operator/ (const MultiArrayView< 2, T, C > &a, T b) |
template<class T, class C> std::ostream & | operator<< (std::ostream &s, const vigra::MultiArrayView< 2, T, C > &m) |
|
add matrices a and b. The result is written into r. All three matrices must have the same shape.
#include "vigra/matrix.hxx" or |
|
Number of columns of a matrix represented as a
#include "vigra/matrix.hxx" or |
|
Create a column vector view for column d of the matrix m
#include "vigra/matrix.hxx" or |
|
calculate the cross product of two matrices representing vectors. That is, x, and y must have a single column of length 3. The result is returned as a temporary matrix.
#include "vigra/matrix.hxx" or |
|
calculate the cross product of two matrices representing vectors. That is, x, y, and r must have a single column of length 3. The result is written into r.
#include "vigra/matrix.hxx" or |
|
calculate the cross product of two vectors of length 3. The result is written into r.
#include "vigra/matrix.hxx" or |
|
create a diagonal matrix from a vector. The vector is given as matrix v, which must either have a single row or column. The result is returned as a temporary matrix. Usage:
vigra::Matrix<double> v(1, len); v = ...; vigra::Matrix<double> m = diagonalMatrix(v);
#include "vigra/matrix.hxx" or |
|
make a diagonal matrix from a vector. The vector is given as matrix v, which must either have a single row or column. The result is witten into the square matrix r.
#include "vigra/matrix.hxx" or |
|
calculate the inner product of two vectors. The vector lenths must match.
#include "vigra/matrix.hxx" or |
|
calculate the inner product of two matrices representing vectors. That is, matrix x must have a single row, and matrix y must have a single column, and the other dimensions must match.
#include "vigra/matrix.hxx" or |
|
create n identity matrix of the given size. Usage:
vigra::Matrix<double> m = vigra::identityMatrix<double>(size);
#include "vigra/matrix.hxx" or |
|
initialize the given square matrix as an identity matrix.
#include "vigra/matrix.hxx" or |
|
create the inverse of square matrix v. The result is returned as a temporary matrix. The inverse is calculated by means of QR decomposition. If v is not invertible,
vigra::Matrix<double> v(n, n); v = ...; vigra::Matrix<double> m = inverse(v);
#include "vigra/linear__solve.hxx" or |
|
invert square matrix v. The result is written into r which must have the same shape. The inverse is calculated by means of QR decomposition. If v is not invertible,
#include "vigra/linear__solve.hxx" or |
|
Check whether matrix m is symmetric.
#include "vigra/matrix.hxx" or |
|
Solve a linear system.
The square matrix a is the coefficient matrix, and the column vectors in b are the right-hand sides of the equation (so, several equations with the same coefficients can be solved in one go). The result is returned int res, whose columns contain the solutions for the correspoinding columns of b. The number of columns of a must equal the number of rows of both b and res, and the number of columns of b and res must be equal. The algorithm uses QR decomposition of a. The algorithm returns
#include "vigra/linear__solve.hxx" or |
|
perform matrix multiplication of matrices a and b. a and b must have matching shapes. The result is returned as a temporary matrix.
#include "vigra/matrix.hxx" or |
|
perform matrix multiplication of matrices a and b. The result is written into r. The three matrices must have matching shapes.
#include "vigra/matrix.hxx" or |
|
Compute the eigensystem of a square, but not necessarily symmetric matrix.
a is a real square matrix, ew is a single-column matrix holding the possibly complex eigenvalues, and ev is a matrix of the same size as a whose columns are the corresponding eigenvectors. Eigenvalues will be sorted from largest to smallest magnitude. The algorithm returns
#include "vigra/eigensystem.hxx" or |
|
calculate the squared Frobenius norm of a matrix. Equal to the sum of squares of the matrix elements. #include "vigra/matrix.hxx" Namespace: vigra |
|
perform matrix multiplication of matrices a and b. a and b must have matching shapes. The result is returned as a temporary matrix.
#include "vigra/matrix.hxx" or |
|
multiply TinyVector a with matrix b. b must be of size
#include "vigra/matrix.hxx" or |
|
multiply matrix a with TinyVector b. a must be of size
#include "vigra/matrix.hxx" or |
|
multiply scalar a with matrix b. The result is returned as a temporary matrix.
#include "vigra/matrix.hxx" or |
|
multiply matrix a with scalar b. The result is returned as a temporary matrix.
#include "vigra/matrix.hxx" or |
|
add matrices a and b. The two matrices must have the same shape. The result is returned as a temporary matrix.
#include "vigra/matrix.hxx" or |
|
negate matrix a. The result is returned as a temporary matrix.
#include "vigra/matrix.hxx" or |
|
subtract matrix b from a. The two matrices must have the same shape. The result is returned as a temporary matrix.
#include "vigra/matrix.hxx" or |
|
divide matrix a by scalar b. The result is returned as a temporary matrix.
#include "vigra/matrix.hxx" or |
|
print a matrix m to the stream s.
#include "vigra/matrix.hxx" or |
|
calculate the outer product of two matrices representing vectors. That is, matrix x must have a single column, and matrix y must have a single row, and the other dimensions must match. The result is returned as a temporary matrix.
#include "vigra/matrix.hxx" or |
|
calculate the outer product of two matrices representing vectors. That is, matrix x must have a single column, and matrix y must have a single row, and the other dimensions must match. The result is written into r.
#include "vigra/matrix.hxx" or |
|
divide matrices a and b pointwise. a and b must have matching shapes. The result is returned as a temporary matrix.
#include "vigra/matrix.hxx" or |
|
divide two matrices a and b pointwise. The result is written into r. All three matrices must have the same shape.
#include "vigra/matrix.hxx" or |
|
multiply matrices a and b pointwise. a and b must have matching shapes. The result is returned as a temporary matrix.
#include "vigra/matrix.hxx" or |
|
multiply two matrices a and b pointwise. The result is written into r. All three matrices must have the same shape.
#include "vigra/matrix.hxx" or |
|
Compute the real roots of a real polynomial using the eigenvalue method.
poly is a real polynomial (compatible to vigra::PolynomialView), and roots a real valued vector (compatible to
#include "vigra/eigensystem.hxx" or
|
|
Compute the roots of a polynomial using the eigenvalue method.
poly is a real polynomial (compatible to vigra::PolynomialView), and roots a complex valued vector (compatible to
#include "vigra/eigensystem.hxx" or
|
|
QR decomposition. a contains the original matrix, results are returned in q and r, where q is a orthogonal matrix, and r is an uppr triangular matrix, and the following relation holds:
assert(a == q * r);
This implementation uses householder transformations. It can be applied in-place, i.e.
#include "vigra/linear__solve.hxx" or |
|
Solve a linear system with right-triangular defining matrix.
The square matrix a must be a right-triangular coefficient matrix as can, for example, be obtained by means of QR decomposition. The column vectors in b are the right-hand sides of the equation (so, several equations with the same coefficients can be solved in one go). The result is returned int x, whose columns contain the solutions for the correspoinding columns of b. The number of columns of a must equal the number of rows of both b and x, and the number of columns of b and x must be equal. This implementation can be applied in-place, i.e.
#include "vigra/linear__solve.hxx" or |
|
Number of rows of a matrix represented as a
#include "vigra/matrix.hxx" or |
|
Create a row vector view for row d of the matrix m
#include "vigra/matrix.hxx" or |
|
divide matrix a by scalar b. The result is written into r. a and r must have the same shape.
#include "vigra/matrix.hxx" or |
|
multiply scalar a with matrix b. The result is written into r. b and r must have the same shape.
#include "vigra/matrix.hxx" or |
|
multiply matrix a with scalar b. The result is written into r. a and r must have the same shape.
#include "vigra/matrix.hxx" or |
|
calculate the squared Frobenius norm of a matrix. Equal to the sum of squares of the matrix elements. #include "vigra/matrix.hxx" Namespace: vigra |
|
subtract matrix b from a. The result is written into r. All three matrices must have the same shape.
#include "vigra/matrix.hxx" or |
|
Compute the eigensystem of a symmetric matrix.
a is a real symmetric matrix, ew is a single-column matrix holding the eigenvalues, and ev is a matrix of the same size as a whose columns are the corresponding eigenvectors. Eigenvalues will be sorted from largest to smallest magnitude. The algorithm returns
#include "vigra/eigensystem.hxx" or |
|
create the transpose of a matrix v. The result is returned as a temporary matrix. Usage:
vigra::Matrix<double> v(rows, cols); v = ...; vigra::Matrix<double> m = transpose(v);
#include "vigra/matrix.hxx" or |
|
transpose matrix v. The result is written into r which must have the correct (i.e. transposed) shape.
#include "vigra/matrix.hxx" or |
© Ullrich Köthe (koethe@informatik.uni-hamburg.de) |
html generated using doxygen and Python
|