LAPACK 3.3.0
|
00001 SUBROUTINE CLARF( SIDE, M, N, V, INCV, TAU, C, LDC, WORK ) 00002 IMPLICIT NONE 00003 * 00004 * -- LAPACK auxiliary routine (version 3.2) -- 00005 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00006 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00007 * November 2006 00008 * 00009 * .. Scalar Arguments .. 00010 CHARACTER SIDE 00011 INTEGER INCV, LDC, M, N 00012 COMPLEX TAU 00013 * .. 00014 * .. Array Arguments .. 00015 COMPLEX C( LDC, * ), V( * ), WORK( * ) 00016 * .. 00017 * 00018 * Purpose 00019 * ======= 00020 * 00021 * CLARF applies a complex elementary reflector H to a complex M-by-N 00022 * matrix C, from either the left or the right. H is represented in the 00023 * form 00024 * 00025 * H = I - tau * v * v' 00026 * 00027 * where tau is a complex scalar and v is a complex vector. 00028 * 00029 * If tau = 0, then H is taken to be the unit matrix. 00030 * 00031 * To apply H' (the conjugate transpose of H), supply conjg(tau) instead 00032 * tau. 00033 * 00034 * Arguments 00035 * ========= 00036 * 00037 * SIDE (input) CHARACTER*1 00038 * = 'L': form H * C 00039 * = 'R': form C * H 00040 * 00041 * M (input) INTEGER 00042 * The number of rows of the matrix C. 00043 * 00044 * N (input) INTEGER 00045 * The number of columns of the matrix C. 00046 * 00047 * V (input) COMPLEX array, dimension 00048 * (1 + (M-1)*abs(INCV)) if SIDE = 'L' 00049 * or (1 + (N-1)*abs(INCV)) if SIDE = 'R' 00050 * The vector v in the representation of H. V is not used if 00051 * TAU = 0. 00052 * 00053 * INCV (input) INTEGER 00054 * The increment between elements of v. INCV <> 0. 00055 * 00056 * TAU (input) COMPLEX 00057 * The value tau in the representation of H. 00058 * 00059 * C (input/output) COMPLEX array, dimension (LDC,N) 00060 * On entry, the M-by-N matrix C. 00061 * On exit, C is overwritten by the matrix H * C if SIDE = 'L', 00062 * or C * H if SIDE = 'R'. 00063 * 00064 * LDC (input) INTEGER 00065 * The leading dimension of the array C. LDC >= max(1,M). 00066 * 00067 * WORK (workspace) COMPLEX array, dimension 00068 * (N) if SIDE = 'L' 00069 * or (M) if SIDE = 'R' 00070 * 00071 * ===================================================================== 00072 * 00073 * .. Parameters .. 00074 COMPLEX ONE, ZERO 00075 PARAMETER ( ONE = ( 1.0E+0, 0.0E+0 ), 00076 $ ZERO = ( 0.0E+0, 0.0E+0 ) ) 00077 * .. 00078 * .. Local Scalars .. 00079 LOGICAL APPLYLEFT 00080 INTEGER I, LASTV, LASTC 00081 * .. 00082 * .. External Subroutines .. 00083 EXTERNAL CGEMV, CGERC 00084 * .. 00085 * .. External Functions .. 00086 LOGICAL LSAME 00087 INTEGER ILACLR, ILACLC 00088 EXTERNAL LSAME, ILACLR, ILACLC 00089 * .. 00090 * .. Executable Statements .. 00091 * 00092 APPLYLEFT = LSAME( SIDE, 'L' ) 00093 LASTV = 0 00094 LASTC = 0 00095 IF( TAU.NE.ZERO ) THEN 00096 ! Set up variables for scanning V. LASTV begins pointing to the end 00097 ! of V. 00098 IF( APPLYLEFT ) THEN 00099 LASTV = M 00100 ELSE 00101 LASTV = N 00102 END IF 00103 IF( INCV.GT.0 ) THEN 00104 I = 1 + (LASTV-1) * INCV 00105 ELSE 00106 I = 1 00107 END IF 00108 ! Look for the last non-zero row in V. 00109 DO WHILE( LASTV.GT.0 .AND. V( I ).EQ.ZERO ) 00110 LASTV = LASTV - 1 00111 I = I - INCV 00112 END DO 00113 IF( APPLYLEFT ) THEN 00114 ! Scan for the last non-zero column in C(1:lastv,:). 00115 LASTC = ILACLC(LASTV, N, C, LDC) 00116 ELSE 00117 ! Scan for the last non-zero row in C(:,1:lastv). 00118 LASTC = ILACLR(M, LASTV, C, LDC) 00119 END IF 00120 END IF 00121 ! Note that lastc.eq.0 renders the BLAS operations null; no special 00122 ! case is needed at this level. 00123 IF( APPLYLEFT ) THEN 00124 * 00125 * Form H * C 00126 * 00127 IF( LASTV.GT.0 ) THEN 00128 * 00129 * w(1:lastc,1) := C(1:lastv,1:lastc)' * v(1:lastv,1) 00130 * 00131 CALL CGEMV( 'Conjugate transpose', LASTV, LASTC, ONE, 00132 $ C, LDC, V, INCV, ZERO, WORK, 1 ) 00133 * 00134 * C(1:lastv,1:lastc) := C(...) - v(1:lastv,1) * w(1:lastc,1)' 00135 * 00136 CALL CGERC( LASTV, LASTC, -TAU, V, INCV, WORK, 1, C, LDC ) 00137 END IF 00138 ELSE 00139 * 00140 * Form C * H 00141 * 00142 IF( LASTV.GT.0 ) THEN 00143 * 00144 * w(1:lastc,1) := C(1:lastc,1:lastv) * v(1:lastv,1) 00145 * 00146 CALL CGEMV( 'No transpose', LASTC, LASTV, ONE, C, LDC, 00147 $ V, INCV, ZERO, WORK, 1 ) 00148 * 00149 * C(1:lastc,1:lastv) := C(...) - w(1:lastc,1) * v(1:lastv,1)' 00150 * 00151 CALL CGERC( LASTC, LASTV, -TAU, WORK, 1, V, INCV, C, LDC ) 00152 END IF 00153 END IF 00154 RETURN 00155 * 00156 * End of CLARF 00157 * 00158 END