LAPACK  3.10.0
LAPACK: Linear Algebra PACKage

◆ dorm22()

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.

Definition at line 161 of file dorm22.f.

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