LAPACK 3.12.1
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches

◆ slarft()

subroutine slarft ( character direct,
character storev,
integer n,
integer k,
real, dimension( ldv, * ) v,
integer ldv,
real, dimension( * ) tau,
real, dimension( ldt, * ) t,
integer ldt )

SLARFT forms the triangular factor T of a block reflector H = I - vtvH

Download SLARFT + dependencies [TGZ] [ZIP] [TXT]

Purpose:
!>
!> SLARFT forms the triangular factor T of a real block reflector H
!> of order n, which is defined as a product of k elementary reflectors.
!>
!> If DIRECT = 'F', H = H(1) H(2) . . . H(k) and T is upper triangular;
!>
!> If DIRECT = 'B', H = H(k) . . . H(2) H(1) and T is lower triangular.
!>
!> If STOREV = 'C', the vector which defines the elementary reflector
!> H(i) is stored in the i-th column of the array V, and
!>
!>    H  =  I - V * T * V**T
!>
!> If STOREV = 'R', the vector which defines the elementary reflector
!> H(i) is stored in the i-th row of the array V, and
!>
!>    H  =  I - V**T * T * V
!>
Parameters
[in]DIRECT
!>          DIRECT is CHARACTER*1
!>          Specifies the order in which the elementary reflectors are
!>          multiplied to form the block reflector:
!>          = 'F': H = H(1) H(2) . . . H(k) (Forward)
!>          = 'B': H = H(k) . . . H(2) H(1) (Backward)
!>
[in]STOREV
!>          STOREV is CHARACTER*1
!>          Specifies how the vectors which define the elementary
!>          reflectors are stored (see also Further Details):
!>          = 'C': columnwise
!>          = 'R': rowwise
!>
[in]N
!>          N is INTEGER
!>          The order of the block reflector H. N >= 0.
!>
[in]K
!>          K is INTEGER
!>          The order of the triangular factor T (= the number of
!>          elementary reflectors). K >= 1.
!>
[in]V
!>          V is REAL array, dimension
!>                               (LDV,K) if STOREV = 'C'
!>                               (LDV,N) if STOREV = 'R'
!>          The matrix V. See further details.
!>
[in]LDV
!>          LDV is INTEGER
!>          The leading dimension of the array V.
!>          If STOREV = 'C', LDV >= max(1,N); if STOREV = 'R', LDV >= K.
!>
[in]TAU
!>          TAU is REAL array, dimension (K)
!>          TAU(i) must contain the scalar factor of the elementary
!>          reflector H(i).
!>
[out]T
!>          T is REAL array, dimension (LDT,K)
!>          The k by k triangular factor T of the block reflector.
!>          If DIRECT = 'F', T is upper triangular; if DIRECT = 'B', T is
!>          lower triangular. The rest of the array is not used.
!>
[in]LDT
!>          LDT is INTEGER
!>          The leading dimension of the array T. LDT >= K.
!>
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Further Details:
!>
!>  The shape of the matrix V and the storage of the vectors which define
!>  the H(i) is best illustrated by the following example with n = 5 and
!>  k = 3. The elements equal to 1 are not stored.
!>
!>  DIRECT = 'F' and STOREV = 'C':         DIRECT = 'F' and STOREV = 'R':
!>
!>               V = (  1       )                 V = (  1 v1 v1 v1 v1 )
!>                   ( v1  1    )                     (     1 v2 v2 v2 )
!>                   ( v1 v2  1 )                     (        1 v3 v3 )
!>                   ( v1 v2 v3 )
!>                   ( v1 v2 v3 )
!>
!>  DIRECT = 'B' and STOREV = 'C':         DIRECT = 'B' and STOREV = 'R':
!>
!>               V = ( v1 v2 v3 )                 V = ( v1 v1  1       )
!>                   ( v1 v2 v3 )                     ( v2 v2 v2  1    )
!>                   (  1 v2 v3 )                     ( v3 v3 v3 v3  1 )
!>                   (     1 v3 )
!>                   (        1 )
!>

Definition at line 160 of file slarft.f.

161*
162* -- LAPACK auxiliary routine --
163* -- LAPACK is a software package provided by Univ. of Tennessee, --
164* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
165*
166* .. Scalar Arguments ..
167 CHARACTER DIRECT, STOREV
168 INTEGER K, LDT, LDV, N
169* ..
170* .. Array Arguments ..
171 REAL T( LDT, * ), TAU( * ), V( LDV, * )
172* ..
173*
174* =====================================================================
175*
176* .. Parameters ..
177 REAL ONE, ZERO
178 parameter( one = 1.0e+0, zero = 0.0e+0 )
179* ..
180* .. Local Scalars ..
181 INTEGER I, J, PREVLASTV, LASTV
182* ..
183* .. External Subroutines ..
184 EXTERNAL sgemv, strmv
185* ..
186* .. External Functions ..
187 LOGICAL LSAME
188 EXTERNAL lsame
189* ..
190* .. Executable Statements ..
191*
192* Quick return if possible
193*
194 IF( n.EQ.0 )
195 $ RETURN
196*
197 IF( lsame( direct, 'F' ) ) THEN
198 prevlastv = n
199 DO i = 1, k
200 prevlastv = max( i, prevlastv )
201 IF( tau( i ).EQ.zero ) THEN
202*
203* H(i) = I
204*
205 DO j = 1, i
206 t( j, i ) = zero
207 END DO
208 ELSE
209*
210* general case
211*
212 IF( lsame( storev, 'C' ) ) THEN
213* Skip any trailing zeros.
214 DO lastv = n, i+1, -1
215 IF( v( lastv, i ).NE.zero ) EXIT
216 END DO
217 DO j = 1, i-1
218 t( j, i ) = -tau( i ) * v( i , j )
219 END DO
220 j = min( lastv, prevlastv )
221*
222* T(1:i-1,i) := - tau(i) * V(i:j,1:i-1)**T * V(i:j,i)
223*
224 CALL sgemv( 'Transpose', j-i, i-1, -tau( i ),
225 $ v( i+1, 1 ), ldv, v( i+1, i ), 1, one,
226 $ t( 1, i ), 1 )
227 ELSE
228* Skip any trailing zeros.
229 DO lastv = n, i+1, -1
230 IF( v( i, lastv ).NE.zero ) EXIT
231 END DO
232 DO j = 1, i-1
233 t( j, i ) = -tau( i ) * v( j , i )
234 END DO
235 j = min( lastv, prevlastv )
236*
237* T(1:i-1,i) := - tau(i) * V(1:i-1,i:j) * V(i,i:j)**T
238*
239 CALL sgemv( 'No transpose', i-1, j-i, -tau( i ),
240 $ v( 1, i+1 ), ldv, v( i, i+1 ), ldv,
241 $ one, t( 1, i ), 1 )
242 END IF
243*
244* T(1:i-1,i) := T(1:i-1,1:i-1) * T(1:i-1,i)
245*
246 CALL strmv( 'Upper', 'No transpose', 'Non-unit', i-1, t,
247 $ ldt, t( 1, i ), 1 )
248 t( i, i ) = tau( i )
249 IF( i.GT.1 ) THEN
250 prevlastv = max( prevlastv, lastv )
251 ELSE
252 prevlastv = lastv
253 END IF
254 END IF
255 END DO
256 ELSE
257 prevlastv = 1
258 DO i = k, 1, -1
259 IF( tau( i ).EQ.zero ) THEN
260*
261* H(i) = I
262*
263 DO j = i, k
264 t( j, i ) = zero
265 END DO
266 ELSE
267*
268* general case
269*
270 IF( i.LT.k ) THEN
271 IF( lsame( storev, 'C' ) ) THEN
272* Skip any leading zeros.
273 DO lastv = 1, i-1
274 IF( v( lastv, i ).NE.zero ) EXIT
275 END DO
276 DO j = i+1, k
277 t( j, i ) = -tau( i ) * v( n-k+i , j )
278 END DO
279 j = max( lastv, prevlastv )
280*
281* T(i+1:k,i) = -tau(i) * V(j:n-k+i,i+1:k)**T * V(j:n-k+i,i)
282*
283 CALL sgemv( 'Transpose', n-k+i-j, k-i, -tau( i ),
284 $ v( j, i+1 ), ldv, v( j, i ), 1, one,
285 $ t( i+1, i ), 1 )
286 ELSE
287* Skip any leading zeros.
288 DO lastv = 1, i-1
289 IF( v( i, lastv ).NE.zero ) EXIT
290 END DO
291 DO j = i+1, k
292 t( j, i ) = -tau( i ) * v( j, n-k+i )
293 END DO
294 j = max( lastv, prevlastv )
295*
296* T(i+1:k,i) = -tau(i) * V(i+1:k,j:n-k+i) * V(i,j:n-k+i)**T
297*
298 CALL sgemv( 'No transpose', k-i, n-k+i-j,
299 $ -tau( i ), v( i+1, j ), ldv, v( i, j ), ldv,
300 $ one, t( i+1, i ), 1 )
301 END IF
302*
303* T(i+1:k,i) := T(i+1:k,i+1:k) * T(i+1:k,i)
304*
305 CALL strmv( 'Lower', 'No transpose', 'Non-unit', k-i,
306 $ t( i+1, i+1 ), ldt, t( i+1, i ), 1 )
307 IF( i.GT.1 ) THEN
308 prevlastv = min( prevlastv, lastv )
309 ELSE
310 prevlastv = lastv
311 END IF
312 END IF
313 t( i, i ) = tau( i )
314 END IF
315 END DO
316 END IF
317 RETURN
318*
319* End of SLARFT
320*
sgemv
subroutine sgemv(trans, m, n, alpha, a, lda, x, incx, beta, y, incy)
SGEMV
Definition sgemv.f:158
lsame
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48
strmv
subroutine strmv(uplo, trans, diag, n, a, lda, x, incx)
STRMV
Definition strmv.f:147
Here is the call graph for this function:
Here is the caller graph for this function: