209 SUBROUTINE sorbdb4( M, P, Q, X11, LDX11, X21, LDX21, THETA,
211 $ TAUP1, TAUP2, TAUQ1, PHANTOM, WORK, LWORK,
219 INTEGER INFO, LWORK, M, P, Q, LDX11, LDX21
222 REAL PHI(*), THETA(*)
223 REAL PHANTOM(*), TAUP1(*), TAUP2(*), TAUQ1(*),
224 $ WORK(*), X11(LDX11,*), X21(LDX21,*)
231 PARAMETER ( NEGONE = -1.0e0, zero = 0.0e0 )
235 INTEGER CHILDINFO, I, ILARF, IORBDB5, J, LLARF,
236 $ LORBDB5, LWORKMIN, LWORKOPT
248 INTRINSIC atan2, cos, max, sin, sqrt
255 lquery = lwork .EQ. -1
259 ELSE IF( p .LT. m-q .OR. m-p .LT. m-q )
THEN
261 ELSE IF( q .LT. m-q .OR. q .GT. m )
THEN
263 ELSE IF( ldx11 .LT. max( 1, p ) )
THEN
265 ELSE IF( ldx21 .LT. max( 1, m-p ) )
THEN
271 IF( info .EQ. 0 )
THEN
273 llarf = max( q-1, p-1, m-p-1 )
276 lworkopt = ilarf + llarf - 1
277 lworkopt = max( lworkopt, iorbdb5 + lorbdb5 - 1 )
279 work(1) = real( lworkopt )
280 IF( lwork .LT. lworkmin .AND. .NOT.lquery )
THEN
284 IF( info .NE. 0 )
THEN
285 CALL xerbla(
'SORBDB4', -info )
287 ELSE IF( lquery )
THEN
299 CALL sorbdb5( p, m-p, q, phantom(1), 1, phantom(p+1), 1,
300 $ x11, ldx11, x21, ldx21, work(iorbdb5),
301 $ lorbdb5, childinfo )
302 CALL sscal( p, negone, phantom(1), 1 )
303 CALL slarfgp( p, phantom(1), phantom(2), 1, taup1(1) )
304 CALL slarfgp( m-p, phantom(p+1), phantom(p+2), 1,
306 theta(i) = atan2( phantom(1), phantom(p+1) )
309 CALL slarf1f(
'L', p, q, phantom(1), 1, taup1(1), x11,
310 $ ldx11, work(ilarf) )
311 CALL slarf1f(
'L', m-p, q, phantom(p+1), 1, taup2(1),
312 $ x21, ldx21, work(ilarf) )
314 CALL sorbdb5( p-i+1, m-p-i+1, q-i+1, x11(i,i-1), 1,
315 $ x21(i,i-1), 1, x11(i,i), ldx11, x21(i,i),
316 $ ldx21, work(iorbdb5), lorbdb5, childinfo )
317 CALL sscal( p-i+1, negone, x11(i,i-1), 1 )
318 CALL slarfgp( p-i+1, x11(i,i-1), x11(i+1,i-1), 1,
320 CALL slarfgp( m-p-i+1, x21(i,i-1), x21(i+1,i-1), 1,
322 theta(i) = atan2( x11(i,i-1), x21(i,i-1) )
325 CALL slarf1f(
'L', p-i+1, q-i+1, x11(i,i-1), 1, taup1(i),
326 $ x11(i,i), ldx11, work(ilarf) )
327 CALL slarf1f(
'L', m-p-i+1, q-i+1, x21(i,i-1), 1,
328 $ taup2(i), x21(i,i), ldx21, work(ilarf) )
331 CALL srot( q-i+1, x11(i,i), ldx11, x21(i,i), ldx21, s, -c )
332 CALL slarfgp( q-i+1, x21(i,i), x21(i,i+1), ldx21, tauq1(i) )
334 CALL slarf1f(
'R', p-i, q-i+1, x21(i,i), ldx21, tauq1(i),
335 $ x11(i+1,i), ldx11, work(ilarf) )
336 CALL slarf1f(
'R', m-p-i, q-i+1, x21(i,i), ldx21, tauq1(i),
337 $ x21(i+1,i), ldx21, work(ilarf) )
338 IF( i .LT. m-q )
THEN
339 s = sqrt( snrm2( p-i, x11(i+1,i), 1 )**2
340 $ + snrm2( m-p-i, x21(i+1,i), 1 )**2 )
341 phi(i) = atan2( s, c )
349 CALL slarfgp( q-i+1, x11(i,i), x11(i,i+1), ldx11, tauq1(i) )
350 CALL slarf1f(
'R', p-i, q-i+1, x11(i,i), ldx11, tauq1(i),
351 $ x11(i+1,i), ldx11, work(ilarf) )
352 CALL slarf1f(
'R', q-p, q-i+1, x11(i,i), ldx11, tauq1(i),
353 $ x21(m-q+1,i), ldx21, work(ilarf) )
359 CALL slarfgp( q-i+1, x21(m-q+i-p,i), x21(m-q+i-p,i+1),
362 CALL slarf1f(
'R', q-i, q-i+1, x21(m-q+i-p,i), ldx21,
363 $ tauq1(i), x21(m-q+i-p+1,i), ldx21,