LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
subroutine zsytrs_rook ( character  UPLO,
integer  N,
integer  NRHS,
complex*16, dimension( lda, * )  A,
integer  LDA,
integer, dimension( * )  IPIV,
complex*16, dimension( ldb, * )  B,
integer  LDB,
integer  INFO 
)

ZSYTRS_ROOK

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

Purpose: 
 ZSYTRS_ROOK solves a system of linear equations A*X = B with
 a complex symmetric matrix A using the factorization A = U*D*U**T or
 A = L*D*L**T computed by ZSYTRF_ROOK.
Parameters
[in]UPLO
          UPLO is CHARACTER*1
          Specifies whether the details of the factorization are stored
          as an upper or lower triangular matrix.
          = 'U':  Upper triangular, form is A = U*D*U**T;
          = 'L':  Lower triangular, form is A = L*D*L**T.
[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]A
          A is COMPLEX*16 array, dimension (LDA,N)
          The block diagonal matrix D and the multipliers used to
          obtain the factor U or L as computed by ZSYTRF_ROOK.
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,N).
[in]IPIV
          IPIV is INTEGER array, dimension (N)
          Details of the interchanges and the block structure of D
          as determined by ZSYTRF_ROOK.
[in,out]B
          B is COMPLEX*16 array, dimension (LDB,NRHS)
          On entry, the right hand side matrix B.
          On exit, the solution matrix 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 = -i, the i-th argument had an illegal value
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2015
Contributors: 
   November 2015, Igor Kozachenko,
                  Computer Science Division,
                  University of California, Berkeley

  September 2007, Sven Hammarling, Nicholas J. Higham, Craig Lucas,
                  School of Mathematics,
                  University of Manchester

Definition at line 138 of file zsytrs_rook.f.

138 *
139 * -- LAPACK computational routine (version 3.6.0) --
140 * -- LAPACK is a software package provided by Univ. of Tennessee, --
141 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
142 * November 2015
143 *
144 * .. Scalar Arguments ..
145  CHARACTER uplo
146  INTEGER info, lda, ldb, n, nrhs
147 * ..
148 * .. Array Arguments ..
149  INTEGER ipiv( * )
150  COMPLEX*16 a( lda, * ), b( ldb, * )
151 * ..
152 *
153 * =====================================================================
154 *
155 * .. Parameters ..
156  COMPLEX*16 cone
157  parameter ( cone = ( 1.0d+0, 0.0d+0 ) )
158 * ..
159 * .. Local Scalars ..
160  LOGICAL upper
161  INTEGER j, k, kp
162  COMPLEX*16 ak, akm1, akm1k, bk, bkm1, denom
163 * ..
164 * .. External Functions ..
165  LOGICAL lsame
166  EXTERNAL lsame
167 * ..
168 * .. External Subroutines ..
169  EXTERNAL zgemv, zgeru, zscal, zswap, xerbla
170 * ..
171 * .. Intrinsic Functions ..
172  INTRINSIC max
173 * ..
174 * .. Executable Statements ..
175 *
176  info = 0
177  upper = lsame( uplo, 'U' )
178  IF( .NOT.upper .AND. .NOT.lsame( uplo, 'L' ) ) THEN
179  info = -1
180  ELSE IF( n.LT.0 ) THEN
181  info = -2
182  ELSE IF( nrhs.LT.0 ) THEN
183  info = -3
184  ELSE IF( lda.LT.max( 1, n ) ) THEN
185  info = -5
186  ELSE IF( ldb.LT.max( 1, n ) ) THEN
187  info = -8
188  END IF
189  IF( info.NE.0 ) THEN
190  CALL xerbla( 'ZSYTRS_ROOK', -info )
191  RETURN
192  END IF
193 *
194 * Quick return if possible
195 *
196  IF( n.EQ.0 .OR. nrhs.EQ.0 )
197  $ RETURN
198 *
199  IF( upper ) THEN
200 *
201 * Solve A*X = B, where A = U*D*U**T.
202 *
203 * First solve U*D*X = B, overwriting B with X.
204 *
205 * K is the main loop index, decreasing from N to 1 in steps of
206 * 1 or 2, depending on the size of the diagonal blocks.
207 *
208  k = n
209  10 CONTINUE
210 *
211 * If K < 1, exit from loop.
212 *
213  IF( k.LT.1 )
214  $ GO TO 30
215 *
216  IF( ipiv( k ).GT.0 ) THEN
217 *
218 * 1 x 1 diagonal block
219 *
220 * Interchange rows K and IPIV(K).
221 *
222  kp = ipiv( k )
223  IF( kp.NE.k )
224  $ CALL zswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
225 *
226 * Multiply by inv(U(K)), where U(K) is the transformation
227 * stored in column K of A.
228 *
229  CALL zgeru( k-1, nrhs, -cone, a( 1, k ), 1, b( k, 1 ), ldb,
230  $ b( 1, 1 ), ldb )
231 *
232 * Multiply by the inverse of the diagonal block.
233 *
234  CALL zscal( nrhs, cone / a( k, k ), b( k, 1 ), ldb )
235  k = k - 1
236  ELSE
237 *
238 * 2 x 2 diagonal block
239 *
240 * Interchange rows K and -IPIV(K) THEN K-1 and -IPIV(K-1)
241 *
242  kp = -ipiv( k )
243  IF( kp.NE.k )
244  $ CALL zswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
245 *
246  kp = -ipiv( k-1 )
247  IF( kp.NE.k-1 )
248  $ CALL zswap( nrhs, b( k-1, 1 ), ldb, b( kp, 1 ), ldb )
249 *
250 * Multiply by inv(U(K)), where U(K) is the transformation
251 * stored in columns K-1 and K of A.
252 *
253  IF( k.GT.2 ) THEN
254  CALL zgeru( k-2, nrhs,-cone, a( 1, k ), 1, b( k, 1 ),
255  $ ldb, b( 1, 1 ), ldb )
256  CALL zgeru( k-2, nrhs,-cone, a( 1, k-1 ), 1, b( k-1, 1 ),
257  $ ldb, b( 1, 1 ), ldb )
258  END IF
259 *
260 * Multiply by the inverse of the diagonal block.
261 *
262  akm1k = a( k-1, k )
263  akm1 = a( k-1, k-1 ) / akm1k
264  ak = a( k, k ) / akm1k
265  denom = akm1*ak - cone
266  DO 20 j = 1, nrhs
267  bkm1 = b( k-1, j ) / akm1k
268  bk = b( k, j ) / akm1k
269  b( k-1, j ) = ( ak*bkm1-bk ) / denom
270  b( k, j ) = ( akm1*bk-bkm1 ) / denom
271  20 CONTINUE
272  k = k - 2
273  END IF
274 *
275  GO TO 10
276  30 CONTINUE
277 *
278 * Next solve U**T *X = B, overwriting B with X.
279 *
280 * K is the main loop index, increasing from 1 to N in steps of
281 * 1 or 2, depending on the size of the diagonal blocks.
282 *
283  k = 1
284  40 CONTINUE
285 *
286 * If K > N, exit from loop.
287 *
288  IF( k.GT.n )
289  $ GO TO 50
290 *
291  IF( ipiv( k ).GT.0 ) THEN
292 *
293 * 1 x 1 diagonal block
294 *
295 * Multiply by inv(U**T(K)), where U(K) is the transformation
296 * stored in column K of A.
297 *
298  IF( k.GT.1 )
299  $ CALL zgemv( 'Transpose', k-1, nrhs, -cone, b,
300  $ ldb, a( 1, k ), 1, cone, b( k, 1 ), ldb )
301 *
302 * Interchange rows K and IPIV(K).
303 *
304  kp = ipiv( k )
305  IF( kp.NE.k )
306  $ CALL zswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
307  k = k + 1
308  ELSE
309 *
310 * 2 x 2 diagonal block
311 *
312 * Multiply by inv(U**T(K+1)), where U(K+1) is the transformation
313 * stored in columns K and K+1 of A.
314 *
315  IF( k.GT.1 ) THEN
316  CALL zgemv( 'Transpose', k-1, nrhs, -cone, b,
317  $ ldb, a( 1, k ), 1, cone, b( k, 1 ), ldb )
318  CALL zgemv( 'Transpose', k-1, nrhs, -cone, b,
319  $ ldb, a( 1, k+1 ), 1, cone, b( k+1, 1 ), ldb )
320  END IF
321 *
322 * Interchange rows K and -IPIV(K) THEN K+1 and -IPIV(K+1).
323 *
324  kp = -ipiv( k )
325  IF( kp.NE.k )
326  $ CALL zswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
327 *
328  kp = -ipiv( k+1 )
329  IF( kp.NE.k+1 )
330  $ CALL zswap( nrhs, b( k+1, 1 ), ldb, b( kp, 1 ), ldb )
331 *
332  k = k + 2
333  END IF
334 *
335  GO TO 40
336  50 CONTINUE
337 *
338  ELSE
339 *
340 * Solve A*X = B, where A = L*D*L**T.
341 *
342 * First solve L*D*X = B, overwriting B with X.
343 *
344 * K is the main loop index, increasing from 1 to N in steps of
345 * 1 or 2, depending on the size of the diagonal blocks.
346 *
347  k = 1
348  60 CONTINUE
349 *
350 * If K > N, exit from loop.
351 *
352  IF( k.GT.n )
353  $ GO TO 80
354 *
355  IF( ipiv( k ).GT.0 ) THEN
356 *
357 * 1 x 1 diagonal block
358 *
359 * Interchange rows K and IPIV(K).
360 *
361  kp = ipiv( k )
362  IF( kp.NE.k )
363  $ CALL zswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
364 *
365 * Multiply by inv(L(K)), where L(K) is the transformation
366 * stored in column K of A.
367 *
368  IF( k.LT.n )
369  $ CALL zgeru( n-k, nrhs, -cone, a( k+1, k ), 1, b( k, 1 ),
370  $ ldb, b( k+1, 1 ), ldb )
371 *
372 * Multiply by the inverse of the diagonal block.
373 *
374  CALL zscal( nrhs, cone / a( k, k ), b( k, 1 ), ldb )
375  k = k + 1
376  ELSE
377 *
378 * 2 x 2 diagonal block
379 *
380 * Interchange rows K and -IPIV(K) THEN K+1 and -IPIV(K+1)
381 *
382  kp = -ipiv( k )
383  IF( kp.NE.k )
384  $ CALL zswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
385 *
386  kp = -ipiv( k+1 )
387  IF( kp.NE.k+1 )
388  $ CALL zswap( nrhs, b( k+1, 1 ), ldb, b( kp, 1 ), ldb )
389 *
390 * Multiply by inv(L(K)), where L(K) is the transformation
391 * stored in columns K and K+1 of A.
392 *
393  IF( k.LT.n-1 ) THEN
394  CALL zgeru( n-k-1, nrhs,-cone, a( k+2, k ), 1, b( k, 1 ),
395  $ ldb, b( k+2, 1 ), ldb )
396  CALL zgeru( n-k-1, nrhs,-cone, a( k+2, k+1 ), 1,
397  $ b( k+1, 1 ), ldb, b( k+2, 1 ), ldb )
398  END IF
399 *
400 * Multiply by the inverse of the diagonal block.
401 *
402  akm1k = a( k+1, k )
403  akm1 = a( k, k ) / akm1k
404  ak = a( k+1, k+1 ) / akm1k
405  denom = akm1*ak - cone
406  DO 70 j = 1, nrhs
407  bkm1 = b( k, j ) / akm1k
408  bk = b( k+1, j ) / akm1k
409  b( k, j ) = ( ak*bkm1-bk ) / denom
410  b( k+1, j ) = ( akm1*bk-bkm1 ) / denom
411  70 CONTINUE
412  k = k + 2
413  END IF
414 *
415  GO TO 60
416  80 CONTINUE
417 *
418 * Next solve L**T *X = B, overwriting B with X.
419 *
420 * K is the main loop index, decreasing from N to 1 in steps of
421 * 1 or 2, depending on the size of the diagonal blocks.
422 *
423  k = n
424  90 CONTINUE
425 *
426 * If K < 1, exit from loop.
427 *
428  IF( k.LT.1 )
429  $ GO TO 100
430 *
431  IF( ipiv( k ).GT.0 ) THEN
432 *
433 * 1 x 1 diagonal block
434 *
435 * Multiply by inv(L**T(K)), where L(K) is the transformation
436 * stored in column K of A.
437 *
438  IF( k.LT.n )
439  $ CALL zgemv( 'Transpose', n-k, nrhs, -cone, b( k+1, 1 ),
440  $ ldb, a( k+1, k ), 1, cone, b( k, 1 ), ldb )
441 *
442 * Interchange rows K and IPIV(K).
443 *
444  kp = ipiv( k )
445  IF( kp.NE.k )
446  $ CALL zswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
447  k = k - 1
448  ELSE
449 *
450 * 2 x 2 diagonal block
451 *
452 * Multiply by inv(L**T(K-1)), where L(K-1) is the transformation
453 * stored in columns K-1 and K of A.
454 *
455  IF( k.LT.n ) THEN
456  CALL zgemv( 'Transpose', n-k, nrhs, -cone, b( k+1, 1 ),
457  $ ldb, a( k+1, k ), 1, cone, b( k, 1 ), ldb )
458  CALL zgemv( 'Transpose', n-k, nrhs, -cone, b( k+1, 1 ),
459  $ ldb, a( k+1, k-1 ), 1, cone, b( k-1, 1 ),
460  $ ldb )
461  END IF
462 *
463 * Interchange rows K and -IPIV(K) THEN K-1 and -IPIV(K-1)
464 *
465  kp = -ipiv( k )
466  IF( kp.NE.k )
467  $ CALL zswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
468 *
469  kp = -ipiv( k-1 )
470  IF( kp.NE.k-1 )
471  $ CALL zswap( nrhs, b( k-1, 1 ), ldb, b( kp, 1 ), ldb )
472 *
473  k = k - 2
474  END IF
475 *
476  GO TO 90
477  100 CONTINUE
478  END IF
479 *
480  RETURN
481 *
482 * End of ZSYTRS_ROOK
483 *
zgemv
subroutine zgemv(TRANS, M, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
ZGEMV
Definition: zgemv.f:160
zswap
subroutine zswap(N, ZX, INCX, ZY, INCY)
ZSWAP
Definition: zswap.f:52
xerbla
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
zgeru
subroutine zgeru(M, N, ALPHA, X, INCX, Y, INCY, A, LDA)
ZGERU
Definition: zgeru.f:132
lsame
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
zscal
subroutine zscal(N, ZA, ZX, INCX)
ZSCAL
Definition: zscal.f:54

Here is the call graph for this function:

Here is the caller graph for this function: