source: trunk/SOURCES/BLAS/slarzt.f @ 22

Last change on this file since 22 was 22, checked in by roche, 9 years ago

Petites adaptations diverses du code pour compilation en gfortran. Ajout d un Makefile flexible a option pour choisir ifort ou gfortran.
File size: 7.7 KB
Line 
1*> \brief \b SLARZT
2*
3*  =========== DOCUMENTATION ===========
4*
5* Online html documentation available at 
6*            http://www.netlib.org/lapack/explore-html/ 
7*
8*> \htmlonly
9*> Download SLARZT + dependencies 
10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/slarzt.f"> 
11*> [TGZ]</a> 
12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/slarzt.f"> 
13*> [ZIP]</a> 
14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/slarzt.f"> 
15*> [TXT]</a>
16*> \endhtmlonly 
17*
18*  Definition:
19*  ===========
20*
21*       SUBROUTINE SLARZT( DIRECT, STOREV, N, K, V, LDV, TAU, T, LDT )
22* 
23*       .. Scalar Arguments ..
24*       CHARACTER          DIRECT, STOREV
25*       INTEGER            K, LDT, LDV, N
26*       ..
27*       .. Array Arguments ..
28*       REAL               T( LDT, * ), TAU( * ), V( LDV, * )
29*       ..
30* 
31*
32*> \par Purpose:
33*  =============
34*>
35*> \verbatim
36*>
37*> SLARZT forms the triangular factor T of a real block reflector
38*> H of order > n, which is defined as a product of k elementary
39*> reflectors.
40*>
41*> If DIRECT = 'F', H = H(1) H(2) . . . H(k) and T is upper triangular;
42*>
43*> If DIRECT = 'B', H = H(k) . . . H(2) H(1) and T is lower triangular.
44*>
45*> If STOREV = 'C', the vector which defines the elementary reflector
46*> H(i) is stored in the i-th column of the array V, and
47*>
48*>    H  =  I - V * T * V**T
49*>
50*> If STOREV = 'R', the vector which defines the elementary reflector
51*> H(i) is stored in the i-th row of the array V, and
52*>
53*>    H  =  I - V**T * T * V
54*>
55*> Currently, only STOREV = 'R' and DIRECT = 'B' are supported.
56*> \endverbatim
57*
58*  Arguments:
59*  ==========
60*
61*> \param[in] DIRECT
62*> \verbatim
63*>          DIRECT is CHARACTER*1
64*>          Specifies the order in which the elementary reflectors are
65*>          multiplied to form the block reflector:
66*>          = 'F': H = H(1) H(2) . . . H(k) (Forward, not supported yet)
67*>          = 'B': H = H(k) . . . H(2) H(1) (Backward)
68*> \endverbatim
69*>
70*> \param[in] STOREV
71*> \verbatim
72*>          STOREV is CHARACTER*1
73*>          Specifies how the vectors which define the elementary
74*>          reflectors are stored (see also Further Details):
75*>          = 'C': columnwise                        (not supported yet)
76*>          = 'R': rowwise
77*> \endverbatim
78*>
79*> \param[in] N
80*> \verbatim
81*>          N is INTEGER
82*>          The order of the block reflector H. N >= 0.
83*> \endverbatim
84*>
85*> \param[in] K
86*> \verbatim
87*>          K is INTEGER
88*>          The order of the triangular factor T (= the number of
89*>          elementary reflectors). K >= 1.
90*> \endverbatim
91*>
92*> \param[in,out] V
93*> \verbatim
94*>          V is REAL array, dimension
95*>                               (LDV,K) if STOREV = 'C'
96*>                               (LDV,N) if STOREV = 'R'
97*>          The matrix V. See further details.
98*> \endverbatim
99*>
100*> \param[in] LDV
101*> \verbatim
102*>          LDV is INTEGER
103*>          The leading dimension of the array V.
104*>          If STOREV = 'C', LDV >= max(1,N); if STOREV = 'R', LDV >= K.
105*> \endverbatim
106*>
107*> \param[in] TAU
108*> \verbatim
109*>          TAU is REAL array, dimension (K)
110*>          TAU(i) must contain the scalar factor of the elementary
111*>          reflector H(i).
112*> \endverbatim
113*>
114*> \param[out] T
115*> \verbatim
116*>          T is REAL array, dimension (LDT,K)
117*>          The k by k triangular factor T of the block reflector.
118*>          If DIRECT = 'F', T is upper triangular; if DIRECT = 'B', T is
119*>          lower triangular. The rest of the array is not used.
120*> \endverbatim
121*>
122*> \param[in] LDT
123*> \verbatim
124*>          LDT is INTEGER
125*>          The leading dimension of the array T. LDT >= K.
126*> \endverbatim
127*
128*  Authors:
129*  ========
130*
131*> \author Univ. of Tennessee 
132*> \author Univ. of California Berkeley 
133*> \author Univ. of Colorado Denver 
134*> \author NAG Ltd. 
135*
136*> \date November 2011
137*
138*> \ingroup realOTHERcomputational
139*
140*> \par Contributors:
141*  ==================
142*>
143*>    A. Petitet, Computer Science Dept., Univ. of Tenn., Knoxville, USA
144*
145*> \par Further Details:
146*  =====================
147*>
148*> \verbatim
149*>
150*>  The shape of the matrix V and the storage of the vectors which define
151*>  the H(i) is best illustrated by the following example with n = 5 and
152*>  k = 3. The elements equal to 1 are not stored; the corresponding
153*>  array elements are modified but restored on exit. The rest of the
154*>  array is not used.
155*>
156*>  DIRECT = 'F' and STOREV = 'C':         DIRECT = 'F' and STOREV = 'R':
157*>
158*>                                              ______V_____
159*>         ( v1 v2 v3 )                        /            \
160*>         ( v1 v2 v3 )                      ( v1 v1 v1 v1 v1 . . . . 1 )
161*>     V = ( v1 v2 v3 )                      ( v2 v2 v2 v2 v2 . . . 1   )
162*>         ( v1 v2 v3 )                      ( v3 v3 v3 v3 v3 . . 1     )
163*>         ( v1 v2 v3 )
164*>            .  .  .
165*>            .  .  .
166*>            1  .  .
167*>               1  .
168*>                  1
169*>
170*>  DIRECT = 'B' and STOREV = 'C':         DIRECT = 'B' and STOREV = 'R':
171*>
172*>                                                        ______V_____
173*>            1                                          /            \
174*>            .  1                           ( 1 . . . . v1 v1 v1 v1 v1 )
175*>            .  .  1                        ( . 1 . . . v2 v2 v2 v2 v2 )
176*>            .  .  .                        ( . . 1 . . v3 v3 v3 v3 v3 )
177*>            .  .  .
178*>         ( v1 v2 v3 )
179*>         ( v1 v2 v3 )
180*>     V = ( v1 v2 v3 )
181*>         ( v1 v2 v3 )
182*>         ( v1 v2 v3 )
183*> \endverbatim
184*>
185*  =====================================================================
186      SUBROUTINE SLARZT( DIRECT, STOREV, N, K, V, LDV, TAU, T, LDT )
187*
188*  -- LAPACK computational routine (version 3.4.0) --
189*  -- LAPACK is a software package provided by Univ. of Tennessee,    --
190*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
191*     November 2011
192*
193*     .. Scalar Arguments ..
194      CHARACTER          DIRECT, STOREV
195      INTEGER            K, LDT, LDV, N
196*     ..
197*     .. Array Arguments ..
198      REAL               T( LDT, * ), TAU( * ), V( LDV, * )
199*     ..
200*
201*  =====================================================================
202*
203*     .. Parameters ..
204      REAL               ZERO
205      PARAMETER          ( ZERO = 0.0E+0 )
206*     ..
207*     .. Local Scalars ..
208      INTEGER            I, INFO, J
209*     ..
210*     .. External Subroutines ..
211      EXTERNAL           SGEMV, STRMV, XERBLA
212*     ..
213*     .. External Functions ..
214      LOGICAL            LSAME
215      EXTERNAL           LSAME
216*     ..
217*     .. Executable Statements ..
218*
219*     Check for currently supported options
220*
221      INFO = 0
222      IF( .NOT.LSAME( DIRECT, 'B' ) ) THEN
223         INFO = -1
224      ELSE IF( .NOT.LSAME( STOREV, 'R' ) ) THEN
225         INFO = -2
226      END IF
227      IF( INFO.NE.0 ) THEN
228         CALL XERBLA( 'SLARZT', -INFO )
229         RETURN
230      END IF
231*
232      DO 20 I = K, 1, -1
233         IF( TAU( I ).EQ.ZERO ) THEN
234*
235*           H(i)  =  I
236*
237            DO 10 J = I, K
238               T( J, I ) = ZERO
239   10       CONTINUE
240         ELSE
241*
242*           general case
243*
244            IF( I.LT.K ) THEN
245*
246*              T(i+1:k,i) = - tau(i) * V(i+1:k,1:n) * V(i,1:n)**T
247*
248               CALL SGEMV( 'No transpose', K-I, N, -TAU( I ),
249     $                     V( I+1, 1 ), LDV, V( I, 1 ), LDV, ZERO,
250     $                     T( I+1, I ), 1 )
251*
252*              T(i+1:k,i) = T(i+1:k,i+1:k) * T(i+1:k,i)
253*
254               CALL STRMV( 'Lower', 'No transpose', 'Non-unit', K-I,
255     $                     T( I+1, I+1 ), LDT, T( I+1, I ), 1 )
256            END IF
257            T( I, I ) = TAU( I )
258         END IF
259   20 CONTINUE
260      RETURN
261*
262*     End of SLARZT
263*
264      END
Note: See TracBrowser for help on using the repository browser.