source: branches/GRISLIv3/SOURCES/BLAS/slarzb.f @ 474

Last change on this file since 474 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: 8.8 KB
Line 
1*> \brief \b SLARZB
2*
3*  =========== DOCUMENTATION ===========
4*
5* Online html documentation available at 
6*            http://www.netlib.org/lapack/explore-html/ 
7*
8*> \htmlonly
9*> Download SLARZB + dependencies 
10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/slarzb.f"> 
11*> [TGZ]</a> 
12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/slarzb.f"> 
13*> [ZIP]</a> 
14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/slarzb.f"> 
15*> [TXT]</a>
16*> \endhtmlonly 
17*
18*  Definition:
19*  ===========
20*
21*       SUBROUTINE SLARZB( SIDE, TRANS, DIRECT, STOREV, M, N, K, L, V,
22*                          LDV, T, LDT, C, LDC, WORK, LDWORK )
23* 
24*       .. Scalar Arguments ..
25*       CHARACTER          DIRECT, SIDE, STOREV, TRANS
26*       INTEGER            K, L, LDC, LDT, LDV, LDWORK, M, N
27*       ..
28*       .. Array Arguments ..
29*       REAL               C( LDC, * ), T( LDT, * ), V( LDV, * ),
30*      $                   WORK( LDWORK, * )
31*       ..
32* 
33*
34*> \par Purpose:
35*  =============
36*>
37*> \verbatim
38*>
39*> SLARZB applies a real block reflector H or its transpose H**T to
40*> a real distributed M-by-N  C from the left or the right.
41*>
42*> Currently, only STOREV = 'R' and DIRECT = 'B' are supported.
43*> \endverbatim
44*
45*  Arguments:
46*  ==========
47*
48*> \param[in] SIDE
49*> \verbatim
50*>          SIDE is CHARACTER*1
51*>          = 'L': apply H or H**T from the Left
52*>          = 'R': apply H or H**T from the Right
53*> \endverbatim
54*>
55*> \param[in] TRANS
56*> \verbatim
57*>          TRANS is CHARACTER*1
58*>          = 'N': apply H (No transpose)
59*>          = 'C': apply H**T (Transpose)
60*> \endverbatim
61*>
62*> \param[in] DIRECT
63*> \verbatim
64*>          DIRECT is CHARACTER*1
65*>          Indicates how H is formed from a product of elementary
66*>          reflectors
67*>          = 'F': H = H(1) H(2) . . . H(k) (Forward, not supported yet)
68*>          = 'B': H = H(k) . . . H(2) H(1) (Backward)
69*> \endverbatim
70*>
71*> \param[in] STOREV
72*> \verbatim
73*>          STOREV is CHARACTER*1
74*>          Indicates how the vectors which define the elementary
75*>          reflectors are stored:
76*>          = 'C': Columnwise                        (not supported yet)
77*>          = 'R': Rowwise
78*> \endverbatim
79*>
80*> \param[in] M
81*> \verbatim
82*>          M is INTEGER
83*>          The number of rows of the matrix C.
84*> \endverbatim
85*>
86*> \param[in] N
87*> \verbatim
88*>          N is INTEGER
89*>          The number of columns of the matrix C.
90*> \endverbatim
91*>
92*> \param[in] K
93*> \verbatim
94*>          K is INTEGER
95*>          The order of the matrix T (= the number of elementary
96*>          reflectors whose product defines the block reflector).
97*> \endverbatim
98*>
99*> \param[in] L
100*> \verbatim
101*>          L is INTEGER
102*>          The number of columns of the matrix V containing the
103*>          meaningful part of the Householder reflectors.
104*>          If SIDE = 'L', M >= L >= 0, if SIDE = 'R', N >= L >= 0.
105*> \endverbatim
106*>
107*> \param[in] V
108*> \verbatim
109*>          V is REAL array, dimension (LDV,NV).
110*>          If STOREV = 'C', NV = K; if STOREV = 'R', NV = L.
111*> \endverbatim
112*>
113*> \param[in] LDV
114*> \verbatim
115*>          LDV is INTEGER
116*>          The leading dimension of the array V.
117*>          If STOREV = 'C', LDV >= L; if STOREV = 'R', LDV >= K.
118*> \endverbatim
119*>
120*> \param[in] T
121*> \verbatim
122*>          T is REAL array, dimension (LDT,K)
123*>          The triangular K-by-K matrix T in the representation of the
124*>          block reflector.
125*> \endverbatim
126*>
127*> \param[in] LDT
128*> \verbatim
129*>          LDT is INTEGER
130*>          The leading dimension of the array T. LDT >= K.
131*> \endverbatim
132*>
133*> \param[in,out] C
134*> \verbatim
135*>          C is REAL array, dimension (LDC,N)
136*>          On entry, the M-by-N matrix C.
137*>          On exit, C is overwritten by H*C or H**T*C or C*H or C*H**T.
138*> \endverbatim
139*>
140*> \param[in] LDC
141*> \verbatim
142*>          LDC is INTEGER
143*>          The leading dimension of the array C. LDC >= max(1,M).
144*> \endverbatim
145*>
146*> \param[out] WORK
147*> \verbatim
148*>          WORK is REAL array, dimension (LDWORK,K)
149*> \endverbatim
150*>
151*> \param[in] LDWORK
152*> \verbatim
153*>          LDWORK is INTEGER
154*>          The leading dimension of the array WORK.
155*>          If SIDE = 'L', LDWORK >= max(1,N);
156*>          if SIDE = 'R', LDWORK >= max(1,M).
157*> \endverbatim
158*
159*  Authors:
160*  ========
161*
162*> \author Univ. of Tennessee 
163*> \author Univ. of California Berkeley 
164*> \author Univ. of Colorado Denver 
165*> \author NAG Ltd. 
166*
167*> \date November 2011
168*
169*> \ingroup realOTHERcomputational
170*
171*> \par Contributors:
172*  ==================
173*>
174*>    A. Petitet, Computer Science Dept., Univ. of Tenn., Knoxville, USA
175*
176*> \par Further Details:
177*  =====================
178*>
179*> \verbatim
180*> \endverbatim
181*>
182*  =====================================================================
183      SUBROUTINE SLARZB( SIDE, TRANS, DIRECT, STOREV, M, N, K, L, V,
184     $                   LDV, T, LDT, C, LDC, WORK, LDWORK )
185*
186*  -- LAPACK computational routine (version 3.4.0) --
187*  -- LAPACK is a software package provided by Univ. of Tennessee,    --
188*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
189*     November 2011
190*
191*     .. Scalar Arguments ..
192      CHARACTER          DIRECT, SIDE, STOREV, TRANS
193      INTEGER            K, L, LDC, LDT, LDV, LDWORK, M, N
194*     ..
195*     .. Array Arguments ..
196      REAL               C( LDC, * ), T( LDT, * ), V( LDV, * ),
197     $                   WORK( LDWORK, * )
198*     ..
199*
200*  =====================================================================
201*
202*     .. Parameters ..
203      REAL               ONE
204      PARAMETER          ( ONE = 1.0E+0 )
205*     ..
206*     .. Local Scalars ..
207      CHARACTER          TRANST
208      INTEGER            I, INFO, J
209*     ..
210*     .. External Functions ..
211      LOGICAL            LSAME
212      EXTERNAL           LSAME
213*     ..
214*     .. External Subroutines ..
215      EXTERNAL           SCOPY, SGEMM, STRMM, XERBLA
216*     ..
217*     .. Executable Statements ..
218*
219*     Quick return if possible
220*
221      IF( M.LE.0 .OR. N.LE.0 )
222     $   RETURN
223*
224*     Check for currently supported options
225*
226      INFO = 0
227      IF( .NOT.LSAME( DIRECT, 'B' ) ) THEN
228         INFO = -3
229      ELSE IF( .NOT.LSAME( STOREV, 'R' ) ) THEN
230         INFO = -4
231      END IF
232      IF( INFO.NE.0 ) THEN
233         CALL XERBLA( 'SLARZB', -INFO )
234         RETURN
235      END IF
236*
237      IF( LSAME( TRANS, 'N' ) ) THEN
238         TRANST = 'T'
239      ELSE
240         TRANST = 'N'
241      END IF
242*
243      IF( LSAME( SIDE, 'L' ) ) THEN
244*
245*        Form  H * C  or  H**T * C
246*
247*        W( 1:n, 1:k ) = C( 1:k, 1:n )**T
248*
249         DO 10 J = 1, K
250            CALL SCOPY( N, C( J, 1 ), LDC, WORK( 1, J ), 1 )
251   10    CONTINUE
252*
253*        W( 1:n, 1:k ) = W( 1:n, 1:k ) + ...
254*                        C( m-l+1:m, 1:n )**T * V( 1:k, 1:l )**T
255*
256         IF( L.GT.0 )
257     $      CALL SGEMM( 'Transpose', 'Transpose', N, K, L, ONE,
258     $                  C( M-L+1, 1 ), LDC, V, LDV, ONE, WORK, LDWORK )
259*
260*        W( 1:n, 1:k ) = W( 1:n, 1:k ) * T**T  or  W( 1:m, 1:k ) * T
261*
262         CALL STRMM( 'Right', 'Lower', TRANST, 'Non-unit', N, K, ONE, T,
263     $               LDT, WORK, LDWORK )
264*
265*        C( 1:k, 1:n ) = C( 1:k, 1:n ) - W( 1:n, 1:k )**T
266*
267         DO 30 J = 1, N
268            DO 20 I = 1, K
269               C( I, J ) = C( I, J ) - WORK( J, I )
270   20       CONTINUE
271   30    CONTINUE
272*
273*        C( m-l+1:m, 1:n ) = C( m-l+1:m, 1:n ) - ...
274*                            V( 1:k, 1:l )**T * W( 1:n, 1:k )**T
275*
276         IF( L.GT.0 )
277     $      CALL SGEMM( 'Transpose', 'Transpose', L, N, K, -ONE, V, LDV,
278     $                  WORK, LDWORK, ONE, C( M-L+1, 1 ), LDC )
279*
280      ELSE IF( LSAME( SIDE, 'R' ) ) THEN
281*
282*        Form  C * H  or  C * H**T
283*
284*        W( 1:m, 1:k ) = C( 1:m, 1:k )
285*
286         DO 40 J = 1, K
287            CALL SCOPY( M, C( 1, J ), 1, WORK( 1, J ), 1 )
288   40    CONTINUE
289*
290*        W( 1:m, 1:k ) = W( 1:m, 1:k ) + ...
291*                        C( 1:m, n-l+1:n ) * V( 1:k, 1:l )**T
292*
293         IF( L.GT.0 )
294     $      CALL SGEMM( 'No transpose', 'Transpose', M, K, L, ONE,
295     $                  C( 1, N-L+1 ), LDC, V, LDV, ONE, WORK, LDWORK )
296*
297*        W( 1:m, 1:k ) = W( 1:m, 1:k ) * T  or  W( 1:m, 1:k ) * T**T
298*
299         CALL STRMM( 'Right', 'Lower', TRANS, 'Non-unit', M, K, ONE, T,
300     $               LDT, WORK, LDWORK )
301*
302*        C( 1:m, 1:k ) = C( 1:m, 1:k ) - W( 1:m, 1:k )
303*
304         DO 60 J = 1, K
305            DO 50 I = 1, M
306               C( I, J ) = C( I, J ) - WORK( I, J )
307   50       CONTINUE
308   60    CONTINUE
309*
310*        C( 1:m, n-l+1:n ) = C( 1:m, n-l+1:n ) - ...
311*                            W( 1:m, 1:k ) * V( 1:k, 1:l )
312*
313         IF( L.GT.0 )
314     $      CALL SGEMM( 'No transpose', 'No transpose', M, L, K, -ONE,
315     $                  WORK, LDWORK, V, LDV, ONE, C( 1, N-L+1 ), LDC )
316*
317      END IF
318*
319      RETURN
320*
321*     End of SLARZB
322*
323      END
Note: See TracBrowser for help on using the repository browser.