LAPACK 3.12.1
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
zgtsvx.f
Go to the documentation of this file.
1*> \brief <b> ZGTSVX computes the solution to system of linear equations A * X = B for GT matrices </b>
2*
3* =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6* http://www.netlib.org/lapack/explore-html/
7*
8*> Download ZGTSVX + dependencies
9*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zgtsvx.f">
10*> [TGZ]</a>
11*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zgtsvx.f">
12*> [ZIP]</a>
13*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zgtsvx.f">
14*> [TXT]</a>
15*
16* Definition:
17* ===========
18*
19* SUBROUTINE ZGTSVX( FACT, TRANS, N, NRHS, DL, D, DU, DLF, DF, DUF,
20* DU2, IPIV, B, LDB, X, LDX, RCOND, FERR, BERR,
21* WORK, RWORK, INFO )
22*
23* .. Scalar Arguments ..
24* CHARACTER FACT, TRANS
25* INTEGER INFO, LDB, LDX, N, NRHS
26* DOUBLE PRECISION RCOND
27* ..
28* .. Array Arguments ..
29* INTEGER IPIV( * )
30* DOUBLE PRECISION BERR( * ), FERR( * ), RWORK( * )
31* COMPLEX*16 B( LDB, * ), D( * ), DF( * ), DL( * ),
32* $ DLF( * ), DU( * ), DU2( * ), DUF( * ),
33* $ WORK( * ), X( LDX, * )
34* ..
35*
36*
37*> \par Purpose:
38* =============
39*>
40*> \verbatim
41*>
42*> ZGTSVX uses the LU factorization to compute the solution to a complex
43*> system of linear equations A * X = B, A**T * X = B, or A**H * X = B,
44*> where A is a tridiagonal matrix of order N and X and B are N-by-NRHS
45*> matrices.
46*>
47*> Error bounds on the solution and a condition estimate are also
48*> provided.
49*> \endverbatim
50*
51*> \par Description:
52* =================
53*>
54*> \verbatim
55*>
56*> The following steps are performed:
57*>
58*> 1. If FACT = 'N', the LU decomposition is used to factor the matrix A
59*> as A = L * U, where L is a product of permutation and unit lower
60*> bidiagonal matrices and U is upper triangular with nonzeros in
61*> only the main diagonal and first two superdiagonals.
62*>
63*> 2. If some U(i,i)=0, so that U is exactly singular, then the routine
64*> returns with INFO = i. Otherwise, the factored form of A is used
65*> to estimate the condition number of the matrix A. If the
66*> reciprocal of the condition number is less than machine precision,
67*> INFO = N+1 is returned as a warning, but the routine still goes on
68*> to solve for X and compute error bounds as described below.
69*>
70*> 3. The system of equations is solved for X using the factored form
71*> of A.
72*>
73*> 4. Iterative refinement is applied to improve the computed solution
74*> matrix and calculate error bounds and backward error estimates
75*> for it.
76*> \endverbatim
77*
78* Arguments:
79* ==========
80*
81*> \param[in] FACT
82*> \verbatim
83*> FACT is CHARACTER*1
84*> Specifies whether or not the factored form of A has been
85*> supplied on entry.
86*> = 'F': DLF, DF, DUF, DU2, and IPIV contain the factored form
87*> of A; DL, D, DU, DLF, DF, DUF, DU2 and IPIV will not
88*> be modified.
89*> = 'N': The matrix will be copied to DLF, DF, and DUF
90*> and factored.
91*> \endverbatim
92*>
93*> \param[in] TRANS
94*> \verbatim
95*> TRANS is CHARACTER*1
96*> Specifies the form of the system of equations:
97*> = 'N': A * X = B (No transpose)
98*> = 'T': A**T * X = B (Transpose)
99*> = 'C': A**H * X = B (Conjugate transpose)
100*> \endverbatim
101*>
102*> \param[in] N
103*> \verbatim
104*> N is INTEGER
105*> The order of the matrix A. N >= 0.
106*> \endverbatim
107*>
108*> \param[in] NRHS
109*> \verbatim
110*> NRHS is INTEGER
111*> The number of right hand sides, i.e., the number of columns
112*> of the matrix B. NRHS >= 0.
113*> \endverbatim
114*>
115*> \param[in] DL
116*> \verbatim
117*> DL is COMPLEX*16 array, dimension (N-1)
118*> The (n-1) subdiagonal elements of A.
119*> \endverbatim
120*>
121*> \param[in] D
122*> \verbatim
123*> D is COMPLEX*16 array, dimension (N)
124*> The n diagonal elements of A.
125*> \endverbatim
126*>
127*> \param[in] DU
128*> \verbatim
129*> DU is COMPLEX*16 array, dimension (N-1)
130*> The (n-1) superdiagonal elements of A.
131*> \endverbatim
132*>
133*> \param[in,out] DLF
134*> \verbatim
135*> DLF is COMPLEX*16 array, dimension (N-1)
136*> If FACT = 'F', then DLF is an input argument and on entry
137*> contains the (n-1) multipliers that define the matrix L from
138*> the LU factorization of A as computed by ZGTTRF.
139*>
140*> If FACT = 'N', then DLF is an output argument and on exit
141*> contains the (n-1) multipliers that define the matrix L from
142*> the LU factorization of A.
143*> \endverbatim
144*>
145*> \param[in,out] DF
146*> \verbatim
147*> DF is COMPLEX*16 array, dimension (N)
148*> If FACT = 'F', then DF is an input argument and on entry
149*> contains the n diagonal elements of the upper triangular
150*> matrix U from the LU factorization of A.
151*>
152*> If FACT = 'N', then DF is an output argument and on exit
153*> contains the n diagonal elements of the upper triangular
154*> matrix U from the LU factorization of A.
155*> \endverbatim
156*>
157*> \param[in,out] DUF
158*> \verbatim
159*> DUF is COMPLEX*16 array, dimension (N-1)
160*> If FACT = 'F', then DUF is an input argument and on entry
161*> contains the (n-1) elements of the first superdiagonal of U.
162*>
163*> If FACT = 'N', then DUF is an output argument and on exit
164*> contains the (n-1) elements of the first superdiagonal of U.
165*> \endverbatim
166*>
167*> \param[in,out] DU2
168*> \verbatim
169*> DU2 is COMPLEX*16 array, dimension (N-2)
170*> If FACT = 'F', then DU2 is an input argument and on entry
171*> contains the (n-2) elements of the second superdiagonal of
172*> U.
173*>
174*> If FACT = 'N', then DU2 is an output argument and on exit
175*> contains the (n-2) elements of the second superdiagonal of
176*> U.
177*> \endverbatim
178*>
179*> \param[in,out] IPIV
180*> \verbatim
181*> IPIV is INTEGER array, dimension (N)
182*> If FACT = 'F', then IPIV is an input argument and on entry
183*> contains the pivot indices from the LU factorization of A as
184*> computed by ZGTTRF.
185*>
186*> If FACT = 'N', then IPIV is an output argument and on exit
187*> contains the pivot indices from the LU factorization of A;
188*> row i of the matrix was interchanged with row IPIV(i).
189*> IPIV(i) will always be either i or i+1; IPIV(i) = i indicates
190*> a row interchange was not required.
191*> \endverbatim
192*>
193*> \param[in] B
194*> \verbatim
195*> B is COMPLEX*16 array, dimension (LDB,NRHS)
196*> The N-by-NRHS right hand side matrix B.
197*> \endverbatim
198*>
199*> \param[in] LDB
200*> \verbatim
201*> LDB is INTEGER
202*> The leading dimension of the array B. LDB >= max(1,N).
203*> \endverbatim
204*>
205*> \param[out] X
206*> \verbatim
207*> X is COMPLEX*16 array, dimension (LDX,NRHS)
208*> If INFO = 0 or INFO = N+1, the N-by-NRHS solution matrix X.
209*> \endverbatim
210*>
211*> \param[in] LDX
212*> \verbatim
213*> LDX is INTEGER
214*> The leading dimension of the array X. LDX >= max(1,N).
215*> \endverbatim
216*>
217*> \param[out] RCOND
218*> \verbatim
219*> RCOND is DOUBLE PRECISION
220*> The estimate of the reciprocal condition number of the matrix
221*> A. If RCOND is less than the machine precision (in
222*> particular, if RCOND = 0), the matrix is singular to working
223*> precision. This condition is indicated by a return code of
224*> INFO > 0.
225*> \endverbatim
226*>
227*> \param[out] FERR
228*> \verbatim
229*> FERR is DOUBLE PRECISION array, dimension (NRHS)
230*> The estimated forward error bound for each solution vector
231*> X(j) (the j-th column of the solution matrix X).
232*> If XTRUE is the true solution corresponding to X(j), FERR(j)
233*> is an estimated upper bound for the magnitude of the largest
234*> element in (X(j) - XTRUE) divided by the magnitude of the
235*> largest element in X(j). The estimate is as reliable as
236*> the estimate for RCOND, and is almost always a slight
237*> overestimate of the true error.
238*> \endverbatim
239*>
240*> \param[out] BERR
241*> \verbatim
242*> BERR is DOUBLE PRECISION array, dimension (NRHS)
243*> The componentwise relative backward error of each solution
244*> vector X(j) (i.e., the smallest relative change in
245*> any element of A or B that makes X(j) an exact solution).
246*> \endverbatim
247*>
248*> \param[out] WORK
249*> \verbatim
250*> WORK is COMPLEX*16 array, dimension (2*N)
251*> \endverbatim
252*>
253*> \param[out] RWORK
254*> \verbatim
255*> RWORK is DOUBLE PRECISION array, dimension (N)
256*> \endverbatim
257*>
258*> \param[out] INFO
259*> \verbatim
260*> INFO is INTEGER
261*> = 0: successful exit
262*> < 0: if INFO = -i, the i-th argument had an illegal value
263*> > 0: if INFO = i, and i is
264*> <= N: U(i,i) is exactly zero. The factorization
265*> has not been completed unless i = N, but the
266*> factor U is exactly singular, so the solution
267*> and error bounds could not be computed.
268*> RCOND = 0 is returned.
269*> = N+1: U is nonsingular, but RCOND is less than machine
270*> precision, meaning that the matrix is singular
271*> to working precision. Nevertheless, the
272*> solution and error bounds are computed because
273*> there are a number of situations where the
274*> computed solution can be more accurate than the
275*> value of RCOND would suggest.
276*> \endverbatim
277*
278* Authors:
279* ========
280*
281*> \author Univ. of Tennessee
282*> \author Univ. of California Berkeley
283*> \author Univ. of Colorado Denver
284*> \author NAG Ltd.
285*
286*> \ingroup gtsvx
287*
288* =====================================================================
289 SUBROUTINE zgtsvx( FACT, TRANS, N, NRHS, DL, D, DU, DLF, DF,
290 $ DUF,
291 $ DU2, IPIV, B, LDB, X, LDX, RCOND, FERR, BERR,
292 $ WORK, RWORK, INFO )
293*
294* -- LAPACK driver routine --
295* -- LAPACK is a software package provided by Univ. of Tennessee, --
296* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
297*
298* .. Scalar Arguments ..
299 CHARACTER FACT, TRANS
300 INTEGER INFO, LDB, LDX, N, NRHS
301 DOUBLE PRECISION RCOND
302* ..
303* .. Array Arguments ..
304 INTEGER IPIV( * )
305 DOUBLE PRECISION BERR( * ), FERR( * ), RWORK( * )
306 COMPLEX*16 B( LDB, * ), D( * ), DF( * ), DL( * ),
307 $ dlf( * ), du( * ), du2( * ), duf( * ),
308 $ work( * ), x( ldx, * )
309* ..
310*
311* =====================================================================
312*
313* .. Parameters ..
314 DOUBLE PRECISION ZERO
315 PARAMETER ( ZERO = 0.0d+0 )
316* ..
317* .. Local Scalars ..
318 LOGICAL NOFACT, NOTRAN
319 CHARACTER NORM
320 DOUBLE PRECISION ANORM
321* ..
322* .. External Functions ..
323 LOGICAL LSAME
324 DOUBLE PRECISION DLAMCH, ZLANGT
325 EXTERNAL LSAME, DLAMCH, ZLANGT
326* ..
327* .. External Subroutines ..
328 EXTERNAL xerbla, zcopy, zgtcon, zgtrfs, zgttrf,
329 $ zgttrs,
330 $ zlacpy
331* ..
332* .. Intrinsic Functions ..
333 INTRINSIC max
334* ..
335* .. Executable Statements ..
336*
337 info = 0
338 nofact = lsame( fact, 'N' )
339 notran = lsame( trans, 'N' )
340 IF( .NOT.nofact .AND. .NOT.lsame( fact, 'F' ) ) THEN
341 info = -1
342 ELSE IF( .NOT.notran .AND. .NOT.lsame( trans, 'T' ) .AND. .NOT.
343 $ lsame( trans, 'C' ) ) THEN
344 info = -2
345 ELSE IF( n.LT.0 ) THEN
346 info = -3
347 ELSE IF( nrhs.LT.0 ) THEN
348 info = -4
349 ELSE IF( ldb.LT.max( 1, n ) ) THEN
350 info = -14
351 ELSE IF( ldx.LT.max( 1, n ) ) THEN
352 info = -16
353 END IF
354 IF( info.NE.0 ) THEN
355 CALL xerbla( 'ZGTSVX', -info )
356 RETURN
357 END IF
358*
359 IF( nofact ) THEN
360*
361* Compute the LU factorization of A.
362*
363 CALL zcopy( n, d, 1, df, 1 )
364 IF( n.GT.1 ) THEN
365 CALL zcopy( n-1, dl, 1, dlf, 1 )
366 CALL zcopy( n-1, du, 1, duf, 1 )
367 END IF
368 CALL zgttrf( n, dlf, df, duf, du2, ipiv, info )
369*
370* Return if INFO is non-zero.
371*
372 IF( info.GT.0 )THEN
373 rcond = zero
374 RETURN
375 END IF
376 END IF
377*
378* Compute the norm of the matrix A.
379*
380 IF( notran ) THEN
381 norm = '1'
382 ELSE
383 norm = 'I'
384 END IF
385 anorm = zlangt( norm, n, dl, d, du )
386*
387* Compute the reciprocal of the condition number of A.
388*
389 CALL zgtcon( norm, n, dlf, df, duf, du2, ipiv, anorm, rcond,
390 $ work,
391 $ info )
392*
393* Compute the solution vectors X.
394*
395 CALL zlacpy( 'Full', n, nrhs, b, ldb, x, ldx )
396 CALL zgttrs( trans, n, nrhs, dlf, df, duf, du2, ipiv, x, ldx,
397 $ info )
398*
399* Use iterative refinement to improve the computed solutions and
400* compute error bounds and backward error estimates for them.
401*
402 CALL zgtrfs( trans, n, nrhs, dl, d, du, dlf, df, duf, du2,
403 $ ipiv,
404 $ b, ldb, x, ldx, ferr, berr, work, rwork, info )
405*
406* Set INFO = N+1 if the matrix is singular to working precision.
407*
408 IF( rcond.LT.dlamch( 'Epsilon' ) )
409 $ info = n + 1
410*
411 RETURN
412*
413* End of ZGTSVX
414*
415 END
subroutine xerbla(srname, info)
Definition cblat2.f:3285
subroutine zcopy(n, zx, incx, zy, incy)
ZCOPY
Definition zcopy.f:81
subroutine zgtcon(norm, n, dl, d, du, du2, ipiv, anorm, rcond, work, info)
ZGTCON
Definition zgtcon.f:139
subroutine zgtrfs(trans, n, nrhs, dl, d, du, dlf, df, duf, du2, ipiv, b, ldb, x, ldx, ferr, berr, work, rwork, info)
ZGTRFS
Definition zgtrfs.f:209
subroutine zgtsvx(fact, trans, n, nrhs, dl, d, du, dlf, df, duf, du2, ipiv, b, ldb, x, ldx, rcond, ferr, berr, work, rwork, info)
ZGTSVX computes the solution to system of linear equations A * X = B for GT matrices
Definition zgtsvx.f:293
subroutine zgttrf(n, dl, d, du, du2, ipiv, info)
ZGTTRF
Definition zgttrf.f:122
subroutine zgttrs(trans, n, nrhs, dl, d, du, du2, ipiv, b, ldb, info)
ZGTTRS
Definition zgttrs.f:137
subroutine zlacpy(uplo, m, n, a, lda, b, ldb)
ZLACPY copies all or part of one two-dimensional array to another.
Definition zlacpy.f:101