LAPACK++
2022.05.00
LAPACK C++ API
|
Functions | |
int64_t | lapack::ptsv (int64_t n, int64_t nrhs, double *D, double *E, double *B, int64_t ldb) |
int64_t | lapack::ptsv (int64_t n, int64_t nrhs, double *D, std::complex< double > *E, std::complex< double > *B, int64_t ldb) |
Computes the solution to a system of linear equations \(A X = B,\) where A is an n-by-n Hermitian positive definite tridiagonal matrix, and X and B are n-by-nrhs matrices. More... | |
int64_t | lapack::ptsv (int64_t n, int64_t nrhs, float *D, float *E, float *B, int64_t ldb) |
int64_t | lapack::ptsv (int64_t n, int64_t nrhs, float *D, std::complex< float > *E, std::complex< float > *B, int64_t ldb) |
int64_t | lapack::ptsvx (lapack::Factored fact, int64_t n, int64_t nrhs, double const *D, std::complex< double > const *E, double *DF, std::complex< double > *EF, std::complex< double > const *B, int64_t ldb, std::complex< double > *X, int64_t ldx, double *rcond, double *ferr, double *berr) |
Uses the factorization \(A = L D L^H\) to compute the solution to a system of linear equations \(A X = B,\) where A is an n-by-n Hermitian positive definite tridiagonal matrix and X and B are n-by-nrhs matrices. More... | |
int64_t | lapack::ptsvx (lapack::Factored fact, int64_t n, int64_t nrhs, double const *D, double const *E, double *DF, double *EF, double const *B, int64_t ldb, double *X, int64_t ldx, double *rcond, double *ferr, double *berr) |
int64_t | lapack::ptsvx (lapack::Factored fact, int64_t n, int64_t nrhs, float const *D, float const *E, float *DF, float *EF, float const *B, int64_t ldb, float *X, int64_t ldx, float *rcond, float *ferr, float *berr) |
int64_t | lapack::ptsvx (lapack::Factored fact, int64_t n, int64_t nrhs, float const *D, std::complex< float > const *E, float *DF, std::complex< float > *EF, std::complex< float > const *B, int64_t ldb, std::complex< float > *X, int64_t ldx, float *rcond, float *ferr, float *berr) |
int64_t lapack::ptsv | ( | int64_t | n, |
int64_t | nrhs, | ||
double * | D, | ||
std::complex< double > * | E, | ||
std::complex< double > * | B, | ||
int64_t | ldb | ||
) |
Computes the solution to a system of linear equations \(A X = B,\) where A is an n-by-n Hermitian positive definite tridiagonal matrix, and X and B are n-by-nrhs matrices.
A is factored as \(A = L D L^H,\) and the factored form of A is then used to solve the system of equations.
Overloaded versions are available for float
, double
, std::complex<float>
, and std::complex<double>
.
[in] | n | The order of the matrix A. n >= 0. |
[in] | nrhs | The number of right hand sides, i.e., the number of columns of the matrix B. nrhs >= 0. |
[in,out] | D | The vector D of length n. On entry, the n diagonal elements of the tridiagonal matrix A. On exit, the n diagonal elements of the diagonal matrix D from the factorization \(A = L D L^H.\) |
[in,out] | E | The vector E of length n-1. On entry, the (n-1) subdiagonal elements of the tridiagonal matrix A. On exit, the (n-1) subdiagonal elements of the unit bidiagonal factor L from the \(L D L^H\) factorization of A. E can also be regarded as the superdiagonal of the unit bidiagonal factor U from the \(U^H D U\) factorization of A. |
[in,out] | B | The n-by-nrhs matrix B, stored in an ldb-by-nrhs array. On entry, the n-by-nrhs right hand side matrix B. On successful exit, the n-by-nrhs solution matrix X. |
[in] | ldb | The leading dimension of the array B. ldb >= max(1,n). |
int64_t lapack::ptsvx | ( | lapack::Factored | fact, |
int64_t | n, | ||
int64_t | nrhs, | ||
double const * | D, | ||
std::complex< double > const * | E, | ||
double * | DF, | ||
std::complex< double > * | EF, | ||
std::complex< double > const * | B, | ||
int64_t | ldb, | ||
std::complex< double > * | X, | ||
int64_t | ldx, | ||
double * | rcond, | ||
double * | ferr, | ||
double * | berr | ||
) |
Uses the factorization \(A = L D L^H\) to compute the solution to a system of linear equations \(A X = B,\) where A is an n-by-n Hermitian positive definite tridiagonal matrix and X and B are n-by-nrhs matrices.
Error bounds on the solution and a condition estimate are also provided.
Overloaded versions are available for float
, double
, std::complex<float>
, and std::complex<double>
.
[in] | fact | Specifies whether or not the factored form of the matrix A is supplied on entry.
|
[in] | n | The order of the matrix A. n >= 0. |
[in] | nrhs | The number of right hand sides, i.e., the number of columns of the matrices B and X. nrhs >= 0. |
[in] | D | The vector D of length n. The n diagonal elements of the tridiagonal matrix A. |
[in] | E | The vector E of length n-1. The (n-1) subdiagonal elements of the tridiagonal matrix A. |
[in,out] | DF | The vector DF of length n.
|
[in,out] | EF | The vector EF of length n-1.
|
[in] | B | The n-by-nrhs matrix B, stored in an ldb-by-nrhs array. The n-by-nrhs right hand side matrix B. |
[in] | ldb | The leading dimension of the array B. ldb >= max(1,n). |
[out] | X | The n-by-nrhs matrix X, stored in an ldx-by-nrhs array. If successful or return value = n+1, the n-by-nrhs solution matrix X. |
[in] | ldx | The leading dimension of the array X. ldx >= max(1,n). |
[out] | rcond | The reciprocal condition number of the matrix A. If rcond is less than the machine precision (in particular, if rcond = 0), the matrix is singular to working precision. This condition is indicated by a return code of return value > 0. |
[out] | ferr | The vector ferr of length nrhs. The forward error bound for each solution vector X(j) (the j-th column of the solution matrix X). If XTRUE is the true solution corresponding to X(j), ferr(j) is an estimated upper bound for the magnitude of the largest element in (X(j) - XTRUE) divided by the magnitude of the largest element in X(j). |
[out] | berr | The vector berr of length nrhs. The componentwise relative backward error of each solution vector X(j) (i.e., the smallest relative change in any element of A or B that makes X(j) an exact solution). |