LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches

◆ slavsy_rook()

subroutine slavsy_rook ( character  uplo,
character  trans,
character  diag,
integer  n,
integer  nrhs,
real, dimension( lda, * )  a,
integer  lda,
integer, dimension( * )  ipiv,
real, dimension( ldb, * )  b,
integer  ldb,
integer  info 
)

SLAVSY_ROOK

Purpose:
 SLAVSY_ROOK  performs one of the matrix-vector operations
    x := A*x  or  x := A'*x,
 where x is an N element vector and A is one of the factors
 from the block U*D*U' or L*D*L' factorization computed by SSYTRF_ROOK.

 If TRANS = 'N', multiplies by U  or U * D  (or L  or L * D)
 If TRANS = 'T', multiplies by U' or D * U' (or L' or D * L')
 If TRANS = 'C', multiplies by U' or D * U' (or L' or D * L')
Parameters
[in]UPLO
          UPLO is CHARACTER*1
          Specifies whether the factor stored in A is upper or lower
          triangular.
          = 'U':  Upper triangular
          = 'L':  Lower triangular
[in]TRANS
          TRANS is CHARACTER*1
          Specifies the operation to be performed:
          = 'N':  x := A*x
          = 'T':  x := A'*x
          = 'C':  x := A'*x
[in]DIAG
          DIAG is CHARACTER*1
          Specifies whether or not the diagonal blocks are unit
          matrices.  If the diagonal blocks are assumed to be unit,
          then A = U or A = L, otherwise A = U*D or A = L*D.
          = 'U':  Diagonal blocks are assumed to be unit matrices.
          = 'N':  Diagonal blocks are assumed to be non-unit matrices.
[in]N
          N is INTEGER
          The number of rows and columns of the matrix A.  N >= 0.
[in]NRHS
          NRHS is INTEGER
          The number of right hand sides, i.e., the number of vectors
          x to be multiplied by A.  NRHS >= 0.
[in]A
          A is REAL array, dimension (LDA,N)
          The block diagonal matrix D and the multipliers used to
          obtain the factor U or L as computed by SSYTRF_ROOK.
          Stored as a 2-D triangular matrix.
[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 SSYTRF_ROOK.

          If UPLO = 'U':
               If IPIV(k) > 0, then rows and columns k and IPIV(k)
               were interchanged and D(k,k) is a 1-by-1 diagonal block.
               (If IPIV( k ) = k, no interchange was done).

               If IPIV(k) < 0 and IPIV(k-1) < 0, then rows and
               columns k and -IPIV(k) were interchanged and rows and
               columns k-1 and -IPIV(k-1) were inerchaged,
               D(k-1:k,k-1:k) is a 2-by-2 diagonal block.

          If UPLO = 'L':
               If IPIV(k) > 0, then rows and columns k and IPIV(k)
               were interchanged and D(k,k) is a 1-by-1 diagonal block.
               (If IPIV( k ) = k, no interchange was done).

               If IPIV(k) < 0 and IPIV(k+1) < 0, then rows and
               columns k and -IPIV(k) were interchanged and rows and
               columns k+1 and -IPIV(k+1) were inerchaged,
               D(k:k+1,k:k+1) is a 2-by-2 diagonal block.
[in,out]B
          B is REAL array, dimension (LDB,NRHS)
          On entry, B contains NRHS vectors of length N.
          On exit, B is overwritten with the product A * B.
[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.

Definition at line 155 of file slavsy_rook.f.

157*
158* -- LAPACK test routine --
159* -- LAPACK is a software package provided by Univ. of Tennessee, --
160* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
161*
162* .. Scalar Arguments ..
163 CHARACTER DIAG, TRANS, UPLO
164 INTEGER INFO, LDA, LDB, N, NRHS
165* ..
166* .. Array Arguments ..
167 INTEGER IPIV( * )
168 REAL A( LDA, * ), B( LDB, * )
169* ..
170*
171* =====================================================================
172*
173* .. Parameters ..
174 REAL ONE
175 parameter( one = 1.0e+0 )
176* ..
177* .. Local Scalars ..
178 LOGICAL NOUNIT
179 INTEGER J, K, KP
180 REAL D11, D12, D21, D22, T1, T2
181* ..
182* .. External Functions ..
183 LOGICAL LSAME
184 EXTERNAL lsame
185* ..
186* .. External Subroutines ..
187 EXTERNAL sgemv, sger, sscal, sswap, xerbla
188* ..
189* .. Intrinsic Functions ..
190 INTRINSIC abs, max
191* ..
192* .. Executable Statements ..
193*
194* Test the input parameters.
195*
196 info = 0
197 IF( .NOT.lsame( uplo, 'U' ) .AND. .NOT.lsame( uplo, 'L' ) ) THEN
198 info = -1
199 ELSE IF( .NOT.lsame( trans, 'N' ) .AND. .NOT.
200 $ lsame( trans, 'T' ) .AND. .NOT.lsame( trans, 'C' ) ) THEN
201 info = -2
202 ELSE IF( .NOT.lsame( diag, 'U' ) .AND. .NOT.lsame( diag, 'N' ) )
203 $ THEN
204 info = -3
205 ELSE IF( n.LT.0 ) THEN
206 info = -4
207 ELSE IF( lda.LT.max( 1, n ) ) THEN
208 info = -6
209 ELSE IF( ldb.LT.max( 1, n ) ) THEN
210 info = -9
211 END IF
212 IF( info.NE.0 ) THEN
213 CALL xerbla( 'SLAVSY_ROOK ', -info )
214 RETURN
215 END IF
216*
217* Quick return if possible.
218*
219 IF( n.EQ.0 )
220 $ RETURN
221*
222 nounit = lsame( diag, 'N' )
223*------------------------------------------
224*
225* Compute B := A * B (No transpose)
226*
227*------------------------------------------
228 IF( lsame( trans, 'N' ) ) THEN
229*
230* Compute B := U*B
231* where U = P(m)*inv(U(m))* ... *P(1)*inv(U(1))
232*
233 IF( lsame( uplo, 'U' ) ) THEN
234*
235* Loop forward applying the transformations.
236*
237 k = 1
238 10 CONTINUE
239 IF( k.GT.n )
240 $ GO TO 30
241 IF( ipiv( k ).GT.0 ) THEN
242*
243* 1 x 1 pivot block
244*
245* Multiply by the diagonal element if forming U * D.
246*
247 IF( nounit )
248 $ CALL sscal( nrhs, a( k, k ), b( k, 1 ), ldb )
249*
250* Multiply by P(K) * inv(U(K)) if K > 1.
251*
252 IF( k.GT.1 ) THEN
253*
254* Apply the transformation.
255*
256 CALL sger( k-1, nrhs, one, a( 1, k ), 1, b( k, 1 ),
257 $ ldb, b( 1, 1 ), ldb )
258*
259* Interchange if P(K) .ne. I.
260*
261 kp = ipiv( k )
262 IF( kp.NE.k )
263 $ CALL sswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
264 END IF
265 k = k + 1
266 ELSE
267*
268* 2 x 2 pivot block
269*
270* Multiply by the diagonal block if forming U * D.
271*
272 IF( nounit ) THEN
273 d11 = a( k, k )
274 d22 = a( k+1, k+1 )
275 d12 = a( k, k+1 )
276 d21 = d12
277 DO 20 j = 1, nrhs
278 t1 = b( k, j )
279 t2 = b( k+1, j )
280 b( k, j ) = d11*t1 + d12*t2
281 b( k+1, j ) = d21*t1 + d22*t2
282 20 CONTINUE
283 END IF
284*
285* Multiply by P(K) * inv(U(K)) if K > 1.
286*
287 IF( k.GT.1 ) THEN
288*
289* Apply the transformations.
290*
291 CALL sger( k-1, nrhs, one, a( 1, k ), 1, b( k, 1 ),
292 $ ldb, b( 1, 1 ), ldb )
293 CALL sger( k-1, nrhs, one, a( 1, k+1 ), 1,
294 $ b( k+1, 1 ), ldb, b( 1, 1 ), ldb )
295*
296* Interchange if a permutation was applied at the
297* K-th step of the factorization.
298*
299* Swap the first of pair with IMAXth
300*
301 kp = abs( ipiv( k ) )
302 IF( kp.NE.k )
303 $ CALL sswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
304*
305* NOW swap the first of pair with Pth
306*
307 kp = abs( ipiv( k+1 ) )
308 IF( kp.NE.k+1 )
309 $ CALL sswap( nrhs, b( k+1, 1 ), ldb, b( kp, 1 ),
310 $ ldb )
311 END IF
312 k = k + 2
313 END IF
314 GO TO 10
315 30 CONTINUE
316*
317* Compute B := L*B
318* where L = P(1)*inv(L(1))* ... *P(m)*inv(L(m)) .
319*
320 ELSE
321*
322* Loop backward applying the transformations to B.
323*
324 k = n
325 40 CONTINUE
326 IF( k.LT.1 )
327 $ GO TO 60
328*
329* Test the pivot index. If greater than zero, a 1 x 1
330* pivot was used, otherwise a 2 x 2 pivot was used.
331*
332 IF( ipiv( k ).GT.0 ) THEN
333*
334* 1 x 1 pivot block:
335*
336* Multiply by the diagonal element if forming L * D.
337*
338 IF( nounit )
339 $ CALL sscal( nrhs, a( k, k ), b( k, 1 ), ldb )
340*
341* Multiply by P(K) * inv(L(K)) if K < N.
342*
343 IF( k.NE.n ) THEN
344 kp = ipiv( k )
345*
346* Apply the transformation.
347*
348 CALL sger( n-k, nrhs, one, a( k+1, k ), 1, b( k, 1 ),
349 $ ldb, b( k+1, 1 ), ldb )
350*
351* Interchange if a permutation was applied at the
352* K-th step of the factorization.
353*
354 IF( kp.NE.k )
355 $ CALL sswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
356 END IF
357 k = k - 1
358*
359 ELSE
360*
361* 2 x 2 pivot block:
362*
363* Multiply by the diagonal block if forming L * D.
364*
365 IF( nounit ) THEN
366 d11 = a( k-1, k-1 )
367 d22 = a( k, k )
368 d21 = a( k, k-1 )
369 d12 = d21
370 DO 50 j = 1, nrhs
371 t1 = b( k-1, j )
372 t2 = b( k, j )
373 b( k-1, j ) = d11*t1 + d12*t2
374 b( k, j ) = d21*t1 + d22*t2
375 50 CONTINUE
376 END IF
377*
378* Multiply by P(K) * inv(L(K)) if K < N.
379*
380 IF( k.NE.n ) THEN
381*
382* Apply the transformation.
383*
384 CALL sger( n-k, nrhs, one, a( k+1, k ), 1, b( k, 1 ),
385 $ ldb, b( k+1, 1 ), ldb )
386 CALL sger( n-k, nrhs, one, a( k+1, k-1 ), 1,
387 $ b( k-1, 1 ), ldb, b( k+1, 1 ), ldb )
388*
389* Interchange if a permutation was applied at the
390* K-th step of the factorization.
391*
392* Swap the second of pair with IMAXth
393*
394 kp = abs( ipiv( k ) )
395 IF( kp.NE.k )
396 $ CALL sswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
397*
398* NOW swap the first of pair with Pth
399*
400 kp = abs( ipiv( k-1 ) )
401 IF( kp.NE.k-1 )
402 $ CALL sswap( nrhs, b( k-1, 1 ), ldb, b( kp, 1 ),
403 $ ldb )
404 END IF
405 k = k - 2
406 END IF
407 GO TO 40
408 60 CONTINUE
409 END IF
410*----------------------------------------
411*
412* Compute B := A' * B (transpose)
413*
414*----------------------------------------
415 ELSE
416*
417* Form B := U'*B
418* where U = P(m)*inv(U(m))* ... *P(1)*inv(U(1))
419* and U' = inv(U'(1))*P(1)* ... *inv(U'(m))*P(m)
420*
421 IF( lsame( uplo, 'U' ) ) THEN
422*
423* Loop backward applying the transformations.
424*
425 k = n
426 70 CONTINUE
427 IF( k.LT.1 )
428 $ GO TO 90
429*
430* 1 x 1 pivot block.
431*
432 IF( ipiv( k ).GT.0 ) THEN
433 IF( k.GT.1 ) THEN
434*
435* Interchange if P(K) .ne. I.
436*
437 kp = ipiv( k )
438 IF( kp.NE.k )
439 $ CALL sswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
440*
441* Apply the transformation
442*
443 CALL sgemv( 'Transpose', k-1, nrhs, one, b, ldb,
444 $ a( 1, k ), 1, one, b( k, 1 ), ldb )
445 END IF
446 IF( nounit )
447 $ CALL sscal( nrhs, a( k, k ), b( k, 1 ), ldb )
448 k = k - 1
449*
450* 2 x 2 pivot block.
451*
452 ELSE
453 IF( k.GT.2 ) THEN
454*
455* Swap the second of pair with Pth
456*
457 kp = abs( ipiv( k ) )
458 IF( kp.NE.k )
459 $ CALL sswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
460*
461* Now swap the first of pair with IMAX(r)th
462*
463 kp = abs( ipiv( k-1 ) )
464 IF( kp.NE.k-1 )
465 $ CALL sswap( nrhs, b( k-1, 1 ), ldb, b( kp, 1 ),
466 $ ldb )
467*
468* Apply the transformations
469*
470 CALL sgemv( 'Transpose', k-2, nrhs, one, b, ldb,
471 $ a( 1, k ), 1, one, b( k, 1 ), ldb )
472 CALL sgemv( 'Transpose', k-2, nrhs, one, b, ldb,
473 $ a( 1, k-1 ), 1, one, b( k-1, 1 ), ldb )
474 END IF
475*
476* Multiply by the diagonal block if non-unit.
477*
478 IF( nounit ) THEN
479 d11 = a( k-1, k-1 )
480 d22 = a( k, k )
481 d12 = a( k-1, k )
482 d21 = d12
483 DO 80 j = 1, nrhs
484 t1 = b( k-1, j )
485 t2 = b( k, j )
486 b( k-1, j ) = d11*t1 + d12*t2
487 b( k, j ) = d21*t1 + d22*t2
488 80 CONTINUE
489 END IF
490 k = k - 2
491 END IF
492 GO TO 70
493 90 CONTINUE
494*
495* Form B := L'*B
496* where L = P(1)*inv(L(1))* ... *P(m)*inv(L(m))
497* and L' = inv(L'(m))*P(m)* ... *inv(L'(1))*P(1)
498*
499 ELSE
500*
501* Loop forward applying the L-transformations.
502*
503 k = 1
504 100 CONTINUE
505 IF( k.GT.n )
506 $ GO TO 120
507*
508* 1 x 1 pivot block
509*
510 IF( ipiv( k ).GT.0 ) THEN
511 IF( k.LT.n ) THEN
512*
513* Interchange if P(K) .ne. I.
514*
515 kp = ipiv( k )
516 IF( kp.NE.k )
517 $ CALL sswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
518*
519* Apply the transformation
520*
521 CALL sgemv( 'Transpose', n-k, nrhs, one, b( k+1, 1 ),
522 $ ldb, a( k+1, k ), 1, one, b( k, 1 ), ldb )
523 END IF
524 IF( nounit )
525 $ CALL sscal( nrhs, a( k, k ), b( k, 1 ), ldb )
526 k = k + 1
527*
528* 2 x 2 pivot block.
529*
530 ELSE
531 IF( k.LT.n-1 ) THEN
532*
533* Swap the first of pair with Pth
534*
535 kp = abs( ipiv( k ) )
536 IF( kp.NE.k )
537 $ CALL sswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
538*
539* Now swap the second of pair with IMAX(r)th
540*
541 kp = abs( ipiv( k+1 ) )
542 IF( kp.NE.k+1 )
543 $ CALL sswap( nrhs, b( k+1, 1 ), ldb, b( kp, 1 ),
544 $ ldb )
545*
546* Apply the transformation
547*
548 CALL sgemv( 'Transpose', n-k-1, nrhs, one,
549 $ b( k+2, 1 ), ldb, a( k+2, k+1 ), 1, one,
550 $ b( k+1, 1 ), ldb )
551 CALL sgemv( 'Transpose', n-k-1, nrhs, one,
552 $ b( k+2, 1 ), ldb, a( k+2, k ), 1, one,
553 $ b( k, 1 ), ldb )
554 END IF
555*
556* Multiply by the diagonal block if non-unit.
557*
558 IF( nounit ) THEN
559 d11 = a( k, k )
560 d22 = a( k+1, k+1 )
561 d21 = a( k+1, k )
562 d12 = d21
563 DO 110 j = 1, nrhs
564 t1 = b( k, j )
565 t2 = b( k+1, j )
566 b( k, j ) = d11*t1 + d12*t2
567 b( k+1, j ) = d21*t1 + d22*t2
568 110 CONTINUE
569 END IF
570 k = k + 2
571 END IF
572 GO TO 100
573 120 CONTINUE
574 END IF
575*
576 END IF
577 RETURN
578*
579* End of SLAVSY_ROOK
580*
xerbla
subroutine xerbla(srname, info)
Definition cblat2.f:3285
sgemv
subroutine sgemv(trans, m, n, alpha, a, lda, x, incx, beta, y, incy)
SGEMV
Definition sgemv.f:158
sger
subroutine sger(m, n, alpha, x, incx, y, incy, a, lda)
SGER
Definition sger.f:130
lsame
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48
sscal
subroutine sscal(n, sa, sx, incx)
SSCAL
Definition sscal.f:79
sswap
subroutine sswap(n, sx, incx, sy, incy)
SSWAP
Definition sswap.f:82
Here is the call graph for this function:
Here is the caller graph for this function: