1: SUBROUTINE DLARF( SIDE, M, N, V, INCV, TAU, C, LDC, WORK )
2: IMPLICIT NONE
3: *
4: * -- LAPACK auxiliary routine (version 3.2) --
5: * -- LAPACK is a software package provided by Univ. of Tennessee, --
6: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
7: * November 2006
8: *
9: * .. Scalar Arguments ..
10: CHARACTER SIDE
11: INTEGER INCV, LDC, M, N
12: DOUBLE PRECISION TAU
13: * ..
14: * .. Array Arguments ..
15: DOUBLE PRECISION C( LDC, * ), V( * ), WORK( * )
16: * ..
17: *
18: * Purpose
19: * =======
20: *
21: * DLARF applies a real elementary reflector H to a real m by n matrix
22: * C, from either the left or the right. H is represented in the form
23: *
24: * H = I - tau * v * v'
25: *
26: * where tau is a real scalar and v is a real vector.
27: *
28: * If tau = 0, then H is taken to be the unit matrix.
29: *
30: * Arguments
31: * =========
32: *
33: * SIDE (input) CHARACTER*1
34: * = 'L': form H * C
35: * = 'R': form C * H
36: *
37: * M (input) INTEGER
38: * The number of rows of the matrix C.
39: *
40: * N (input) INTEGER
41: * The number of columns of the matrix C.
42: *
43: * V (input) DOUBLE PRECISION array, dimension
44: * (1 + (M-1)*abs(INCV)) if SIDE = 'L'
45: * or (1 + (N-1)*abs(INCV)) if SIDE = 'R'
46: * The vector v in the representation of H. V is not used if
47: * TAU = 0.
48: *
49: * INCV (input) INTEGER
50: * The increment between elements of v. INCV <> 0.
51: *
52: * TAU (input) DOUBLE PRECISION
53: * The value tau in the representation of H.
54: *
55: * C (input/output) DOUBLE PRECISION array, dimension (LDC,N)
56: * On entry, the m by n matrix C.
57: * On exit, C is overwritten by the matrix H * C if SIDE = 'L',
58: * or C * H if SIDE = 'R'.
59: *
60: * LDC (input) INTEGER
61: * The leading dimension of the array C. LDC >= max(1,M).
62: *
63: * WORK (workspace) DOUBLE PRECISION array, dimension
64: * (N) if SIDE = 'L'
65: * or (M) if SIDE = 'R'
66: *
67: * =====================================================================
68: *
69: * .. Parameters ..
70: DOUBLE PRECISION ONE, ZERO
71: PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 )
72: * ..
73: * .. Local Scalars ..
74: LOGICAL APPLYLEFT
75: INTEGER I, LASTV, LASTC
76: * ..
77: * .. External Subroutines ..
78: EXTERNAL DGEMV, DGER
79: * ..
80: * .. External Functions ..
81: LOGICAL LSAME
82: INTEGER ILADLR, ILADLC
83: EXTERNAL LSAME, ILADLR, ILADLC
84: * ..
85: * .. Executable Statements ..
86: *
87: APPLYLEFT = LSAME( SIDE, 'L' )
88: LASTV = 0
89: LASTC = 0
90: IF( TAU.NE.ZERO ) THEN
91: ! Set up variables for scanning V. LASTV begins pointing to the end
92: ! of V.
93: IF( APPLYLEFT ) THEN
94: LASTV = M
95: ELSE
96: LASTV = N
97: END IF
98: IF( INCV.GT.0 ) THEN
99: I = 1 + (LASTV-1) * INCV
100: ELSE
101: I = 1
102: END IF
103: ! Look for the last non-zero row in V.
104: DO WHILE( LASTV.GT.0 .AND. V( I ).EQ.ZERO )
105: LASTV = LASTV - 1
106: I = I - INCV
107: END DO
108: IF( APPLYLEFT ) THEN
109: ! Scan for the last non-zero column in C(1:lastv,:).
110: LASTC = ILADLC(LASTV, N, C, LDC)
111: ELSE
112: ! Scan for the last non-zero row in C(:,1:lastv).
113: LASTC = ILADLR(M, LASTV, C, LDC)
114: END IF
115: END IF
116: ! Note that lastc.eq.0 renders the BLAS operations null; no special
117: ! case is needed at this level.
118: IF( APPLYLEFT ) THEN
119: *
120: * Form H * C
121: *
122: IF( LASTV.GT.0 ) THEN
123: *
124: * w(1:lastc,1) := C(1:lastv,1:lastc)' * v(1:lastv,1)
125: *
126: CALL DGEMV( 'Transpose', LASTV, LASTC, ONE, C, LDC, V, INCV,
127: $ ZERO, WORK, 1 )
128: *
129: * C(1:lastv,1:lastc) := C(...) - v(1:lastv,1) * w(1:lastc,1)'
130: *
131: CALL DGER( LASTV, LASTC, -TAU, V, INCV, WORK, 1, C, LDC )
132: END IF
133: ELSE
134: *
135: * Form C * H
136: *
137: IF( LASTV.GT.0 ) THEN
138: *
139: * w(1:lastc,1) := C(1:lastc,1:lastv) * v(1:lastv,1)
140: *
141: CALL DGEMV( 'No transpose', LASTC, LASTV, ONE, C, LDC,
142: $ V, INCV, ZERO, WORK, 1 )
143: *
144: * C(1:lastc,1:lastv) := C(...) - w(1:lastc,1) * v(1:lastv,1)'
145: *
146: CALL DGER( LASTC, LASTV, -TAU, WORK, 1, V, INCV, C, LDC )
147: END IF
148: END IF
149: RETURN
150: *
151: * End of DLARF
152: *
153: END
CVSweb interface <joel.bertrand@systella.fr>