?hesvxx
?hesvxx
Uses extra precise iterative refinement to compute the solution to the system of linear equations with a Hermitian indefinite coefficient matrix A applying the diagonal pivoting factorization.
Syntax
lapack_int LAPACKE_chesvxx
(
int
matrix_layout
,
char
fact
,
char
uplo
,
lapack_int
n
,
lapack_int
nrhs
,
lapack_complex_float*
a
,
lapack_int
lda
,
lapack_complex_float*
af
,
lapack_int
ldaf
,
lapack_int*
ipiv
,
char*
equed
,
float*
s
,
lapack_complex_float*
b
,
lapack_int
ldb
,
lapack_complex_float*
x
,
lapack_int
ldx
,
float*
rcond
,
float*
rpvgrw
,
float*
berr
,
lapack_int
n_err_bnds
,
float*
err_bnds_norm
,
float*
err_bnds_comp
,
lapack_int
nparams
,
const float*
params
);
lapack_int LAPACKE_zhesvxx
(
int
matrix_layout
,
char
fact
,
char
uplo
,
lapack_int
n
,
lapack_int
nrhs
,
lapack_complex_double*
a
,
lapack_int
lda
,
lapack_complex_double*
af
,
lapack_int
ldaf
,
lapack_int*
ipiv
,
char*
equed
,
double*
s
,
lapack_complex_double*
b
,
lapack_int
ldb
,
lapack_complex_double*
x
,
lapack_int
ldx
,
double*
rcond
,
double*
rpvgrw
,
double*
berr
,
lapack_int
n_err_bnds
,
double*
err_bnds_norm
,
double*
err_bnds_comp
,
lapack_int
nparams
,
const double*
params
);
Include Files
- mkl.h
Description
The routine uses the diagonal pivoting factorization to compute the solution to a complex/double complex system of linear equations , where
A*X
= B
A
is an n
-by-n
Hermitian matrix, the columns of matrix B
are individual right-hand sides, and the columns of X
are the corresponding solutions.Both normwise and maximum componentwise error bounds are also provided on request. The routine returns a solution with a small guaranteed error (
O(eps)
, where eps
is the working machine precision) unless the matrix is very ill-conditioned, in which case a warning is returned. Relevant condition numbers are also calculated and returned.The routine accepts user-provided factorizations and equilibration factors; see definitions of the
fact
and equed
options. Solving with refinement and using a factorization from a previous call of the routine also produces a solution with O(eps)
errors or warnings but that may not be true for general user-provided factorizations and equilibration factors if they differ from what the routine would itself produce.The routine
?hesvxx
performs the following steps:- If, scaling factors are computed to equilibrate the system:fact='E'diag(s)*A*diag(s) *inv(diag(s))*X=diag(s)*BWhether or not the system will be equilibrated depends on the scaling of the matrixA, but if equilibration is used,Ais overwritten byanddiag(s)*A*diag(s)Bby.diag(s)*B
- Iforfact='N', the LU decomposition is used to factor the matrix'E'A(after equilibration if) asfact='E', ifA=U*D*UT,uplo='U'or, ifA=L*D*LT,uplo='L'whereUorLis a product of permutation and unit upper (lower) triangular matrices, andDis a symmetric and block diagonal with 1-by-1 and 2-by-2 diagonal blocks.
- If some, so thatD(i,i)=0Dis exactly singular, the routine returns with. Otherwise, the factored form ofinfo=iAis used to estimate the condition number of the matrixA(see thercondparameter). If the reciprocal of the condition number is less than machine precision, the routine still goes on to solve forXand compute error bounds.
- The system of equations is solved forXusing the factored form ofA.
- By default, unlessparams[0]is set to zero, the routine applies iterative refinement to get a small error and error bounds. Refinement calculates the residual to at least twice the working precision.
- If equilibration was used, the matrixXis premultiplied byso that it solves the original system before equilibration.diag(r)
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','N', or'E'.Specifies whether or not the factored form of the matrixAis supplied on entry, and if not, whether the matrixAshould be equilibrated before it is factored.If, on entry,fact='F'afandipivcontain the factored form ofA. Ifequedis not'N', the matrixAhas been equilibrated with scaling factors given bys. Parametersa,af, andipivare not modified.If, the matrixfact='N'Awill be copied toafand factored.If, the matrixfact='E'Awill be equilibrated, if necessary, copied toafand factored.
- uplo
- Must be'U'or'L'.Indicates whether the upper or lower triangular part ofAis stored:If, the upper triangle ofuplo='U'Ais stored.If, the lower triangle ofuplo='L'Ais stored.
- 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.
- a,af,b
- Arrays:a(size max(,lda*n))af(size max(,ldaf*n))b, (size max(.ldb*nrhs) for column major layout and max(ldb*n) for row major layout),The arrayacontains the Hermitian matrixAas specified byuplo. Ifuplo='U', the leadingn-by-nupper triangular part ofacontains the upper triangular part of the matrixAand the strictly lower triangular part ofais not referenced. Ifuplo='L', the leadingn-by-nlower triangular part ofacontains the lower triangular part of the matrixAand the strictly upper triangular part ofais not referenced.The arrayafis an input argument if. It contains the block diagonal matrixfact='F'Dand the multipliers used to obtain the factorUandLfrom the factorizationorA=U*D*UTas computed byA=L*D*LT?hetrf.The arraybcontains the matrixBwhose columns are the right-hand sides for the systems of equations.
- lda
- The leading dimension of the arraya;.lda≥max(1,n)
- ldaf
- The leading dimension of the arrayaf;.ldaf≥max(1,n)
- ipiv
- Array, size at leastmax(1,. The arrayn)ipivis an input argument if. It contains details of the interchanges and the block structure offact='F'Das determined by?sytrf.If> 0, rows and columnsipiv[k-1]kandwere interchanged andipiv[k-1]D) is a 1-by-1 diagonal block.k,kIfanduplo='U',ipiv[i] =ipiv[i- 1] =m< 0Dhas a 2-by-2 diagonal block in rows and columnsiandi+ 1, and thei-th row and column ofAwere interchanged with them-th row and column.Ifanduplo='L',ipiv[i] =ipiv[i- 1] =m< 0Dhas a 2-by-2 diagonal block in rows and columnsiandi+ 1, and the (i+ 1)-st row and column ofAwere interchanged with them-th row and column.
- equed
- Must be'N'or'Y'.equedis an input argument if. It specifies the form of equilibration that was done:fact='F'If, no equilibration was done (always true ifequed='N').fact='N'if, both row and column equilibration was done, that is,equed='Y'Ahas been replaced by.diag(s)*A*diag(s)
- s
- Array, size (n). The arrayscontains the scale factors forA. If,equed='Y'Ais multiplied on the left and right bydiag(.s)This array is an input argument ifonly; otherwise it is an output argument.fact='F'Ifandfact='F'equed='Y', each element ofsmust be positive.Each element ofsshould be a power of the radix to ensure a reliable solution and error estimates. Scaling by powers of the radix does not cause rounding errors unless the result underflows or overflows. Rounding errors during scaling lead to refining with a matrix that is not equivalent to the input matrix, producing error estimates that may not be reliable.
- ldb
- The leading dimension of the arrayb;.ldb≥max(1,n) for column major layout andldb≥nrhsfor row major layout
- ldx
- The leading dimension of the output arrayx;.ldx≥max(1,n) for column major layout andldx≥nrhsfor row major layout
- n_err_bnds
- Number of error bounds to return for each right hand side and each type (normwise or componentwise). Seeerr_bnds_normanderr_bnds_compdescriptions in theOutput Argumentssection below.
- nparams
- Specifies the number of parameters set inparams. If≤0, theparamsarray is never referenced and default values are used.
- params
- Array, sizemax(1,. Specifies algorithm parameters. If an entry is less than 0.0, that entry is filled with the default value used for that parameter. Only positions up tonparams)nparamsare accessed; defaults are used for higher-numbered parameters. If defaults are acceptable, you can passnparams= 0, which prevents the source code from accessing theparamsargument.: Whether to perform iterative refinement or not. Default: 1.0 (for single precision flavors), 1.0D+0 (for double precision flavors).params[0]
- =0.0
- No refinement is performed and no error bounds are computed.
- =1.0
- Use the extra-precise refinement algorithm.
(Other values are reserved for future use.): Maximum number of residual computations allowed for refinement.params[1]- Default
- 10
- Aggressive
- Set to 100 to permit convergence using approximate factorizations or factorizations other thanLU. If the factorization uses a technique other than Gaussian elimination, the guarantees inerr_bnds_normanderr_bnds_compmay no longer be trustworthy.
: Flag determining if the code will attempt to find a solution with a small componentwise relative error in the double-precision algorithm. Positive is true, 0.0 is false. Default: 1.0 (attempt componentwise convergence).params[2]
Output Parameters
- x
- Array, sizemax(1,.ldx*nrhs) for column major layout and max(1,ldx*n) for row major layoutIf, the arrayinfo= 0xcontains the solutionn-by-nrhsmatrixXto the original system of equations. Note thatAandBare modified on exit if, and the solution to the equilibrated system is:equed≠'N'inv(.diag(s))*X
- a
- Ifandfact='E'equed='Y', overwritten by.diag(s)*A*diag(s)
- af
- If,fact='N'afis an output argument and on exit returns the block diagonal matrixDand the multipliers used to obtain the factorUorLfrom the factorizationorA=U*D*UT.A=L*D*LT
- b
- Ifequed='N',Bis not modified.Ifequed='Y',Bis overwritten by.diag(s)*B
- s
- This array is an output argument if. Each element of this array is a power of the radix. See the description offact≠'F'sinInput Argumentssection.
- rcond
- Reciprocal scaled condition number. An estimate of the reciprocal Skeel condition number of the matrixAafter equilibration (if done). Ifrcondis less than the machine precision, in particular, ifrcond= 0, the matrix is singular to working precision. Note that the error may still be small even if this number is very small and the matrix appears ill-conditioned.
- rpvgrw
- Contains the reciprocal pivot growth factor:If this is much less than 1, the stability of theLUfactorization of the (equlibrated) matrixAcould be poor. This also means that the solutionX, estimated condition numbers, and error bounds could be unreliable. If factorization fails with0 <, this parameter contains the reciprocal pivot growth factor for the leadinginfo≤ninfocolumns ofA.
- berr
- Array, size at leastmax(1,.nrhs)Contains the component-wise relative backward error for each solution vectorx, that is, the smallest relative change in any element ofjAorBthat makesxan exact solution.j
- err_bnds_norm
- Array of size. For each right-hand side, contains information about various error bounds and condition numbers corresponding to the normwise relative errornrhs*n_err_bnds, which is defined as follows:Normwise relative error in thei-th solution vectorThe array is indexed by the type of error information as described below. There are currently up to three pieces of information returned.
- err=1
- "Trust/don't trust" boolean. Trust the answer if the reciprocal condition number is less than the thresholdsqrt(n)*