LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
subroutine dorm22 ( character  SIDE,
character  TRANS,
integer  M,
integer  N,
integer  N1,
integer  N2,
double precision, dimension( ldq, * )  Q,
integer  LDQ,
double precision, dimension( ldc, * )  C,
integer  LDC,
double precision, dimension( * )  WORK,
integer  LWORK,
integer  INFO 
)

DORM22 multiplies a general matrix by a banded orthogonal matrix.

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

Purpose 
  DORM22 overwrites the general real M-by-N matrix C with

                  SIDE = 'L'     SIDE = 'R'
  TRANS = 'N':      Q * C          C * Q
  TRANS = 'T':      Q**T * C       C * Q**T

  where Q is a real orthogonal matrix of order NQ, with NQ = M if
  SIDE = 'L' and NQ = N if SIDE = 'R'.
  The orthogonal matrix Q processes a 2-by-2 block structure

         [  Q11  Q12  ]
     Q = [            ]
         [  Q21  Q22  ],

  where Q12 is an N1-by-N1 lower triangular matrix and Q21 is an
  N2-by-N2 upper triangular matrix.
Parameters
[in]SIDE
          SIDE is CHARACTER*1
          = 'L': apply Q or Q**T from the Left;
          = 'R': apply Q or Q**T from the Right.
[in]TRANS
          TRANS is CHARACTER*1
          = 'N':  apply Q (No transpose);
          = 'C':  apply Q**T (Conjugate transpose).
[in]M
          M is INTEGER
          The number of rows of the matrix C. M >= 0.
[in]N
          N is INTEGER
          The number of columns of the matrix C. N >= 0.
[in]N1
[in]N2
          N1 is INTEGER
          N2 is INTEGER
          The dimension of Q12 and Q21, respectively. N1, N2 >= 0.
          The following requirement must be satisfied:
          N1 + N2 = M if SIDE = 'L' and N1 + N2 = N if SIDE = 'R'.
[in]Q
          Q is DOUBLE PRECISION array, dimension
                                       (LDQ,M) if SIDE = 'L'
                                       (LDQ,N) if SIDE = 'R'
[in]LDQ
          LDQ is INTEGER
          The leading dimension of the array Q.
          LDQ >= max(1,M) if SIDE = 'L'; LDQ >= max(1,N) if SIDE = 'R'.
[in,out]C
          C is DOUBLE PRECISION array, dimension (LDC,N)
          On entry, the M-by-N matrix C.
          On exit, C is overwritten by Q*C or Q**T*C or C*Q**T or C*Q.
[in]LDC
          LDC is INTEGER
          The leading dimension of the array C. LDC >= max(1,M).
[out]WORK
          WORK is DOUBLE PRECISION array, dimension (MAX(1,LWORK))
          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
[in]LWORK
          LWORK is INTEGER
          The dimension of the array WORK.
          If SIDE = 'L', LWORK >= max(1,N);
          if SIDE = 'R', LWORK >= max(1,M).
          For optimum performance LWORK >= M*N.

          If LWORK = -1, then a workspace query is assumed; the routine
          only calculates the optimal size of the WORK array, returns
          this value as the first entry of the WORK array, and no error
          message related to LWORK is issued by XERBLA.
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
January 2015

Definition at line 165 of file dorm22.f.

165 *
166 * -- LAPACK computational routine (version 3.6.0) --
167 * -- LAPACK is a software package provided by Univ. of Tennessee, --
168 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
169 * January 2015
170 *
171  IMPLICIT NONE
172 *
173 * .. Scalar Arguments ..
174  CHARACTER side, trans
175  INTEGER m, n, n1, n2, ldq, ldc, lwork, info
176 * ..
177 * .. Array Arguments ..
178  DOUBLE PRECISION q( ldq, * ), c( ldc, * ), work( * )
179 * ..
180 *
181 * =====================================================================
182 *
183 * .. Parameters ..
184  DOUBLE PRECISION one
185  parameter ( one = 1.0d+0 )
186 *
187 * .. Local Scalars ..
188  LOGICAL left, lquery, notran
189  INTEGER i, ldwork, len, lwkopt, nb, nq, nw
190 * ..
191 * .. External Functions ..
192  LOGICAL lsame
193  EXTERNAL lsame
194 * ..
195 * .. External Subroutines ..
196  EXTERNAL dgemm, dlacpy, dtrmm, xerbla
197 * ..
198 * .. Intrinsic Functions ..
199  INTRINSIC dble, max, min
200 * ..
201 * .. Executable Statements ..
202 *
203 * Test the input arguments
204 *
205  info = 0
206  left = lsame( side, 'L' )
207  notran = lsame( trans, 'N' )
208  lquery = ( lwork.EQ.-1 )
209 *
210 * NQ is the order of Q;
211 * NW is the minimum dimension of WORK.
212 *
213  IF( left ) THEN
214  nq = m
215  ELSE
216  nq = n
217  END IF
218  nw = nq
219  IF( n1.EQ.0 .OR. n2.EQ.0 ) nw = 1
220  IF( .NOT.left .AND. .NOT.lsame( side, 'R' ) ) THEN
221  info = -1
222  ELSE IF( .NOT.lsame( trans, 'N' ) .AND. .NOT.lsame( trans, 'T' ) )
223  $ THEN
224  info = -2
225  ELSE IF( m.LT.0 ) THEN
226  info = -3
227  ELSE IF( n.LT.0 ) THEN
228  info = -4
229  ELSE IF( n1.LT.0 .OR. n1+n2.NE.nq ) THEN
230  info = -5
231  ELSE IF( n2.LT.0 ) THEN
232  info = -6
233  ELSE IF( ldq.LT.max( 1, nq ) ) THEN
234  info = -8
235  ELSE IF( ldc.LT.max( 1, m ) ) THEN
236  info = -10
237  ELSE IF( lwork.LT.nw .AND. .NOT.lquery ) THEN
238  info = -12
239  END IF
240 *
241  IF( info.EQ.0 ) THEN
242  lwkopt = m*n
243  work( 1 ) = dble( lwkopt )
244  END IF
245 *
246  IF( info.NE.0 ) THEN
247  CALL xerbla( 'DORM22', -info )
248  RETURN
249  ELSE IF( lquery ) THEN
250  RETURN
251  END IF
252 *
253 * Quick return if possible
254 *
255  IF( m.EQ.0 .OR. n.EQ.0 ) THEN
256  work( 1 ) = 1
257  RETURN
258  END IF
259 *
260 * Degenerate cases (N1 = 0 or N2 = 0) are handled using DTRMM.
261 *
262  IF( n1.EQ.0 ) THEN
263  CALL dtrmm( side, 'Upper', trans, 'Non-Unit', m, n, one,
264  $ q, ldq, c, ldc )
265  work( 1 ) = one
266  RETURN
267  ELSE IF( n2.EQ.0 ) THEN
268  CALL dtrmm( side, 'Lower', trans, 'Non-Unit', m, n, one,
269  $ q, ldq, c, ldc )
270  work( 1 ) = one
271  RETURN
272  END IF
273 *
274 * Compute the largest chunk size available from the workspace.
275 *
276  nb = max( 1, min( lwork, lwkopt ) / nq )
277 *
278  IF( left ) THEN
279  IF( notran ) THEN
280  DO i = 1, n, nb
281  len = min( nb, n-i+1 )
282  ldwork = m
283 *
284 * Multiply bottom part of C by Q12.
285 *
286  CALL dlacpy( 'All', n1, len, c( n2+1, i ), ldc, work,
287  $ ldwork )
288  CALL dtrmm( 'Left', 'Lower', 'No Transpose', 'Non-Unit',
289  $ n1, len, one, q( 1, n2+1 ), ldq, work,
290  $ ldwork )
291 *
292 * Multiply top part of C by Q11.
293 *
294  CALL dgemm( 'No Transpose', 'No Transpose', n1, len, n2,
295  $ one, q, ldq, c( 1, i ), ldc, one, work,
296  $ ldwork )
297 *
298 * Multiply top part of C by Q21.
299 *
300  CALL dlacpy( 'All', n2, len, c( 1, i ), ldc,
301  $ work( n1+1 ), ldwork )
302  CALL dtrmm( 'Left', 'Upper', 'No Transpose', 'Non-Unit',
303  $ n2, len, one, q( n1+1, 1 ), ldq,
304  $ work( n1+1 ), ldwork )
305 *
306 * Multiply bottom part of C by Q22.
307 *
308  CALL dgemm( 'No Transpose', 'No Transpose', n2, len, n1,
309  $ one, q( n1+1, n2+1 ), ldq, c( n2+1, i ), ldc,
310  $ one, work( n1+1 ), ldwork )
311 *
312 * Copy everything back.
313 *
314  CALL dlacpy( 'All', m, len, work, ldwork, c( 1, i ),
315  $ ldc )
316  END DO
317  ELSE
318  DO i = 1, n, nb
319  len = min( nb, n-i+1 )
320  ldwork = m
321 *
322 * Multiply bottom part of C by Q21**T.
323 *
324  CALL dlacpy( 'All', n2, len, c( n1+1, i ), ldc, work,
325  $ ldwork )
326  CALL dtrmm( 'Left', 'Upper', 'Transpose', 'Non-Unit',
327  $ n2, len, one, q( n1+1, 1 ), ldq, work,
328  $ ldwork )
329 *
330 * Multiply top part of C by Q11**T.
331 *
332  CALL dgemm( 'Transpose', 'No Transpose', n2, len, n1,
333  $ one, q, ldq, c( 1, i ), ldc, one, work,
334  $ ldwork )
335 *
336 * Multiply top part of C by Q12**T.
337 *
338  CALL dlacpy( 'All', n1, len, c( 1, i ), ldc,
339  $ work( n2+1 ), ldwork )
340  CALL dtrmm( 'Left', 'Lower', 'Transpose', 'Non-Unit',
341  $ n1, len, one, q( 1, n2+1 ), ldq,
342  $ work( n2+1 ), ldwork )
343 *
344 * Multiply bottom part of C by Q22**T.
345 *
346  CALL dgemm( 'Transpose', 'No Transpose', n1, len, n2,
347  $ one, q( n1+1, n2+1 ), ldq, c( n1+1, i ), ldc,
348  $ one, work( n2+1 ), ldwork )
349 *
350 * Copy everything back.
351 *
352  CALL dlacpy( 'All', m, len, work, ldwork, c( 1, i ),
353  $ ldc )
354  END DO
355  END IF
356  ELSE
357  IF( notran ) THEN
358  DO i = 1, m, nb
359  len = min( nb, m-i+1 )
360  ldwork = len
361 *
362 * Multiply right part of C by Q21.
363 *
364  CALL dlacpy( 'All', len, n2, c( i, n1+1 ), ldc, work,
365  $ ldwork )
366  CALL dtrmm( 'Right', 'Upper', 'No Transpose', 'Non-Unit',
367  $ len, n2, one, q( n1+1, 1 ), ldq, work,
368  $ ldwork )
369 *
370 * Multiply left part of C by Q11.
371 *
372  CALL dgemm( 'No Transpose', 'No Transpose', len, n2, n1,
373  $ one, c( i, 1 ), ldc, q, ldq, one, work,
374  $ ldwork )
375 *
376 * Multiply left part of C by Q12.
377 *
378  CALL dlacpy( 'All', len, n1, c( i, 1 ), ldc,
379  $ work( 1 + n2*ldwork ), ldwork )
380  CALL dtrmm( 'Right', 'Lower', 'No Transpose', 'Non-Unit',
381  $ len, n1, one, q( 1, n2+1 ), ldq,
382  $ work( 1 + n2*ldwork ), ldwork )
383 *
384 * Multiply right part of C by Q22.
385 *
386  CALL dgemm( 'No Transpose', 'No Transpose', len, n1, n2,
387  $ one, c( i, n1+1 ), ldc, q( n1+1, n2+1 ), ldq,
388  $ one, work( 1 + n2*ldwork ), ldwork )
389 *
390 * Copy everything back.
391 *
392  CALL dlacpy( 'All', len, n, work, ldwork, c( i, 1 ),
393  $ ldc )
394  END DO
395  ELSE
396  DO i = 1, m, nb
397  len = min( nb, m-i+1 )
398  ldwork = len
399 *
400 * Multiply right part of C by Q12**T.
401 *
402  CALL dlacpy( 'All', len, n1, c( i, n2+1 ), ldc, work,
403  $ ldwork )
404  CALL dtrmm( 'Right', 'Lower', 'Transpose', 'Non-Unit',
405  $ len, n1, one, q( 1, n2+1 ), ldq, work,
406  $ ldwork )
407 *
408 * Multiply left part of C by Q11**T.
409 *
410  CALL dgemm( 'No Transpose', 'Transpose', len, n1, n2,
411  $ one, c( i, 1 ), ldc, q, ldq, one, work,
412  $ ldwork )
413 *
414 * Multiply left part of C by Q21**T.
415 *
416  CALL dlacpy( 'All', len, n2, c( i, 1 ), ldc,
417  $ work( 1 + n1*ldwork ), ldwork )
418  CALL dtrmm( 'Right', 'Upper', 'Transpose', 'Non-Unit',
419  $ len, n2, one, q( n1+1, 1 ), ldq,
420  $ work( 1 + n1*ldwork ), ldwork )
421 *
422 * Multiply right part of C by Q22**T.
423 *
424  CALL dgemm( 'No Transpose', 'Transpose', len, n2, n1,
425  $ one, c( i, n2+1 ), ldc, q( n1+1, n2+1 ), ldq,
426  $ one, work( 1 + n1*ldwork ), ldwork )
427 *
428 * Copy everything back.
429 *
430  CALL dlacpy( 'All', len, n, work, ldwork, c( i, 1 ),
431  $ ldc )
432  END DO
433  END IF
434  END IF
435 *
436  work( 1 ) = dble( lwkopt )
437  RETURN
438 *
439 * End of DORM22
440 *
dtrmm
subroutine dtrmm(SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB)
DTRMM
Definition: dtrmm.f:179
dlacpy
subroutine dlacpy(UPLO, M, N, A, LDA, B, LDB)
DLACPY copies all or part of one two-dimensional array to another.
Definition: dlacpy.f:105
dgemm
subroutine dgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
DGEMM
Definition: dgemm.f:189
xerbla
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
lsame
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55

Here is the call graph for this function:

Here is the caller graph for this function: