LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
subroutine clavsp ( character  UPLO,
character  TRANS,
character  DIAG,
integer  N,
integer  NRHS,
complex, dimension( * )  A,
integer, dimension( * )  IPIV,
complex, dimension( ldb, * )  B,
integer  LDB,
integer  INFO 
)

CLAVSP

Purpose: 
    CLAVSP  performs one of the matrix-vector operations
       x := A*x  or  x := A^T*x,
    where x is an N element vector and  A is one of the factors
    from the symmetric factorization computed by CSPTRF.
    CSPTRF produces a factorization of the form
         U * D * U^T     or     L * D * L^T,
    where U (or L) is a product of permutation and unit upper (lower)
    triangular matrices, U^T (or L^T) is the transpose of
    U (or L), and D is symmetric and block diagonal with 1 x 1 and
    2 x 2 diagonal blocks.  The multipliers for the transformations
    and the upper or lower triangular parts of the diagonal blocks
    are stored columnwise in packed format in the linear array A.

    If TRANS = 'N' or 'n', CLAVSP multiplies either by U or U * D
    (or L or L * D).
    If TRANS = 'C' or 'c', CLAVSP multiplies either by U^T or D * U^T
    (or L^T or D * L^T ).
  UPLO   - CHARACTER*1
           On entry, UPLO specifies whether the triangular matrix
           stored in A is upper or lower triangular.
              UPLO = 'U' or 'u'   The matrix is upper triangular.
              UPLO = 'L' or 'l'   The matrix is lower triangular.
           Unchanged on exit.

  TRANS  - CHARACTER*1
           On entry, TRANS specifies the operation to be performed as
           follows:
              TRANS = 'N' or 'n'   x := A*x.
              TRANS = 'T' or 't'   x := A^T*x.
           Unchanged on exit.

  DIAG   - CHARACTER*1
           On entry, DIAG specifies whether the diagonal blocks are
           assumed to be unit matrices, as follows:
              DIAG = 'U' or 'u'   Diagonal blocks are unit matrices.
              DIAG = 'N' or 'n'   Diagonal blocks are non-unit.
           Unchanged on exit.

  N      - INTEGER
           On entry, N specifies the order of the matrix A.
           N must be at least zero.
           Unchanged on exit.

  NRHS   - INTEGER
           On entry, NRHS specifies the number of right hand sides,
           i.e., the number of vectors x to be multiplied by A.
           NRHS must be at least zero.
           Unchanged on exit.

  A      - COMPLEX array, dimension( N*(N+1)/2 )
           On entry, A contains a block diagonal matrix and the
           multipliers of the transformations used to obtain it,
           stored as a packed triangular matrix.
           Unchanged on exit.

  IPIV   - INTEGER array, dimension( N )
           On entry, IPIV contains the vector of pivot indices as
           determined by CSPTRF.
           If IPIV( K ) = K, no interchange was done.
           If IPIV( K ) <> K but IPIV( K ) > 0, then row K was inter-
           changed with row IPIV( K ) and a 1 x 1 pivot block was used.
           If IPIV( K ) < 0 and UPLO = 'U', then row K-1 was exchanged
           with row | IPIV( K ) | and a 2 x 2 pivot block was used.
           If IPIV( K ) < 0 and UPLO = 'L', then row K+1 was exchanged
           with row | IPIV( K ) | and a 2 x 2 pivot block was used.

  B      - COMPLEX array, dimension( LDB, NRHS )
           On entry, B contains NRHS vectors of length N.
           On exit, B is overwritten with the product A * B.

  LDB    - INTEGER
           On entry, LDB contains the leading dimension of B as
           declared in the calling program.  LDB must be at least
           max( 1, N ).
           Unchanged on exit.

  INFO   - INTEGER
           INFO is the error flag.
           On exit, a value of 0 indicates a successful exit.
           A negative value, say -K, indicates that the K-th argument
           has an illegal value.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2011

Definition at line 133 of file clavsp.f.

133 *
134 * -- LAPACK test routine (version 3.4.0) --
135 * -- LAPACK is a software package provided by Univ. of Tennessee, --
136 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
137 * November 2011
138 *
139 * .. Scalar Arguments ..
140  CHARACTER diag, trans, uplo
141  INTEGER info, ldb, n, nrhs
142 * ..
143 * .. Array Arguments ..
144  INTEGER ipiv( * )
145  COMPLEX a( * ), b( ldb, * )
146 * ..
147 *
148 * =====================================================================
149 *
150 * .. Parameters ..
151  COMPLEX one
152  parameter ( one = ( 1.0e+0, 0.0e+0 ) )
153 * ..
154 * .. Local Scalars ..
155  LOGICAL nounit
156  INTEGER j, k, kc, kcnext, kp
157  COMPLEX d11, d12, d21, d22, t1, t2
158 * ..
159 * .. External Functions ..
160  LOGICAL lsame
161  EXTERNAL lsame
162 * ..
163 * .. External Subroutines ..
164  EXTERNAL cgemv, cgeru, cscal, cswap, xerbla
165 * ..
166 * .. Intrinsic Functions ..
167  INTRINSIC abs, max
168 * ..
169 * .. Executable Statements ..
170 *
171 * Test the input parameters.
172 *
173  info = 0
174  IF( .NOT.lsame( uplo, 'U' ) .AND. .NOT.lsame( uplo, 'L' ) ) THEN
175  info = -1
176  ELSE IF( .NOT.lsame( trans, 'N' ) .AND. .NOT.lsame( trans, 'T' ) )
177  $ THEN
178  info = -2
179  ELSE IF( .NOT.lsame( diag, 'U' ) .AND. .NOT.lsame( diag, 'N' ) )
180  $ THEN
181  info = -3
182  ELSE IF( n.LT.0 ) THEN
183  info = -4
184  ELSE IF( ldb.LT.max( 1, n ) ) THEN
185  info = -8
186  END IF
187  IF( info.NE.0 ) THEN
188  CALL xerbla( 'CLAVSP ', -info )
189  RETURN
190  END IF
191 *
192 * Quick return if possible.
193 *
194  IF( n.EQ.0 )
195  $ RETURN
196 *
197  nounit = lsame( diag, 'N' )
198 *------------------------------------------
199 *
200 * Compute B := A * B (No transpose)
201 *
202 *------------------------------------------
203  IF( lsame( trans, 'N' ) ) THEN
204 *
205 * Compute B := U*B
206 * where U = P(m)*inv(U(m))* ... *P(1)*inv(U(1))
207 *
208  IF( lsame( uplo, 'U' ) ) THEN
209 *
210 * Loop forward applying the transformations.
211 *
212  k = 1
213  kc = 1
214  10 CONTINUE
215  IF( k.GT.n )
216  $ GO TO 30
217 *
218 * 1 x 1 pivot block
219 *
220  IF( ipiv( k ).GT.0 ) THEN
221 *
222 * Multiply by the diagonal element if forming U * D.
223 *
224  IF( nounit )
225  $ CALL cscal( nrhs, a( kc+k-1 ), b( k, 1 ), ldb )
226 *
227 * Multiply by P(K) * inv(U(K)) if K > 1.
228 *
229  IF( k.GT.1 ) THEN
230 *
231 * Apply the transformation.
232 *
233  CALL cgeru( k-1, nrhs, one, a( kc ), 1, b( k, 1 ),
234  $ ldb, b( 1, 1 ), ldb )
235 *
236 * Interchange if P(K) != I.
237 *
238  kp = ipiv( k )
239  IF( kp.NE.k )
240  $ CALL cswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
241  END IF
242  kc = kc + k
243  k = k + 1
244  ELSE
245 *
246 * 2 x 2 pivot block
247 *
248  kcnext = kc + k
249 *
250 * Multiply by the diagonal block if forming U * D.
251 *
252  IF( nounit ) THEN
253  d11 = a( kcnext-1 )
254  d22 = a( kcnext+k )
255  d12 = a( kcnext+k-1 )
256  d21 = d12
257  DO 20 j = 1, nrhs
258  t1 = b( k, j )
259  t2 = b( k+1, j )
260  b( k, j ) = d11*t1 + d12*t2
261  b( k+1, j ) = d21*t1 + d22*t2
262  20 CONTINUE
263  END IF
264 *
265 * Multiply by P(K) * inv(U(K)) if K > 1.
266 *
267  IF( k.GT.1 ) THEN
268 *
269 * Apply the transformations.
270 *
271  CALL cgeru( k-1, nrhs, one, a( kc ), 1, b( k, 1 ),
272  $ ldb, b( 1, 1 ), ldb )
273  CALL cgeru( k-1, nrhs, one, a( kcnext ), 1,
274  $ b( k+1, 1 ), ldb, b( 1, 1 ), ldb )
275 *
276 * Interchange if P(K) != I.
277 *
278  kp = abs( ipiv( k ) )
279  IF( kp.NE.k )
280  $ CALL cswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
281  END IF
282  kc = kcnext + k + 1
283  k = k + 2
284  END IF
285  GO TO 10
286  30 CONTINUE
287 *
288 * Compute B := L*B
289 * where L = P(1)*inv(L(1))* ... *P(m)*inv(L(m)) .
290 *
291  ELSE
292 *
293 * Loop backward applying the transformations to B.
294 *
295  k = n
296  kc = n*( n+1 ) / 2 + 1
297  40 CONTINUE
298  IF( k.LT.1 )
299  $ GO TO 60
300  kc = kc - ( n-k+1 )
301 *
302 * Test the pivot index. If greater than zero, a 1 x 1
303 * pivot was used, otherwise a 2 x 2 pivot was used.
304 *
305  IF( ipiv( k ).GT.0 ) THEN
306 *
307 * 1 x 1 pivot block:
308 *
309 * Multiply by the diagonal element if forming L * D.
310 *
311  IF( nounit )
312  $ CALL cscal( nrhs, a( kc ), b( k, 1 ), ldb )
313 *
314 * Multiply by P(K) * inv(L(K)) if K < N.
315 *
316  IF( k.NE.n ) THEN
317  kp = ipiv( k )
318 *
319 * Apply the transformation.
320 *
321  CALL cgeru( n-k, nrhs, one, a( kc+1 ), 1, b( k, 1 ),
322  $ ldb, b( k+1, 1 ), ldb )
323 *
324 * Interchange if a permutation was applied at the
325 * K-th step of the factorization.
326 *
327  IF( kp.NE.k )
328  $ CALL cswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
329  END IF
330  k = k - 1
331 *
332  ELSE
333 *
334 * 2 x 2 pivot block:
335 *
336  kcnext = kc - ( n-k+2 )
337 *
338 * Multiply by the diagonal block if forming L * D.
339 *
340  IF( nounit ) THEN
341  d11 = a( kcnext )
342  d22 = a( kc )
343  d21 = a( kcnext+1 )
344  d12 = d21
345  DO 50 j = 1, nrhs
346  t1 = b( k-1, j )
347  t2 = b( k, j )
348  b( k-1, j ) = d11*t1 + d12*t2
349  b( k, j ) = d21*t1 + d22*t2
350  50 CONTINUE
351  END IF
352 *
353 * Multiply by P(K) * inv(L(K)) if K < N.
354 *
355  IF( k.NE.n ) THEN
356 *
357 * Apply the transformation.
358 *
359  CALL cgeru( n-k, nrhs, one, a( kc+1 ), 1, b( k, 1 ),
360  $ ldb, b( k+1, 1 ), ldb )
361  CALL cgeru( n-k, nrhs, one, a( kcnext+2 ), 1,
362  $ b( k-1, 1 ), ldb, b( k+1, 1 ), ldb )
363 *
364 * Interchange if a permutation was applied at the
365 * K-th step of the factorization.
366 *
367  kp = abs( ipiv( k ) )
368  IF( kp.NE.k )
369  $ CALL cswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
370  END IF
371  kc = kcnext
372  k = k - 2
373  END IF
374  GO TO 40
375  60 CONTINUE
376  END IF
377 *-------------------------------------------------
378 *
379 * Compute B := A^T * B (transpose)
380 *
381 *-------------------------------------------------
382  ELSE
383 *
384 * Form B := U^T*B
385 * where U = P(m)*inv(U(m))* ... *P(1)*inv(U(1))
386 * and U^T = inv(U^T(1))*P(1)* ... *inv(U^T(m))*P(m)
387 *
388  IF( lsame( uplo, 'U' ) ) THEN
389 *
390 * Loop backward applying the transformations.
391 *
392  k = n
393  kc = n*( n+1 ) / 2 + 1
394  70 IF( k.LT.1 )
395  $ GO TO 90
396  kc = kc - k
397 *
398 * 1 x 1 pivot block.
399 *
400  IF( ipiv( k ).GT.0 ) THEN
401  IF( k.GT.1 ) THEN
402 *
403 * Interchange if P(K) != I.
404 *
405  kp = ipiv( k )
406  IF( kp.NE.k )
407  $ CALL cswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
408 *
409 * Apply the transformation:
410 * y := y - B' * conjg(x)
411 * where x is a column of A and y is a row of B.
412 *
413  CALL cgemv( 'Transpose', k-1, nrhs, one, b, ldb,
414  $ a( kc ), 1, one, b( k, 1 ), ldb )
415  END IF
416  IF( nounit )
417  $ CALL cscal( nrhs, a( kc+k-1 ), b( k, 1 ), ldb )
418  k = k - 1
419 *
420 * 2 x 2 pivot block.
421 *
422  ELSE
423  kcnext = kc - ( k-1 )
424  IF( k.GT.2 ) THEN
425 *
426 * Interchange if P(K) != I.
427 *
428  kp = abs( ipiv( k ) )
429  IF( kp.NE.k-1 )
430  $ CALL cswap( nrhs, b( k-1, 1 ), ldb, b( kp, 1 ),
431  $ ldb )
432 *
433 * Apply the transformations.
434 *
435  CALL cgemv( 'Transpose', k-2, nrhs, one, b, ldb,
436  $ a( kc ), 1, one, b( k, 1 ), ldb )
437 *
438  CALL cgemv( 'Transpose', k-2, nrhs, one, b, ldb,
439  $ a( kcnext ), 1, one, b( k-1, 1 ), ldb )
440  END IF
441 *
442 * Multiply by the diagonal block if non-unit.
443 *
444  IF( nounit ) THEN
445  d11 = a( kc-1 )
446  d22 = a( kc+k-1 )
447  d12 = a( kc+k-2 )
448  d21 = d12
449  DO 80 j = 1, nrhs
450  t1 = b( k-1, j )
451  t2 = b( k, j )
452  b( k-1, j ) = d11*t1 + d12*t2
453  b( k, j ) = d21*t1 + d22*t2
454  80 CONTINUE
455  END IF
456  kc = kcnext
457  k = k - 2
458  END IF
459  GO TO 70
460  90 CONTINUE
461 *
462 * Form B := L^T*B
463 * where L = P(1)*inv(L(1))* ... *P(m)*inv(L(m))
464 * and L^T = inv(L(m))*P(m)* ... *inv(L(1))*P(1)
465 *
466  ELSE
467 *
468 * Loop forward applying the L-transformations.
469 *
470  k = 1
471  kc = 1
472  100 CONTINUE
473  IF( k.GT.n )
474  $ GO TO 120
475 *
476 * 1 x 1 pivot block
477 *
478  IF( ipiv( k ).GT.0 ) THEN
479  IF( k.LT.n ) THEN
480 *
481 * Interchange if P(K) != I.
482 *
483  kp = ipiv( k )
484  IF( kp.NE.k )
485  $ CALL cswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
486 *
487 * Apply the transformation
488 *
489  CALL cgemv( 'Transpose', n-k, nrhs, one, b( k+1, 1 ),
490  $ ldb, a( kc+1 ), 1, one, b( k, 1 ), ldb )
491  END IF
492  IF( nounit )
493  $ CALL cscal( nrhs, a( kc ), b( k, 1 ), ldb )
494  kc = kc + n - k + 1
495  k = k + 1
496 *
497 * 2 x 2 pivot block.
498 *
499  ELSE
500  kcnext = kc + n - k + 1
501  IF( k.LT.n-1 ) THEN
502 *
503 * Interchange if P(K) != I.
504 *
505  kp = abs( ipiv( k ) )
506  IF( kp.NE.k+1 )
507  $ CALL cswap( nrhs, b( k+1, 1 ), ldb, b( kp, 1 ),
508  $ ldb )
509 *
510 * Apply the transformation
511 *
512  CALL cgemv( 'Transpose', n-k-1, nrhs, one,
513  $ b( k+2, 1 ), ldb, a( kcnext+1 ), 1, one,
514  $ b( k+1, 1 ), ldb )
515 *
516  CALL cgemv( 'Transpose', n-k-1, nrhs, one,
517  $ b( k+2, 1 ), ldb, a( kc+2 ), 1, one,
518  $ b( k, 1 ), ldb )
519  END IF
520 *
521 * Multiply by the diagonal block if non-unit.
522 *
523  IF( nounit ) THEN
524  d11 = a( kc )
525  d22 = a( kcnext )
526  d21 = a( kc+1 )
527  d12 = d21
528  DO 110 j = 1, nrhs
529  t1 = b( k, j )
530  t2 = b( k+1, j )
531  b( k, j ) = d11*t1 + d12*t2
532  b( k+1, j ) = d21*t1 + d22*t2
533  110 CONTINUE
534  END IF
535  kc = kcnext + ( n-k )
536  k = k + 2
537  END IF
538  GO TO 100
539  120 CONTINUE
540  END IF
541 *
542  END IF
543  RETURN
544 *
545 * End of CLAVSP
546 *
xerbla
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
cscal
subroutine cscal(N, CA, CX, INCX)
CSCAL
Definition: cscal.f:54
cgemv
subroutine cgemv(TRANS, M, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
CGEMV
Definition: cgemv.f:160
cswap
subroutine cswap(N, CX, INCX, CY, INCY)
CSWAP
Definition: cswap.f:52
cgeru
subroutine cgeru(M, N, ALPHA, X, INCX, Y, INCY, A, LDA)
CGERU
Definition: cgeru.f:132
lsame
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55

Here is the call graph for this function:

Here is the caller graph for this function: