SRC/dormtr.f

   1 *> \brief \b DORMTR
   2 *
   3 *  =========== DOCUMENTATION ===========
   4 *
   5 * Online html documentation available at
   6 *            http://www.netlib.org/lapack/explore-html/
   7 *
   8 *> \htmlonly
   9 *> Download DORMTR + dependencies
  10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dormtr.f">
  11 *> [TGZ]</a>
  12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dormtr.f">
  13 *> [ZIP]</a>
  14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dormtr.f">
  15 *> [TXT]</a>
  16 *> \endhtmlonly
  17 *
  18 *  Definition:
  19 *  ===========
  20 *
  21 *       SUBROUTINE DORMTR( SIDE, UPLO, TRANS, M, N, A, LDA, TAU, C, LDC,
  22 *                          WORK, LWORK, INFO )
  23 *
  24 *       .. Scalar Arguments ..
  25 *       CHARACTER          SIDE, TRANS, UPLO
  26 *       INTEGER            INFO, LDA, LDC, LWORK, M, N
  27 *       ..
  28 *       .. Array Arguments ..
  29 *       DOUBLE PRECISION   A( LDA, * ), C( LDC, * ), TAU( * ), WORK( * )
  30 *       ..
  31 *
  32 *
  33 *> \par Purpose:
  34 *  =============
  35 *>
  36 *> \verbatim
  37 *>
  38 *> DORMTR overwrites the general real M-by-N matrix C with
  39 *>
  40 *>                 SIDE = 'L'     SIDE = 'R'
  41 *> TRANS = 'N':      Q * C          C * Q
  42 *> TRANS = 'T':      Q**T * C       C * Q**T
  43 *>
  44 *> where Q is a real orthogonal matrix of order nq, with nq = m if
  45 *> SIDE = 'L' and nq = n if SIDE = 'R'. Q is defined as the product of
  46 *> nq-1 elementary reflectors, as returned by DSYTRD:
  47 *>
  48 *> if UPLO = 'U', Q = H(nq-1) . . . H(2) H(1);
  49 *>
  50 *> if UPLO = 'L', Q = H(1) H(2) . . . H(nq-1).
  51 *> \endverbatim
  52 *
  53 *  Arguments:
  54 *  ==========
  55 *
  56 *> \param[in] SIDE
  57 *> \verbatim
  58 *>          SIDE is CHARACTER*1
  59 *>          = 'L': apply Q or Q**T from the Left;
  60 *>          = 'R': apply Q or Q**T from the Right.
  61 *> \endverbatim
  62 *>
  63 *> \param[in] UPLO
  64 *> \verbatim
  65 *>          UPLO is CHARACTER*1
  66 *>          = 'U': Upper triangle of A contains elementary reflectors
  67 *>                 from DSYTRD;
  68 *>          = 'L': Lower triangle of A contains elementary reflectors
  69 *>                 from DSYTRD.
  70 *> \endverbatim
  71 *>
  72 *> \param[in] TRANS
  73 *> \verbatim
  74 *>          TRANS is CHARACTER*1
  75 *>          = 'N':  No transpose, apply Q;
  76 *>          = 'T':  Transpose, apply Q**T.
  77 *> \endverbatim
  78 *>
  79 *> \param[in] M
  80 *> \verbatim
  81 *>          M is INTEGER
  82 *>          The number of rows of the matrix C. M >= 0.
  83 *> \endverbatim
  84 *>
  85 *> \param[in] N
  86 *> \verbatim
  87 *>          N is INTEGER
  88 *>          The number of columns of the matrix C. N >= 0.
  89 *> \endverbatim
  90 *>
  91 *> \param[in] A
  92 *> \verbatim
  93 *>          A is DOUBLE PRECISION array, dimension
  94 *>                               (LDA,M) if SIDE = 'L'
  95 *>                               (LDA,N) if SIDE = 'R'
  96 *>          The vectors which define the elementary reflectors, as
  97 *>          returned by DSYTRD.
  98 *> \endverbatim
  99 *>
 100 *> \param[in] LDA
 101 *> \verbatim
 102 *>          LDA is INTEGER
 103 *>          The leading dimension of the array A.
 104 *>          LDA >= max(1,M) if SIDE = 'L'; LDA >= max(1,N) if SIDE = 'R'.
 105 *> \endverbatim
 106 *>
 107 *> \param[in] TAU
 108 *> \verbatim
 109 *>          TAU is DOUBLE PRECISION array, dimension
 110 *>                               (M-1) if SIDE = 'L'
 111 *>                               (N-1) if SIDE = 'R'
 112 *>          TAU(i) must contain the scalar factor of the elementary
 113 *>          reflector H(i), as returned by DSYTRD.
 114 *> \endverbatim
 115 *>
 116 *> \param[in,out] C
 117 *> \verbatim
 118 *>          C is DOUBLE PRECISION array, dimension (LDC,N)
 119 *>          On entry, the M-by-N matrix C.
 120 *>          On exit, C is overwritten by Q*C or Q**T*C or C*Q**T or C*Q.
 121 *> \endverbatim
 122 *>
 123 *> \param[in] LDC
 124 *> \verbatim
 125 *>          LDC is INTEGER
 126 *>          The leading dimension of the array C. LDC >= max(1,M).
 127 *> \endverbatim
 128 *>
 129 *> \param[out] WORK
 130 *> \verbatim
 131 *>          WORK is DOUBLE PRECISION array, dimension (MAX(1,LWORK))
 132 *>          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
 133 *> \endverbatim
 134 *>
 135 *> \param[in] LWORK
 136 *> \verbatim
 137 *>          LWORK is INTEGER
 138 *>          The dimension of the array WORK.
 139 *>          If SIDE = 'L', LWORK >= max(1,N);
 140 *>          if SIDE = 'R', LWORK >= max(1,M).
 141 *>          For optimum performance LWORK >= N*NB if SIDE = 'L', and
 142 *>          LWORK >= M*NB if SIDE = 'R', where NB is the optimal
 143 *>          blocksize.
 144 *>
 145 *>          If LWORK = -1, then a workspace query is assumed; the routine
 146 *>          only calculates the optimal size of the WORK array, returns
 147 *>          this value as the first entry of the WORK array, and no error
 148 *>          message related to LWORK is issued by XERBLA.
 149 *> \endverbatim
 150 *>
 151 *> \param[out] INFO
 152 *> \verbatim
 153 *>          INFO is INTEGER
 154 *>          = 0:  successful exit
 155 *>          < 0:  if INFO = -i, the i-th argument had an illegal value
 156 *> \endverbatim
 157 *
 158 *  Authors:
 159 *  ========
 160 *
 161 *> \author Univ. of Tennessee
 162 *> \author Univ. of California Berkeley
 163 *> \author Univ. of Colorado Denver
 164 *> \author NAG Ltd.
 165 *
 166 *> \date November 2011
 167 *
 168 *> \ingroup doubleOTHERcomputational
 169 *
 170 *  =====================================================================
 171       SUBROUTINE DORMTR( SIDE, UPLO, TRANS, M, N, A, LDA, TAU, C, LDC,
 172      $                   WORK, LWORK, INFO )
 173 *
 174 *  -- LAPACK computational routine (version 3.4.0) --
 175 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
 176 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
 177 *     November 2011
 178 *
 179 *     .. Scalar Arguments ..
 180       CHARACTER          SIDE, TRANS, UPLO
 181       INTEGER            INFO, LDA, LDC, LWORK, M, N
 182 *     ..
 183 *     .. Array Arguments ..
 184       DOUBLE PRECISION   A( LDA, * ), C( LDC, * ), TAU( * ), WORK( * )
 185 *     ..
 186 *
 187 *  =====================================================================
 188 *
 189 *     .. Local Scalars ..
 190       LOGICAL            LEFT, LQUERY, UPPER
 191       INTEGER            I1, I2, IINFO, LWKOPT, MI, NB, NI, NQ, NW
 192 *     ..
 193 *     .. External Functions ..
 194       LOGICAL            LSAME
 195       INTEGER            ILAENV
 196       EXTERNAL           LSAME, ILAENV
 197 *     ..
 198 *     .. External Subroutines ..
 199       EXTERNAL           DORMQL, DORMQR, XERBLA
 200 *     ..
 201 *     .. Intrinsic Functions ..
 202       INTRINSIC          MAX
 203 *     ..
 204 *     .. Executable Statements ..
 205 *
 206 *     Test the input arguments
 207 *
 208       INFO = 0
 209       LEFT = LSAME( SIDE, 'L' )
 210       UPPER = LSAME( UPLO, 'U' )
 211       LQUERY = ( LWORK.EQ.-1 )
 212 *
 213 *     NQ is the order of Q and NW is the minimum dimension of WORK
 214 *
 215       IF( LEFT ) THEN
 216          NQ = M
 217          NW = N
 218       ELSE
 219          NQ = N
 220          NW = M
 221       END IF
 222       IF( .NOT.LEFT .AND. .NOT.LSAME( SIDE, 'R' ) ) THEN
 223          INFO = -1
 224       ELSE IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN
 225          INFO = -2
 226       ELSE IF( .NOT.LSAME( TRANS, 'N' ) .AND. .NOT.LSAME( TRANS, 'T' ) )
 227      $          THEN
 228          INFO = -3
 229       ELSE IF( M.LT.0 ) THEN
 230          INFO = -4
 231       ELSE IF( N.LT.0 ) THEN
 232          INFO = -5
 233       ELSE IF( LDA.LT.MAX( 1, NQ ) ) THEN
 234          INFO = -7
 235       ELSE IF( LDC.LT.MAX( 1, M ) ) THEN
 236          INFO = -10
 237       ELSE IF( LWORK.LT.MAX( 1, NW ) .AND. .NOT.LQUERY ) THEN
 238          INFO = -12
 239       END IF
 240 *
 241       IF( INFO.EQ.0 ) THEN
 242          IF( UPPER ) THEN
 243             IF( LEFT ) THEN
 244                NB = ILAENV( 1, 'DORMQL', SIDE // TRANS, M-1, N, M-1,
 245      $              -1 )
 246             ELSE
 247                NB = ILAENV( 1, 'DORMQL', SIDE // TRANS, M, N-1, N-1,
 248      $              -1 )
 249             END IF
 250          ELSE
 251             IF( LEFT ) THEN
 252                NB = ILAENV( 1, 'DORMQR', SIDE // TRANS, M-1, N, M-1,
 253      $              -1 )
 254             ELSE
 255                NB = ILAENV( 1, 'DORMQR', SIDE // TRANS, M, N-1, N-1,
 256      $              -1 )
 257             END IF
 258          END IF
 259          LWKOPT = MAX( 1, NW )*NB
 260          WORK( 1 ) = LWKOPT
 261       END IF
 262 *
 263       IF( INFO.NE.0 ) THEN
 264          CALL XERBLA( 'DORMTR', -INFO )
 265          RETURN
 266       ELSE IF( LQUERY ) THEN
 267          RETURN
 268       END IF
 269 *
 270 *     Quick return if possible
 271 *
 272       IF( M.EQ.0 .OR. N.EQ.0 .OR. NQ.EQ.1 ) THEN
 273          WORK( 1 ) = 1
 274          RETURN
 275       END IF
 276 *
 277       IF( LEFT ) THEN
 278          MI = M - 1
 279          NI = N
 280       ELSE
 281          MI = M
 282          NI = N - 1
 283       END IF
 284 *
 285       IF( UPPER ) THEN
 286 *
 287 *        Q was determined by a call to DSYTRD with UPLO = 'U'
 288 *
 289          CALL DORMQL( SIDE, TRANS, MI, NI, NQ-1, A( 1, 2 ), LDA, TAU, C,
 290      $                LDC, WORK, LWORK, IINFO )
 291       ELSE
 292 *
 293 *        Q was determined by a call to DSYTRD with UPLO = 'L'
 294 *
 295          IF( LEFT ) THEN
 296             I1 = 2
 297             I2 = 1
 298          ELSE
 299             I1 = 1
 300             I2 = 2
 301          END IF
 302          CALL DORMQR( SIDE, TRANS, MI, NI, NQ-1, A( 2, 1 ), LDA, TAU,
 303      $                C( I1, I2 ), LDC, WORK, LWORK, IINFO )
 304       END IF
 305       WORK( 1 ) = LWKOPT
 306       RETURN
 307 *
 308 *     End of DORMTR
 309 *
 310       END