LAPACK++
LAPACK C++ API
General matrix: LU: tridiagonal

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)
 

Detailed Description

Function Documentation

◆ gtsv()

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>.

Parameters
[in]nThe order of the matrix A. n >= 0.
[in]nrhsThe number of right hand sides, i.e., the number of columns of the matrix B. nrhs >= 0.
[in,out]DLThe 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]DThe 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]DUThe 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]BThe 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]ldbThe leading dimension of the array B. ldb >= max(1,n).
Returns
= 0: successful exit
> 0: if return value = i, U(i,i) is exactly zero, and the solution has not been computed. The factorization has not been completed unless i = n.

◆ gtsvx()

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>.

Parameters
[in]factSpecifies whether or not the factored form of A has been supplied on entry.
  • lapack::Factored::Factored: DLF, DF, DUF, DU2, and ipiv contain the factored form of A; DL, D, DU, DLF, DF, DUF, DU2 and ipiv will not be modified.
  • lapack::Factored::NotFactored: The matrix will be copied to DLF, DF, and DUF and factored.
[in]transSpecifies the form of the system of equations:
  • lapack::Op::NoTrans: \(A X = B\) (No transpose)
  • lapack::Op::Trans: \(A^T X = B\) (Transpose)
  • lapack::Op::ConjTrans: \(A^H X = B\) (Conjugate transpose)
[in]nThe order of the matrix A. n >= 0.
[in]nrhsThe number of right hand sides, i.e., the number of columns of the matrix B. nrhs >= 0.
[in]DLThe vector DL of length n-1. The (n-1) subdiagonal elements of A.
[in]DThe vector D of length n. The n diagonal elements of A.
[in]DUThe vector DU of length n-1. The (n-1) superdiagonal elements of A.
[in,out]DLFThe vector DLF of length n-1.
  • If fact = Factored, then DLF is an input argument and on entry contains the (n-1) multipliers that define the matrix L from the LU factorization of A as computed by lapack::gttrf.
  • If fact = NotFactored, then DLF is an output argument and on exit contains the (n-1) multipliers that define the matrix L from the LU factorization of A.
[in,out]DFThe vector DF of length n.
  • If fact = Factored, then DF is an input argument and on entry contains the n diagonal elements of the upper triangular matrix U from the LU factorization of A.
  • If fact = NotFactored, then DF is an output argument and on exit contains the n diagonal elements of the upper triangular matrix U from the LU factorization of A.
[in,out]DUFThe vector DUF of length n-1.
  • If fact = Factored, then DUF is an input argument and on entry contains the (n-1) elements of the first superdiagonal of U.
  • If fact = NotFactored, then DUF is an output argument and on exit contains the (n-1) elements of the first superdiagonal of U.
[in,out]DU2The vector DU2 of length n-2.
  • If fact = Factored, then DU2 is an input argument and on entry contains the (n-2) elements of the second superdiagonal of U.
  • If fact = NotFactored, then DU2 is an output argument and on exit contains the (n-2) elements of the second superdiagonal of U.
[in,out]ipivThe vector ipiv of length n.
  • If fact = Factored, then ipiv is an input argument and on entry contains the pivot indices from the LU factorization of A as computed by lapack::gttrf.
  • If fact = NotFactored, then ipiv is an output argument and on exit contains the pivot indices from the LU factorization of A; row i of the matrix was interchanged with row ipiv(i). ipiv(i) will always be either i or i+1; ipiv(i) = i indicates a row interchange was not required.
[in]BThe n-by-nrhs matrix B, stored in an ldb-by-nrhs array. The n-by-nrhs right hand side matrix B.
[in]ldbThe leading dimension of the array B. ldb >= max(1,n).
[out]XThe 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]ldxThe leading dimension of the array X. ldx >= max(1,n).
[out]rcondThe 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]ferrThe 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]berrThe 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).
Returns
= 0: successful exit
> 0 and <= n: if return value = i, U(i,i) is exactly zero. The factorization has not been completed unless i = n, but the factor U is exactly singular, so the solution and error bounds could not be computed. rcond = 0 is returned.
= n+1: U is nonsingular, but rcond is less than machine precision, meaning that the matrix is singular to working precision. Nevertheless, the solution and error bounds are computed because there are a number of situations where the computed solution can be more accurate than the value of rcond would suggest.