bdlaapp.f

Go to the documentation of this file.

      SUBROUTINE bdlaapp( ISIDE, M, N, NB, A, LDA, NITRAF, ITRAF,
     $                    DTRAF, WORK )
      IMPLICIT NONE
*
*     .. Scalar Arguments ..
      INTEGER            ISIDE, LDA, M, N, NB, NITRAF
*     ..
*     .. Array Arguments ..
      INTEGER            ITRAF( * )
      DOUBLE PRECISION   A( LDA, * ), DTRAF( * ), WORK( * )
*
*
*  Purpose
*  =======
*
*  BDLAAPP computes
*
*          B = Q**T * A       or       B = A * Q,
*
*  where A is an M-by-N matrix and Q is an orthogonal matrix represented
*  by the parameters in the arrays ITRAF and DTRAF as described in
*  BDTREXC.
*
*  This is an auxiliary routine called by BDTRSEN.
*
*  Arguments
*  =========
*
*  ISIDE   (input) INTEGER
*          Specifies whether Q multiplies A from the left or right as
*          follows:
*          = 0: compute B = Q**T * A;
*          = 1: compute B = A * Q.
*
*  M       (input) INTEGER
*          The number of rows of A.
*
*  N       (input) INTEGER
*          The number of columns of A.
*
*  NB      (input) INTEGER
*          If ISIDE = 0, the Q is applied block column-wise to the rows
*          of A and NB specifies the maximal width of the block columns.
*          If ISIDE = 1, this variable is not referenced.
*
*  A       (input/output) DOUBLE PRECISION array, dimension (LDA,N)
*          On entry, the matrix A.
*          On exit, A is overwritten by B.
*
*  LDA     (input) INTEGER
*          The leading dimension of the array A. LDA >= max(1,N).
*
*  NITRAF  (input) INTEGER
*          Length of the array ITRAF. NITRAF >= 0.
*
*  ITRAF   (input) INTEGER array, length NITRAF
*          List of parameters for representing the transformation
*          matrix Q, see BDTREXC.
*
*  DTRAF   (output) DOUBLE PRECISION array, length k, where
*          List of parameters for representing the transformation
*          matrix Q, see BDTREXC.
*
*  WORK    (workspace) DOUBLE PRECISION array, dimension (N)
*
*  =====================================================================
*
 
*     .. Parameters ..
      DOUBLE PRECISION   ZERO, ONE
      parameter( zero = 0.0d+0, one = 1.0d+0 )
*     ..
*     .. Local Scalars ..
      INTEGER            I, IT, J, NNB, PD
      DOUBLE PRECISION   TAU
*     ..
*     .. External Subroutines ..
      EXTERNAL           dlarfx, drot
*     .. Intrinsic Functions ..
      INTRINSIC          min
*     ..
*     .. Executable Statements ..
*
*     Quick return if possible.
*
      IF( m.LE.0 .OR. n.LE.0 )
     $   RETURN
*
      IF( iside.EQ.0 ) THEN
*
*        Apply Q from left.
*
         DO 20 j = 1, n, nb
            pd = 1
            nnb = min( nb, n - j + 1 )
            DO 10 i = 1, nitraf
               it = itraf(i)
               IF( it.LE.m ) THEN
*
*                 Apply Givens rotation.
*
                  CALL drot( nnb, a(it,j), lda, a(it+1,j), lda,
     $                       dtraf(pd), dtraf(pd+1) )
                  pd = pd + 2
               ELSE IF( it.LE.2*m ) THEN
*
*                 Apply Householder reflector of first kind.
*
                  tau = dtraf(pd)
                  dtraf(pd) = one
                  CALL dlarfx( 'Left', 3, nnb, dtraf(pd), tau,
     $                         a(it-m,j), lda, work )
                  dtraf(pd) = tau
                  pd = pd + 3
               ELSE
*
*                 Apply Householder reflector of second kind.
*
                  tau = dtraf(pd+2)
                  dtraf(pd+2) = one
                  CALL dlarfx( 'Left', 3, nnb, dtraf(pd), tau,
     $                         a(it-2*m,j), lda, work )
                  dtraf(pd+2) = tau
                  pd = pd + 3
               END IF
   10       CONTINUE
   20    CONTINUE
      ELSE
         pd = 1
         DO 30 i = 1, nitraf
            it = itraf(i)
            IF( it.LE.n ) THEN
*
*              Apply Givens rotation.
*
               CALL drot( m, a(1,it), 1, a(1,it+1), 1, dtraf(pd),
     $                    dtraf(pd+1) )
               pd = pd + 2
            ELSE IF( it.LE.2*n ) THEN
*
*              Apply Householder reflector of first kind.
*
               tau = dtraf(pd)
               dtraf(pd) = one
               CALL dlarfx( 'Right', m, 3, dtraf(pd), tau, a(1,it-n),
     $                      lda, work )
               dtraf(pd) = tau
               pd = pd + 3
            ELSE
*
*              Apply Householder reflector of second kind.
*
               tau = dtraf(pd+2)
               dtraf(pd+2) = one
               CALL dlarfx( 'Right', m, 3, dtraf(pd), tau, a(1,it-2*n),
     $                      lda, work )
               dtraf(pd+2) = tau
               pd = pd + 3
            END IF
   30    CONTINUE
      END IF
*
      RETURN
*
*     End of BDLAAPP
*
      END