The namespace for linear algebra classes and functions. More...
Classes | |
class | lanczos |
Lanczos diagonalization. More... | |
class | linear_solver |
A generic solver for the linear system ![]() | |
class | linear_solver_arma |
Armadillo linear solver. More... | |
class | linear_solver_eigen_colQR |
Eigen linear solver using QR decomposition with column pivoting. More... | |
class | linear_solver_eigen_fullLU |
Eigen linear solver using LU decomposition with full pivoting. More... | |
class | linear_solver_eigen_fullQR |
Eigen linear solver using QR decomposition with full pivoting. More... | |
class | linear_solver_eigen_houseQR |
Eigen linear solver using QR decomposition with column pivoting. More... | |
class | linear_solver_eigen_LDLT |
Eigen linear solver using LDLT decomposition with full pivoting. More... | |
class | linear_solver_eigen_LLT |
Eigen linear solver using LLT decomposition with full pivoting. More... | |
class | linear_solver_eigen_partLU |
Eigen linear solver using LU decomposition with partial pivoting. More... | |
class | linear_solver_HH |
Generic Householder linear solver. More... | |
class | linear_solver_LU |
Generic linear solver using LU decomposition. More... | |
class | linear_solver_QR |
Generic linear solver using QR decomposition. More... | |
class | ubvector_2_mem |
Allocation object for 2 arrays of equal size. More... | |
class | ubvector_4_mem |
Allocation object for 4 arrays of equal size. More... | |
class | ubvector_5_mem |
Allocation object for 5 arrays of equal size. More... | |
Functions | |
template<class mat_t , class vec_t , class vec2_t > | |
int | bidiag_decomp (size_t M, size_t N, mat_t &A, vec_t &tau_U, vec2_t &tau_V) |
Factor a matrix into bidiagonal form. More... | |
template<class mat_t , class vec_t , class mat2_t , class vec2_t , class mat3_t , class vec3_t , class vec4_t > | |
int | bidiag_unpack (size_t M, size_t N, const mat_t &A, const vec_t &tau_U, mat2_t &U, const vec2_t &tau_V, mat3_t &V, vec3_t &diag, vec4_t &superdiag) |
Unpack a matrix A with the bidiagonal decomposition and create matrices U , V , diagonal diag and superdiagonal superdiag . More... | |
template<class mat_t , class vec_t , class vec2_t , class mat2_t > | |
int | bidiag_unpack2 (size_t M, size_t N, mat_t &A, vec_t &tau_U, vec2_t &tau_V, mat2_t &V) |
Unpack a matrix A with the bidiagonal decomposition and create matrix V . | |
template<class mat_t , class vec_t , class vec2_t > | |
int | bidiag_unpack_B (size_t M, size_t N, const mat_t &A, vec_t &diag, vec2_t &superdiag) |
Unpack the diagonal and superdiagonal of the bidiagonal decomposition of A into diag and superdiag . | |
template<> | |
int | cholesky_decomp< Eigen::MatrixXd > (const size_t M, Eigen::MatrixXd &A, bool err_on_fail) |
Eigen specialization of cholesky_decomp() | |
template<class mat_t > | |
int | cholesky_decomp (const size_t M, mat_t &A, bool err_on_fail=true) |
Compute the in-place Cholesky decomposition of a symmetric positive-definite square matrix. More... | |
template<class mat_t , class vec_t , class vec2_t > | |
void | cholesky_solve (const size_t N, const mat_t &LLT, const vec_t &b, vec2_t &x) |
Solve a symmetric positive-definite linear system after a Cholesky decomposition. More... | |
template<class mat_t , class vec_t > | |
void | cholesky_svx (const size_t N, const mat_t &LLT, vec_t &x) |
Solve a linear system in place using a Cholesky decomposition. | |
template<class mat_t > | |
void | cholesky_invert (const size_t N, mat_t &LLT) |
Compute the inverse of a symmetric positive definite matrix given the Cholesky decomposition. More... | |
template<class mat_t , class vec_t > | |
int | cholesky_decomp_unit (const size_t N, mat_t &A, vec_t &D) |
Cholesky decomposition with unit-diagonal triangular parts. More... | |
void | create_givens (const double a, const double b, double &c, double &s) |
Create a Givens rotation matrix. More... | |
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. More... | |
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. More... | |
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_svx (size_t N, size_t M, mat_t &A, vec_t &x) |
Solve a linear system after Householder decomposition in place. More... | |
template<class mat_t , class vec_t , class vec2_t > | |
int | HH_solve (size_t n, mat_t &A, const vec_t &b, vec2_t &x) |
Solve linear system after Householder decomposition. | |
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 v[1] through v[n-1] (inclusive) More... | |
template<class mat_t > | |
double | householder_transform_subcol (mat_t &A, const size_t ir, const size_t ic, const size_t M) |
Compute the Householder transform of a vector formed with n rows of a column of a matrix. More... | |
template<class mat_t > | |
double | householder_transform_subrow (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 columns of a row of a matrix. More... | |
template<class vec_t , class mat_t > | |
void | householder_hm (const size_t M, const size_t N, double tau, const vec_t &v, mat_t &A) |
Apply a Householder transformation ![]() ![]() (M,N) More... | |
template<class mat_t > | |
void | householder_hm_subcol (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 to the lower-right part of M when the transformation is stored in a column of M2 . More... | |
template<class mat_t > | |
void | householder_hm_subrow (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 to the lower-right part of M when the transformation is stored in a row of M2 . More... | |
template<class vec_t , class vec2_t > | |
void | householder_hv (const size_t N, double tau, const vec_t &v, vec2_t &w) |
Apply a Householder transformation v to vector w . | |
template<class mat_t , class vec_t > | |
void | householder_hv_subcol (const mat_t &A, vec_t &w, double tau, const size_t ie, const size_t N) |
Apply a Householder transformation v to vector w where v is stored as a column in a matrix A . More... | |
template<class mat_t > | |
void | householder_hm1 (const size_t M, const size_t N, double tau, mat_t &A) |
Apply a Householder transformation ![]() | |
template<class mat_t > | |
void | householder_hm1_sub (const size_t M, const size_t N, double tau, mat_t &A, size_t ir, size_t ic) |
Apply a Householder transformation ![]() | |
template<class vec_t , class mat_t > | |
void | householder_mh (const size_t M, const size_t N, double tau, const vec_t &v, mat_t &A) |
Apply the Householder transformation (v,tau) to the right-hand side of the matrix A . | |
template<class mat_t , class mat2_t > | |
void | householder_mh_subrow (mat_t &M, const size_t ir, const size_t ic, const size_t nr, const size_t nc, const mat2_t &M2, const size_t ir2, const size_t ic2, double tau) |
Apply the Householder transformation (v,tau) to the right-hand side of the matrix A . More... | |
template<class mat_t > | |
int | diagonal_has_zero (const size_t N, mat_t &A) |
Return 1 if at least one of the elements in the diagonal is zero. | |
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 . More... | |
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. More... | |
template<class mat_t , class mat_row_t > | |
int | LU_decomp_alt (const size_t N, mat_t &A, o2scl::permutation &p, int &signum) |
An alternate form of LU decomposition with matrix row objects. | |
template<class mat_t , class vec_t , class vec2_t > | |
int | LU_solve (const size_t N, const mat_t &LU, const o2scl::permutation &p, const vec_t &b, vec2_t &x) |
Solve a linear system after LU decomposition. More... | |
template<class mat_t , class mat2_t , class vec_t , class vec2_t , class vec3_t > | |
int | LU_refine (const size_t N, const mat_t &A, const mat2_t &LU, const o2scl::permutation &p, const vec_t &b, vec2_t &x, vec3_t &residual) |
Refine the solution of a linear system. More... | |
template<class mat_t , class mat2_t , class mat_col_t > | |
int | LU_invert (const size_t N, const mat_t &LU, const o2scl::permutation &p, mat2_t &inverse) |
Compute the inverse of a matrix from its LU decomposition. More... | |
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. More... | |
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. More... | |
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. More... | |
template<> | |
void | QR_decomp_unpack< arma::mat, arma::mat, arma::mat > (const size_t M, const size_t N, arma::mat &A, arma::mat &Q, arma::mat &R) |
Armadillo specialization of QR_decomp_unpack(). | |
template<> | |
void | QR_decomp_unpack< Eigen::MatrixXd, Eigen::MatrixXd, Eigen::MatrixXd > (const size_t M, const size_t N, Eigen::MatrixXd &A, Eigen::MatrixXd &Q, Eigen::MatrixXd &R) |
Eigen specialization of QR_decomp_unpack(). | |
template<class mat_t , class vec_t > | |
void | 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 > | |
void | 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 > | |
void | 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 mat_t , class vec_t , class vec2_t > | |
void | 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 , class vec3_t > | |
void | QR_solve (size_t N, const mat_t &QR, const vec_t &tau, const vec2_t &b, vec3_t &x) |
Solve the system A x = b using the QR factorization. | |
template<class mat1_t , class mat2_t , class vec1_t , class vec2_t > | |
void | 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,. More... | |
template<class mat_t , class mat2_t , class mat3_t > | |
void | QR_decomp_unpack (const size_t M, const size_t N, mat_t &A, mat2_t &Q, mat3_t &R) |
Compute the unpacked QR decomposition of matrix A . More... | |
template<class mat_t , class vec_t , class vec2_t > | |
void | QRPT_decomp (size_t M, size_t N, mat_t &A, vec_t &tau, o2scl::permutation &p, int &signum, vec2_t &norm) |
Compute the QR decomposition of matrix A . | |
template<class mat_t , class mat2_t , class vec_t , class vec2_t > | |
void | SV_decomp (size_t M, size_t N, mat_t &A, mat2_t &V, vec_t &S, vec2_t &work) |
Factorise a general matrix into its SV decomposition using the Golub-Reinsch algorithm. More... | |
template<class mat_t , class mat2_t , class mat3_t , class vec_t , class vec2_t > | |
void | SV_decomp_mod (size_t M, size_t N, mat_t &A, mat2_t &X, mat3_t &V, vec_t &S, vec2_t &work) |
SV decomposition by the modified Golub-Reinsch algorithm which is better for ![]() | |
template<class mat_t , class mat2_t , class vec_t , class vec2_t , class vec3_t > | |
void | SV_solve (size_t M, size_t N, mat_t &U, mat2_t &V, vec_t &S, vec2_t &b, vec3_t &x) |
Solve the system A x = b using the SV decomposition. More... | |
template<class mat_t , class mat2_t , class vec_t > | |
void | SV_decomp_jacobi (size_t M, size_t N, mat_t &A, mat2_t &Q, vec_t &S) |
SV decomposition using one-sided Jacobi orthogonalization. More... | |
template<class mat_t , class vec_t > | |
void | balance_columns (size_t M, size_t N, mat_t &A, vec_t &D) |
Balance a general matrix A by scaling the columns by the diagonal matrix D. More... | |
template<class vec_t , class vec2_t > | |
void | chop_small_elements (size_t N, vec_t &d, vec2_t &f) |
Zero out small elements in f according to the scales set in d . More... | |
template<class vec_t , class vec2_t > | |
double | trailing_eigenvalue (size_t n, const vec_t &d, const vec2_t &f) |
Desc. More... | |
void | create_schur (double d0, double f0, double d1, double &c, double &s) |
Desc. More... | |
template<class vec_t , class vec2_t , class mat_t , class mat2_t > | |
void | svd2 (size_t M, size_t N, vec_t &d, vec2_t &f, mat_t &U, mat2_t &V) |
2-variable SVD More... | |
template<class vec_t , class vec2_t , class mat_t , class mat2_t > | |
void | svd2_sub (size_t M, size_t N, vec_t &d, vec2_t &f, mat_t &U, mat2_t &V, size_t a) |
Shifted 2-variable SVD. More... | |
template<class vec_t , class vec2_t , class mat_t > | |
void | chase_out_intermediate_zero (size_t M, size_t n, vec_t &d, vec2_t &f, mat_t &U, size_t k0) |
Desc. More... | |
template<class vec_t , class vec2_t , class mat_t > | |
void | chase_out_trailing_zero (size_t N, size_t n, vec_t &d, vec2_t &f, mat_t &V) |
Desc. More... | |
template<class vec_t , class vec2_t , class mat_t > | |
void | chase_out_trailing_zero_sub (size_t N, size_t n, vec_t &d, vec2_t &f, mat_t &V, size_t a) |
Desc. More... | |
template<class vec_t , class vec2_t , class mat_t , class mat2_t > | |
void | qrstep (size_t M, size_t N, size_t n, vec_t &d, vec2_t &f, mat_t &U, mat2_t &V) |
Desc. More... | |
template<class vec_t , class vec2_t , class mat_t , class mat2_t > | |
void | qrstep_sub (size_t M, size_t N, size_t n, vec_t &d, vec2_t &f, mat_t &U, mat2_t &V, size_t a) |
A special form of qrstep() for SV_decomp() More... | |
template<class vec_t , class vec2_t , class vec3_t , class vec4_t , class mem_t , class mem_vec_t > | |
void | solve_tridiag_sym (const vec_t &diag, const vec2_t &offdiag, const vec3_t &b, vec4_t &x, size_t N, mem_t &m) |
Solve a symmetric tridiagonal linear system with user-specified memory. More... | |
template<class vec_t , class vec2_t , class vec3_t , class vec4_t , class vec5_t , class mem_t , class mem_vec_t > | |
void | solve_tridiag_nonsym (const vec_t &diag, const vec2_t &abovediag, const vec3_t &belowdiag, const vec4_t &rhs, vec5_t &x, size_t N, mem_t &m) |
Solve an asymmetric tridiagonal linear system with user-specified memory. More... | |
template<class vec_t , class vec2_t , class vec3_t , class vec4_t , class mem_t , class mem_vec_t > | |
void | solve_cyc_tridiag_sym (const vec_t &diag, const vec2_t &offdiag, const vec3_t &b, vec4_t &x, size_t N, mem_t &m) |
Solve a symmetric cyclic tridiagonal linear system with user specified memory. More... | |
template<class vec_t , class vec2_t , class vec3_t , class vec4_t , class vec5_t , class mem_t , class mem_vec_t > | |
void | solve_cyc_tridiag_nonsym (const vec_t &diag, const vec2_t &abovediag, const vec3_t &belowdiag, const vec4_t &rhs, vec5_t &x, size_t N, mem_t &m) |
Solve an asymmetric cyclic tridiagonal linear system with user-specified memory. More... | |
template<class vec_t , class vec2_t , class vec3_t , class vec4_t > | |
void | solve_tridiag_sym (const vec_t &diag, const vec2_t &offdiag, const vec3_t &b, vec4_t &x, size_t N) |
Solve a symmetric tridiagonal linear system. | |
template<class vec_t , class vec2_t , class vec3_t , class vec4_t , class vec5_t > | |
void | solve_tridiag_nonsym (const vec_t &diag, const vec2_t &abovediag, const vec3_t &belowdiag, const vec4_t &rhs, vec5_t &x, size_t N) |
Solve an asymmetric tridiagonal linear system. | |
template<class vec_t , class vec2_t , class vec3_t , class vec4_t > | |
void | solve_cyc_tridiag_sym (const vec_t &diag, const vec2_t &offdiag, const vec3_t &b, vec4_t &x, size_t N) |
Solve a symmetric cyclic tridiagonal linear system. | |
template<class vec_t , class vec2_t , class vec3_t , class vec4_t , class vec5_t > | |
void | solve_cyc_tridiag_nonsym (const vec_t &diag, const vec2_t &abovediag, const vec3_t &belowdiag, const vec4_t &rhs, vec5_t &x, size_t N) |
Solve an asymmetric cyclic tridiagonal linear system. | |
See Linear Algebra for more complete information about linear algebra in O2scl .
This namespace documentation is in the file src/base/cblas.h
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 | ||
) |
This performs and
.
Definition at line 86 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 | ||
) |
This performs and
.
Definition at line 58 of file givens_base.h.
void o2scl_linalg::balance_columns | ( | size_t | M, |
size_t | N, | ||
mat_t & | A, | ||
vec_t & | D | ||
) |
The function computes where
is a diagonal matrix.
Definition at line 615 of file svd_base.h.
int o2scl_linalg::bidiag_decomp | ( | size_t | M, |
size_t | N, | ||
mat_t & | A, | ||
vec_t & | tau_U, | ||
vec2_t & | tau_V | ||
) |
Factor matrix A of size (M,N)
with into
where U and V are orthogonal and B is upper bidiagonal.
After the function call, the matrix is stored the diagonal and first superdiagonal of
A
. The matrices and
are stored as packed sets of Householder transformations in the lower and upper triangular parts of
A
, respectively.
Adapted from the GSL version which was based on algorithm 5.4.2 in Golub96.
Definition at line 65 of file bidiag_base.h.
int o2scl_linalg::bidiag_unpack | ( | size_t | M, |
size_t | N, | ||
const mat_t & | A, | ||
const vec_t & | tau_U, | ||
mat2_t & | U, | ||
const vec2_t & | tau_V, | ||
mat3_t & | V, | ||
vec3_t & | diag, | ||
vec4_t & | superdiag | ||
) |
Given a matrix A
of size (M,N)
with created by bidiag_decomp(), this function creates the matrix
U
of size (M,N)
, the matrix V
of size (N,N)
, the diagonal diag
of size N
and the super-diagonal superdiag
of size N-1
.
Definition at line 113 of file bidiag_base.h.
void o2scl_linalg::chase_out_intermediate_zero | ( | size_t | M, |
size_t | n, | ||
vec_t & | d, | ||
vec2_t & | f, | ||
mat_t & | U, | ||
size_t | k0 | ||
) |
The vector d
should be of size n
, the vector f
should be of size n
, and the matrix U should be of size (M,n)
Used in qrstep() and qrstep_sub().
Definition at line 457 of file svdstep_base.h.
void o2scl_linalg::chase_out_trailing_zero | ( | size_t | N, |
size_t | n, | ||
vec_t & | d, | ||
vec2_t & | f, | ||
mat_t & | V | ||
) |
The vector d
should be of size n
, the vector f
should be of size n
, and the matrix V should be of size (N,n)
Used in qrstep().
Definition at line 507 of file svdstep_base.h.
void o2scl_linalg::chase_out_trailing_zero_sub | ( | size_t | N, |
size_t | n, | ||
vec_t & | d, | ||
vec2_t & | f, | ||
mat_t & | V, | ||
size_t | a | ||
) |
The vector d
should be of size n
, the vector f
should be of size n
, and the matrix V should be of size (N,n)
Used in qrstep_sub().
Definition at line 556 of file svdstep_base.h.
int o2scl_linalg::cholesky_decomp | ( | const size_t | M, |
mat_t & | A, | ||
bool | err_on_fail = true |
||
) |
On input, the upper triangular part of A is ignored (only the lower triangular part and diagonal are used). On output, the diagonal and lower triangular part contain the matrix L and the upper triangular part contains L^T.
If the matrix is not positive-definite, the error handler will be called, unless err_on_fail
is false, in which case a non-zero value will be returned.
Definition at line 62 of file cholesky_base.h.
int o2scl_linalg::cholesky_decomp_unit | ( | const size_t | N, |
mat_t & | A, | ||
vec_t & | D | ||
) |
A = L D L^T, where diag(L) = (1,1,...,1). Upon exit, A contains L and L^T as for Cholesky, and the diagonal of A is (1,1,...,1). The vector Dis set to the diagonal elements of the diagonal matrix D.
Definition at line 313 of file cholesky_base.h.
void o2scl_linalg::cholesky_invert | ( | const size_t | N, |
mat_t & | LLT | ||
) |
Given a Cholesky decomposition produced by cholesky_decomp(), this function returns the inverse of that matrix in LLT
.
Definition at line 218 of file cholesky_base.h.
void o2scl_linalg::cholesky_solve | ( | const size_t | N, |
const mat_t & | LLT, | ||
const vec_t & | b, | ||
vec2_t & | x | ||
) |
Given the Cholesky decomposition of a matrix A in LLT
, this function solves the system A*x=b
.
Definition at line 171 of file cholesky_base.h.
void o2scl_linalg::chop_small_elements | ( | size_t | N, |
vec_t & | d, | ||
vec2_t & | f | ||
) |
The parameter N
is the size of d
.
Used in SV_decomp().
Definition at line 58 of file svdstep_base.h.
void o2scl_linalg::create_givens | ( | const double | a, |
const double | b, | ||
double & | c, | ||
double & | s | ||
) |
Given values a
and b
, create entries c
and s
of a matrix for which
void o2scl_linalg::create_schur | ( | double | d0, |
double | f0, | ||
double | d1, | ||
double & | c, | ||
double & | s | ||
) |
Used in svd2() and svd2_sub().
Definition at line 120 of file svdstep_base.h.
int o2scl_linalg::HH_svx | ( | size_t | N, |
size_t | M, | ||
mat_t & | A, | ||
vec_t & | x | ||
) |
void o2scl_linalg::householder_hm | ( | const size_t | M, |
const size_t | N, | ||
double | tau, | ||
const vec_t & | v, | ||
mat_t & | A | ||
) |
The vector v
must have at least N
entries, with the exception that the vector element v[0]
is never referenced by this function.
Definition at line 210 of file householder_base.h.
void o2scl_linalg::householder_hm1_sub | ( | const size_t | M, |
const size_t | N, | ||
double | tau, | ||
mat_t & | A, | ||
size_t | ir, | ||
size_t | ic | ||
) |
Used in SV_decomp_mod() and bidiag_unpack2().
Definition at line 436 of file householder_base.h.
void o2scl_linalg::householder_hm_subcol | ( | 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 | ||
) |
This applies a householder transformation (v,tau)
to a lower-right submatrix of M
. The submatrix has nr-ir
rows and nc-ic
columns and starts at row ir
of column ic
of the original matrix M
. The vector containing the transformation is taken from a column of M2
starting at row ir2
and column ic2
. The matrix M2
must have at least ic2+1
columns and at least nr-ir+ir2
rows.
This function is used in QR_decomp() and QR_unpack() .
Definition at line 255 of file householder_base.h.
void o2scl_linalg::householder_hm_subrow | ( | 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 | ||
) |
This applies a householder transformation (v,tau)
to a lower-right submatrix of M
. The submatrix has nr-ir
rows and nc-ic
columns and starts at row ir
of column ic
of the original matrix M
. The vector containing the transformation is taken from a row of M2
starting at row ir2
and column ic2
. The matrix M2
must have at least ir2+1
rows and nr-ir+ic2
columns.
Used in bidiag_unpack().
Definition at line 301 of file householder_base.h.
void o2scl_linalg::householder_hv_subcol | ( | const mat_t & | A, |
vec_t & | w, | ||
double | tau, | ||
const size_t | ie, | ||
const size_t | N | ||
) |
Used in QR_QTvec().
Definition at line 364 of file householder_base.h.
void o2scl_linalg::householder_mh_subrow | ( | mat_t & | M, |
const size_t | ir, | ||
const size_t | ic, | ||
const size_t | nr, | ||
const size_t | nc, | ||
const mat2_t & | M2, | ||
const size_t | ir2, | ||
const size_t | ic2, | ||
double | tau | ||
) |
Used in bidiag_decomp().
Definition at line 519 of file householder_base.h.
double o2scl_linalg::householder_transform | ( | const size_t | n, |
vec_t & | v | ||
) |
On exit, this function returns the value of . If
n
is less than or equal to 1 then this function returns zero without calling the error handler.
Definition at line 68 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 | M | ||
) |
This performs a Householder transform of a vector defined by a column of a matrix A
which starts at element A(ir,ic)
and ends at element A(M-1,ic)
. If M-1
is equal to ir+1
, this function quietly does nothing.
Used in QR_decomp() and SV_decomp_mod().
Definition at line 118 of file householder_base.h.
double o2scl_linalg::householder_transform_subrow | ( | mat_t & | A, |
const size_t | ir, | ||
const size_t | ic, | ||
const size_t | N | ||
) |
This performs a Householder transform of a vector defined by a row of a matrix A
which starts at element A(ir,ic)
and ends at element A(ir,N-1)
If N-1
is equal to ic
, this function quietly does nothing.
Definition at line 165 of file householder_base.h.
int o2scl_linalg::LU_decomp | ( | const size_t | N, |
mat_t & | A, | ||
o2scl::permutation & | p, | ||
int & | signum | ||
) |
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 | ||
) |
int o2scl_linalg::LU_invert | ( | const size_t | N, |
const mat_t & | LU, | ||
const o2scl::permutation & | p, | ||
mat2_t & | inverse | ||
) |
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.
double o2scl_linalg::LU_lndet | ( | const size_t | N, |
const mat_t & | LU | ||
) |
int o2scl_linalg::LU_refine | ( | const size_t | N, |
const mat_t & | A, | ||
const mat2_t & | LU, | ||
const o2scl::permutation & | p, | ||
const vec_t & | b, | ||
vec2_t & | x, | ||
vec3_t & | residual | ||
) |
int o2scl_linalg::LU_sgndet | ( | const size_t | N, |
const mat_t & | LU, | ||
int | signum | ||
) |
int o2scl_linalg::LU_solve | ( | const size_t | N, |
const mat_t & | LU, | ||
const o2scl::permutation & | p, | ||
const vec_t & | b, | ||
vec2_t & | x | ||
) |
int o2scl_linalg::LU_svx | ( | const size_t | N, |
const mat_t & | LU, | ||
const o2scl::permutation & | p, | ||
vec_t & | x | ||
) |
void o2scl_linalg::QR_decomp_unpack | ( | const size_t | M, |
const size_t | N, | ||
mat_t & | A, | ||
mat2_t & | Q, | ||
mat3_t & | R | ||
) |
void o2scl_linalg::QR_update | ( | size_t | M, |
size_t | N, | ||
mat1_t & | Q, | ||
mat2_t & | R, | ||
vec1_t & | w, | ||
vec2_t & | v | ||
) |
The parameters M
and N
are the number of rows and columns of the matrix 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)
void o2scl_linalg::qrstep | ( | size_t | M, |
size_t | N, | ||
size_t | n, | ||
vec_t & | d, | ||
vec2_t & | f, | ||
mat_t & | U, | ||
mat2_t & | V | ||
) |
The vector d
should be of size n
, the vector f
should be of size n
, the matrix U should be of size (M,N)
, and the matrix V
should be of size (N,N)
.
Definition at line 603 of file svdstep_base.h.
void o2scl_linalg::qrstep_sub | ( | size_t | M, |
size_t | N, | ||
size_t | n, | ||
vec_t & | d, | ||
vec2_t & | f, | ||
mat_t & | U, | ||
mat2_t & | V, | ||
size_t | a | ||
) |
The vector d
should be of size n
, the vector f
should be of size n
, the matrix U should be of size (M,n)
, and the matrix V
should be of size (N,n)
.
A version of qrstep(), but starting at the a'th column of U, the a'th column of V, and the a'th entries of d
and f
.
This function is used by SV_decomp().
Definition at line 762 of file svdstep_base.h.
void o2scl_linalg::solve_cyc_tridiag_nonsym | ( | const vec_t & | diag, |
const vec2_t & | abovediag, | ||
const vec3_t & | belowdiag, | ||
const vec4_t & | rhs, | ||
vec5_t & | x, | ||
size_t | N, | ||
mem_t & | m | ||
) |
This function solves the system where
is a matrix of the form
* * diag[0] abovediag[0] 0 ..... belowdiag[N-1] * belowdiag[0] diag[1] abovediag[1] ..... * 0 belowdiag[1] diag[2] * 0 0 belowdiag[2] ..... * ... ... * abovediag[N-1] ...
This function solves the following system without the corner elements and then use Sherman-Morrison formula to compensate for them.
Definition at line 338 of file tridiag_base.h.
void o2scl_linalg::solve_cyc_tridiag_sym | ( | const vec_t & | diag, |
const vec2_t & | offdiag, | ||
const vec3_t & | b, | ||
vec4_t & | x, | ||
size_t | N, | ||
mem_t & | m | ||
) |
This function solves the system where
is a matrix of the form
* * diag[0] offdiag[0] 0 ..... offdiag[N-1] * offdiag[0] diag[1] offdiag[1] ..... * 0 offdiag[1] diag[2] * 0 0 offdiag[2] ..... * ... ... * offdiag[N-1] ...
See EngelnMullges96 .
Definition at line 241 of file tridiag_base.h.
void o2scl_linalg::solve_tridiag_nonsym | ( | const vec_t & | diag, |
const vec2_t & | abovediag, | ||
const vec3_t & | belowdiag, | ||
const vec4_t & | rhs, | ||
vec5_t & | x, | ||
size_t | N, | ||
mem_t & | m | ||
) |
This function solves the system where
is a matrix of the form
* * diag[0] abovediag[0] 0 ..... * belowdiag[0] diag[1] abovediag[1] ..... * 0 belowdiag[1] diag[2] * 0 0 belowdiag[2] .....
This function uses plain Gauss elimination, not bothering with the zeroes.
Definition at line 183 of file tridiag_base.h.
void o2scl_linalg::solve_tridiag_sym | ( | const vec_t & | diag, |
const vec2_t & | offdiag, | ||
const vec3_t & | b, | ||
vec4_t & | x, | ||
size_t | N, | ||
mem_t & | m | ||
) |
This function solves the system where
is a matrix of the form
* * diag[0] offdiag[0] 0 ..... * offdiag[0] diag[1] offdiag[1] ..... * 0 offdiag[1] diag[2] * 0 0 offdiag[2] .....
given the N
diagonal elements in diag
, N-1
diagonal elements in offdiag
, and the N
elements b
from the RHS.
See EngelnMullges96 .
Definition at line 118 of file tridiag_base.h.
void o2scl_linalg::SV_decomp | ( | size_t | M, |
size_t | N, | ||
mat_t & | A, | ||
mat2_t & | V, | ||
vec_t & | S, | ||
vec2_t & | work | ||
) |
This factors matrix A
of size (M,N)
into
where U
is a column-orthogonal matrix of size (M,N)
(stored in A
), D
is a diagonal matrix of size (N,N)
(stored in the vector S
of size N
), and V
is a orthogonal matrix of size (N,N)
. The vector work
is a workspace vector of size N
. The matrices U
and V
are constructed so that
This algorithm requres .
Definition at line 73 of file svd_base.h.
void o2scl_linalg::SV_decomp_jacobi | ( | size_t | M, |
size_t | N, | ||
mat_t & | A, | ||
mat2_t & | Q, | ||
vec_t & | S | ||
) |
This factors matrix A
of size (M,N)
into
where U
is a column-orthogonal matrix of size (M,N)
(stored in A
), D
is a diagonal matrix of size (N,N)
(stored in the vector S
of size N
), and V
is a orthogonal matrix of size (N,N)
.
This function computes singular values to higher relative accuracy than SV_decomp() and SV_decomp_mod().
Definition at line 438 of file svd_base.h.
void o2scl_linalg::SV_decomp_mod | ( | size_t | M, |
size_t | N, | ||
mat_t & | A, | ||
mat2_t & | X, | ||
mat3_t & | V, | ||
vec_t & | S, | ||
vec2_t & | work | ||
) |
This factors matrix A
of size (M,N)
into
where U
is a column-orthogonal matrix of size (M,N)
(stored in A
), D
is a diagonal matrix of size (N,N)
(stored in the vector S
of size N
), and V
is a orthogonal matrix of size (N,N)
. The vector work
is a workspace vector of size N
and the matrix X
is a workspace of size (N,N)
.
Definition at line 290 of file svd_base.h.
void o2scl_linalg::SV_solve | ( | size_t | M, |
size_t | N, | ||
mat_t & | U, | ||
mat2_t & | V, | ||
vec_t & | S, | ||
vec2_t & | b, | ||
vec3_t & | x | ||
) |
Solves a linear system using the output of SV_decomp(). Only non-zero singular values are used in computing solution. In the over-determined case, , the system is solved in the least-squares sense.
Definition at line 373 of file svd_base.h.
void o2scl_linalg::svd2 | ( | size_t | M, |
size_t | N, | ||
vec_t & | d, | ||
vec2_t & | f, | ||
mat_t & | U, | ||
mat2_t & | V | ||
) |
The parameter M
is the number of rows in U
and N
is the number of rows in V
. Both U and V assumed to have two columns.
Used in qrstep().
Definition at line 176 of file svdstep_base.h.
void o2scl_linalg::svd2_sub | ( | size_t | M, |
size_t | N, | ||
vec_t & | d, | ||
vec2_t & | f, | ||
mat_t & | U, | ||
mat2_t & | V, | ||
size_t | a | ||
) |
The parameter M
is the number of rows in U
and N
is the number of rows in V
. Both U and V assumed to have two columns.
Used in qrstep_sub().
Definition at line 316 of file svdstep_base.h.
double o2scl_linalg::trailing_eigenvalue | ( | size_t | n, |
const vec_t & | d, | ||
const vec2_t & | f | ||
) |
The parameter n
is the size of the vector d
.
Used in qrstep() and qrstep_sub().
Definition at line 85 of file svdstep_base.h.
Documentation generated with Doxygen. Provided under the
GNU Free Documentation License (see License Information).