1 SUBROUTINE pclarzb( SIDE, TRANS, DIRECT, STOREV, M, N, K, L, V,
2 $ IV, JV, DESCV, T, C, IC, JC, DESCC, WORK )
9 CHARACTER DIRECT, SIDE, , TRANS
10 INTEGER IC, IV, JC, JV, K, L, M
13 INTEGER DESCC( * ), ( * )
14 COMPLEX C( * ), T( * ), V( * ), WORK( * )
221 INTEGER BLOCK_CYCLIC_2D, CSRC_, CTXT_, DLEN_, DTYPE_,
222 $ lld_, mb_, m_, nb_, n_, rsrc_
223 parameter( block_cyclic_2d = 1, dlen_ = 9, dtype_ = 1,
224 $ ctxt_ = 2, m_ = 3, n_ = 4, mb_
225 $ rsrc_ = 7, csrc_ = 8, lld_ = 9 )
227 parameter( one = ( 1.0e+0, 0.0e+0 ),
228 $ zero = ( 0.0e+0, 0.0e+0 ) )
232 CHARACTER COLBTOP, TRANST
233 INTEGER ICCOL1, ICCOL2, ICOFFC1, ICOFFC2, ICOFFV,
234 $ icrow1, icrow2, ictxt, iibeg, iic1, iic2,
235 $ iiend, iinxt, iiv, ileft, info, ioffc2, ioffv,
236 $ ipt, ipv, ipw, iroffc1, iroffc2, itop, ivcol,
237 $ ivrow, j, jjbeg, jjend, jjnxt, jjc1, jjc2, jjv,
238 $ ldc, ldv, lv, lw, mbc, mbv, mpc1, mpc2, mpc20,
239 $ mqv, mqv0, mycol, mydist, myrow, nbc, nbv,
240 $ npcol, nprow, nqc1, nqc2, nqcall, nqv
245 $ clamov,
claset, ctrbr2d, ctrbs2d,
254 INTEGER ICEIL, NUMROC
255 EXTERNAL iceil, lsame, numroc
261 IF( m.LE.0 .OR. n.LE.0 .OR. k.LE.0 )
266 ictxt = descc( ctxt_ )
272 IF( .NOT.lsame( direct,
'B' ) )
THEN
274 ELSE IF( .NOT.lsame( storev,
'R' ) )
THEN
278 CALL pxerbla( ictxt,
'PCLARZB', -info
279 CALL blacs_abort( ictxt, 1 )
283 left = lsame( side,
'L' )
284 IF( lsame( trans,
'N' ) )
THEN
290 CALL infog2l( iv, jv, descv, nprow, npcol, myrow, mycol, iiv, jjv,
294 icoffv = mod( jv-1, nbv )
295 nqv = numroc( l+icoffv, nbv, mycol, ivcol, npcol )
299 iiv =
min( iiv, ldv )
300 jjv =
min( jjv,
max( 1, numroc( descv( n_ ), nbv, mycol,
301 $ descv( csrc_ ), npcol ) ) )
302 ioffv = iiv + ( jjv-1 ) * ldv
305 nqcall = numroc( descc( n_ ), nbc, mycol, descc( csrc_ ), npcol )
306 CALL infog2l( ic, jc, descc, nprow, npcol, myrow, mycol, iic1,
307 $ jjc1, icrow1, iccol1 )
309 iic1 =
min( iic1, ldc )
310 jjc1 =
min( jjc1,
max( 1, nqcall ) )
313 iroffc1 = mod( ic-1, mbc )
314 mpc1 = numroc( k+iroffc1, mbc, myrow, icrow1, nprow )
315 IF( myrow.EQ.icrow1 )
316 $ mpc1 = mpc1 - iroffc1
317 icoffc1 = mod( jc-1, nbc )
318 nqc1 = numroc( n+icoffc1, nbc, mycol, iccol1, npcol )
319 IF( mycol.EQ.iccol1 )
320 $ nqc1 = nqc1 - icoffc1
321 CALL infog2l( ic+m-l, jc, descc, nprow, npcol, myrow, mycol,
322 $ iic2, jjc2, icrow2, iccol2
323 iroffc2 = mod( ic+m-l-1, mbc )
324 mpc2 = numroc( l+iroffc2, mbc, myrow, icrow2, nprow )
325 IF( myrow.EQ.icrow2 )
326 $ mpc2 = mpc2 - iroffc2
330 iroffc1 = mod( ic-1, mbc )
331 mpc1 = numroc( m+iroffc1, mbc, myrow, icrow1, nprow )
332 IF( myrow.EQ.icrow1 )
333 $ mpc1 = mpc1 - iroffc1
334 icoffc1 = mod( jc-1, nbc )
335 nqc1 = numroc( k+icoffc1, nbc, mycol, iccol1, npcol )
336 IF( mycol.EQ.iccol1 )
337 $ nqc1 = nqc1 - icoffc1
338 CALL infog2l( ic, jc+n-l, descc, nprow, npcol, myrow, mycol,
339 $ iic2, jjc2, icrow2, iccol2 )
342 icoffc2 = mod( jc+n-l-1, nbc )
343 nqc2 = numroc( l+icoffc2, nbc, mycol, iccol2, npcol )
344 IF( mycol.EQ.iccol2 )
345 $ nqc2 = nqc2 - icoffc2
347 iic2 =
min( iic2, ldc )
348 jjc2 =
min( jjc2, nqcall )
349 ioffc2 = iic2 + ( jjc2-1 ) * ldc
351 IF( lsame( side,
'L' ) )
THEN
358 mqv0 = numroc( m+icoffv, nbv, mycol, ivcol, npcol )
359 IF( mycol.EQ.ivcol )
THEN
364 IF( myrow.EQ.icrow2 )
THEN
365 mpc20 = mpc2 + iroffc2
376 ipw = ipv + mpc20 * k
381 IF( myrow.EQ.ivrow )
THEN
382 IF( mycol.EQ.ivcol )
THEN
383 CALL clamov(
'All', k, mqv, v( ioffv ), ldv,
384 $ work( ipw+icoffv*lw ), lw )
386 CALL clamov(
'All', k, mqv, v( ioffv ), ldv,
393 CALL pbctran( ictxt,
'Rowwise',
'Conjugate transpose', k,
394 $ m+icoffv, descv( nb_ ), work( ipw ), lw, zero,
395 $ work( ipv ), lv, ivrow, ivcol, icrow2, -1,
400 IF( myrow.EQ.icrow2 )
401 $ ipv = ipv + iroffc2
409 CALL cgemm(
'Transpose''No transpose', nqc2, k, mpc2,
410 $ one, c( ioffc2 ), ldc, work( ipv ), lv, zero,
413 CALL claset(
'All', nqc2, k, zero, zero, work( ipw ), lw
419 mydist = mod( myrow-icrow1+nprow, nprow )
420 itop =
max( 0, mydist * mbc - iroffc1 )
422 iiend = iic1 + mpc1 - 1
423 iinxt =
min( iceil( iibeg, mbc ) * mbc, iiend )
426 IF( iibeg.LE.iinxt )
THEN
428 $ one, c( iibeg+(jjc1-1)*ldc ), ldc, one,
429 $ work( ipw+itop ), lw )
430 mydist = mydist + nprow
431 itop = mydist * mbc - iroffc1
433 iinxt =
min( iinxt+mbc, iiend )
438 CALL cgsum2d( ictxt,
'Columnwise',
' ', nqc2, k, work( ipw ),
443 IF( myrow.EQ.ivrow )
THEN
444 IF( mycol.EQ.ivcol )
THEN
448 CALL ctrbs2d( ictxt,
'Rowwise',
' ',
'Lower',
'Non unit',
451 CALL ctrbr2d( ictxt,
'Rowwise',
' ',
'Lower',
'Non unit',
452 $ k, k, t, mbv, myrow, ivcol )
454 CALL ctrmm(
'Right',
'Lower', transt,
'Non unit', nqc2, k,
455 $ one, t, mbv, work( ipw ), lw )
457 CALL cgebs2d( ictxt,
'Columnwise',
' ', nqc2, k,
460 CALL cgebr2d( ictxt,
'Columnwise',
' ', nqc2, k,
461 $ work( ipw ), lw, ivrow, mycol )
467 mydist = mod( myrow-icrow1+nprow, nprow )
468 itop =
max( 0, mydist * mbc - iroffc1 )
470 iiend = iic1 + mpc1 - 1
471 iinxt =
min( iceil( iibeg, mbc ) * mbc, iiend )
474 IF( iibeg.LE.iinxt )
THEN
475 CALL pbcmatadd( ictxt,
'Transpose', iinxt-iibeg+1, nqc2,
476 $ -one, work( ipw+itop ), lw, one,
477 $ c( iibeg+(jjc1-1)*ldc ), ldc )
478 mydist = mydist + nprow
479 itop = mydist * mbc - iroffc1
481 iinxt =
min( iinxt+mbc, iiend )
491 CALL clacgv( mpc2, work( ipv+(j-1)*lv ), 1 )
493 CALL cgemm(
'No transpose', 'transpose
', MPC2, NQC2, K, -ONE,
494 $ WORK( IPV ), LV, WORK( IPW ), LW, ONE,
512 CALL PB_TOPGET( ICTXT, 'broadcast
', 'columnwise
', COLBTOP )
513.EQ.
IF( MYROWIVROW ) THEN
514 CALL CGEBS2D( ICTXT, 'columnwise
', COLBTOP, K, NQC2,
517 $ CALL CTRBS2D( ICTXT, 'columnwise
', COLBTOP, 'lower
',
518 $ 'non unit
', K, K, T, MBV )
519 CALL CLAMOV( 'all
', K, NQC2, V( IOFFV ), LDV, WORK( IPV ),
522 CALL CGEBR2D( ICTXT, 'columnwise
', COLBTOP, K, NQC2,
523 $ WORK( IPV ), LV, IVROW, MYCOL )
525 $ CALL CTRBR2D( ICTXT, 'columnwise
', COLBTOP, 'lower
',
526 $ 'non unit
', K, K, T, MBV, IVROW, MYCOL )
533 CALL CGEMM( 'no transpose
', 'transpose
', MPC2, K, NQC2,
534 $ ONE, C( IOFFC2 ), LDC, WORK( IPV ), LV, ZERO,
537 CALL CLASET( 'all
', MPC2, K, ZERO, ZERO, WORK( IPW ), LW )
543 MYDIST = MOD( MYCOL-ICCOL1+NPCOL, NPCOL )
544 ILEFT = MAX( 0, MYDIST * NBC - ICOFFC1 )
546 JJEND = JJC1 + NQC1 - 1
547 JJNXT = MIN( ICEIL( JJBEG, NBC ) * NBC, JJEND )
550.LE.
IF( JJBEGJJNXT ) THEN
551 CALL PBCMATADD( ICTXT, 'no transpose
', MPC2,
552 $ JJNXT-JJBEG+1, ONE,
553 $ C( IIC1+(JJBEG-1)*LDC ), LDC, ONE,
554 $ WORK( IPW+ILEFT*LW ), LW )
555 MYDIST = MYDIST + NPCOL
556 ILEFT = MYDIST * NBC - ICOFFC1
558 JJNXT = MIN( JJNXT+NBC, JJEND )
563 CALL CGSUM2D( ICTXT, 'rowwise
', ' ', MPC2, K, WORK( IPW ),
568.EQ.
IF( MYCOLIVCOL ) THEN
570 CALL CLACGV( K-J+1, T( J+(J-1)*MBV ), 1 )
572 CALL CTRMM( 'right
', 'lower
', TRANS, 'non unit
', MPC2, K,
573 $ ONE, T, MBV, WORK( IPW ), LW )
574 CALL CGEBS2D( ICTXT, 'rowwise
', ' ', MPC2, K, WORK( IPW ),
577 CALL CLACGV( K-J+1, T( J+(J-1)*MBV ), 1 )
580 CALL CGEBR2D( ICTXT, 'rowwise
', ' ', MPC2, K, WORK( IPW ),
587 MYDIST = MOD( MYCOL-ICCOL1+NPCOL, NPCOL )
588 ILEFT = MAX( 0, MYDIST * NBC - ICOFFC1 )
590 JJEND = JJC1 + NQC1 - 1
591 JJNXT = MIN( ICEIL( JJBEG, NBC ) * NBC, JJEND )
594.LE.
IF( JJBEGJJNXT ) THEN
595 CALL PBCMATADD( ICTXT, 'no transpose
', MPC2,
596 $ JJNXT-JJBEG+1, -ONE,
597 $ WORK( IPW+ILEFT*LW ), LW, ONE,
598 $ C( IIC1+(JJBEG-1)*LDC ), LDC )
599 MYDIST = MYDIST + NPCOL
600 ILEFT = MYDIST * NBC - ICOFFC1
602 JJNXT = MIN( JJNXT+NBC, JJEND )
612 CALL CLACGV( K, WORK( IPV+(J-1)*LV ), 1 )
615 $ CALL CGEMM( 'no transpose',
'No transpose', mpc2, nqc2, k,
616 $ -one, work( ipw ), lw, work( ipv ), lv, one,