SRC/dlasr.f

   1 *> \brief \b DLASR applies a sequence of plane rotations to a general rectangular matrix.
   2 *
   3 *  =========== DOCUMENTATION ===========
   4 *
   5 * Online html documentation available at
   6 *            http://www.netlib.org/lapack/explore-html/
   7 *
   8 *> \htmlonly
   9 *> Download DLASR + dependencies
  10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlasr.f">
  11 *> [TGZ]</a>
  12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlasr.f">
  13 *> [ZIP]</a>
  14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlasr.f">
  15 *> [TXT]</a>
  16 *> \endhtmlonly
  17 *
  18 *  Definition:
  19 *  ===========
  20 *
  21 *       SUBROUTINE DLASR( SIDE, PIVOT, DIRECT, M, N, C, S, A, LDA )
  22 *
  23 *       .. Scalar Arguments ..
  24 *       CHARACTER          DIRECT, PIVOT, SIDE
  25 *       INTEGER            LDA, M, N
  26 *       ..
  27 *       .. Array Arguments ..
  28 *       DOUBLE PRECISION   A( LDA, * ), C( * ), S( * )
  29 *       ..
  30 *
  31 *
  32 *> \par Purpose:
  33 *  =============
  34 *>
  35 *> \verbatim
  36 *>
  37 *> DLASR applies a sequence of plane rotations to a real matrix A,
  38 *> from either the left or the right.
  39 *>
  40 *> When SIDE = 'L', the transformation takes the form
  41 *>
  42 *>    A := P*A
  43 *>
  44 *> and when SIDE = 'R', the transformation takes the form
  45 *>
  46 *>    A := A*P**T
  47 *>
  48 *> where P is an orthogonal matrix consisting of a sequence of z plane
  49 *> rotations, with z = M when SIDE = 'L' and z = N when SIDE = 'R',
  50 *> and P**T is the transpose of P.
  51 *>
  52 *> When DIRECT = 'F' (Forward sequence), then
  53 *>
  54 *>    P = P(z-1) * ... * P(2) * P(1)
  55 *>
  56 *> and when DIRECT = 'B' (Backward sequence), then
  57 *>
  58 *>    P = P(1) * P(2) * ... * P(z-1)
  59 *>
  60 *> where P(k) is a plane rotation matrix defined by the 2-by-2 rotation
  61 *>
  62 *>    R(k) = (  c(k)  s(k) )
  63 *>         = ( -s(k)  c(k) ).
  64 *>
  65 *> When PIVOT = 'V' (Variable pivot), the rotation is performed
  66 *> for the plane (k,k+1), i.e., P(k) has the form
  67 *>
  68 *>    P(k) = (  1                                            )
  69 *>           (       ...                                     )
  70 *>           (              1                                )
  71 *>           (                   c(k)  s(k)                  )
  72 *>           (                  -s(k)  c(k)                  )
  73 *>           (                                1              )
  74 *>           (                                     ...       )
  75 *>           (                                            1  )
  76 *>
  77 *> where R(k) appears as a rank-2 modification to the identity matrix in
  78 *> rows and columns k and k+1.
  79 *>
  80 *> When PIVOT = 'T' (Top pivot), the rotation is performed for the
  81 *> plane (1,k+1), so P(k) has the form
  82 *>
  83 *>    P(k) = (  c(k)                    s(k)                 )
  84 *>           (         1                                     )
  85 *>           (              ...                              )
  86 *>           (                     1                         )
  87 *>           ( -s(k)                    c(k)                 )
  88 *>           (                                 1             )
  89 *>           (                                      ...      )
  90 *>           (                                             1 )
  91 *>
  92 *> where R(k) appears in rows and columns 1 and k+1.
  93 *>
  94 *> Similarly, when PIVOT = 'B' (Bottom pivot), the rotation is
  95 *> performed for the plane (k,z), giving P(k) the form
  96 *>
  97 *>    P(k) = ( 1                                             )
  98 *>           (      ...                                      )
  99 *>           (             1                                 )
 100 *>           (                  c(k)                    s(k) )
 101 *>           (                         1                     )
 102 *>           (                              ...              )
 103 *>           (                                     1         )
 104 *>           (                 -s(k)                    c(k) )
 105 *>
 106 *> where R(k) appears in rows and columns k and z.  The rotations are
 107 *> performed without ever forming P(k) explicitly.
 108 *> \endverbatim
 109 *
 110 *  Arguments:
 111 *  ==========
 112 *
 113 *> \param[in] SIDE
 114 *> \verbatim
 115 *>          SIDE is CHARACTER*1
 116 *>          Specifies whether the plane rotation matrix P is applied to
 117 *>          A on the left or the right.
 118 *>          = 'L':  Left, compute A := P*A
 119 *>          = 'R':  Right, compute A:= A*P**T
 120 *> \endverbatim
 121 *>
 122 *> \param[in] PIVOT
 123 *> \verbatim
 124 *>          PIVOT is CHARACTER*1
 125 *>          Specifies the plane for which P(k) is a plane rotation
 126 *>          matrix.
 127 *>          = 'V':  Variable pivot, the plane (k,k+1)
 128 *>          = 'T':  Top pivot, the plane (1,k+1)
 129 *>          = 'B':  Bottom pivot, the plane (k,z)
 130 *> \endverbatim
 131 *>
 132 *> \param[in] DIRECT
 133 *> \verbatim
 134 *>          DIRECT is CHARACTER*1
 135 *>          Specifies whether P is a forward or backward sequence of
 136 *>          plane rotations.
 137 *>          = 'F':  Forward, P = P(z-1)*...*P(2)*P(1)
 138 *>          = 'B':  Backward, P = P(1)*P(2)*...*P(z-1)
 139 *> \endverbatim
 140 *>
 141 *> \param[in] M
 142 *> \verbatim
 143 *>          M is INTEGER
 144 *>          The number of rows of the matrix A.  If m <= 1, an immediate
 145 *>          return is effected.
 146 *> \endverbatim
 147 *>
 148 *> \param[in] N
 149 *> \verbatim
 150 *>          N is INTEGER
 151 *>          The number of columns of the matrix A.  If n <= 1, an
 152 *>          immediate return is effected.
 153 *> \endverbatim
 154 *>
 155 *> \param[in] C
 156 *> \verbatim
 157 *>          C is DOUBLE PRECISION array, dimension
 158 *>                  (M-1) if SIDE = 'L'
 159 *>                  (N-1) if SIDE = 'R'
 160 *>          The cosines c(k) of the plane rotations.
 161 *> \endverbatim
 162 *>
 163 *> \param[in] S
 164 *> \verbatim
 165 *>          S is DOUBLE PRECISION array, dimension
 166 *>                  (M-1) if SIDE = 'L'
 167 *>                  (N-1) if SIDE = 'R'
 168 *>          The sines s(k) of the plane rotations.  The 2-by-2 plane
 169 *>          rotation part of the matrix P(k), R(k), has the form
 170 *>          R(k) = (  c(k)  s(k) )
 171 *>                 ( -s(k)  c(k) ).
 172 *> \endverbatim
 173 *>
 174 *> \param[in,out] A
 175 *> \verbatim
 176 *>          A is DOUBLE PRECISION array, dimension (LDA,N)
 177 *>          The M-by-N matrix A.  On exit, A is overwritten by P*A if
 178 *>          SIDE = 'R' or by A*P**T if SIDE = 'L'.
 179 *> \endverbatim
 180 *>
 181 *> \param[in] LDA
 182 *> \verbatim
 183 *>          LDA is INTEGER
 184 *>          The leading dimension of the array A.  LDA >= max(1,M).
 185 *> \endverbatim
 186 *
 187 *  Authors:
 188 *  ========
 189 *
 190 *> \author Univ. of Tennessee
 191 *> \author Univ. of California Berkeley
 192 *> \author Univ. of Colorado Denver
 193 *> \author NAG Ltd.
 194 *
 195 *> \date September 2012
 196 *
 197 *> \ingroup auxOTHERauxiliary
 198 *
 199 *  =====================================================================
 200       SUBROUTINE DLASR( SIDE, PIVOT, DIRECT, M, N, C, S, A, LDA )
 201 *
 202 *  -- LAPACK auxiliary routine (version 3.4.2) --
 203 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
 204 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
 205 *     September 2012
 206 *
 207 *     .. Scalar Arguments ..
 208       CHARACTER          DIRECT, PIVOT, SIDE
 209       INTEGER            LDA, M, N
 210 *     ..
 211 *     .. Array Arguments ..
 212       DOUBLE PRECISION   A( LDA, * ), C( * ), S( * )
 213 *     ..
 214 *
 215 *  =====================================================================
 216 *
 217 *     .. Parameters ..
 218       DOUBLE PRECISION   ONE, ZERO
 219       PARAMETER          ( ONE = 1.0D+0, ZERO = 0.0D+0 )
 220 *     ..
 221 *     .. Local Scalars ..
 222       INTEGER            I, INFO, J
 223       DOUBLE PRECISION   CTEMP, STEMP, TEMP
 224 *     ..
 225 *     .. External Functions ..
 226       LOGICAL            LSAME
 227       EXTERNAL           LSAME
 228 *     ..
 229 *     .. External Subroutines ..
 230       EXTERNAL           XERBLA
 231 *     ..
 232 *     .. Intrinsic Functions ..
 233       INTRINSIC          MAX
 234 *     ..
 235 *     .. Executable Statements ..
 236 *
 237 *     Test the input parameters
 238 *
 239       INFO = 0
 240       IF( .NOT.( LSAME( SIDE, 'L' ) .OR. LSAME( SIDE, 'R' ) ) ) THEN
 241          INFO = 1
 242       ELSE IF( .NOT.( LSAME( PIVOT, 'V' ) .OR. LSAME( PIVOT,
 243      $         'T' ) .OR. LSAME( PIVOT, 'B' ) ) ) THEN
 244          INFO = 2
 245       ELSE IF( .NOT.( LSAME( DIRECT, 'F' ) .OR. LSAME( DIRECT, 'B' ) ) )
 246      $          THEN
 247          INFO = 3
 248       ELSE IF( M.LT.0 ) THEN
 249          INFO = 4
 250       ELSE IF( N.LT.0 ) THEN
 251          INFO = 5
 252       ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
 253          INFO = 9
 254       END IF
 255       IF( INFO.NE.0 ) THEN
 256          CALL XERBLA( 'DLASR ', INFO )
 257          RETURN
 258       END IF
 259 *
 260 *     Quick return if possible
 261 *
 262       IF( ( M.EQ.0 ) .OR. ( N.EQ.0 ) )
 263      $   RETURN
 264       IF( LSAME( SIDE, 'L' ) ) THEN
 265 *
 266 *        Form  P * A
 267 *
 268          IF( LSAME( PIVOT, 'V' ) ) THEN
 269             IF( LSAME( DIRECT, 'F' ) ) THEN
 270                DO 20 J = 1, M - 1
 271                   CTEMP = C( J )
 272                   STEMP = S( J )
 273                   IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
 274                      DO 10 I = 1, N
 275                         TEMP = A( J+1, I )
 276                         A( J+1, I ) = CTEMP*TEMP - STEMP*A( J, I )
 277                         A( J, I ) = STEMP*TEMP + CTEMP*A( J, I )
 278    10                CONTINUE
 279                   END IF
 280    20          CONTINUE
 281             ELSE IF( LSAME( DIRECT, 'B' ) ) THEN
 282                DO 40 J = M - 1, 1, -1
 283                   CTEMP = C( J )
 284                   STEMP = S( J )
 285                   IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
 286                      DO 30 I = 1, N
 287                         TEMP = A( J+1, I )
 288                         A( J+1, I ) = CTEMP*TEMP - STEMP*A( J, I )
 289                         A( J, I ) = STEMP*TEMP + CTEMP*A( J, I )
 290    30                CONTINUE
 291                   END IF
 292    40          CONTINUE
 293             END IF
 294          ELSE IF( LSAME( PIVOT, 'T' ) ) THEN
 295             IF( LSAME( DIRECT, 'F' ) ) THEN
 296                DO 60 J = 2, M
 297                   CTEMP = C( J-1 )
 298                   STEMP = S( J-1 )
 299                   IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
 300                      DO 50 I = 1, N
 301                         TEMP = A( J, I )
 302                         A( J, I ) = CTEMP*TEMP - STEMP*A( 1, I )
 303                         A( 1, I ) = STEMP*TEMP + CTEMP*A( 1, I )
 304    50                CONTINUE
 305                   END IF
 306    60          CONTINUE
 307             ELSE IF( LSAME( DIRECT, 'B' ) ) THEN
 308                DO 80 J = M, 2, -1
 309                   CTEMP = C( J-1 )
 310                   STEMP = S( J-1 )
 311                   IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
 312                      DO 70 I = 1, N
 313                         TEMP = A( J, I )
 314                         A( J, I ) = CTEMP*TEMP - STEMP*A( 1, I )
 315                         A( 1, I ) = STEMP*TEMP + CTEMP*A( 1, I )
 316    70                CONTINUE
 317                   END IF
 318    80          CONTINUE
 319             END IF
 320          ELSE IF( LSAME( PIVOT, 'B' ) ) THEN
 321             IF( LSAME( DIRECT, 'F' ) ) THEN
 322                DO 100 J = 1, M - 1
 323                   CTEMP = C( J )
 324                   STEMP = S( J )
 325                   IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
 326                      DO 90 I = 1, N
 327                         TEMP = A( J, I )
 328                         A( J, I ) = STEMP*A( M, I ) + CTEMP*TEMP
 329                         A( M, I ) = CTEMP*A( M, I ) - STEMP*TEMP
 330    90                CONTINUE
 331                   END IF
 332   100          CONTINUE
 333             ELSE IF( LSAME( DIRECT, 'B' ) ) THEN
 334                DO 120 J = M - 1, 1, -1
 335                   CTEMP = C( J )
 336                   STEMP = S( J )
 337                   IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
 338                      DO 110 I = 1, N
 339                         TEMP = A( J, I )
 340                         A( J, I ) = STEMP*A( M, I ) + CTEMP*TEMP
 341                         A( M, I ) = CTEMP*A( M, I ) - STEMP*TEMP
 342   110                CONTINUE
 343                   END IF
 344   120          CONTINUE
 345             END IF
 346          END IF
 347       ELSE IF( LSAME( SIDE, 'R' ) ) THEN
 348 *
 349 *        Form A * P**T
 350 *
 351          IF( LSAME( PIVOT, 'V' ) ) THEN
 352             IF( LSAME( DIRECT, 'F' ) ) THEN
 353                DO 140 J = 1, N - 1
 354                   CTEMP = C( J )
 355                   STEMP = S( J )
 356                   IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
 357                      DO 130 I = 1, M
 358                         TEMP = A( I, J+1 )
 359                         A( I, J+1 ) = CTEMP*TEMP - STEMP*A( I, J )
 360                         A( I, J ) = STEMP*TEMP + CTEMP*A( I, J )
 361   130                CONTINUE
 362                   END IF
 363   140          CONTINUE
 364             ELSE IF( LSAME( DIRECT, 'B' ) ) THEN
 365                DO 160 J = N - 1, 1, -1
 366                   CTEMP = C( J )
 367                   STEMP = S( J )
 368                   IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
 369                      DO 150 I = 1, M
 370                         TEMP = A( I, J+1 )
 371                         A( I, J+1 ) = CTEMP*TEMP - STEMP*A( I, J )
 372                         A( I, J ) = STEMP*TEMP + CTEMP*A( I, J )
 373   150                CONTINUE
 374                   END IF
 375   160          CONTINUE
 376             END IF
 377          ELSE IF( LSAME( PIVOT, 'T' ) ) THEN
 378             IF( LSAME( DIRECT, 'F' ) ) THEN
 379                DO 180 J = 2, N
 380                   CTEMP = C( J-1 )
 381                   STEMP = S( J-1 )
 382                   IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
 383                      DO 170 I = 1, M
 384                         TEMP = A( I, J )
 385                         A( I, J ) = CTEMP*TEMP - STEMP*A( I, 1 )
 386                         A( I, 1 ) = STEMP*TEMP + CTEMP*A( I, 1 )
 387   170                CONTINUE
 388                   END IF
 389   180          CONTINUE
 390             ELSE IF( LSAME( DIRECT, 'B' ) ) THEN
 391                DO 200 J = N, 2, -1
 392                   CTEMP = C( J-1 )
 393                   STEMP = S( J-1 )
 394                   IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
 395                      DO 190 I = 1, M
 396                         TEMP = A( I, J )
 397                         A( I, J ) = CTEMP*TEMP - STEMP*A( I, 1 )
 398                         A( I, 1 ) = STEMP*TEMP + CTEMP*A( I, 1 )
 399   190                CONTINUE
 400                   END IF
 401   200          CONTINUE
 402             END IF
 403          ELSE IF( LSAME( PIVOT, 'B' ) ) THEN
 404             IF( LSAME( DIRECT, 'F' ) ) THEN
 405                DO 220 J = 1, N - 1
 406                   CTEMP = C( J )
 407                   STEMP = S( J )
 408                   IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
 409                      DO 210 I = 1, M
 410                         TEMP = A( I, J )
 411                         A( I, J ) = STEMP*A( I, N ) + CTEMP*TEMP
 412                         A( I, N ) = CTEMP*A( I, N ) - STEMP*TEMP
 413   210                CONTINUE
 414                   END IF
 415   220          CONTINUE
 416             ELSE IF( LSAME( DIRECT, 'B' ) ) THEN
 417                DO 240 J = N - 1, 1, -1
 418                   CTEMP = C( J )
 419                   STEMP = S( J )
 420                   IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
 421                      DO 230 I = 1, M
 422                         TEMP = A( I, J )
 423                         A( I, J ) = STEMP*A( I, N ) + CTEMP*TEMP
 424                         A( I, N ) = CTEMP*A( I, N ) - STEMP*TEMP
 425   230                CONTINUE
 426                   END IF
 427   240          CONTINUE
 428             END IF
 429          END IF
 430       END IF
 431 *
 432       RETURN
 433 *
 434 *     End of DLASR
 435 *
 436       END