?gtsvx
?gtsvx
Computes the solution to the real or complex system of linear equations with a tridiagonal coefficient matrix A and multiple right-hand sides, and provides error bounds on the solution.
Syntax
lapack_int LAPACKE_sgtsvx
(
int
matrix_layout
,
char
fact
,
char
trans
,
lapack_int
n
,
lapack_int
nrhs
,
const float*
dl
,
const float*
d
,
const float*
du
,
float*
dlf
,
float*
df
,
float*
duf
,
float*
du2
,
lapack_int*
ipiv
,
const float*
b
,
lapack_int
ldb
,
float*
x
,
lapack_int
ldx
,
float*
rcond
,
float*
ferr
,
float*
berr
);
lapack_int LAPACKE_dgtsvx
(
int
matrix_layout
,
char
fact
,
char
trans
,
lapack_int
n
,
lapack_int
nrhs
,
const double*
dl
,
const double*
d
,
const double*
du
,
double*
dlf
,
double*
df
,
double*
duf
,
double*
du2
,
lapack_int*
ipiv
,
const double*
b
,
lapack_int
ldb
,
double*
x
,
lapack_int
ldx
,
double*
rcond
,
double*
ferr
,
double*
berr
);
lapack_int LAPACKE_cgtsvx
(
int
matrix_layout
,
char
fact
,
char
trans
,
lapack_int
n
,
lapack_int
nrhs
,
const lapack_complex_float*
dl
,
const lapack_complex_float*
d
,
const lapack_complex_float*
du
,
lapack_complex_float*
dlf
,
lapack_complex_float*
df
,
lapack_complex_float*
duf
,
lapack_complex_float*
du2
,
lapack_int*
ipiv
,
const lapack_complex_float*
b
,
lapack_int
ldb
,
lapack_complex_float*
x
,
lapack_int
ldx
,
float*
rcond
,
float*
ferr
,
float*
berr
);
lapack_int LAPACKE_zgtsvx
(
int
matrix_layout
,
char
fact
,
char
trans
,
lapack_int
n
,
lapack_int
nrhs
,
const lapack_complex_double*
dl
,
const lapack_complex_double*
d
,
const lapack_complex_double*
du
,
lapack_complex_double*
dlf
,
lapack_complex_double*
df
,
lapack_complex_double*
duf
,
lapack_complex_double*
du2
,
lapack_int*
ipiv
,
const lapack_complex_double*
b
,
lapack_int
ldb
,
lapack_complex_double*
x
,
lapack_int
ldx
,
double*
rcond
,
double*
ferr
,
double*
berr
);
Include Files
- mkl.h
Description
The routine uses the ,
LU
factorization to compute the solution to a real or 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
, the columns of matrix B
are individual right-hand sides, and the columns of X
are the corresponding solutions.Error bounds on the solution and a condition estimate are also provided.
The routine
?gtsvx
performs the following steps:- If, thefact='N'LUdecomposition is used to factor the matrixAas, whereA=L*ULis a product of permutation and unit lower bidiagonal matrices andUis an upper triangular matrix with nonzeroes in only the main diagonal and first two superdiagonals.
- If some= 0, so thatUi,iUis exactly singular, then the routine returns with. Otherwise, the factored form ofinfo=iAis used to estimate the condition number of the matrixA. If the reciprocal of the condition number is less than machine precision,is returned as a warning, but the routine still goes on to solve forinfo=n+ 1Xand compute error bounds as described below.
- The system of equations is solved forXusing the factored form ofA.
- Iterative refinement is applied to improve the computed solution matrix and calculate error bounds and backward error estimates for it.
Input Parameters
- matrix_layout
- Specifies whether matrix storage layout is row major (LAPACK_ROW_MAJOR) or column major (LAPACK_COL_MAJOR).
- fact
- Must be'F'or'N'.Specifies whether or not the factored form of the matrixAhas been supplied on entry.If: on entry,fact='F'dlf,df,duf,du2, andipivcontain the factored form ofA; arraysdl,d,du,dlf,df,duf,du2, andipivwill not be modified.If, the matrixfact='N'Awill be copied todlf,df, anddufand factored.
- trans
- Must be'N','T', or'C'.Specifies the form of the system of equations:If, the system has the formtrans='N'(No transpose).A*X=BIf, the system has the formtrans='T'AT*X=B(Transpose).If, the system has the formtrans='C'AH*X=B(Conjugate transpose).
- n
- The number of linear equations, the order of the matrixA;n≥0.
- nrhs
- The number of right hand sides, the number of columns of the matricesBandX;nrhs≥0.
- dl,d,du,dlf,df,duf,du2,b
- Arrays:dl, size (n-1), contains the subdiagonal elements ofA.d, size (n), contains the diagonal elements ofA.du, size (n-1), contains the superdiagonal elements ofA.dlf, size (n-1). If, thenfact='F'dlfis an input argument and on entry contains the (n-1) multipliers that define the matrixLfrom theLUfactorization ofAas computed by?gttrf.df, size (n). If, thenfact='F'dfis an input argument and on entry contains thendiagonal elements of the upper triangular matrixUfrom theLUfactorization ofA.duf, size (n-1). If, thenfact='F'dufis an input argument and on entry contains the (n-1) elements of the first superdiagonal ofU.du2, size (n-2). If, thenfact='F'du2is an input argument and on entry contains the (n-2) elements of the second superdiagonal ofU.b, size max(contains the right-hand side matrixldb*nrhs) for column major layout and max(ldb*n) for row major layout,B.
- ldb
- The leading dimension ofb;.ldb≥max(1,n) for column major layout andldb≥nrhsfor row major layout
- ldx
- The leading dimension ofx;.ldx≥max(1,n) for column major layout andldx≥nrhsfor row major layout
- ipiv
- Array, size at leastmax(1,. Ifn), thenfact='F'ipivis an input argument and on entry contains the pivot indices, as returned by?gttrf.
Output Parameters
- x
- Array, sizemax(1,.ldx*nrhs) for column major layout and max(1,ldx*n) for row major layoutIforinfo= 0, the arrayinfo=n+1xcontains the solution matrixX.
- dlf
- If, thenfact='N'dlfis an output argument and on exit contains the(multipliers that define the matrixn-1)Lfrom theLUfactorization of A.
- df
- If, thenfact='N'dfis an output argument and on exit contains thendiagonal elements of the upper triangular matrixUfrom theLUfactorization ofA.
- duf
- If, thenfact='N'dufis an output argument and on exit contains the(elements of the first superdiagonal ofn-1)U.
- du2
- If, thenfact='N'du2is an output argument and on exit contains the(elements of the second superdiagonal ofn-2)U.
- ipiv
- The arrayipivis an output argument ifand, on exit, contains the pivot indices from the factorizationfact='N'; rowA=L*Uiof the matrix was interchanged with rowipiv[i-1]. The value ofipiv[i-1] will always beiori+1;ipiv[i-1]=iindicates a row interchange was not required.
- rcond
- An estimate of the reciprocal condition number of the matrixA. Ifrcondis less than the machine precision (in particular, ifrcond=0), the matrix is singular to working precision. This condition is indicated by a return code ofinfo>0.
- ferr
- Array, size at leastmax(1,. Contains the estimated forward error bound for each solution vectornrhs)(thexjj-th column of the solution matrixX). Ifxtrueis the true solution corresponding to,xjis an estimated upper bound for the magnitude of the largest element inferr[j-1]divided by the magnitude of the largest element inx-jxtrue. The estimate is as reliable as the estimate forxjrcond, and is almost always a slight overestimate of the true error.
- berr
- Array, size at leastmax(1,. Contains the component-wise relative backward error for each solution vectornrhs), that is, the smallest relative change in any element ofxjAorBthat makesan exact solution.xj
Return Values
This function returns a value
info
.If , the execution is successful.
info
= 0If , parameter
info
= -i
i
had an illegal value. If , and , then is exactly zero. The factorization has not been completed unless , but the factor , and , then
info
= i
i
≤
n
U
i
, i
i
= n
U
is exactly singular, so the solution and error bounds could not be computed; rcond
= 0 is returned. If info
= i
i
= n
+ 1U
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.