4919
 4920
 4921
 4922
 4923
 4924
 4925
 4926      CHARACTER*1        UPLO
 4927      INTEGER            IA, ICTXT, INCX, INCY, INFO, IX, IY, JA, JX,
 4928     $                   JY, M, N
 4929      REAL               ALPHA, ERR
 4930
 4931
 4932      INTEGER            DESCA( * ), DESCX( * ), DESCY( * )
 4933      REAL               A( * ), G( * ), PA( * ), X( * ), Y( * )
 4934
 4935
 4936
 4937
 4938
 4939
 4940
 4941
 4942
 4943
 4944
 4945
 4946
 4947
 4948
 4949
 4950
 4951
 4952
 4953
 4954
 4955
 4956
 4957
 4958
 4959
 4960
 4961
 4962
 4963
 4964
 4965
 4966
 4967
 4968
 4969
 4970
 4971
 4972
 4973
 4974
 4975
 4976
 4977
 4978
 4979
 4980
 4981
 4982
 4983
 4984
 4985
 4986
 4987
 4988
 4989
 4990
 4991
 4992
 4993
 4994
 4995
 4996
 4997
 4998
 4999
 5000
 5001
 5002
 5003
 5004
 5005
 5006
 5007
 5008
 5009
 5010
 5011
 5012
 5013
 5014
 5015
 5016
 5017
 5018
 5019
 5020
 5021
 5022
 5023
 5024
 5025
 5026
 5027
 5028
 5029
 5030
 5031
 5032
 5033
 5034
 5035
 5036
 5037
 5038
 5039
 5040
 5041
 5042
 5043
 5044
 5045
 5046
 5047
 5048
 5049
 5050
 5051
 5052
 5053
 5054
 5055
 5056
 5057
 5058
 5059
 5060
 5061
 5062
 5063
 5064
 5065
 5066
 5067
 5068
 5069
 5070
 5071
 5072
 5073
 5074
 5075
 5076
 5077
 5078
 5079
 5080
 5081
 5082
 5083
 5084
 5085
 5086
 5087
 5088
 5089
 5090
 5091
 5092
 5093
 5094
 5095
 5096
 5097
 5098
 5099
 5100
 5101
 5102
 5103
 5104
 5105
 5106
 5107
 5108      INTEGER            BLOCK_CYCLIC_2D_INB, CSRC_, CTXT_, DLEN_,
 5109     $                   DTYPE_, IMB_, INB_, LLD_, MB_, M_, NB_, N_,
 5110     $                   RSRC_
 5111      parameter( block_cyclic_2d_inb = 2, dlen_ = 11,
 5112     $                   dtype_ = 1, ctxt_ = 2, m_ = 3, n_ = 4,
 5113     $                   imb_ = 5, inb_ = 6, mb_ = 7, nb_ = 8,
 5114     $                   rsrc_ = 9, csrc_ = 10, lld_ = 11 )
 5115      REAL               ZERO, ONE
 5116      parameter( zero = 0.0e+0, one = 1.0e+0 )
 5117
 5118
 5119      LOGICAL            COLREP, LOWER, ROWREP, UPPER
 5120      INTEGER            I, IACOL, IAROW, IB, IBEG, ICURROW, IEND, IIA,
 5121     $                   IN, IOFFA, IOFFXI, IOFFXJ, IOFFYI, IOFFYJ, J,
 5122     $                   JJA, KK, LDA, LDPA, LDX, LDY, MYCOL, MYROW,
 5123     $                   NPCOL, NPROW
 5124      REAL               EPS, ERRI, GTMP, ATMP
 5125
 5126
 5127      EXTERNAL           blacs_gridinfo, igsum2d, 
pb_infog2l, sgamx2d
 
 5128
 5129
 5130      LOGICAL            LSAME
 5131      REAL               PSLAMCH
 5133
 5134
 5135      INTRINSIC          abs, 
max, 
min, mod, sqrt
 
 5136
 5137
 5138
 5139      CALL blacs_gridinfo( ictxt, nprow, npcol, myrow, mycol )
 5140
 5142
 5143      upper = 
lsame( uplo, 
'U' )
 
 5144      lower = 
lsame( uplo, 
'L' )
 
 5145
 5146      lda = 
max( 1, desca( m_ ) )
 
 5147      ldx = 
max( 1, descx( m_ ) )
 
 5148      ldy = 
max( 1, descy( m_ ) )
 
 5149
 5150
 5151
 5152
 5153
 5154      DO 70 j = 1, n
 5155
 5156         ioffxj = ix + ( jx - 1 ) * ldx + ( j - 1 ) * incx
 5157         ioffyj = iy + ( jy - 1 ) * ldy + ( j - 1 ) * incy
 5158
 5159         IF( lower ) THEN
 5160            ibeg = j
 5161            iend = m
 5162            DO 10 i = 1, j-1
 5163               g( i ) = zero
 5164   10       CONTINUE
 5165         ELSE IF( upper ) THEN
 5166            ibeg = 1
 5167            iend = j
 5168            DO 20 i = j+1, m
 5169               g( i ) = zero
 5170   20       CONTINUE
 5171         ELSE
 5172            ibeg = 1
 5173            iend = m
 5174         END IF
 5175
 5176         DO 30 i = ibeg, iend
 5177            ioffa = ia + i - 1 + ( ja + j - 2 ) * lda
 5178            ioffxi = ix + ( jx - 1 ) * ldx + ( i - 1 ) * incx
 5179            ioffyi = iy + ( jy - 1 ) * ldy + ( i - 1 ) * incy
 5180            atmp = x( ioffxi ) * y( ioffyj )
 5181            atmp = atmp + y( ioffyi ) * x( ioffxj )
 5182            gtmp = abs( x( ioffxi ) * y( ioffyj ) )
 5183            gtmp = gtmp + abs( y( ioffyi ) * x( ioffxj ) )
 5184            g( i ) = abs( alpha ) * gtmp + abs( a( ioffa ) )
 5185            a( ioffa ) = alpha*atmp + a( ioffa )
 5186
 5187   30    CONTINUE
 5188
 5189
 5190
 5191         info = 0
 5192         err  = zero
 5193         ldpa = desca( lld_ )
 5194         ioffa = ia + ( ja + j - 2 ) * lda
 5195         CALL pb_infog2l( ia, ja+j-1, desca, nprow, npcol, myrow, mycol,
 
 5196     $                    iia, jja, iarow, iacol )
 5197         rowrep = ( iarow.EQ.-1 )
 5198         colrep = ( iacol.EQ.-1 )
 5199
 5200         IF( mycol.EQ.iacol .OR. colrep ) THEN
 5201
 5202            icurrow = iarow
 5203            ib = desca( imb_ ) - ia + 1
 5204            IF( ib.LE.0 )
 5205     $         ib = ( ( -ib ) / desca( mb_ ) + 1 ) * desca( mb_ ) + ib
 5207            in = ia + ib - 1
 5208
 5209            DO 40 i = ia, in
 5210
 5211               IF( myrow.EQ.icurrow .OR. rowrep ) THEN
 5212                  erri = abs( pa( iia+(jja-1)*ldpa ) - a( ioffa ) )/eps
 5213                  IF( g( i-ia+1 ).NE.zero )
 5214     $               erri = erri / g( i-ia+1 )
 5215                  err = 
max( err, erri )
 
 5216                  IF( err*sqrt( eps ).GE.one )
 5217     $               info = 1
 5218                  iia = iia + 1
 5219               END IF
 5220
 5221               ioffa = ioffa + 1
 5222
 5223   40       CONTINUE
 5224
 5225            icurrow = mod( icurrow+1, nprow )
 5226
 5227            DO 60 i = in+1, ia+m-1, desca( mb_ )
 5228               ib = 
min( ia+m-i, desca( mb_ ) )
 
 5229
 5230               DO 50 kk = 0, ib-1
 5231
 5232                  IF( myrow.EQ.icurrow .OR. rowrep ) THEN
 5233                     erri = abs( pa( iia+(jja-1)*ldpa )-a( ioffa ) )/eps
 5234                     IF( g( i+kk-ia+1 ).NE.zero )
 5235     $                  erri = erri / g( i+kk-ia+1 )
 5236                     err = 
max( err, erri )
 
 5237                     IF( err*sqrt( eps ).GE.one )
 5238     $                  info = 1
 5239                     iia = iia + 1
 5240                  END IF
 5241
 5242                  ioffa = ioffa + 1
 5243
 5244   50          CONTINUE
 5245
 5246               icurrow = mod( icurrow+1, nprow )
 5247
 5248   60       CONTINUE
 5249
 5250         END IF
 5251
 5252
 5253
 5254         CALL igsum2d( ictxt, 'All', ' ', 1, 1, info, 1, -1, mycol )
 5255         CALL sgamx2d( ictxt, 'All', ' ', 1, 1, err, 1, i, j, -1, -1,
 5256     $                 mycol )
 5257         IF( info.NE.0 )
 5258     $      GO TO 80
 5259
 5260   70 CONTINUE
 5261
 5262   80 CONTINUE
 5263
 5264      RETURN
 5265
 5266
 5267
subroutine pb_infog2l(i, j, desc, nprow, npcol, myrow, mycol, ii, jj, prow, pcol)
 
real function pslamch(ictxt, cmach)