LTL  2.0.x
Classes
lapack.h File Reference

Classes

struct  lapack_syev_dispatch< T >
 
struct  lapack_syev_dispatch< double >
 
struct  lapack_syev_dispatch< float >
 
struct  lapack_sbev_dispatch< T >
 
struct  lapack_sbev_dispatch< double >
 
struct  lapack_sbev_dispatch< float >
 

Functions

template<typename T >
bool lapack_gesv (MArray< T, 2 > &A, MArray< T, 1 > &b, MArray< int, 1 > &ipiv)
 
template<typename T >
bool lapack_getrf (MArray< T, 2 > &A, MArray< int, 1 > &ipiv)
 
template<typename T >
bool lapack_getri (MArray< T, 2 > &A, MArray< int, 1 > &ipiv)
 
template<typename T >
bool lapack_gesv (MArray< T, 2 > &A, MArray< T, 2 > &B, MArray< int, 1 > ipiv)
 
template<typename T >
bool lapack_syev (const MArray< T, 2 > &A, MArray< T, 1 > &val, MArray< T, 2 > &vec)
 
template<typename T >
bool lapack_sbev (const MArray< T, 2 > &AB, MArray< T, 1 > &val, MArray< T, 2 > &vec)
 

Function Documentation

template<typename T >
bool lapack_gesv ( MArray< T, 2 > &  A,
MArray< T, 1 > &  b,
MArray< int, 1 > &  ipiv 
)

GESV computes the solution to a system of linear equations A * x = b, where A is an N-by-N matrix and x and b are N-vectors.

LU decomposition with partial pivoting and row interchanges is used to factor A as A = P * L * U, where P is a permutation matrix, L is unit lower triangular, and U is upper triangular. The factored form of A is then used to solve the system of equations A * x = b.

The matrix A is replaced by L and U (or A^-1, see below). The vector b is overwritten with the solution of Ax=b.

The permutation matrix is returned in ipiv. Call getri() to compute the inverse of A.

References LTL_ASSERT.

template<typename T >
bool lapack_getrf ( MArray< T, 2 > &  A,
MArray< int, 1 > &  ipiv 
)

DGETRF - compute an LU factorization of a general M-by-N matrix A using partial pivoting with row interchanges

DGETRF computes an LU factorization of a general M-by-N matrix A using partial pivoting with row interchanges.

The factorization has the form A = P * L * U where P is a permutation matrix, L is lower triangular with unit diagonal elements (lower trapezoidal if m > n), and U is upper triangular (upper trapezoidal if m < n).

This is the right-looking Level 3 BLAS version of the algorithm.

template<typename T >
bool lapack_getri ( MArray< T, 2 > &  A,
MArray< int, 1 > &  ipiv 
)

GETRI computes the inverse of a matrix using the LU factorization computed by DGETRF.

This method inverts U and then computes inv(A) by solving the system inv(A)*L = inv(U) for inv(A).

A (input/output) DOUBLE PRECISION array, dimension (LDA,N) On entry, the factors L and U from the factorization A = P*L*U as computed by DGETRF. On exit, if INFO = 0, the inverse of the original matrix A.

IPIV (input) INTEGER array, dimension (N) The pivot indices from DGETRF; for 1<=i<=N, row i of the matrix was interchanged with row IPIV(i).

References LTL_ASSERT.

template<typename T >
bool lapack_gesv ( MArray< T, 2 > &  A,
MArray< T, 2 > &  B,
MArray< int, 1 >  ipiv 
)

GESV computes the solution to Nrhs systems of linear equations A * X = B, where A is an N-by-N matrix and X and B are N x Nrhs-matrices.

LU decomposition with partial pivoting and row interchanges is used to factor A as A = P * L * U, where P is a permutation matrix, L is unit lower triangular, and U is upper triangular. The factored form of A is then used to solve the system of equations A * x = b.

The matrix A is replaced by L and U. The Matrix B is overwritten with the solutions of A x_i = b_i, where x_i and b_i are the column vectors of X and B

References LTL_ASSERT.

template<typename T >
bool lapack_syev ( const MArray< T, 2 > &  A,
MArray< T, 1 > &  val,
MArray< T, 2 > &  vec 
)

Purpose

SYEV computes all eigenvalues and, optionally, eigenvectors of a real symmetric matrix A.

References lapack_syev_dispatch< T >::call(), and LTL_ASSERT.

template<typename T >
bool lapack_sbev ( const MArray< T, 2 > &  AB,
MArray< T, 1 > &  val,
MArray< T, 2 > &  vec 
)

Purpose

SYEV computes all eigenvalues and, optionally, eigenvectors of a real symmetric matrix A in banded storage (assuming superdiagonals are stored.

References lapack_sbev_dispatch< T >::call(), and LTL_ASSERT.