LAPACK  3.6.0
LAPACK: Linear Algebra PACKage
Collaboration diagram for complex16:

Functions

subroutine zgtcon (NORM, N, DL, D, DU, DU2, IPIV, ANORM, RCOND, WORK, INFO)
 ZGTCON More...
 
subroutine zgtrfs (TRANS, N, NRHS, DL, D, DU, DLF, DF, DUF, DU2, IPIV, B, LDB, X, LDX, FERR, BERR, WORK, RWORK, INFO)
 ZGTRFS More...
 
subroutine zgttrf (N, DL, D, DU, DU2, IPIV, INFO)
 ZGTTRF More...
 
subroutine zgttrs (TRANS, N, NRHS, DL, D, DU, DU2, IPIV, B, LDB, INFO)
 ZGTTRS More...
 
subroutine zgtts2 (ITRANS, N, NRHS, DL, D, DU, DU2, IPIV, B, LDB)
 ZGTTS2 solves a system of linear equations with a tridiagonal matrix using the LU factorization computed by sgttrf. More...
 

Detailed Description

This is the group of complex16 computational functions for GT matrices

Function Documentation

subroutine zgtcon ( character  NORM,
integer  N,
complex*16, dimension( * )  DL,
complex*16, dimension( * )  D,
complex*16, dimension( * )  DU,
complex*16, dimension( * )  DU2,
integer, dimension( * )  IPIV,
double precision  ANORM,
double precision  RCOND,
complex*16, dimension( * )  WORK,
integer  INFO 
)

ZGTCON

Download ZGTCON + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 ZGTCON estimates the reciprocal of the condition number of a complex
 tridiagonal matrix A using the LU factorization as computed by
 ZGTTRF.

 An estimate is obtained for norm(inv(A)), and the reciprocal of the
 condition number is computed as RCOND = 1 / (ANORM * norm(inv(A))).
Parameters
[in]NORM
          NORM is CHARACTER*1
          Specifies whether the 1-norm condition number or the
          infinity-norm condition number is required:
          = '1' or 'O':  1-norm;
          = 'I':         Infinity-norm.
[in]N
          N is INTEGER
          The order of the matrix A.  N >= 0.
[in]DL
          DL is COMPLEX*16 array, dimension (N-1)
          The (n-1) multipliers that define the matrix L from the
          LU factorization of A as computed by ZGTTRF.
[in]D
          D is COMPLEX*16 array, dimension (N)
          The n diagonal elements of the upper triangular matrix U from
          the LU factorization of A.
[in]DU
          DU is COMPLEX*16 array, dimension (N-1)
          The (n-1) elements of the first superdiagonal of U.
[in]DU2
          DU2 is COMPLEX*16 array, dimension (N-2)
          The (n-2) elements of the second superdiagonal of U.
[in]IPIV
          IPIV is INTEGER array, dimension (N)
          The pivot indices; for 1 <= i <= n, 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]ANORM
          ANORM is DOUBLE PRECISION
          If NORM = '1' or 'O', the 1-norm of the original matrix A.
          If NORM = 'I', the infinity-norm of the original matrix A.
[out]RCOND
          RCOND is DOUBLE PRECISION
          The reciprocal of the condition number of the matrix A,
          computed as RCOND = 1/(ANORM * AINVNM), where AINVNM is an
          estimate of the 1-norm of inv(A) computed in this routine.
[out]WORK
          WORK is COMPLEX*16 array, dimension (2*N)
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
September 2012

Definition at line 143 of file zgtcon.f.

143 *
144 * -- LAPACK computational routine (version 3.4.2) --
145 * -- LAPACK is a software package provided by Univ. of Tennessee, --
146 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
147 * September 2012
148 *
149 * .. Scalar Arguments ..
150  CHARACTER norm
151  INTEGER info, n
152  DOUBLE PRECISION anorm, rcond
153 * ..
154 * .. Array Arguments ..
155  INTEGER ipiv( * )
156  COMPLEX*16 d( * ), dl( * ), du( * ), du2( * ), work( * )
157 * ..
158 *
159 * =====================================================================
160 *
161 * .. Parameters ..
162  DOUBLE PRECISION one, zero
163  parameter( one = 1.0d+0, zero = 0.0d+0 )
164 * ..
165 * .. Local Scalars ..
166  LOGICAL onenrm
167  INTEGER i, kase, kase1
168  DOUBLE PRECISION ainvnm
169 * ..
170 * .. Local Arrays ..
171  INTEGER isave( 3 )
172 * ..
173 * .. External Functions ..
174  LOGICAL lsame
175  EXTERNAL lsame
176 * ..
177 * .. External Subroutines ..
178  EXTERNAL xerbla, zgttrs, zlacn2
179 * ..
180 * .. Intrinsic Functions ..
181  INTRINSIC dcmplx
182 * ..
183 * .. Executable Statements ..
184 *
185 * Test the input arguments.
186 *
187  info = 0
188  onenrm = norm.EQ.'1' .OR. lsame( norm, 'O' )
189  IF( .NOT.onenrm .AND. .NOT.lsame( norm, 'I' ) ) THEN
190  info = -1
191  ELSE IF( n.LT.0 ) THEN
192  info = -2
193  ELSE IF( anorm.LT.zero ) THEN
194  info = -8
195  END IF
196  IF( info.NE.0 ) THEN
197  CALL xerbla( 'ZGTCON', -info )
198  RETURN
199  END IF
200 *
201 * Quick return if possible
202 *
203  rcond = zero
204  IF( n.EQ.0 ) THEN
205  rcond = one
206  RETURN
207  ELSE IF( anorm.EQ.zero ) THEN
208  RETURN
209  END IF
210 *
211 * Check that D(1:N) is non-zero.
212 *
213  DO 10 i = 1, n
214  IF( d( i ).EQ.dcmplx( zero ) )
215  $ RETURN
216  10 CONTINUE
217 *
218  ainvnm = zero
219  IF( onenrm ) THEN
220  kase1 = 1
221  ELSE
222  kase1 = 2
223  END IF
224  kase = 0
225  20 CONTINUE
226  CALL zlacn2( n, work( n+1 ), work, ainvnm, kase, isave )
227  IF( kase.NE.0 ) THEN
228  IF( kase.EQ.kase1 ) THEN
229 *
230 * Multiply by inv(U)*inv(L).
231 *
232  CALL zgttrs( 'No transpose', n, 1, dl, d, du, du2, ipiv,
233  $ work, n, info )
234  ELSE
235 *
236 * Multiply by inv(L**H)*inv(U**H).
237 *
238  CALL zgttrs( 'Conjugate transpose', n, 1, dl, d, du, du2,
239  $ ipiv, work, n, info )
240  END IF
241  GO TO 20
242  END IF
243 *
244 * Compute the estimate of the reciprocal condition number.
245 *
246  IF( ainvnm.NE.zero )
247  $ rcond = ( one / ainvnm ) / anorm
248 *
249  RETURN
250 *
251 * End of ZGTCON
252 *
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine zlacn2(N, V, X, EST, KASE, ISAVE)
ZLACN2 estimates the 1-norm of a square matrix, using reverse communication for evaluating matrix-vec...
Definition: zlacn2.f:135
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine zgttrs(TRANS, N, NRHS, DL, D, DU, DU2, IPIV, B, LDB, INFO)
ZGTTRS
Definition: zgttrs.f:140

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine zgtrfs ( character  TRANS,
integer  N,
integer  NRHS,
complex*16, dimension( * )  DL,
complex*16, dimension( * )  D,
complex*16, dimension( * )  DU,
complex*16, dimension( * )  DLF,
complex*16, dimension( * )  DF,
complex*16, dimension( * )  DUF,
complex*16, dimension( * )  DU2,
integer, dimension( * )  IPIV,
complex*16, dimension( ldb, * )  B,
integer  LDB,
complex*16, dimension( ldx, * )  X,
integer  LDX,
double precision, dimension( * )  FERR,
double precision, dimension( * )  BERR,
complex*16, dimension( * )  WORK,
double precision, dimension( * )  RWORK,
integer  INFO 
)

ZGTRFS

Download ZGTRFS + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 ZGTRFS improves the computed solution to a system of linear
 equations when the coefficient matrix is tridiagonal, and provides
 error bounds and backward error estimates for the solution.
Parameters
[in]TRANS
          TRANS is CHARACTER*1
          Specifies the form of the system of equations:
          = 'N':  A * X = B     (No transpose)
          = 'T':  A**T * X = B  (Transpose)
          = 'C':  A**H * X = B  (Conjugate transpose)
[in]N
          N is INTEGER
          The order of the matrix A.  N >= 0.
[in]NRHS
          NRHS is INTEGER
          The number of right hand sides, i.e., the number of columns
          of the matrix B.  NRHS >= 0.
[in]DL
          DL is COMPLEX*16 array, dimension (N-1)
          The (n-1) subdiagonal elements of A.
[in]D
          D is COMPLEX*16 array, dimension (N)
          The diagonal elements of A.
[in]DU
          DU is COMPLEX*16 array, dimension (N-1)
          The (n-1) superdiagonal elements of A.
[in]DLF
          DLF is COMPLEX*16 array, dimension (N-1)
          The (n-1) multipliers that define the matrix L from the
          LU factorization of A as computed by ZGTTRF.
[in]DF
          DF is COMPLEX*16 array, dimension (N)
          The n diagonal elements of the upper triangular matrix U from
          the LU factorization of A.
[in]DUF
          DUF is COMPLEX*16 array, dimension (N-1)
          The (n-1) elements of the first superdiagonal of U.
[in]DU2
          DU2 is COMPLEX*16 array, dimension (N-2)
          The (n-2) elements of the second superdiagonal of U.
[in]IPIV
          IPIV is INTEGER array, dimension (N)
          The pivot indices; for 1 <= i <= n, 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]B
          B is COMPLEX*16 array, dimension (LDB,NRHS)
          The right hand side matrix B.
[in]LDB
          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).
[in,out]X
          X is COMPLEX*16 array, dimension (LDX,NRHS)
          On entry, the solution matrix X, as computed by ZGTTRS.
          On exit, the improved solution matrix X.
[in]LDX
          LDX is INTEGER
          The leading dimension of the array X.  LDX >= max(1,N).
[out]FERR
          FERR is DOUBLE PRECISION array, dimension (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
          BERR is DOUBLE PRECISION array, dimension (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).
[out]WORK
          WORK is COMPLEX*16 array, dimension (2*N)
[out]RWORK
          RWORK is DOUBLE PRECISION array, dimension (N)
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
Internal Parameters:
  ITMAX is the maximum number of steps of iterative refinement.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
September 2012

Definition at line 212 of file zgtrfs.f.

212 *
213 * -- LAPACK computational routine (version 3.4.2) --
214 * -- LAPACK is a software package provided by Univ. of Tennessee, --
215 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
216 * September 2012
217 *
218 * .. Scalar Arguments ..
219  CHARACTER trans
220  INTEGER info, ldb, ldx, n, nrhs
221 * ..
222 * .. Array Arguments ..
223  INTEGER ipiv( * )
224  DOUBLE PRECISION berr( * ), ferr( * ), rwork( * )
225  COMPLEX*16 b( ldb, * ), d( * ), df( * ), dl( * ),
226  $ dlf( * ), du( * ), du2( * ), duf( * ),
227  $ work( * ), x( ldx, * )
228 * ..
229 *
230 * =====================================================================
231 *
232 * .. Parameters ..
233  INTEGER itmax
234  parameter( itmax = 5 )
235  DOUBLE PRECISION zero, one
236  parameter( zero = 0.0d+0, one = 1.0d+0 )
237  DOUBLE PRECISION two
238  parameter( two = 2.0d+0 )
239  DOUBLE PRECISION three
240  parameter( three = 3.0d+0 )
241 * ..
242 * .. Local Scalars ..
243  LOGICAL notran
244  CHARACTER transn, transt
245  INTEGER count, i, j, kase, nz
246  DOUBLE PRECISION eps, lstres, s, safe1, safe2, safmin
247  COMPLEX*16 zdum
248 * ..
249 * .. Local Arrays ..
250  INTEGER isave( 3 )
251 * ..
252 * .. External Subroutines ..
253  EXTERNAL xerbla, zaxpy, zcopy, zgttrs, zlacn2, zlagtm
254 * ..
255 * .. Intrinsic Functions ..
256  INTRINSIC abs, dble, dcmplx, dimag, max
257 * ..
258 * .. External Functions ..
259  LOGICAL lsame
260  DOUBLE PRECISION dlamch
261  EXTERNAL lsame, dlamch
262 * ..
263 * .. Statement Functions ..
264  DOUBLE PRECISION cabs1
265 * ..
266 * .. Statement Function definitions ..
267  cabs1( zdum ) = abs( dble( zdum ) ) + abs( dimag( zdum ) )
268 * ..
269 * .. Executable Statements ..
270 *
271 * Test the input parameters.
272 *
273  info = 0
274  notran = lsame( trans, 'N' )
275  IF( .NOT.notran .AND. .NOT.lsame( trans, 'T' ) .AND. .NOT.
276  $ lsame( trans, 'C' ) ) THEN
277  info = -1
278  ELSE IF( n.LT.0 ) THEN
279  info = -2
280  ELSE IF( nrhs.LT.0 ) THEN
281  info = -3
282  ELSE IF( ldb.LT.max( 1, n ) ) THEN
283  info = -13
284  ELSE IF( ldx.LT.max( 1, n ) ) THEN
285  info = -15
286  END IF
287  IF( info.NE.0 ) THEN
288  CALL xerbla( 'ZGTRFS', -info )
289  RETURN
290  END IF
291 *
292 * Quick return if possible
293 *
294  IF( n.EQ.0 .OR. nrhs.EQ.0 ) THEN
295  DO 10 j = 1, nrhs
296  ferr( j ) = zero
297  berr( j ) = zero
298  10 CONTINUE
299  RETURN
300  END IF
301 *
302  IF( notran ) THEN
303  transn = 'N'
304  transt = 'C'
305  ELSE
306  transn = 'C'
307  transt = 'N'
308  END IF
309 *
310 * NZ = maximum number of nonzero elements in each row of A, plus 1
311 *
312  nz = 4
313  eps = dlamch( 'Epsilon' )
314  safmin = dlamch( 'Safe minimum' )
315  safe1 = nz*safmin
316  safe2 = safe1 / eps
317 *
318 * Do for each right hand side
319 *
320  DO 110 j = 1, nrhs
321 *
322  count = 1
323  lstres = three
324  20 CONTINUE
325 *
326 * Loop until stopping criterion is satisfied.
327 *
328 * Compute residual R = B - op(A) * X,
329 * where op(A) = A, A**T, or A**H, depending on TRANS.
330 *
331  CALL zcopy( n, b( 1, j ), 1, work, 1 )
332  CALL zlagtm( trans, n, 1, -one, dl, d, du, x( 1, j ), ldx, one,
333  $ work, n )
334 *
335 * Compute abs(op(A))*abs(x) + abs(b) for use in the backward
336 * error bound.
337 *
338  IF( notran ) THEN
339  IF( n.EQ.1 ) THEN
340  rwork( 1 ) = cabs1( b( 1, j ) ) +
341  $ cabs1( d( 1 ) )*cabs1( x( 1, j ) )
342  ELSE
343  rwork( 1 ) = cabs1( b( 1, j ) ) +
344  $ cabs1( d( 1 ) )*cabs1( x( 1, j ) ) +
345  $ cabs1( du( 1 ) )*cabs1( x( 2, j ) )
346  DO 30 i = 2, n - 1
347  rwork( i ) = cabs1( b( i, j ) ) +
348  $ cabs1( dl( i-1 ) )*cabs1( x( i-1, j ) ) +
349  $ cabs1( d( i ) )*cabs1( x( i, j ) ) +
350  $ cabs1( du( i ) )*cabs1( x( i+1, j ) )
351  30 CONTINUE
352  rwork( n ) = cabs1( b( n, j ) ) +
353  $ cabs1( dl( n-1 ) )*cabs1( x( n-1, j ) ) +
354  $ cabs1( d( n ) )*cabs1( x( n, j ) )
355  END IF
356  ELSE
357  IF( n.EQ.1 ) THEN
358  rwork( 1 ) = cabs1( b( 1, j ) ) +
359  $ cabs1( d( 1 ) )*cabs1( x( 1, j ) )
360  ELSE
361  rwork( 1 ) = cabs1( b( 1, j ) ) +
362  $ cabs1( d( 1 ) )*cabs1( x( 1, j ) ) +
363  $ cabs1( dl( 1 ) )*cabs1( x( 2, j ) )
364  DO 40 i = 2, n - 1
365  rwork( i ) = cabs1( b( i, j ) ) +
366  $ cabs1( du( i-1 ) )*cabs1( x( i-1, j ) ) +
367  $ cabs1( d( i ) )*cabs1( x( i, j ) ) +
368  $ cabs1( dl( i ) )*cabs1( x( i+1, j ) )
369  40 CONTINUE
370  rwork( n ) = cabs1( b( n, j ) ) +
371  $ cabs1( du( n-1 ) )*cabs1( x( n-1, j ) ) +
372  $ cabs1( d( n ) )*cabs1( x( n, j ) )
373  END IF
374  END IF
375 *
376 * Compute componentwise relative backward error from formula
377 *
378 * max(i) ( abs(R(i)) / ( abs(op(A))*abs(X) + abs(B) )(i) )
379 *
380 * where abs(Z) is the componentwise absolute value of the matrix
381 * or vector Z. If the i-th component of the denominator is less
382 * than SAFE2, then SAFE1 is added to the i-th components of the
383 * numerator and denominator before dividing.
384 *
385  s = zero
386  DO 50 i = 1, n
387  IF( rwork( i ).GT.safe2 ) THEN
388  s = max( s, cabs1( work( i ) ) / rwork( i ) )
389  ELSE
390  s = max( s, ( cabs1( work( i ) )+safe1 ) /
391  $ ( rwork( i )+safe1 ) )
392  END IF
393  50 CONTINUE
394  berr( j ) = s
395 *
396 * Test stopping criterion. Continue iterating if
397 * 1) The residual BERR(J) is larger than machine epsilon, and
398 * 2) BERR(J) decreased by at least a factor of 2 during the
399 * last iteration, and
400 * 3) At most ITMAX iterations tried.
401 *
402  IF( berr( j ).GT.eps .AND. two*berr( j ).LE.lstres .AND.
403  $ count.LE.itmax ) THEN
404 *
405 * Update solution and try again.
406 *
407  CALL zgttrs( trans, n, 1, dlf, df, duf, du2, ipiv, work, n,
408  $ info )
409  CALL zaxpy( n, dcmplx( one ), work, 1, x( 1, j ), 1 )
410  lstres = berr( j )
411  count = count + 1
412  GO TO 20
413  END IF
414 *
415 * Bound error from formula
416 *
417 * norm(X - XTRUE) / norm(X) .le. FERR =
418 * norm( abs(inv(op(A)))*
419 * ( abs(R) + NZ*EPS*( abs(op(A))*abs(X)+abs(B) ))) / norm(X)
420 *
421 * where
422 * norm(Z) is the magnitude of the largest component of Z
423 * inv(op(A)) is the inverse of op(A)
424 * abs(Z) is the componentwise absolute value of the matrix or
425 * vector Z
426 * NZ is the maximum number of nonzeros in any row of A, plus 1
427 * EPS is machine epsilon
428 *
429 * The i-th component of abs(R)+NZ*EPS*(abs(op(A))*abs(X)+abs(B))
430 * is incremented by SAFE1 if the i-th component of
431 * abs(op(A))*abs(X) + abs(B) is less than SAFE2.
432 *
433 * Use ZLACN2 to estimate the infinity-norm of the matrix
434 * inv(op(A)) * diag(W),
435 * where W = abs(R) + NZ*EPS*( abs(op(A))*abs(X)+abs(B) )))
436 *
437  DO 60 i = 1, n
438  IF( rwork( i ).GT.safe2 ) THEN
439  rwork( i ) = cabs1( work( i ) ) + nz*eps*rwork( i )
440  ELSE
441  rwork( i ) = cabs1( work( i ) ) + nz*eps*rwork( i ) +
442  $ safe1
443  END IF
444  60 CONTINUE
445 *
446  kase = 0
447  70 CONTINUE
448  CALL zlacn2( n, work( n+1 ), work, ferr( j ), kase, isave )
449  IF( kase.NE.0 ) THEN
450  IF( kase.EQ.1 ) THEN
451 *
452 * Multiply by diag(W)*inv(op(A)**H).
453 *
454  CALL zgttrs( transt, n, 1, dlf, df, duf, du2, ipiv, work,
455  $ n, info )
456  DO 80 i = 1, n
457  work( i ) = rwork( i )*work( i )
458  80 CONTINUE
459  ELSE
460 *
461 * Multiply by inv(op(A))*diag(W).
462 *
463  DO 90 i = 1, n
464  work( i ) = rwork( i )*work( i )
465  90 CONTINUE
466  CALL zgttrs( transn, n, 1, dlf, df, duf, du2, ipiv, work,
467  $ n, info )
468  END IF
469  GO TO 70
470  END IF
471 *
472 * Normalize error.
473 *
474  lstres = zero
475  DO 100 i = 1, n
476  lstres = max( lstres, cabs1( x( i, j ) ) )
477  100 CONTINUE
478  IF( lstres.NE.zero )
479  $ ferr( j ) = ferr( j ) / lstres
480 *
481  110 CONTINUE
482 *
483  RETURN
484 *
485 * End of ZGTRFS
486 *
subroutine zlagtm(TRANS, N, NRHS, ALPHA, DL, D, DU, X, LDX, BETA, B, LDB)
ZLAGTM performs a matrix-matrix product of the form C = αAB+βC, where A is a tridiagonal matrix...
Definition: zlagtm.f:147
subroutine zcopy(N, ZX, INCX, ZY, INCY)
ZCOPY
Definition: zcopy.f:52
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine zlacn2(N, V, X, EST, KASE, ISAVE)
ZLACN2 estimates the 1-norm of a square matrix, using reverse communication for evaluating matrix-vec...
Definition: zlacn2.f:135
double precision function dlamch(CMACH)
DLAMCH
Definition: dlamch.f:65
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine zgttrs(TRANS, N, NRHS, DL, D, DU, DU2, IPIV, B, LDB, INFO)
ZGTTRS
Definition: zgttrs.f:140
subroutine zaxpy(N, ZA, ZX, INCX, ZY, INCY)
ZAXPY
Definition: zaxpy.f:53

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine zgttrf ( integer  N,
complex*16, dimension( * )  DL,
complex*16, dimension( * )  D,
complex*16, dimension( * )  DU,
complex*16, dimension( * )  DU2,
integer, dimension( * )  IPIV,
integer  INFO 
)

ZGTTRF

Download ZGTTRF + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 ZGTTRF computes an LU factorization of a complex tridiagonal matrix A
 using elimination with partial pivoting and row interchanges.

 The factorization has the form
    A = L * U
 where L is a product of permutation and unit lower bidiagonal
 matrices and U is upper triangular with nonzeros in only the main
 diagonal and first two superdiagonals.
Parameters
[in]N
          N is INTEGER
          The order of the matrix A.
[in,out]DL
          DL is COMPLEX*16 array, dimension (N-1)
          On entry, DL must contain the (n-1) sub-diagonal elements of
          A.

          On exit, DL is overwritten by the (n-1) multipliers that
          define the matrix L from the LU factorization of A.
[in,out]D
          D is COMPLEX*16 array, dimension (N)
          On entry, D must contain the diagonal elements of A.

          On exit, D is overwritten by the n diagonal elements of the
          upper triangular matrix U from the LU factorization of A.
[in,out]DU
          DU is COMPLEX*16 array, dimension (N-1)
          On entry, DU must contain the (n-1) super-diagonal elements
          of A.

          On exit, DU is overwritten by the (n-1) elements of the first
          super-diagonal of U.
[out]DU2
          DU2 is COMPLEX*16 array, dimension (N-2)
          On exit, DU2 is overwritten by the (n-2) elements of the
          second super-diagonal of U.
[out]IPIV
          IPIV is INTEGER array, dimension (N)
          The pivot indices; for 1 <= i <= n, 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.
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -k, the k-th argument had an illegal value
          > 0:  if INFO = k, U(k,k) is exactly zero. The factorization
                has been completed, but the factor U is exactly
                singular, and division by zero will occur if it is used
                to solve a system of equations.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
September 2012

Definition at line 126 of file zgttrf.f.

126 *
127 * -- LAPACK computational routine (version 3.4.2) --
128 * -- LAPACK is a software package provided by Univ. of Tennessee, --
129 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
130 * September 2012
131 *
132 * .. Scalar Arguments ..
133  INTEGER info, n
134 * ..
135 * .. Array Arguments ..
136  INTEGER ipiv( * )
137  COMPLEX*16 d( * ), dl( * ), du( * ), du2( * )
138 * ..
139 *
140 * =====================================================================
141 *
142 * .. Parameters ..
143  DOUBLE PRECISION zero
144  parameter( zero = 0.0d+0 )
145 * ..
146 * .. Local Scalars ..
147  INTEGER i
148  COMPLEX*16 fact, temp, zdum
149 * ..
150 * .. External Subroutines ..
151  EXTERNAL xerbla
152 * ..
153 * .. Intrinsic Functions ..
154  INTRINSIC abs, dble, dimag
155 * ..
156 * .. Statement Functions ..
157  DOUBLE PRECISION cabs1
158 * ..
159 * .. Statement Function definitions ..
160  cabs1( zdum ) = abs( dble( zdum ) ) + abs( dimag( zdum ) )
161 * ..
162 * .. Executable Statements ..
163 *
164  info = 0
165  IF( n.LT.0 ) THEN
166  info = -1
167  CALL xerbla( 'ZGTTRF', -info )
168  RETURN
169  END IF
170 *
171 * Quick return if possible
172 *
173  IF( n.EQ.0 )
174  $ RETURN
175 *
176 * Initialize IPIV(i) = i and DU2(i) = 0
177 *
178  DO 10 i = 1, n
179  ipiv( i ) = i
180  10 CONTINUE
181  DO 20 i = 1, n - 2
182  du2( i ) = zero
183  20 CONTINUE
184 *
185  DO 30 i = 1, n - 2
186  IF( cabs1( d( i ) ).GE.cabs1( dl( i ) ) ) THEN
187 *
188 * No row interchange required, eliminate DL(I)
189 *
190  IF( cabs1( d( i ) ).NE.zero ) THEN
191  fact = dl( i ) / d( i )
192  dl( i ) = fact
193  d( i+1 ) = d( i+1 ) - fact*du( i )
194  END IF
195  ELSE
196 *
197 * Interchange rows I and I+1, eliminate DL(I)
198 *
199  fact = d( i ) / dl( i )
200  d( i ) = dl( i )
201  dl( i ) = fact
202  temp = du( i )
203  du( i ) = d( i+1 )
204  d( i+1 ) = temp - fact*d( i+1 )
205  du2( i ) = du( i+1 )
206  du( i+1 ) = -fact*du( i+1 )
207  ipiv( i ) = i + 1
208  END IF
209  30 CONTINUE
210  IF( n.GT.1 ) THEN
211  i = n - 1
212  IF( cabs1( d( i ) ).GE.cabs1( dl( i ) ) ) THEN
213  IF( cabs1( d( i ) ).NE.zero ) THEN
214  fact = dl( i ) / d( i )
215  dl( i ) = fact
216  d( i+1 ) = d( i+1 ) - fact*du( i )
217  END IF
218  ELSE
219  fact = d( i ) / dl( i )
220  d( i ) = dl( i )
221  dl( i ) = fact
222  temp = du( i )
223  du( i ) = d( i+1 )
224  d( i+1 ) = temp - fact*d( i+1 )
225  ipiv( i ) = i + 1
226  END IF
227  END IF
228 *
229 * Check for a zero on the diagonal of U.
230 *
231  DO 40 i = 1, n
232  IF( cabs1( d( i ) ).EQ.zero ) THEN
233  info = i
234  GO TO 50
235  END IF
236  40 CONTINUE
237  50 CONTINUE
238 *
239  RETURN
240 *
241 * End of ZGTTRF
242 *
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine zgttrs ( character  TRANS,
integer  N,
integer  NRHS,
complex*16, dimension( * )  DL,
complex*16, dimension( * )  D,
complex*16, dimension( * )  DU,
complex*16, dimension( * )  DU2,
integer, dimension( * )  IPIV,
complex*16, dimension( ldb, * )  B,
integer  LDB,
integer  INFO 
)

ZGTTRS

Download ZGTTRS + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 ZGTTRS solves one of the systems of equations
    A * X = B,  A**T * X = B,  or  A**H * X = B,
 with a tridiagonal matrix A using the LU factorization computed
 by ZGTTRF.
Parameters
[in]TRANS
          TRANS is CHARACTER*1
          Specifies the form of the system of equations.
          = 'N':  A * X = B     (No transpose)
          = 'T':  A**T * X = B  (Transpose)
          = 'C':  A**H * X = B  (Conjugate transpose)
[in]N
          N is INTEGER
          The order of the matrix A.
[in]NRHS
          NRHS is INTEGER
          The number of right hand sides, i.e., the number of columns
          of the matrix B.  NRHS >= 0.
[in]DL
          DL is COMPLEX*16 array, dimension (N-1)
          The (n-1) multipliers that define the matrix L from the
          LU factorization of A.
[in]D
          D is COMPLEX*16 array, dimension (N)
          The n diagonal elements of the upper triangular matrix U from
          the LU factorization of A.
[in]DU
          DU is COMPLEX*16 array, dimension (N-1)
          The (n-1) elements of the first super-diagonal of U.
[in]DU2
          DU2 is COMPLEX*16 array, dimension (N-2)
          The (n-2) elements of the second super-diagonal of U.
[in]IPIV
          IPIV is INTEGER array, dimension (N)
          The pivot indices; for 1 <= i <= n, 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,out]B
          B is COMPLEX*16 array, dimension (LDB,NRHS)
          On entry, the matrix of right hand side vectors B.
          On exit, B is overwritten by the solution vectors X.
[in]LDB
          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -k, the k-th argument had an illegal value
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
September 2012

Definition at line 140 of file zgttrs.f.

140 *
141 * -- LAPACK computational routine (version 3.4.2) --
142 * -- LAPACK is a software package provided by Univ. of Tennessee, --
143 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
144 * September 2012
145 *
146 * .. Scalar Arguments ..
147  CHARACTER trans
148  INTEGER info, ldb, n, nrhs
149 * ..
150 * .. Array Arguments ..
151  INTEGER ipiv( * )
152  COMPLEX*16 b( ldb, * ), d( * ), dl( * ), du( * ), du2( * )
153 * ..
154 *
155 * =====================================================================
156 *
157 * .. Local Scalars ..
158  LOGICAL notran
159  INTEGER itrans, j, jb, nb
160 * ..
161 * .. External Functions ..
162  INTEGER ilaenv
163  EXTERNAL ilaenv
164 * ..
165 * .. External Subroutines ..
166  EXTERNAL xerbla, zgtts2
167 * ..
168 * .. Intrinsic Functions ..
169  INTRINSIC max, min
170 * ..
171 * .. Executable Statements ..
172 *
173  info = 0
174  notran = ( trans.EQ.'N' .OR. trans.EQ.'n' )
175  IF( .NOT.notran .AND. .NOT.( trans.EQ.'T' .OR. trans.EQ.
176  $ 't' ) .AND. .NOT.( trans.EQ.'C' .OR. trans.EQ.'c' ) ) THEN
177  info = -1
178  ELSE IF( n.LT.0 ) THEN
179  info = -2
180  ELSE IF( nrhs.LT.0 ) THEN
181  info = -3
182  ELSE IF( ldb.LT.max( n, 1 ) ) THEN
183  info = -10
184  END IF
185  IF( info.NE.0 ) THEN
186  CALL xerbla( 'ZGTTRS', -info )
187  RETURN
188  END IF
189 *
190 * Quick return if possible
191 *
192  IF( n.EQ.0 .OR. nrhs.EQ.0 )
193  $ RETURN
194 *
195 * Decode TRANS
196 *
197  IF( notran ) THEN
198  itrans = 0
199  ELSE IF( trans.EQ.'T' .OR. trans.EQ.'t' ) THEN
200  itrans = 1
201  ELSE
202  itrans = 2
203  END IF
204 *
205 * Determine the number of right-hand sides to solve at a time.
206 *
207  IF( nrhs.EQ.1 ) THEN
208  nb = 1
209  ELSE
210  nb = max( 1, ilaenv( 1, 'ZGTTRS', trans, n, nrhs, -1, -1 ) )
211  END IF
212 *
213  IF( nb.GE.nrhs ) THEN
214  CALL zgtts2( itrans, n, nrhs, dl, d, du, du2, ipiv, b, ldb )
215  ELSE
216  DO 10 j = 1, nrhs, nb
217  jb = min( nrhs-j+1, nb )
218  CALL zgtts2( itrans, n, jb, dl, d, du, du2, ipiv, b( 1, j ),
219  $ ldb )
220  10 CONTINUE
221  END IF
222 *
223 * End of ZGTTRS
224 *
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine zgtts2(ITRANS, N, NRHS, DL, D, DU, DU2, IPIV, B, LDB)
ZGTTS2 solves a system of linear equations with a tridiagonal matrix using the LU factorization compu...
Definition: zgtts2.f:130
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
Definition: tstiee.f:83

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine zgtts2 ( integer  ITRANS,
integer  N,
integer  NRHS,
complex*16, dimension( * )  DL,
complex*16, dimension( * )  D,
complex*16, dimension( * )  DU,
complex*16, dimension( * )  DU2,
integer, dimension( * )  IPIV,
complex*16, dimension( ldb, * )  B,
integer  LDB 
)

ZGTTS2 solves a system of linear equations with a tridiagonal matrix using the LU factorization computed by sgttrf.

Download ZGTTS2 + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 ZGTTS2 solves one of the systems of equations
    A * X = B,  A**T * X = B,  or  A**H * X = B,
 with a tridiagonal matrix A using the LU factorization computed
 by ZGTTRF.
Parameters
[in]ITRANS
          ITRANS is INTEGER
          Specifies the form of the system of equations.
          = 0:  A * X = B     (No transpose)
          = 1:  A**T * X = B  (Transpose)
          = 2:  A**H * X = B  (Conjugate transpose)
[in]N
          N is INTEGER
          The order of the matrix A.
[in]NRHS
          NRHS is INTEGER
          The number of right hand sides, i.e., the number of columns
          of the matrix B.  NRHS >= 0.
[in]DL
          DL is COMPLEX*16 array, dimension (N-1)
          The (n-1) multipliers that define the matrix L from the
          LU factorization of A.
[in]D
          D is COMPLEX*16 array, dimension (N)
          The n diagonal elements of the upper triangular matrix U from
          the LU factorization of A.
[in]DU
          DU is COMPLEX*16 array, dimension (N-1)
          The (n-1) elements of the first super-diagonal of U.
[in]DU2
          DU2 is COMPLEX*16 array, dimension (N-2)
          The (n-2) elements of the second super-diagonal of U.
[in]IPIV
          IPIV is INTEGER array, dimension (N)
          The pivot indices; for 1 <= i <= n, 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,out]B
          B is COMPLEX*16 array, dimension (LDB,NRHS)
          On entry, the matrix of right hand side vectors B.
          On exit, B is overwritten by the solution vectors X.
[in]LDB
          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
September 2012

Definition at line 130 of file zgtts2.f.

130 *
131 * -- LAPACK computational routine (version 3.4.2) --
132 * -- LAPACK is a software package provided by Univ. of Tennessee, --
133 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
134 * September 2012
135 *
136 * .. Scalar Arguments ..
137  INTEGER itrans, ldb, n, nrhs
138 * ..
139 * .. Array Arguments ..
140  INTEGER ipiv( * )
141  COMPLEX*16 b( ldb, * ), d( * ), dl( * ), du( * ), du2( * )
142 * ..
143 *
144 * =====================================================================
145 *
146 * .. Local Scalars ..
147  INTEGER i, j
148  COMPLEX*16 temp
149 * ..
150 * .. Intrinsic Functions ..
151  INTRINSIC dconjg
152 * ..
153 * .. Executable Statements ..
154 *
155 * Quick return if possible
156 *
157  IF( n.EQ.0 .OR. nrhs.EQ.0 )
158  $ RETURN
159 *
160  IF( itrans.EQ.0 ) THEN
161 *
162 * Solve A*X = B using the LU factorization of A,
163 * overwriting each right hand side vector with its solution.
164 *
165  IF( nrhs.LE.1 ) THEN
166  j = 1
167  10 CONTINUE
168 *
169 * Solve L*x = b.
170 *
171  DO 20 i = 1, n - 1
172  IF( ipiv( i ).EQ.i ) THEN
173  b( i+1, j ) = b( i+1, j ) - dl( i )*b( i, j )
174  ELSE
175  temp = b( i, j )
176  b( i, j ) = b( i+1, j )
177  b( i+1, j ) = temp - dl( i )*b( i, j )
178  END IF
179  20 CONTINUE
180 *
181 * Solve U*x = b.
182 *
183  b( n, j ) = b( n, j ) / d( n )
184  IF( n.GT.1 )
185  $ b( n-1, j ) = ( b( n-1, j )-du( n-1 )*b( n, j ) ) /
186  $ d( n-1 )
187  DO 30 i = n - 2, 1, -1
188  b( i, j ) = ( b( i, j )-du( i )*b( i+1, j )-du2( i )*
189  $ b( i+2, j ) ) / d( i )
190  30 CONTINUE
191  IF( j.LT.nrhs ) THEN
192  j = j + 1
193  GO TO 10
194  END IF
195  ELSE
196  DO 60 j = 1, nrhs
197 *
198 * Solve L*x = b.
199 *
200  DO 40 i = 1, n - 1
201  IF( ipiv( i ).EQ.i ) THEN
202  b( i+1, j ) = b( i+1, j ) - dl( i )*b( i, j )
203  ELSE
204  temp = b( i, j )
205  b( i, j ) = b( i+1, j )
206  b( i+1, j ) = temp - dl( i )*b( i, j )
207  END IF
208  40 CONTINUE
209 *
210 * Solve U*x = b.
211 *
212  b( n, j ) = b( n, j ) / d( n )
213  IF( n.GT.1 )
214  $ b( n-1, j ) = ( b( n-1, j )-du( n-1 )*b( n, j ) ) /
215  $ d( n-1 )
216  DO 50 i = n - 2, 1, -1
217  b( i, j ) = ( b( i, j )-du( i )*b( i+1, j )-du2( i )*
218  $ b( i+2, j ) ) / d( i )
219  50 CONTINUE
220  60 CONTINUE
221  END IF
222  ELSE IF( itrans.EQ.1 ) THEN
223 *
224 * Solve A**T * X = B.
225 *
226  IF( nrhs.LE.1 ) THEN
227  j = 1
228  70 CONTINUE
229 *
230 * Solve U**T * x = b.
231 *
232  b( 1, j ) = b( 1, j ) / d( 1 )
233  IF( n.GT.1 )
234  $ b( 2, j ) = ( b( 2, j )-du( 1 )*b( 1, j ) ) / d( 2 )
235  DO 80 i = 3, n
236  b( i, j ) = ( b( i, j )-du( i-1 )*b( i-1, j )-du2( i-2 )*
237  $ b( i-2, j ) ) / d( i )
238  80 CONTINUE
239 *
240 * Solve L**T * x = b.
241 *
242  DO 90 i = n - 1, 1, -1
243  IF( ipiv( i ).EQ.i ) THEN
244  b( i, j ) = b( i, j ) - dl( i )*b( i+1, j )
245  ELSE
246  temp = b( i+1, j )
247  b( i+1, j ) = b( i, j ) - dl( i )*temp
248  b( i, j ) = temp
249  END IF
250  90 CONTINUE
251  IF( j.LT.nrhs ) THEN
252  j = j + 1
253  GO TO 70
254  END IF
255  ELSE
256  DO 120 j = 1, nrhs
257 *
258 * Solve U**T * x = b.
259 *
260  b( 1, j ) = b( 1, j ) / d( 1 )
261  IF( n.GT.1 )
262  $ b( 2, j ) = ( b( 2, j )-du( 1 )*b( 1, j ) ) / d( 2 )
263  DO 100 i = 3, n
264  b( i, j ) = ( b( i, j )-du( i-1 )*b( i-1, j )-
265  $ du2( i-2 )*b( i-2, j ) ) / d( i )
266  100 CONTINUE
267 *
268 * Solve L**T * x = b.
269 *
270  DO 110 i = n - 1, 1, -1
271  IF( ipiv( i ).EQ.i ) THEN
272  b( i, j ) = b( i, j ) - dl( i )*b( i+1, j )
273  ELSE
274  temp = b( i+1, j )
275  b( i+1, j ) = b( i, j ) - dl( i )*temp
276  b( i, j ) = temp
277  END IF
278  110 CONTINUE
279  120 CONTINUE
280  END IF
281  ELSE
282 *
283 * Solve A**H * X = B.
284 *
285  IF( nrhs.LE.1 ) THEN
286  j = 1
287  130 CONTINUE
288 *
289 * Solve U**H * x = b.
290 *
291  b( 1, j ) = b( 1, j ) / dconjg( d( 1 ) )
292  IF( n.GT.1 )
293  $ b( 2, j ) = ( b( 2, j )-dconjg( du( 1 ) )*b( 1, j ) ) /
294  $ dconjg( d( 2 ) )
295  DO 140 i = 3, n
296  b( i, j ) = ( b( i, j )-dconjg( du( i-1 ) )*b( i-1, j )-
297  $ dconjg( du2( i-2 ) )*b( i-2, j ) ) /
298  $ dconjg( d( i ) )
299  140 CONTINUE
300 *
301 * Solve L**H * x = b.
302 *
303  DO 150 i = n - 1, 1, -1
304  IF( ipiv( i ).EQ.i ) THEN
305  b( i, j ) = b( i, j ) - dconjg( dl( i ) )*b( i+1, j )
306  ELSE
307  temp = b( i+1, j )
308  b( i+1, j ) = b( i, j ) - dconjg( dl( i ) )*temp
309  b( i, j ) = temp
310  END IF
311  150 CONTINUE
312  IF( j.LT.nrhs ) THEN
313  j = j + 1
314  GO TO 130
315  END IF
316  ELSE
317  DO 180 j = 1, nrhs
318 *
319 * Solve U**H * x = b.
320 *
321  b( 1, j ) = b( 1, j ) / dconjg( d( 1 ) )
322  IF( n.GT.1 )
323  $ b( 2, j ) = ( b( 2, j )-dconjg( du( 1 ) )*b( 1, j ) )
324  $ / dconjg( d( 2 ) )
325  DO 160 i = 3, n
326  b( i, j ) = ( b( i, j )-dconjg( du( i-1 ) )*
327  $ b( i-1, j )-dconjg( du2( i-2 ) )*
328  $ b( i-2, j ) ) / dconjg( d( i ) )
329  160 CONTINUE
330 *
331 * Solve L**H * x = b.
332 *
333  DO 170 i = n - 1, 1, -1
334  IF( ipiv( i ).EQ.i ) THEN
335  b( i, j ) = b( i, j ) - dconjg( dl( i ) )*
336  $ b( i+1, j )
337  ELSE
338  temp = b( i+1, j )
339  b( i+1, j ) = b( i, j ) - dconjg( dl( i ) )*temp
340  b( i, j ) = temp
341  END IF
342  170 CONTINUE
343  180 CONTINUE
344  END IF
345  END IF
346 *
347 * End of ZGTTS2
348 *

Here is the caller graph for this function: