← Back to branch summary

~ubuntu-branches/debian/sid/lammps/sid

« back to all changes in this revision

Viewing changes to lib/linalg/dorm2r.f

  • Committer: Package Import Robot
  • Author(s): Anton Gladky
  • Date: 2015-04-29 23:44:49 UTC
  • mfrom: (5.1.3 experimental)
  • Revision ID: package-import@ubuntu.com-20150429234449-mbhy9utku6hp6oq8
Tags: 0~20150313.gitfa668e1-1
Upload into unstable.

Show diffs side-by-side

added added

removed removed

Lines of Context:
lib/linalg/dorm2r.f
 
1
*> \brief \b DORM2R multiplies a general matrix by the orthogonal matrix from a QR factorization determined by sgeqrf (unblocked algorithm).
 
2
*
 
3
*  =========== DOCUMENTATION ===========
 
4
*
 
5
* Online html documentation available at 
 
6
*            http://www.netlib.org/lapack/explore-html/ 
 
7
*
 
8
*> \htmlonly
 
9
*> Download DORM2R + dependencies 
 
10
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dorm2r.f"> 
 
11
*> [TGZ]</a> 
 
12
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dorm2r.f"> 
 
13
*> [ZIP]</a> 
 
14
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dorm2r.f"> 
 
15
*> [TXT]</a>
 
16
*> \endhtmlonly 
 
17
*
 
18
*  Definition:
 
19
*  ===========
 
20
*
 
21
*       SUBROUTINE DORM2R( SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC,
 
22
*                          WORK, INFO )
 
23
* 
 
24
*       .. Scalar Arguments ..
 
25
*       CHARACTER          SIDE, TRANS
 
26
*       INTEGER            INFO, K, LDA, LDC, M, N
 
27
*       ..
 
28
*       .. Array Arguments ..
 
29
*       DOUBLE PRECISION   A( LDA, * ), C( LDC, * ), TAU( * ), WORK( * )
 
30
*       ..
 
31
*  
 
32
*
 
33
*> \par Purpose:
 
34
*  =============
 
35
*>
 
36
*> \verbatim
 
37
*>
 
38
*> DORM2R overwrites the general real m by n matrix C with
 
39
*>
 
40
*>       Q * C  if SIDE = 'L' and TRANS = 'N', or
 
41
*>
 
42
*>       Q**T* C  if SIDE = 'L' and TRANS = 'T', or
 
43
*>
 
44
*>       C * Q  if SIDE = 'R' and TRANS = 'N', or
 
45
*>
 
46
*>       C * Q**T if SIDE = 'R' and TRANS = 'T',
 
47
*>
 
48
*> where Q is a real orthogonal matrix defined as the product of k
 
49
*> elementary reflectors
 
50
*>
 
51
*>       Q = H(1) H(2) . . . H(k)
 
52
*>
 
53
*> as returned by DGEQRF. Q is of order m if SIDE = 'L' and of order n
 
54
*> if SIDE = 'R'.
 
55
*> \endverbatim
 
56
*
 
57
*  Arguments:
 
58
*  ==========
 
59
*
 
60
*> \param[in] SIDE
 
61
*> \verbatim
 
62
*>          SIDE is CHARACTER*1
 
63
*>          = 'L': apply Q or Q**T from the Left
 
64
*>          = 'R': apply Q or Q**T from the Right
 
65
*> \endverbatim
 
66
*>
 
67
*> \param[in] TRANS
 
68
*> \verbatim
 
69
*>          TRANS is CHARACTER*1
 
70
*>          = 'N': apply Q  (No transpose)
 
71
*>          = 'T': apply Q**T (Transpose)
 
72
*> \endverbatim
 
73
*>
 
74
*> \param[in] M
 
75
*> \verbatim
 
76
*>          M is INTEGER
 
77
*>          The number of rows of the matrix C. M >= 0.
 
78
*> \endverbatim
 
79
*>
 
80
*> \param[in] N
 
81
*> \verbatim
 
82
*>          N is INTEGER
 
83
*>          The number of columns of the matrix C. N >= 0.
 
84
*> \endverbatim
 
85
*>
 
86
*> \param[in] K
 
87
*> \verbatim
 
88
*>          K is INTEGER
 
89
*>          The number of elementary reflectors whose product defines
 
90
*>          the matrix Q.
 
91
*>          If SIDE = 'L', M >= K >= 0;
 
92
*>          if SIDE = 'R', N >= K >= 0.
 
93
*> \endverbatim
 
94
*>
 
95
*> \param[in] A
 
96
*> \verbatim
 
97
*>          A is DOUBLE PRECISION array, dimension (LDA,K)
 
98
*>          The i-th column must contain the vector which defines the
 
99
*>          elementary reflector H(i), for i = 1,2,...,k, as returned by
 
100
*>          DGEQRF in the first k columns of its array argument A.
 
101
*>          A is modified by the routine but restored on exit.
 
102
*> \endverbatim
 
103
*>
 
104
*> \param[in] LDA
 
105
*> \verbatim
 
106
*>          LDA is INTEGER
 
107
*>          The leading dimension of the array A.
 
108
*>          If SIDE = 'L', LDA >= max(1,M);
 
109
*>          if SIDE = 'R', LDA >= max(1,N).
 
110
*> \endverbatim
 
111
*>
 
112
*> \param[in] TAU
 
113
*> \verbatim
 
114
*>          TAU is DOUBLE PRECISION array, dimension (K)
 
115
*>          TAU(i) must contain the scalar factor of the elementary
 
116
*>          reflector H(i), as returned by DGEQRF.
 
117
*> \endverbatim
 
118
*>
 
119
*> \param[in,out] C
 
120
*> \verbatim
 
121
*>          C is DOUBLE PRECISION array, dimension (LDC,N)
 
122
*>          On entry, the m by n matrix C.
 
123
*>          On exit, C is overwritten by Q*C or Q**T*C or C*Q**T or C*Q.
 
124
*> \endverbatim
 
125
*>
 
126
*> \param[in] LDC
 
127
*> \verbatim
 
128
*>          LDC is INTEGER
 
129
*>          The leading dimension of the array C. LDC >= max(1,M).
 
130
*> \endverbatim
 
131
*>
 
132
*> \param[out] WORK
 
133
*> \verbatim
 
134
*>          WORK is DOUBLE PRECISION array, dimension
 
135
*>                                   (N) if SIDE = 'L',
 
136
*>                                   (M) if SIDE = 'R'
 
137
*> \endverbatim
 
138
*>
 
139
*> \param[out] INFO
 
140
*> \verbatim
 
141
*>          INFO is INTEGER
 
142
*>          = 0: successful exit
 
143
*>          < 0: if INFO = -i, the i-th argument had an illegal value
 
144
*> \endverbatim
 
145
*
 
146
*  Authors:
 
147
*  ========
 
148
*
 
149
*> \author Univ. of Tennessee 
 
150
*> \author Univ. of California Berkeley 
 
151
*> \author Univ. of Colorado Denver 
 
152
*> \author NAG Ltd. 
 
153
*
 
154
*> \date September 2012
 
155
*
 
156
*> \ingroup doubleOTHERcomputational
 
157
*
 
158
*  =====================================================================
 
159
      SUBROUTINE DORM2R( SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC,
 
160
     $                   WORK, INFO )
 
161
*
 
162
*  -- LAPACK computational routine (version 3.4.2) --
 
163
*  -- LAPACK is a software package provided by Univ. of Tennessee,    --
 
164
*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
 
165
*     September 2012
 
166
*
 
167
*     .. Scalar Arguments ..
 
168
      CHARACTER          SIDE, TRANS
 
169
      INTEGER            INFO, K, LDA, LDC, M, N
 
170
*     ..
 
171
*     .. Array Arguments ..
 
172
      DOUBLE PRECISION   A( LDA, * ), C( LDC, * ), TAU( * ), WORK( * )
 
173
*     ..
 
174
*
 
175
*  =====================================================================
 
176
*
 
177
*     .. Parameters ..
 
178
      DOUBLE PRECISION   ONE
 
179
      PARAMETER          ( ONE = 1.0D+0 )
 
180
*     ..
 
181
*     .. Local Scalars ..
 
182
      LOGICAL            LEFT, NOTRAN
 
183
      INTEGER            I, I1, I2, I3, IC, JC, MI, NI, NQ
 
184
      DOUBLE PRECISION   AII
 
185
*     ..
 
186
*     .. External Functions ..
 
187
      LOGICAL            LSAME
 
188
      EXTERNAL           LSAME
 
189
*     ..
 
190
*     .. External Subroutines ..
 
191
      EXTERNAL           DLARF, XERBLA
 
192
*     ..
 
193
*     .. Intrinsic Functions ..
 
194
      INTRINSIC          MAX
 
195
*     ..
 
196
*     .. Executable Statements ..
 
197
*
 
198
*     Test the input arguments
 
199
*
 
200
      INFO = 0
 
201
      LEFT = LSAME( SIDE, 'L' )
 
202
      NOTRAN = LSAME( TRANS, 'N' )
 
203
*
 
204
*     NQ is the order of Q
 
205
*
 
206
      IF( LEFT ) THEN
 
207
         NQ = M
 
208
      ELSE
 
209
         NQ = N
 
210
      END IF
 
211
      IF( .NOT.LEFT .AND. .NOT.LSAME( SIDE, 'R' ) ) THEN
 
212
         INFO = -1
 
213
      ELSE IF( .NOT.NOTRAN .AND. .NOT.LSAME( TRANS, 'T' ) ) THEN
 
214
         INFO = -2
 
215
      ELSE IF( M.LT.0 ) THEN
 
216
         INFO = -3
 
217
      ELSE IF( N.LT.0 ) THEN
 
218
         INFO = -4
 
219
      ELSE IF( K.LT.0 .OR. K.GT.NQ ) THEN
 
220
         INFO = -5
 
221
      ELSE IF( LDA.LT.MAX( 1, NQ ) ) THEN
 
222
         INFO = -7
 
223
      ELSE IF( LDC.LT.MAX( 1, M ) ) THEN
 
224
         INFO = -10
 
225
      END IF
 
226
      IF( INFO.NE.0 ) THEN
 
227
         CALL XERBLA( 'DORM2R', -INFO )
 
228
         RETURN
 
229
      END IF
 
230
*
 
231
*     Quick return if possible
 
232
*
 
233
      IF( M.EQ.0 .OR. N.EQ.0 .OR. K.EQ.0 )
 
234
     $   RETURN
 
235
*
 
236
      IF( ( LEFT .AND. .NOT.NOTRAN ) .OR. ( .NOT.LEFT .AND. NOTRAN ) )
 
237
     $     THEN
 
238
         I1 = 1
 
239
         I2 = K
 
240
         I3 = 1
 
241
      ELSE
 
242
         I1 = K
 
243
         I2 = 1
 
244
         I3 = -1
 
245
      END IF
 
246
*
 
247
      IF( LEFT ) THEN
 
248
         NI = N
 
249
         JC = 1
 
250
      ELSE
 
251
         MI = M
 
252
         IC = 1
 
253
      END IF
 
254
*
 
255
      DO 10 I = I1, I2, I3
 
256
         IF( LEFT ) THEN
 
257
*
 
258
*           H(i) is applied to C(i:m,1:n)
 
259
*
 
260
            MI = M - I + 1
 
261
            IC = I
 
262
         ELSE
 
263
*
 
264
*           H(i) is applied to C(1:m,i:n)
 
265
*
 
266
            NI = N - I + 1
 
267
            JC = I
 
268
         END IF
 
269
*
 
270
*        Apply H(i)
 
271
*
 
272
         AII = A( I, I )
 
273
         A( I, I ) = ONE
 
274
         CALL DLARF( SIDE, MI, NI, A( I, I ), 1, TAU( I ), C( IC, JC ),
 
275
     $               LDC, WORK )
 
276
         A( I, I ) = AII
 
277
   10 CONTINUE
 
278
      RETURN
 
279
*
 
280
*     End of DORM2R
 
281
*
 
282
      END