LAPACK 3.12.0 LAPACK: Linear Algebra PACKage
Searching...
No Matches

## ◆ zlaqz3()

 subroutine zlaqz3 ( logical, intent(in) ilschur, logical, intent(in) ilq, logical, intent(in) ilz, integer, intent(in) n, integer, intent(in) ilo, integer, intent(in) ihi, integer, intent(in) nshifts, integer, intent(in) nblock_desired, complex*16, dimension( * ), intent(inout) alpha, complex*16, dimension( * ), intent(inout) beta, complex*16, dimension( lda, * ), intent(inout) a, integer, intent(in) lda, complex*16, dimension( ldb, * ), intent(inout) b, integer, intent(in) ldb, complex*16, dimension( ldq, * ), intent(inout) q, integer, intent(in) ldq, complex*16, dimension( ldz, * ), intent(inout) z, integer, intent(in) ldz, complex*16, dimension( ldqc, * ), intent(inout) qc, integer, intent(in) ldqc, complex*16, dimension( ldzc, * ), intent(inout) zc, integer, intent(in) ldzc, complex*16, dimension( * ), intent(inout) work, integer, intent(in) lwork, integer, intent(out) info )

ZLAQZ3

Purpose:
` ZLAQZ3 Executes a single multishift QZ sweep`
Parameters
 [in] ILSCHUR ``` ILSCHUR is LOGICAL Determines whether or not to update the full Schur form``` [in] ILQ ``` ILQ is LOGICAL Determines whether or not to update the matrix Q``` [in] ILZ ``` ILZ is LOGICAL Determines whether or not to update the matrix Z``` [in] N ``` N is INTEGER The order of the matrices A, B, Q, and Z. N >= 0.``` [in] ILO ` ILO is INTEGER` [in] IHI ` IHI is INTEGER` [in] NSHIFTS ``` NSHIFTS is INTEGER The desired number of shifts to use``` [in] NBLOCK_DESIRED ``` NBLOCK_DESIRED is INTEGER The desired size of the computational windows``` [in] ALPHA ``` ALPHA is COMPLEX*16 array. SR contains the alpha parts of the shifts to use.``` [in] BETA ``` BETA is COMPLEX*16 array. SS contains the scale of the shifts to use.``` [in,out] A ` A is COMPLEX*16 array, dimension (LDA, N)` [in] LDA ``` LDA is INTEGER The leading dimension of the array A. LDA >= max( 1, N ).``` [in,out] B ` B is COMPLEX*16 array, dimension (LDB, N)` [in] LDB ``` LDB is INTEGER The leading dimension of the array B. LDB >= max( 1, N ).``` [in,out] Q ` Q is COMPLEX*16 array, dimension (LDQ, N)` [in] LDQ ` LDQ is INTEGER` [in,out] Z ` Z is COMPLEX*16 array, dimension (LDZ, N)` [in] LDZ ` LDZ is INTEGER` [in,out] QC ` QC is COMPLEX*16 array, dimension (LDQC, NBLOCK_DESIRED)` [in] LDQC ` LDQC is INTEGER` [in,out] ZC ` ZC is COMPLEX*16 array, dimension (LDZC, NBLOCK_DESIRED)` [in] LDZC ` LDZ is INTEGER` [out] WORK ``` WORK is COMPLEX*16 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. LWORK >= max(1,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```
Date
May 2020

Definition at line 204 of file zlaqz3.f.

208 IMPLICIT NONE
209
210* Function arguments
211 LOGICAL, INTENT( IN ) :: ILSCHUR, ILQ, ILZ
212 INTEGER, INTENT( IN ) :: N, ILO, IHI, LDA, LDB, LDQ, LDZ, LWORK,
213 \$ NSHIFTS, NBLOCK_DESIRED, LDQC, LDZC
214
215 COMPLEX*16, INTENT( INOUT ) :: A( LDA, * ), B( LDB, * ), Q( LDQ,
216 \$ * ), Z( LDZ, * ), QC( LDQC, * ), ZC( LDZC, * ), WORK( * ),
217 \$ ALPHA( * ), BETA( * )
218
219 INTEGER, INTENT( OUT ) :: INFO
220
221* Parameters
222 COMPLEX*16 CZERO, CONE
223 parameter( czero = ( 0.0d+0, 0.0d+0 ), cone = ( 1.0d+0,
224 \$ 0.0d+0 ) )
225 DOUBLE PRECISION :: ZERO, ONE, HALF
226 parameter( zero = 0.0d0, one = 1.0d0, half = 0.5d0 )
227
228* Local scalars
229 INTEGER :: I, J, NS, ISTARTM, ISTOPM, SHEIGHT, SWIDTH, K, NP,
230 \$ ISTARTB, ISTOPB, ISHIFT, NBLOCK, NPOS
231 DOUBLE PRECISION :: SAFMIN, SAFMAX, C, SCALE
232 COMPLEX*16 :: TEMP, TEMP2, TEMP3, S
233
234* External Functions
235 EXTERNAL :: xerbla, zlaset, zlartg, zrot, zlaqz1, zgemm, zlacpy
236 DOUBLE PRECISION, EXTERNAL :: DLAMCH
237
238 info = 0
239 IF ( nblock_desired .LT. nshifts+1 ) THEN
240 info = -8
241 END IF
242 IF ( lwork .EQ.-1 ) THEN
243* workspace query, quick return
244 work( 1 ) = n*nblock_desired
245 RETURN
246 ELSE IF ( lwork .LT. n*nblock_desired ) THEN
247 info = -25
248 END IF
249
250 IF( info.NE.0 ) THEN
251 CALL xerbla( 'ZLAQZ3', -info )
252 RETURN
253 END IF
254
255*
256* Executable statements
257*
258
259* Get machine constants
260 safmin = dlamch( 'SAFE MINIMUM' )
261 safmax = one/safmin
262
263 IF ( ilo .GE. ihi ) THEN
264 RETURN
265 END IF
266
267 IF ( ilschur ) THEN
268 istartm = 1
269 istopm = n
270 ELSE
271 istartm = ilo
272 istopm = ihi
273 END IF
274
275 ns = nshifts
276 npos = max( nblock_desired-ns, 1 )
277
278
279* The following block introduces the shifts and chases
280* them down one by one just enough to make space for
281* the other shifts. The near-the-diagonal block is
282* of size (ns+1) x ns.
283
284 CALL zlaset( 'FULL', ns+1, ns+1, czero, cone, qc, ldqc )
285 CALL zlaset( 'FULL', ns, ns, czero, cone, zc, ldzc )
286
287 DO i = 1, ns
288* Introduce the shift
289 scale = sqrt( abs( alpha( i ) ) ) * sqrt( abs( beta( i ) ) )
290 IF( scale .GE. safmin .AND. scale .LE. safmax ) THEN
291 alpha( i ) = alpha( i )/scale
292 beta( i ) = beta( i )/scale
293 END IF
294
295 temp2 = beta( i )*a( ilo, ilo )-alpha( i )*b( ilo, ilo )
296 temp3 = beta( i )*a( ilo+1, ilo )
297
298 IF ( abs( temp2 ) .GT. safmax .OR.
299 \$ abs( temp3 ) .GT. safmax ) THEN
300 temp2 = cone
301 temp3 = czero
302 END IF
303
304 CALL zlartg( temp2, temp3, c, s, temp )
305 CALL zrot( ns, a( ilo, ilo ), lda, a( ilo+1, ilo ), lda, c,
306 \$ s )
307 CALL zrot( ns, b( ilo, ilo ), ldb, b( ilo+1, ilo ), ldb, c,
308 \$ s )
309 CALL zrot( ns+1, qc( 1, 1 ), 1, qc( 1, 2 ), 1, c,
310 \$ dconjg( s ) )
311
312* Chase the shift down
313 DO j = 1, ns-i
314
315 CALL zlaqz1( .true., .true., j, 1, ns, ihi-ilo+1, a( ilo,
316 \$ ilo ), lda, b( ilo, ilo ), ldb, ns+1, 1, qc,
317 \$ ldqc, ns, 1, zc, ldzc )
318
319 END DO
320
321 END DO
322
323* Update the rest of the pencil
324
325* Update A(ilo:ilo+ns,ilo+ns:istopm) and B(ilo:ilo+ns,ilo+ns:istopm)
326* from the left with Qc(1:ns+1,1:ns+1)'
327 sheight = ns+1
328 swidth = istopm-( ilo+ns )+1
329 IF ( swidth > 0 ) THEN
330 CALL zgemm( 'C', 'N', sheight, swidth, sheight, cone, qc, ldqc,
331 \$ a( ilo, ilo+ns ), lda, czero, work, sheight )
332 CALL zlacpy( 'ALL', sheight, swidth, work, sheight, a( ilo,
333 \$ ilo+ns ), lda )
334 CALL zgemm( 'C', 'N', sheight, swidth, sheight, cone, qc, ldqc,
335 \$ b( ilo, ilo+ns ), ldb, czero, work, sheight )
336 CALL zlacpy( 'ALL', sheight, swidth, work, sheight, b( ilo,
337 \$ ilo+ns ), ldb )
338 END IF
339 IF ( ilq ) THEN
340 CALL zgemm( 'N', 'N', n, sheight, sheight, cone, q( 1, ilo ),
341 \$ ldq, qc, ldqc, czero, work, n )
342 CALL zlacpy( 'ALL', n, sheight, work, n, q( 1, ilo ), ldq )
343 END IF
344
345* Update A(istartm:ilo-1,ilo:ilo+ns-1) and B(istartm:ilo-1,ilo:ilo+ns-1)
346* from the right with Zc(1:ns,1:ns)
347 sheight = ilo-1-istartm+1
348 swidth = ns
349 IF ( sheight > 0 ) THEN
350 CALL zgemm( 'N', 'N', sheight, swidth, swidth, cone,
351 \$ a( istartm, ilo ), lda, zc, ldzc, czero, work,
352 \$ sheight )
353 CALL zlacpy( 'ALL', sheight, swidth, work, sheight, a( istartm,
354 \$ ilo ), lda )
355 CALL zgemm( 'N', 'N', sheight, swidth, swidth, cone,
356 \$ b( istartm, ilo ), ldb, zc, ldzc, czero, work,
357 \$ sheight )
358 CALL zlacpy( 'ALL', sheight, swidth, work, sheight, b( istartm,
359 \$ ilo ), ldb )
360 END IF
361 IF ( ilz ) THEN
362 CALL zgemm( 'N', 'N', n, swidth, swidth, cone, z( 1, ilo ),
363 \$ ldz, zc, ldzc, czero, work, n )
364 CALL zlacpy( 'ALL', n, swidth, work, n, z( 1, ilo ), ldz )
365 END IF
366
367* The following block chases the shifts down to the bottom
368* right block. If possible, a shift is moved down npos
369* positions at a time
370
371 k = ilo
372 DO WHILE ( k < ihi-ns )
373 np = min( ihi-ns-k, npos )
374* Size of the near-the-diagonal block
375 nblock = ns+np
376* istartb points to the first row we will be updating
377 istartb = k+1
378* istopb points to the last column we will be updating
379 istopb = k+nblock-1
380
381 CALL zlaset( 'FULL', ns+np, ns+np, czero, cone, qc, ldqc )
382 CALL zlaset( 'FULL', ns+np, ns+np, czero, cone, zc, ldzc )
383
384* Near the diagonal shift chase
385 DO i = ns-1, 0, -1
386 DO j = 0, np-1
387* Move down the block with index k+i+j, updating
388* the (ns+np x ns+np) block:
389* (k:k+ns+np,k:k+ns+np-1)
390 CALL zlaqz1( .true., .true., k+i+j, istartb, istopb, ihi,
391 \$ a, lda, b, ldb, nblock, k+1, qc, ldqc,
392 \$ nblock, k, zc, ldzc )
393 END DO
394 END DO
395
396* Update rest of the pencil
397
398* Update A(k+1:k+ns+np, k+ns+np:istopm) and
399* B(k+1:k+ns+np, k+ns+np:istopm)
400* from the left with Qc(1:ns+np,1:ns+np)'
401 sheight = ns+np
402 swidth = istopm-( k+ns+np )+1
403 IF ( swidth > 0 ) THEN
404 CALL zgemm( 'C', 'N', sheight, swidth, sheight, cone, qc,
405 \$ ldqc, a( k+1, k+ns+np ), lda, czero, work,
406 \$ sheight )
407 CALL zlacpy( 'ALL', sheight, swidth, work, sheight, a( k+1,
408 \$ k+ns+np ), lda )
409 CALL zgemm( 'C', 'N', sheight, swidth, sheight, cone, qc,
410 \$ ldqc, b( k+1, k+ns+np ), ldb, czero, work,
411 \$ sheight )
412 CALL zlacpy( 'ALL', sheight, swidth, work, sheight, b( k+1,
413 \$ k+ns+np ), ldb )
414 END IF
415 IF ( ilq ) THEN
416 CALL zgemm( 'N', 'N', n, nblock, nblock, cone, q( 1, k+1 ),
417 \$ ldq, qc, ldqc, czero, work, n )
418 CALL zlacpy( 'ALL', n, nblock, work, n, q( 1, k+1 ), ldq )
419 END IF
420
421* Update A(istartm:k,k:k+ns+npos-1) and B(istartm:k,k:k+ns+npos-1)
422* from the right with Zc(1:ns+np,1:ns+np)
423 sheight = k-istartm+1
424 swidth = nblock
425 IF ( sheight > 0 ) THEN
426 CALL zgemm( 'N', 'N', sheight, swidth, swidth, cone,
427 \$ a( istartm, k ), lda, zc, ldzc, czero, work,
428 \$ sheight )
429 CALL zlacpy( 'ALL', sheight, swidth, work, sheight,
430 \$ a( istartm, k ), lda )
431 CALL zgemm( 'N', 'N', sheight, swidth, swidth, cone,
432 \$ b( istartm, k ), ldb, zc, ldzc, czero, work,
433 \$ sheight )
434 CALL zlacpy( 'ALL', sheight, swidth, work, sheight,
435 \$ b( istartm, k ), ldb )
436 END IF
437 IF ( ilz ) THEN
438 CALL zgemm( 'N', 'N', n, nblock, nblock, cone, z( 1, k ),
439 \$ ldz, zc, ldzc, czero, work, n )
440 CALL zlacpy( 'ALL', n, nblock, work, n, z( 1, k ), ldz )
441 END IF
442
443 k = k+np
444
445 END DO
446
447* The following block removes the shifts from the bottom right corner
448* one by one. Updates are initially applied to A(ihi-ns+1:ihi,ihi-ns:ihi).
449
450 CALL zlaset( 'FULL', ns, ns, czero, cone, qc, ldqc )
451 CALL zlaset( 'FULL', ns+1, ns+1, czero, cone, zc, ldzc )
452
453* istartb points to the first row we will be updating
454 istartb = ihi-ns+1
455* istopb points to the last column we will be updating
456 istopb = ihi
457
458 DO i = 1, ns
459* Chase the shift down to the bottom right corner
460 DO ishift = ihi-i, ihi-1
461 CALL zlaqz1( .true., .true., ishift, istartb, istopb, ihi,
462 \$ a, lda, b, ldb, ns, ihi-ns+1, qc, ldqc, ns+1,
463 \$ ihi-ns, zc, ldzc )
464 END DO
465
466 END DO
467
468* Update rest of the pencil
469
470* Update A(ihi-ns+1:ihi, ihi+1:istopm)
471* from the left with Qc(1:ns,1:ns)'
472 sheight = ns
473 swidth = istopm-( ihi+1 )+1
474 IF ( swidth > 0 ) THEN
475 CALL zgemm( 'C', 'N', sheight, swidth, sheight, cone, qc, ldqc,
476 \$ a( ihi-ns+1, ihi+1 ), lda, czero, work, sheight )
477 CALL zlacpy( 'ALL', sheight, swidth, work, sheight,
478 \$ a( ihi-ns+1, ihi+1 ), lda )
479 CALL zgemm( 'C', 'N', sheight, swidth, sheight, cone, qc, ldqc,
480 \$ b( ihi-ns+1, ihi+1 ), ldb, czero, work, sheight )
481 CALL zlacpy( 'ALL', sheight, swidth, work, sheight,
482 \$ b( ihi-ns+1, ihi+1 ), ldb )
483 END IF
484 IF ( ilq ) THEN
485 CALL zgemm( 'N', 'N', n, ns, ns, cone, q( 1, ihi-ns+1 ), ldq,
486 \$ qc, ldqc, czero, work, n )
487 CALL zlacpy( 'ALL', n, ns, work, n, q( 1, ihi-ns+1 ), ldq )
488 END IF
489
490* Update A(istartm:ihi-ns,ihi-ns:ihi)
491* from the right with Zc(1:ns+1,1:ns+1)
492 sheight = ihi-ns-istartm+1
493 swidth = ns+1
494 IF ( sheight > 0 ) THEN
495 CALL zgemm( 'N', 'N', sheight, swidth, swidth, cone,
496 \$ a( istartm, ihi-ns ), lda, zc, ldzc, czero, work,
497 \$ sheight )
498 CALL zlacpy( 'ALL', sheight, swidth, work, sheight, a( istartm,
499 \$ ihi-ns ), lda )
500 CALL zgemm( 'N', 'N', sheight, swidth, swidth, cone,
501 \$ b( istartm, ihi-ns ), ldb, zc, ldzc, czero, work,
502 \$ sheight )
503 CALL zlacpy( 'ALL', sheight, swidth, work, sheight, b( istartm,
504 \$ ihi-ns ), ldb )
505 END IF
506 IF ( ilz ) THEN
507 CALL zgemm( 'N', 'N', n, ns+1, ns+1, cone, z( 1, ihi-ns ), ldz,
508 \$ zc, ldzc, czero, work, n )
509 CALL zlacpy( 'ALL', n, ns+1, work, n, z( 1, ihi-ns ), ldz )
510 END IF
511
subroutine xerbla(srname, info)
Definition cblat2.f:3285
subroutine zgemm(transa, transb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc)
ZGEMM
Definition zgemm.f:188
subroutine zlacpy(uplo, m, n, a, lda, b, ldb)
ZLACPY copies all or part of one two-dimensional array to another.
Definition zlacpy.f:103
double precision function dlamch(cmach)
DLAMCH
Definition dlamch.f:69
subroutine zlaqz1(ilq, ilz, k, istartm, istopm, ihi, a, lda, b, ldb, nq, qstart, q, ldq, nz, zstart, z, ldz)
ZLAQZ1
Definition zlaqz1.f:173
subroutine zlartg(f, g, c, s, r)
ZLARTG generates a plane rotation with real cosine and complex sine.
Definition zlartg.f90:116
subroutine zlaset(uplo, m, n, alpha, beta, a, lda)
ZLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition zlaset.f:106
subroutine zrot(n, cx, incx, cy, incy, c, s)
ZROT applies a plane rotation with real cosine and complex sine to a pair of complex vectors.
Definition zrot.f:103
Here is the call graph for this function:
Here is the caller graph for this function: