LAPACK++
LAPACK C++ API
|
Functions | |
int64_t | lapack::gtsv (int64_t n, int64_t nrhs, double *DL, double *D, double *DU, double *B, int64_t ldb) |
int64_t | lapack::gtsv (int64_t n, int64_t nrhs, float *DL, float *D, float *DU, float *B, int64_t ldb) |
int64_t | lapack::gtsv (int64_t n, int64_t nrhs, std::complex< double > *DL, std::complex< double > *D, std::complex< double > *DU, std::complex< double > *B, int64_t ldb) |
Solves the equation. More... | |
int64_t | lapack::gtsv (int64_t n, int64_t nrhs, std::complex< float > *DL, std::complex< float > *D, std::complex< float > *DU, std::complex< float > *B, int64_t ldb) |
int64_t | lapack::gtsvx (lapack::Factored fact, lapack::Op trans, int64_t n, int64_t nrhs, double const *DL, double const *D, double const *DU, double *DLF, double *DF, double *DUF, double *DU2, int64_t *ipiv, double const *B, int64_t ldb, double *X, int64_t ldx, double *rcond, double *ferr, double *berr) |
int64_t | lapack::gtsvx (lapack::Factored fact, lapack::Op trans, int64_t n, int64_t nrhs, float const *DL, float const *D, float const *DU, float *DLF, float *DF, float *DUF, float *DU2, int64_t *ipiv, float const *B, int64_t ldb, float *X, int64_t ldx, float *rcond, float *ferr, float *berr) |
int64_t | lapack::gtsvx (lapack::Factored fact, lapack::Op trans, int64_t n, int64_t nrhs, std::complex< double > const *DL, std::complex< double > const *D, std::complex< double > const *DU, std::complex< double > *DLF, std::complex< double > *DF, std::complex< double > *DUF, std::complex< double > *DU2, int64_t *ipiv, std::complex< double > const *B, int64_t ldb, std::complex< double > *X, int64_t ldx, double *rcond, double *ferr, double *berr) |
Uses the LU factorization to compute the solution to a complex system of linear equations. More... | |
int64_t | lapack::gtsvx (lapack::Factored fact, lapack::Op trans, int64_t n, int64_t nrhs, std::complex< float > const *DL, std::complex< float > const *D, std::complex< float > const *DU, std::complex< float > *DLF, std::complex< float > *DF, std::complex< float > *DUF, std::complex< float > *DU2, int64_t *ipiv, std::complex< float > const *B, int64_t ldb, std::complex< float > *X, int64_t ldx, float *rcond, float *ferr, float *berr) |
int64_t lapack::gtsv | ( | int64_t | n, |
int64_t | nrhs, | ||
std::complex< double > * | DL, | ||
std::complex< double > * | D, | ||
std::complex< double > * | DU, | ||
std::complex< double > * | B, | ||
int64_t | ldb | ||
) |
Solves the equation.
\[ A X = B, \]
where A is an n-by-n tridiagonal matrix, by Gaussian elimination with partial pivoting.
Note that the equation \(A^T X = B\) may be solved by interchanging the order of the arguments DU and DL.
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] | DL | The vector DL of length n-1. On entry, DL must contain the (n-1) subdiagonal elements of A. On exit, DL is overwritten by the (n-2) elements of the second superdiagonal of the upper triangular matrix U from the LU factorization of A, in DL(1), ..., DL(n-2). |
[in,out] | D | The vector D of length n. On entry, D must contain the diagonal elements of A. On exit, D is overwritten by the n diagonal elements of U. |
[in,out] | DU | The vector DU of length n-1. On entry, DU must contain the (n-1) superdiagonal elements of A. On exit, DU is overwritten by the (n-1) elements of the first superdiagonal of U. |
[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::gtsvx | ( | lapack::Factored | fact, |
lapack::Op | trans, | ||
int64_t | n, | ||
int64_t | nrhs, | ||
std::complex< double > const * | DL, | ||
std::complex< double > const * | D, | ||
std::complex< double > const * | DU, | ||
std::complex< double > * | DLF, | ||
std::complex< double > * | DF, | ||
std::complex< double > * | DUF, | ||
std::complex< double > * | DU2, | ||
int64_t * | ipiv, | ||
std::complex< double > const * | B, | ||
int64_t | ldb, | ||
std::complex< double > * | X, | ||
int64_t | ldx, | ||
double * | rcond, | ||
double * | ferr, | ||
double * | berr | ||
) |
Uses the LU factorization to compute the solution to a complex system of linear equations.
\[ A X = B, \]
\[ A^T X = B, \]
or
\[ A^H X = B, \]
where A is a tridiagonal matrix of order n 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 A has been supplied on entry.
|
[in] | trans | Specifies the form of the system of equations:
|
[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] | DL | The vector DL of length n-1. The (n-1) subdiagonal elements of A. |
[in] | D | The vector D of length n. The n diagonal elements of A. |
[in] | DU | The vector DU of length n-1. The (n-1) superdiagonal elements of A. |
[in,out] | DLF | The vector DLF of length n-1.
|
[in,out] | DF | The vector DF of length n.
|
[in,out] | DUF | The vector DUF of length n-1.
|
[in,out] | DU2 | The vector DU2 of length n-2.
|
[in,out] | ipiv | The vector ipiv of length n.
|
[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 estimate of 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 estimated 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). The estimate is as reliable as the estimate for rcond, and is almost always a slight overestimate of the true error. |
[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). |