2
    3
    4
    5
    6
    7
    8
    9      INTEGER            II, JJ, M
   10      REAL               H33, H43H34, H44
   11
   12
   13      INTEGER            DESCA( * )
   14      REAL               A( * ), V( * )
   15
   16
   17
   18
   19
   20
   21
   22
   23
   24
   25
   26
   27
   28
   29
   30
   31
   32
   33
   34
   35
   36
   37
   38
   39
   40
   41
   42
   43
   44
   45
   46
   47
   48
   49
   50
   51
   52
   53
   54
   55
   56
   57
   58
   59
   60
   61
   62
   63
   64
   65
   66
   67
   68
   69
   70
   71
   72
   73
   74
   75
   76
   77
   78
   79
   80
   81
   82
   83
   84
   85
   86
   87
   88
   89
   90
   91
   92
   93
   94
   95
   96
   97
   98
   99
  100
  101
  102
  103
  104
  105
  106
  107
  108
  109
  110
  111
  112
  113      INTEGER            BLOCK_CYCLIC_2D, CSRC_, CTXT_, DLEN_, DTYPE_,
  114     $                   LLD_, MB_, M_, NB_, N_, RSRC_
  115      parameter( block_cyclic_2d = 1, dlen_ = 9, dtype_ = 1,
  116     $                     ctxt_ = 2, m_ = 3, n_ = 4, mb_ = 5, nb_ = 6,
  117     $                     rsrc_ = 7, csrc_ = 8, lld_ = 9 )
  118
  119
  120      INTEGER            CONTXT, DOWN, HBL, ICOL, IROW, JSRC, LDA, LEFT,
  121     $                   MODKM1, MYCOL, MYROW, NPCOL, NPROW, NUM, RIGHT,
  122     $                   RSRC, UP
  123      REAL               H22, H33S, H44S, S, V1, V2
  124
  125
  126      REAL               BUF( 4 )
  127      REAL               H11( 1 )
  128      REAL               H12( 1 )
  129      REAL               H21( 1 )
  130      REAL               V3( 1 )
  131
  132
  133      EXTERNAL           blacs_gridinfo, sgerv2d, sgesd2d, 
infog2l 
  134
  135
  136      INTRINSIC          abs, mod
  137
  138
  139
  140      hbl = desca( mb_ )
  141      contxt = desca( ctxt_ )
  142      lda = desca( lld_ )
  143      CALL blacs_gridinfo( contxt, nprow, npcol, myrow, mycol )
  144      left = mod( mycol+npcol-1, npcol )
  145      right = mod( mycol+1, npcol )
  146      up = mod( myrow+nprow-1, nprow )
  147      down = mod( myrow+1, nprow )
  148      num = nprow*npcol
  149
  150
  151
  152      modkm1 = mod( m+1, hbl )
  153      IF( modkm1.EQ.0 ) THEN
  154         IF( ( myrow.EQ.ii ) .AND. ( right.EQ.jj ) .AND.
  155     $       ( npcol.GT.1 ) ) THEN
  156            CALL infog2l( m+2, m+1, desca, nprow, npcol, myrow, mycol,
 
  157     $                    irow, icol, rsrc, jsrc )
  158            buf( 1 ) = a( ( icol-1 )*lda+irow )
  159            CALL sgesd2d( contxt, 1, 1, buf, 1, ii, jj )
  160         END IF
  161         IF( ( down.EQ.ii ) .AND. ( right.EQ.jj ) .AND. ( num.GT.1 ) )
  162     $        THEN
  163            CALL infog2l( m, m, desca, nprow, npcol, myrow, mycol, irow,
 
  164     $                    icol, rsrc, jsrc )
  165            buf( 1 ) = a( ( icol-1 )*lda+irow )
  166            buf( 2 ) = a( ( icol-1 )*lda+irow+1 )
  167            buf( 3 ) = a( icol*lda+irow )
  168            buf( 4 ) = a( icol*lda+irow+1 )
  169            CALL sgesd2d( contxt, 4, 1, buf, 4, ii, jj )
  170         END IF
  171         IF( ( myrow.EQ.ii ) .AND. ( mycol.EQ.jj ) ) THEN
  172            CALL infog2l( m+2, m+2, desca, nprow, npcol, myrow, mycol,
 
  173     $                    irow, icol, rsrc, jsrc )
  174            IF( npcol.GT.1 ) THEN
  175               CALL sgerv2d( contxt, 1, 1, v3, 1, myrow, left )
  176            ELSE
  177               v3( 1 ) = a( ( icol-2 )*lda+irow )
  178            END IF
  179            IF( num.GT.1 ) THEN
  180               CALL sgerv2d( contxt, 4, 1, buf, 4, up, left )
  181               h11( 1 ) = buf( 1 )
  182               h21( 1 ) = buf( 2 )
  183               h12( 1 ) = buf( 3 )
  184               h22 = buf( 4 )
  185            ELSE
  186               h11( 1 ) = a( ( icol-3 )*lda+irow-2 )
  187               h21( 1 ) = a( ( icol-3 )*lda+irow-1 )
  188               h12( 1 ) = a( ( icol-2 )*lda+irow-2 )
  189               h22 = a( ( icol-2 )*lda+irow-1 )
  190            END IF
  191         END IF
  192      END IF
  193      IF( modkm1.EQ.1 ) THEN
  194         IF( ( down.EQ.ii ) .AND. ( right.EQ.jj ) .AND. ( num.GT.1 ) )
  195     $        THEN
  196            CALL infog2l( m, m, desca, nprow, npcol, myrow, mycol, irow,
 
  197     $                    icol, rsrc, jsrc )
  198            CALL sgesd2d( contxt, 1, 1, a( ( icol-1 )*lda+irow ), 1, ii,
  199     $                    jj )
  200         END IF
  201         IF( ( down.EQ.ii ) .AND. ( mycol.EQ.jj ) .AND. ( nprow.GT.1 ) )
  202     $        THEN
  203            CALL infog2l( m, m+1, desca, nprow, npcol, myrow, mycol,
 
  204     $                    irow, icol, rsrc, jsrc )
  205            CALL sgesd2d( contxt, 1, 1, a( ( icol-1 )*lda+irow ), 1, ii,
  206     $                    jj )
  207         END IF
  208         IF( ( myrow.EQ.ii ) .AND. ( right.EQ.jj ) .AND.
  209     $       ( npcol.GT.1 ) ) THEN
  210            CALL infog2l( m+1, m, desca, nprow, npcol, myrow, mycol,
 
  211     $                    irow, icol, rsrc, jsrc )
  212            CALL sgesd2d( contxt, 1, 1, a( ( icol-1 )*lda+irow ), 1, ii,
  213     $                    jj )
  214         END IF
  215         IF( ( myrow.EQ.ii ) .AND. ( mycol.EQ.jj ) ) THEN
  216            CALL infog2l( m+2, m+2, desca, nprow, npcol, myrow, mycol,
 
  217     $                    irow, icol, rsrc, jsrc )
  218            IF( num.GT.1 ) THEN
  219               CALL sgerv2d( contxt, 1, 1, h11, 1, up, left )
  220            ELSE
  221               h11( 1 ) = a( ( icol-3 )*lda+irow-2 )
  222            END IF
  223            IF( nprow.GT.1 ) THEN
  224               CALL sgerv2d( contxt, 1, 1, h12, 1, up, mycol )
  225            ELSE
  226               h12( 1 ) = a( ( icol-2 )*lda+irow-2 )
  227            END IF
  228            IF( npcol.GT.1 ) THEN
  229               CALL sgerv2d( contxt, 1, 1, h21, 1, myrow, left )
  230            ELSE
  231               h21( 1 ) = a( ( icol-3 )*lda+irow-1 )
  232            END IF
  233            h22 = a( ( icol-2 )*lda+irow-1 )
  234            v3( 1 ) = a( ( icol-2 )*lda+irow )
  235         END IF
  236      END IF
  237      IF( ( myrow.NE.ii ) .OR. ( mycol.NE.jj ) )
  238     $   RETURN
  239
  240      IF( modkm1.GT.1 ) THEN
  241         CALL infog2l( m+2, m+2, desca, nprow, npcol, myrow, mycol,
 
  242     $                 irow, icol, rsrc, jsrc )
  243         h11( 1 ) = a( ( icol-3 )*lda+irow-2 )
  244         h21( 1 ) = a( ( icol-3 )*lda+irow-1 )
  245         h12( 1 ) = a( ( icol-2 )*lda+irow-2 )
  246         h22 = a( ( icol-2 )*lda+irow-1 )
  247         v3( 1 ) = a( ( icol-2 )*lda+irow )
  248      END IF
  249
  250      h44s = h44 - h11( 1 )
  251      h33s = h33 - h11( 1 )
  252      v1 = ( h33s*h44s-h43h34 ) / h21( 1 ) + h12( 1 )
  253      v2 = h22 - h11( 1 ) - h33s - h44s
  254      s = abs( v1 ) + abs( v2 ) + abs( v3( 1 ) )
  255      v1 = v1 / s
  256      v2 = v2 / s
  257      v3( 1 ) = v3( 1 ) / s
  258      v( 1 ) = v1
  259      v( 2 ) = v2
  260      v( 3 ) = v3( 1 )
  261
  262      RETURN
  263
  264
  265
subroutine infog2l(grindx, gcindx, desc, nprow, npcol, myrow, mycol, lrindx, lcindx, rsrc, csrc)