LAPACK: SRC/clasr.f Source File

Go to the documentation of this file.
00001       SUBROUTINE CLASR( SIDE, PIVOT, DIRECT, M, N, C, S, A, LDA )
00002 *
00003 *  -- LAPACK auxiliary routine (version 3.2) --
00004 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
00005 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
00006 *     November 2006
00007 *
00008 *     .. Scalar Arguments ..
00009       CHARACTER          DIRECT, PIVOT, SIDE
00010       INTEGER            LDA, M, N
00011 *     ..
00012 *     .. Array Arguments ..
00013       REAL               C( * ), S( * )
00014       COMPLEX            A( LDA, * )
00015 *     ..
00016 *
00017 *  Purpose
00018 *  =======
00019 *
00020 *  CLASR applies a sequence of real plane rotations to a complex matrix
00021 *  A, from either the left or the right.
00022 *
00023 *  When SIDE = 'L', the transformation takes the form
00024 *
00025 *     A := P*A
00026 *
00027 *  and when SIDE = 'R', the transformation takes the form
00028 *
00029 *     A := A*P**T
00030 *
00031 *  where P is an orthogonal matrix consisting of a sequence of z plane
00032 *  rotations, with z = M when SIDE = 'L' and z = N when SIDE = 'R',
00033 *  and P**T is the transpose of P.
00034 *  
00035 *  When DIRECT = 'F' (Forward sequence), then
00036 *  
00037 *     P = P(z-1) * ... * P(2) * P(1)
00038 *  
00039 *  and when DIRECT = 'B' (Backward sequence), then
00040 *  
00041 *     P = P(1) * P(2) * ... * P(z-1)
00042 *  
00043 *  where P(k) is a plane rotation matrix defined by the 2-by-2 rotation
00044 *  
00045 *     R(k) = (  c(k)  s(k) )
00046 *          = ( -s(k)  c(k) ).
00047 *  
00048 *  When PIVOT = 'V' (Variable pivot), the rotation is performed
00049 *  for the plane (k,k+1), i.e., P(k) has the form
00050 *  
00051 *     P(k) = (  1                                            )
00052 *            (       ...                                     )
00053 *            (              1                                )
00054 *            (                   c(k)  s(k)                  )
00055 *            (                  -s(k)  c(k)                  )
00056 *            (                                1              )
00057 *            (                                     ...       )
00058 *            (                                            1  )
00059 *  
00060 *  where R(k) appears as a rank-2 modification to the identity matrix in
00061 *  rows and columns k and k+1.
00062 *  
00063 *  When PIVOT = 'T' (Top pivot), the rotation is performed for the
00064 *  plane (1,k+1), so P(k) has the form
00065 *  
00066 *     P(k) = (  c(k)                    s(k)                 )
00067 *            (         1                                     )
00068 *            (              ...                              )
00069 *            (                     1                         )
00070 *            ( -s(k)                    c(k)                 )
00071 *            (                                 1             )
00072 *            (                                      ...      )
00073 *            (                                             1 )
00074 *  
00075 *  where R(k) appears in rows and columns 1 and k+1.
00076 *  
00077 *  Similarly, when PIVOT = 'B' (Bottom pivot), the rotation is
00078 *  performed for the plane (k,z), giving P(k) the form
00079 *  
00080 *     P(k) = ( 1                                             )
00081 *            (      ...                                      )
00082 *            (             1                                 )
00083 *            (                  c(k)                    s(k) )
00084 *            (                         1                     )
00085 *            (                              ...              )
00086 *            (                                     1         )
00087 *            (                 -s(k)                    c(k) )
00088 *  
00089 *  where R(k) appears in rows and columns k and z.  The rotations are
00090 *  performed without ever forming P(k) explicitly.
00091 *
00092 *  Arguments
00093 *  =========
00094 *
00095 *  SIDE    (input) CHARACTER*1
00096 *          Specifies whether the plane rotation matrix P is applied to
00097 *          A on the left or the right.
00098 *          = 'L':  Left, compute A := P*A
00099 *          = 'R':  Right, compute A:= A*P**T
00100 *
00101 *  PIVOT   (input) CHARACTER*1
00102 *          Specifies the plane for which P(k) is a plane rotation
00103 *          matrix.
00104 *          = 'V':  Variable pivot, the plane (k,k+1)
00105 *          = 'T':  Top pivot, the plane (1,k+1)
00106 *          = 'B':  Bottom pivot, the plane (k,z)
00107 *
00108 *  DIRECT  (input) CHARACTER*1
00109 *          Specifies whether P is a forward or backward sequence of
00110 *          plane rotations.
00111 *          = 'F':  Forward, P = P(z-1)*...*P(2)*P(1)
00112 *          = 'B':  Backward, P = P(1)*P(2)*...*P(z-1)
00113 *
00114 *  M       (input) INTEGER
00115 *          The number of rows of the matrix A.  If m <= 1, an immediate
00116 *          return is effected.
00117 *
00118 *  N       (input) INTEGER
00119 *          The number of columns of the matrix A.  If n <= 1, an
00120 *          immediate return is effected.
00121 *
00122 *  C       (input) REAL array, dimension
00123 *                  (M-1) if SIDE = 'L'
00124 *                  (N-1) if SIDE = 'R'
00125 *          The cosines c(k) of the plane rotations.
00126 *
00127 *  S       (input) REAL array, dimension
00128 *                  (M-1) if SIDE = 'L'
00129 *                  (N-1) if SIDE = 'R'
00130 *          The sines s(k) of the plane rotations.  The 2-by-2 plane
00131 *          rotation part of the matrix P(k), R(k), has the form
00132 *          R(k) = (  c(k)  s(k) )
00133 *                 ( -s(k)  c(k) ).
00134 *
00135 *  A       (input/output) COMPLEX array, dimension (LDA,N)
00136 *          The M-by-N matrix A.  On exit, A is overwritten by P*A if
00137 *          SIDE = 'R' or by A*P**T if SIDE = 'L'.
00138 *
00139 *  LDA     (input) INTEGER
00140 *          The leading dimension of the array A.  LDA >= max(1,M).
00141 *
00142 *  =====================================================================
00143 *
00144 *     .. Parameters ..
00145       REAL               ONE, ZERO
00146       PARAMETER          ( ONE = 1.0E+0, ZERO = 0.0E+0 )
00147 *     ..
00148 *     .. Local Scalars ..
00149       INTEGER            I, INFO, J
00150       REAL               CTEMP, STEMP
00151       COMPLEX            TEMP
00152 *     ..
00153 *     .. Intrinsic Functions ..
00154       INTRINSIC          MAX
00155 *     ..
00156 *     .. External Functions ..
00157       LOGICAL            LSAME
00158       EXTERNAL           LSAME
00159 *     ..
00160 *     .. External Subroutines ..
00161       EXTERNAL           XERBLA
00162 *     ..
00163 *     .. Executable Statements ..
00164 *
00165 *     Test the input parameters
00166 *
00167       INFO = 0
00168       IF( .NOT.( LSAME( SIDE, 'L' ) .OR. LSAME( SIDE, 'R' ) ) ) THEN
00169          INFO = 1
00170       ELSE IF( .NOT.( LSAME( PIVOT, 'V' ) .OR. LSAME( PIVOT,
00171      $         'T' ) .OR. LSAME( PIVOT, 'B' ) ) ) THEN
00172          INFO = 2
00173       ELSE IF( .NOT.( LSAME( DIRECT, 'F' ) .OR. LSAME( DIRECT, 'B' ) ) )
00174      $          THEN
00175          INFO = 3
00176       ELSE IF( M.LT.0 ) THEN
00177          INFO = 4
00178       ELSE IF( N.LT.0 ) THEN
00179          INFO = 5
00180       ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
00181          INFO = 9
00182       END IF
00183       IF( INFO.NE.0 ) THEN
00184          CALL XERBLA( 'CLASR ', INFO )
00185          RETURN
00186       END IF
00187 *
00188 *     Quick return if possible
00189 *
00190       IF( ( M.EQ.0 ) .OR. ( N.EQ.0 ) )
00191      $   RETURN
00192       IF( LSAME( SIDE, 'L' ) ) THEN
00193 *
00194 *        Form  P * A
00195 *
00196          IF( LSAME( PIVOT, 'V' ) ) THEN
00197             IF( LSAME( DIRECT, 'F' ) ) THEN
00198                DO 20 J = 1, M - 1
00199                   CTEMP = C( J )
00200                   STEMP = S( J )
00201                   IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
00202                      DO 10 I = 1, N
00203                         TEMP = A( J+1, I )
00204                         A( J+1, I ) = CTEMP*TEMP - STEMP*A( J, I )
00205                         A( J, I ) = STEMP*TEMP + CTEMP*A( J, I )
00206    10                CONTINUE
00207                   END IF
00208    20          CONTINUE
00209             ELSE IF( LSAME( DIRECT, 'B' ) ) THEN
00210                DO 40 J = M - 1, 1, -1
00211                   CTEMP = C( J )
00212                   STEMP = S( J )
00213                   IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
00214                      DO 30 I = 1, N
00215                         TEMP = A( J+1, I )
00216                         A( J+1, I ) = CTEMP*TEMP - STEMP*A( J, I )
00217                         A( J, I ) = STEMP*TEMP + CTEMP*A( J, I )
00218    30                CONTINUE
00219                   END IF
00220    40          CONTINUE
00221             END IF
00222          ELSE IF( LSAME( PIVOT, 'T' ) ) THEN
00223             IF( LSAME( DIRECT, 'F' ) ) THEN
00224                DO 60 J = 2, M
00225                   CTEMP = C( J-1 )
00226                   STEMP = S( J-1 )
00227                   IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
00228                      DO 50 I = 1, N
00229                         TEMP = A( J, I )
00230                         A( J, I ) = CTEMP*TEMP - STEMP*A( 1, I )
00231                         A( 1, I ) = STEMP*TEMP + CTEMP*A( 1, I )
00232    50                CONTINUE
00233                   END IF
00234    60          CONTINUE
00235             ELSE IF( LSAME( DIRECT, 'B' ) ) THEN
00236                DO 80 J = M, 2, -1
00237                   CTEMP = C( J-1 )
00238                   STEMP = S( J-1 )
00239                   IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
00240                      DO 70 I = 1, N
00241                         TEMP = A( J, I )
00242                         A( J, I ) = CTEMP*TEMP - STEMP*A( 1, I )
00243                         A( 1, I ) = STEMP*TEMP + CTEMP*A( 1, I )
00244    70                CONTINUE
00245                   END IF
00246    80          CONTINUE
00247             END IF
00248          ELSE IF( LSAME( PIVOT, 'B' ) ) THEN
00249             IF( LSAME( DIRECT, 'F' ) ) THEN
00250                DO 100 J = 1, M - 1
00251                   CTEMP = C( J )
00252                   STEMP = S( J )
00253                   IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
00254                      DO 90 I = 1, N
00255                         TEMP = A( J, I )
00256                         A( J, I ) = STEMP*A( M, I ) + CTEMP*TEMP
00257                         A( M, I ) = CTEMP*A( M, I ) - STEMP*TEMP
00258    90                CONTINUE
00259                   END IF
00260   100          CONTINUE
00261             ELSE IF( LSAME( DIRECT, 'B' ) ) THEN
00262                DO 120 J = M - 1, 1, -1
00263                   CTEMP = C( J )
00264                   STEMP = S( J )
00265                   IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
00266                      DO 110 I = 1, N
00267                         TEMP = A( J, I )
00268                         A( J, I ) = STEMP*A( M, I ) + CTEMP*TEMP
00269                         A( M, I ) = CTEMP*A( M, I ) - STEMP*TEMP
00270   110                CONTINUE
00271                   END IF
00272   120          CONTINUE
00273             END IF
00274          END IF
00275       ELSE IF( LSAME( SIDE, 'R' ) ) THEN
00276 *
00277 *        Form A * P'
00278 *
00279          IF( LSAME( PIVOT, 'V' ) ) THEN
00280             IF( LSAME( DIRECT, 'F' ) ) THEN
00281                DO 140 J = 1, N - 1
00282                   CTEMP = C( J )
00283                   STEMP = S( J )
00284                   IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
00285                      DO 130 I = 1, M
00286                         TEMP = A( I, J+1 )
00287                         A( I, J+1 ) = CTEMP*TEMP - STEMP*A( I, J )
00288                         A( I, J ) = STEMP*TEMP + CTEMP*A( I, J )
00289   130                CONTINUE
00290                   END IF
00291   140          CONTINUE
00292             ELSE IF( LSAME( DIRECT, 'B' ) ) THEN
00293                DO 160 J = N - 1, 1, -1
00294                   CTEMP = C( J )
00295                   STEMP = S( J )
00296                   IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
00297                      DO 150 I = 1, M
00298                         TEMP = A( I, J+1 )
00299                         A( I, J+1 ) = CTEMP*TEMP - STEMP*A( I, J )
00300                         A( I, J ) = STEMP*TEMP + CTEMP*A( I, J )
00301   150                CONTINUE
00302                   END IF
00303   160          CONTINUE
00304             END IF
00305          ELSE IF( LSAME( PIVOT, 'T' ) ) THEN
00306             IF( LSAME( DIRECT, 'F' ) ) THEN
00307                DO 180 J = 2, N
00308                   CTEMP = C( J-1 )
00309                   STEMP = S( J-1 )
00310                   IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
00311                      DO 170 I = 1, M
00312                         TEMP = A( I, J )
00313                         A( I, J ) = CTEMP*TEMP - STEMP*A( I, 1 )
00314                         A( I, 1 ) = STEMP*TEMP + CTEMP*A( I, 1 )
00315   170                CONTINUE
00316                   END IF
00317   180          CONTINUE
00318             ELSE IF( LSAME( DIRECT, 'B' ) ) THEN
00319                DO 200 J = N, 2, -1
00320                   CTEMP = C( J-1 )
00321                   STEMP = S( J-1 )
00322                   IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
00323                      DO 190 I = 1, M
00324                         TEMP = A( I, J )
00325                         A( I, J ) = CTEMP*TEMP - STEMP*A( I, 1 )
00326                         A( I, 1 ) = STEMP*TEMP + CTEMP*A( I, 1 )
00327   190                CONTINUE
00328                   END IF
00329   200          CONTINUE
00330             END IF
00331          ELSE IF( LSAME( PIVOT, 'B' ) ) THEN
00332             IF( LSAME( DIRECT, 'F' ) ) THEN
00333                DO 220 J = 1, N - 1
00334                   CTEMP = C( J )
00335                   STEMP = S( J )
00336                   IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
00337                      DO 210 I = 1, M
00338                         TEMP = A( I, J )
00339                         A( I, J ) = STEMP*A( I, N ) + CTEMP*TEMP
00340                         A( I, N ) = CTEMP*A( I, N ) - STEMP*TEMP
00341   210                CONTINUE
00342                   END IF
00343   220          CONTINUE
00344             ELSE IF( LSAME( DIRECT, 'B' ) ) THEN
00345                DO 240 J = N - 1, 1, -1
00346                   CTEMP = C( J )
00347                   STEMP = S( J )
00348                   IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
00349                      DO 230 I = 1, M
00350                         TEMP = A( I, J )
00351                         A( I, J ) = STEMP*A( I, N ) + CTEMP*TEMP
00352                         A( I, N ) = CTEMP*A( I, N ) - STEMP*TEMP
00353   230                CONTINUE
00354                   END IF
00355   240          CONTINUE
00356             END IF
00357          END IF
00358       END IF
00359 *
00360       RETURN
00361 *
00362 *     End of CLASR
00363 *
00364       END