212 SUBROUTINE zunbdb4( M, P, Q, X11, LDX11, X21, LDX21, THETA, PHI,
213 $ taup1, taup2, tauq1, phantom, work, lwork,
222 INTEGER INFO, LWORK, M, P, Q, LDX11, LDX21
225 DOUBLE PRECISION PHI(*), THETA(*)
226 COMPLEX*16 PHANTOM(*), TAUP1(*), TAUP2(*), TAUQ1(*),
227 $ work(*), x11(ldx11,*), x21(ldx21,*)
233 COMPLEX*16 NEGONE, ONE, ZERO
234 parameter( negone = (-1.0d0,0.0d0), one = (1.0d0,0.0d0),
235 $ zero = (0.0d0,0.0d0) )
238 DOUBLE PRECISION C, S
239 INTEGER CHILDINFO, I, ILARF, IORBDB5, J, LLARF,
240 $ lorbdb5, lworkmin, lworkopt
247 DOUBLE PRECISION DZNRM2
251 INTRINSIC atan2, cos, max, sin, sqrt
258 lquery = lwork .EQ. -1
262 ELSE IF( p .LT. m-q .OR. m-p .LT. m-q )
THEN
264 ELSE IF( q .LT. m-q .OR. q .GT. m )
THEN
266 ELSE IF( ldx11 .LT. max( 1, p ) )
THEN
268 ELSE IF( ldx21 .LT. max( 1, m-p ) )
THEN
274 IF( info .EQ. 0 )
THEN
276 llarf = max( q-1, p-1, m-p-1 )
279 lworkopt = ilarf + llarf - 1
280 lworkopt = max( lworkopt, iorbdb5 + lorbdb5 - 1 )
283 IF( lwork .LT. lworkmin .AND. .NOT.lquery )
THEN
287 IF( info .NE. 0 )
THEN
288 CALL
xerbla(
'ZUNBDB4', -info )
290 ELSE IF( lquery )
THEN
302 CALL
zunbdb5( p, m-p, q, phantom(1), 1, phantom(p+1), 1,
303 $ x11, ldx11, x21, ldx21, work(iorbdb5),
304 $ lorbdb5, childinfo )
305 CALL
zscal( p, negone, phantom(1), 1 )
306 CALL
zlarfgp( p, phantom(1), phantom(2), 1, taup1(1) )
307 CALL
zlarfgp( m-p, phantom(p+1), phantom(p+2), 1, taup2(1) )
308 theta(i) = atan2( dble( phantom(1) ), dble( phantom(p+1) ) )
313 CALL
zlarf(
'L', p, q, phantom(1), 1, dconjg(taup1(1)), x11,
314 $ ldx11, work(ilarf) )
315 CALL
zlarf(
'L', m-p, q, phantom(p+1), 1, dconjg(taup2(1)),
316 $ x21, ldx21, work(ilarf) )
318 CALL
zunbdb5( p-i+1, m-p-i+1, q-i+1, x11(i,i-1), 1,
319 $ x21(i,i-1), 1, x11(i,i), ldx11, x21(i,i),
320 $ ldx21, work(iorbdb5), lorbdb5, childinfo )
321 CALL
zscal( p-i+1, negone, x11(i,i-1), 1 )
322 CALL
zlarfgp( p-i+1, x11(i,i-1), x11(i+1,i-1), 1, taup1(i) )
323 CALL
zlarfgp( m-p-i+1, x21(i,i-1), x21(i+1,i-1), 1,
325 theta(i) = atan2( dble( x11(i,i-1) ), dble( x21(i,i-1) ) )
330 CALL
zlarf(
'L', p-i+1, q-i+1, x11(i,i-1), 1,
331 $ dconjg(taup1(i)), x11(i,i), ldx11, work(ilarf) )
332 CALL
zlarf(
'L', m-p-i+1, q-i+1, x21(i,i-1), 1,
333 $ dconjg(taup2(i)), x21(i,i), ldx21, work(ilarf) )
336 CALL
zdrot( q-i+1, x11(i,i), ldx11, x21(i,i), ldx21, s, -c )
337 CALL
zlacgv( q-i+1, x21(i,i), ldx21 )
338 CALL
zlarfgp( q-i+1, x21(i,i), x21(i,i+1), ldx21, tauq1(i) )
341 CALL
zlarf(
'R', p-i, q-i+1, x21(i,i), ldx21, tauq1(i),
342 $ x11(i+1,i), ldx11, work(ilarf) )
343 CALL
zlarf(
'R', m-p-i, q-i+1, x21(i,i), ldx21, tauq1(i),
344 $ x21(i+1,i), ldx21, work(ilarf) )
345 CALL
zlacgv( q-i+1, x21(i,i), ldx21 )
346 IF( i .LT. m-q )
THEN
347 s = sqrt( dznrm2( p-i, x11(i+1,i), 1, x11(i+1,i),
348 $ 1 )**2 + dznrm2( m-p-i, x21(i+1,i), 1, x21(i+1,i),
350 phi(i) = atan2( s, c )
358 CALL
zlacgv( q-i+1, x11(i,i), ldx11 )
359 CALL
zlarfgp( q-i+1, x11(i,i), x11(i,i+1), ldx11, tauq1(i) )
361 CALL
zlarf(
'R', p-i, q-i+1, x11(i,i), ldx11, tauq1(i),
362 $ x11(i+1,i), ldx11, work(ilarf) )
363 CALL
zlarf(
'R', q-p, q-i+1, x11(i,i), ldx11, tauq1(i),
364 $ x21(m-q+1,i), ldx21, work(ilarf) )
365 CALL
zlacgv( q-i+1, x11(i,i), ldx11 )
371 CALL
zlacgv( q-i+1, x21(m-q+i-p,i), ldx21 )
372 CALL
zlarfgp( q-i+1, x21(m-q+i-p,i), x21(m-q+i-p,i+1), ldx21,
375 CALL
zlarf(
'R', q-i, q-i+1, x21(m-q+i-p,i), ldx21, tauq1(i),
376 $ x21(m-q+i-p+1,i), ldx21, work(ilarf) )
377 CALL
zlacgv( q-i+1, x21(m-q+i-p,i), ldx21 )
subroutine zunbdb4(M, P, Q, X11, LDX11, X21, LDX21, THETA, PHI, TAUP1, TAUP2, TAUQ1, PHANTOM, WORK, LWORK, INFO)
ZUNBDB4
subroutine zdrot(N, CX, INCX, CY, INCY, C, S)
ZDROT
subroutine zunbdb5(M1, M2, N, X1, INCX1, X2, INCX2, Q1, LDQ1, Q2, LDQ2, WORK, LWORK, INFO)
ZUNBDB5
subroutine zlarfgp(N, ALPHA, X, INCX, TAU)
ZLARFGP generates an elementary reflector (Householder matrix) with non-negatibe beta.
subroutine zlarf(SIDE, M, N, V, INCV, TAU, C, LDC, WORK)
ZLARF applies an elementary reflector to a general rectangular matrix.
subroutine xerbla(SRNAME, INFO)
XERBLA
subroutine zlacgv(N, X, INCX)
ZLACGV conjugates a complex vector.
subroutine zscal(N, ZA, ZX, INCX)
ZSCAL