LAPACK  3.10.0
LAPACK: Linear Algebra PACKage

◆ zsytrs_rook()

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.
Contributors:
   December 2016, 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 134 of file zsytrs_rook.f.

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