Data Structures | |
class | linear_solver |
A generic solver for the linear system ![]() | |
class | linear_solver_lu |
Generic linear solver using LU decomposition. More... | |
class | linear_solver_qr |
Generic linear solver using QR decomposition. More... | |
class | linear_solver_hh |
Generic Householder linear solver. More... | |
class | gsl_solver_LU |
GSL solver by LU decomposition. More... | |
class | gsl_solver_QR |
GSL solver by QR decomposition. More... | |
class | gsl_solver_HH |
GSL Householder solver. More... | |
Functions | |
void | create_givens (const double a, const double b, double &c, double &s) |
Desc. | |
template<class mat1_t , class mat2_t > | |
void | apply_givens_qr (size_t M, size_t N, mat1_t &Q, mat2_t &R, size_t i, size_t j, double c, double s) |
Apply a rotation to matrices from the QR decomposition. | |
template<class mat1_t , class mat2_t > | |
void | apply_givens_lq (size_t M, size_t N, mat1_t &Q, mat2_t &L, size_t i, size_t j, double c, double s) |
Apply a rotation to matrices from the LQ decomposition. | |
template<class vec_t > | |
void | apply_givens_vec (vec_t &v, size_t i, size_t j, double c, double s) |
Apply a rotation to a vector, ![]() | |
template<class mat_t , class vec_t > | |
int | HH_solve (size_t n, mat_t &A, const vec_t &b, vec_t &x) |
Desc. | |
template<class mat_t , class vec_t > | |
int | HH_svx (size_t N, size_t M, mat_t &A, vec_t &x) |
Desc. | |
template<class vec_t > | |
double | householder_transform (const size_t n, vec_t &v) |
Replace the vector v with a householder vector and a coefficient tau that annihilates the last n-1 elements of v . | |
template<class mat_t > | |
double | householder_transform_subcol (mat_t &A, const size_t ir, const size_t ic, const size_t n) |
Compute the householder transform of a vector formed with the last n rows of a column of a matrix. | |
template<class vec_t , class mat_t > | |
int | householder_hm (const size_t M, const size_t N, double tau, const vec_t &v, mat_t &A) |
Apply a householder transformation v,tau to matrix m. | |
template<class mat_t > | |
int | householder_hm_sub (mat_t &M, const size_t ir, const size_t ic, const size_t nr, const size_t nc, const mat_t &M2, const size_t ir2, const size_t ic2, double tau) |
Apply a householder transformation v , tau to submatrix of m. | |
template<class vec_t > | |
int | householder_hv (const size_t N, double tau, const vec_t &v, vec_t &w) |
Apply a householder transformation v to vector w . | |
template<class mat_t , class vec_t > | |
int | householder_hv_sub (const mat_t &M, vec_t &w, double tau, const size_t ie, const size_t N) |
Apply a householder transformation v to vector w . | |
template<class mat1_t , class mat2_t > | |
int | householder_hm_sub2 (const size_t M, const size_t ic, double tau, const mat1_t &mv, mat2_t &A) |
Special version of householder transformation for QR_unpack(). | |
template<class mat_t > | |
int | LU_decomp (const size_t N, mat_t &A, o2scl::permutation &p, int &signum) |
Compute the LU decomposition of the matrix A . | |
template<class mat_t , class vec_t > | |
int | LU_solve (const size_t N, const mat_t &LU, const o2scl::permutation &p, const vec_t &b, vec_t &x) |
Solve a linear system after LU decomposition. | |
template<class mat_t , class vec_t > | |
int | LU_svx (const size_t N, const mat_t &LU, const o2scl::permutation &p, vec_t &x) |
Solve a linear system after LU decomposition in place. | |
template<class mat_t , class vec_t > | |
int | LU_refine (const size_t N, const mat_t &A, const mat_t &LU, const o2scl::permutation &p, const vec_t &b, vec_t &x, vec_t &residual) |
Refine the solution of a linear system. | |
template<class mat_t , class mat_col_t > | |
int | LU_invert (const size_t N, const mat_t &LU, const o2scl::permutation &p, mat_t &inverse) |
Compute the inverse of a matrix from its LU decomposition. | |
template<class mat_t > | |
double | LU_det (const size_t N, const mat_t &LU, int signum) |
Compute the determinant of a matrix from its LU decomposition. | |
template<class mat_t > | |
double | LU_lndet (const size_t N, const mat_t &LU) |
Compute the logarithm of the absolute value of the determinant of a matrix from its LU decomposition. | |
template<class mat_t > | |
int | LU_sgndet (const size_t N, const mat_t &LU, int signum) |
Compute the sign of the determinant of a matrix from its LU decomposition. | |
template<class mat_t , class vec_t > | |
int | QR_decomp (size_t M, size_t N, mat_t &A, vec_t &tau) |
Compute the QR decomposition of matrix A . | |
template<class mat_t , class vec_t , class vec2_t > | |
int | QR_solve (size_t N, const mat_t &QR, const vec_t &tau, const vec2_t &b, vec2_t &x) |
Solve the system A x = b using the QR factorization. | |
template<class mat_t , class vec_t , class vec2_t > | |
int | QR_svx (size_t M, size_t N, const mat_t &QR, const vec_t &tau, vec2_t &x) |
Solve the system A x = b in place using the QR factorization. | |
template<class mat_t , class vec_t , class vec2_t > | |
int | QR_QTvec (const size_t M, const size_t N, const mat_t &QR, const vec_t &tau, vec2_t &v) |
Form the product Q^T v from a QR factorized matrix. | |
template<class mat1_t , class mat2_t , class mat3_t , class vec_t > | |
int | QR_unpack (const size_t M, const size_t N, const mat1_t &QR, const vec_t &tau, mat2_t &Q, mat3_t &R) |
Unpack the QR matrix to the individual Q and R components. | |
template<class mat1_t , class mat2_t , class vec1_t , class vec2_t > | |
int | QR_update (size_t M, size_t N, mat1_t &Q, mat2_t &R, vec1_t &w, vec2_t &v) |
Update a QR factorisation for A= Q R , A' = A + u v^T,. |
void o2scl_linalg::apply_givens_lq | ( | size_t | M, | |
size_t | N, | |||
mat1_t & | Q, | |||
mat2_t & | L, | |||
size_t | i, | |||
size_t | j, | |||
double | c, | |||
double | s | |||
) | [inline] |
Apply a rotation to matrices from the LQ decomposition.
This performs and
.
Definition at line 83 of file givens_base.h.
void o2scl_linalg::apply_givens_qr | ( | size_t | M, | |
size_t | N, | |||
mat1_t & | Q, | |||
mat2_t & | R, | |||
size_t | i, | |||
size_t | j, | |||
double | c, | |||
double | s | |||
) | [inline] |
Apply a rotation to matrices from the QR decomposition.
This performs and
.
Definition at line 57 of file givens_base.h.
int o2scl_linalg::householder_hm_sub | ( | mat_t & | M, | |
const size_t | ir, | |||
const size_t | ic, | |||
const size_t | nr, | |||
const size_t | nc, | |||
const mat_t & | M2, | |||
const size_t | ir2, | |||
const size_t | ic2, | |||
double | tau | |||
) | [inline] |
Apply a householder transformation v
, tau
to submatrix of m.
Used in QR_decomp().
Definition at line 163 of file householder_base.h.
int o2scl_linalg::householder_hv_sub | ( | const mat_t & | M, | |
vec_t & | w, | |||
double | tau, | |||
const size_t | ie, | |||
const size_t | N | |||
) | [inline] |
Apply a householder transformation v
to vector w
.
Used in QR_QTvec().
Definition at line 231 of file householder_base.h.
double o2scl_linalg::householder_transform_subcol | ( | mat_t & | A, | |
const size_t | ir, | |||
const size_t | ic, | |||
const size_t | n | |||
) | [inline] |
Compute the householder transform of a vector formed with the last n
rows of a column of a matrix.
Used in QR_decomp().
Definition at line 95 of file householder_base.h.
int o2scl_linalg::LU_decomp | ( | const size_t | N, | |
mat_t & | A, | |||
o2scl::permutation & | p, | |||
int & | signum | |||
) | [inline] |
Compute the LU decomposition of the matrix A
.
On output the diagonal and upper triangular part of the input matrix A contain the matrix U. The lower triangular part of the input matrix (excluding the diagonal) contains L. The diagonal elements of L are unity, and are not stored.
The permutation matrix P is encoded in the permutation p. The j-th column of the matrix P is given by the k-th column of the identity matrix, where k = p_j the j-th element of the permutation vector. The sign of the permutation is given by signum. It has the value (-1)^n, where n is the number of interchanges in the permutation.
The algorithm used in the decomposition is Gaussian Elimination with partial pivoting (Golub & Van Loan, Matrix Computations, Algorithm 3.4.1).
double o2scl_linalg::LU_det | ( | const size_t | N, | |
const mat_t & | LU, | |||
int | signum | |||
) | [inline] |
Compute the determinant of a matrix from its LU decomposition.
These functions compute the determinant of a matrix A from its LU decomposition, LU. The determinant is computed as the product of the diagonal elements of U and the sign of the row permutation signum.
int o2scl_linalg::LU_invert | ( | const size_t | N, | |
const mat_t & | LU, | |||
const o2scl::permutation & | p, | |||
mat_t & | inverse | |||
) | [inline] |
Compute the inverse of a matrix from its LU decomposition.
These functions compute the inverse of a matrix A from its LU decomposition (LU,p), storing the result in the matrix inverse. The inverse is computed by solving the system A x = b for each column of the identity matrix. It is preferable to avoid direct use of the inverse whenever possible, as the linear solver functions can obtain the same result more efficiently and reliably (consult any introductory textbook on numerical linear algebra for details).
double o2scl_linalg::LU_lndet | ( | const size_t | N, | |
const mat_t & | LU | |||
) | [inline] |
Compute the logarithm of the absolute value of the determinant of a matrix from its LU decomposition.
These functions compute the logarithm of the absolute value of the determinant of a matrix A, , from its LU decomposition, LU. This function may be useful if the direct computation of the determinant would overflow or underflow.
int o2scl_linalg::LU_refine | ( | const size_t | N, | |
const mat_t & | A, | |||
const mat_t & | LU, | |||
const o2scl::permutation & | p, | |||
const vec_t & | b, | |||
vec_t & | x, | |||
vec_t & | residual | |||
) | [inline] |
int o2scl_linalg::LU_sgndet | ( | const size_t | N, | |
const mat_t & | LU, | |||
int | signum | |||
) | [inline] |
int o2scl_linalg::LU_solve | ( | const size_t | N, | |
const mat_t & | LU, | |||
const o2scl::permutation & | p, | |||
const vec_t & | b, | |||
vec_t & | x | |||
) | [inline] |
int o2scl_linalg::LU_svx | ( | const size_t | N, | |
const mat_t & | LU, | |||
const o2scl::permutation & | p, | |||
vec_t & | x | |||
) | [inline] |
int o2scl_linalg::QR_update | ( | size_t | M, | |
size_t | N, | |||
mat1_t & | Q, | |||
mat2_t & | R, | |||
vec1_t & | w, | |||
vec2_t & | v | |||
) | [inline] |
Update a QR factorisation for A= Q R , A' = A + u v^T,.
M
and N
are the number of rows and columns of R
.
* Q' R' = QR + u v^T * = Q (R + Q^T u v^T) * = Q (R + w v^T) * * where w = Q^T u. * * Algorithm from Golub and Van Loan, "Matrix Computations", Section * 12.5 (Updating Matrix Factorizations, Rank-One Changes)
Documentation generated with Doxygen and provided under the GNU Free Documentation License. See License Information for details.
Project hosting provided by
,
O2scl Sourceforge Project Page