3266
3267
3268
3269
3270
3271
3272
3273
3274
3275
3276
3277
3278 COMPLEX*16 ZERO
3279 parameter( zero = ( 0.0d0, 0.0d0 ) )
3280 DOUBLE PRECISION RZERO, RONE
3281 parameter( rzero = 0.0d0, rone = 1.0d0 )
3282
3283 COMPLEX*16 ALPHA, BETA
3284 DOUBLE PRECISION EPS, ERR
3285 INTEGER KK, LDA, LDB, LDC, LDCC, M, N, NOUT
3286 LOGICAL FATAL, MV
3287 CHARACTER*1 TRANSA, TRANSB
3288
3289 COMPLEX*16 A( LDA, * ), B( LDB, * ), C( LDC, * ),
3290 $ CC( LDCC, * ), CT( * )
3291 DOUBLE PRECISION G( * )
3292
3293 COMPLEX*16 CL
3294 DOUBLE PRECISION ERRI
3295 INTEGER I, J, K
3296 LOGICAL CTRANA, CTRANB, TRANA, TRANB
3297
3298 INTRINSIC abs, dimag, dconjg, max, dble, sqrt
3299
3300 DOUBLE PRECISION ABS1
3301
3302 abs1( cl ) = abs( dble( cl ) ) + abs( dimag( cl ) )
3303
3304 trana = transa.EQ.'T'.OR.transa.EQ.'C'
3305 tranb = transb.EQ.'T'.OR.transb.EQ.'C'
3306 ctrana = transa.EQ.'C'
3307 ctranb = transb.EQ.'C'
3308
3309
3310
3311
3312
3313 DO 220 j = 1, n
3314
3315 DO 10 i = 1, m
3316 ct( i ) = zero
3317 g( i ) = rzero
3318 10 CONTINUE
3319 IF( .NOT.trana.AND..NOT.tranb )THEN
3320 DO 30 k = 1, kk
3321 DO 20 i = 1, m
3322 ct( i ) = ct( i ) + a( i, k )*b( k, j )
3323 g( i ) = g( i ) + abs1( a( i, k ) )*abs1( b( k, j ) )
3324 20 CONTINUE
3325 30 CONTINUE
3326 ELSE IF( trana.AND..NOT.tranb )THEN
3327 IF( ctrana )THEN
3328 DO 50 k = 1, kk
3329 DO 40 i = 1, m
3330 ct( i ) = ct( i ) + dconjg( a( k, i ) )*b( k, j )
3331 g( i ) = g( i ) + abs1( a( k, i ) )*
3332 $ abs1( b( k, j ) )
3333 40 CONTINUE
3334 50 CONTINUE
3335 ELSE
3336 DO 70 k = 1, kk
3337 DO 60 i = 1, m
3338 ct( i ) = ct( i ) + a( k, i )*b( k, j )
3339 g( i ) = g( i ) + abs1( a( k, i ) )*
3340 $ abs1( b( k, j ) )
3341 60 CONTINUE
3342 70 CONTINUE
3343 END IF
3344 ELSE IF( .NOT.trana.AND.tranb )THEN
3345 IF( ctranb )THEN
3346 DO 90 k = 1, kk
3347 DO 80 i = 1, m
3348 ct( i ) = ct( i ) + a( i, k )*dconjg( b( j, k ) )
3349 g( i ) = g( i ) + abs1( a( i, k ) )*
3350 $ abs1( b( j, k ) )
3351 80 CONTINUE
3352 90 CONTINUE
3353 ELSE
3354 DO 110 k = 1, kk
3355 DO 100 i = 1, m
3356 ct( i ) = ct( i ) + a( i, k )*b( j, k )
3357 g( i ) = g( i ) + abs1( a( i, k ) )*
3358 $ abs1( b( j, k ) )
3359 100 CONTINUE
3360 110 CONTINUE
3361 END IF
3362 ELSE IF( trana.AND.tranb )THEN
3363 IF( ctrana )THEN
3364 IF( ctranb )THEN
3365 DO 130 k = 1, kk
3366 DO 120 i = 1, m
3367 ct( i ) = ct( i ) + dconjg( a( k, i ) )*
3368 $ dconjg( b( j, k ) )
3369 g( i ) = g( i ) + abs1( a( k, i ) )*
3370 $ abs1( b( j, k ) )
3371 120 CONTINUE
3372 130 CONTINUE
3373 ELSE
3374 DO 150 k = 1, kk
3375 DO 140 i = 1, m
3376 ct( i ) = ct( i ) + dconjg( a( k, i ) )*
3377 $ b( j, k )
3378 g( i ) = g( i ) + abs1( a( k, i ) )*
3379 $ abs1( b( j, k ) )
3380 140 CONTINUE
3381 150 CONTINUE
3382 END IF
3383 ELSE
3384 IF( ctranb )THEN
3385 DO 170 k = 1, kk
3386 DO 160 i = 1, m
3387 ct( i ) = ct( i ) + a( k, i )*
3388 $ dconjg( b( j, k ) )
3389 g( i ) = g( i ) + abs1( a( k, i ) )*
3390 $ abs1( b( j, k ) )
3391 160 CONTINUE
3392 170 CONTINUE
3393 ELSE
3394 DO 190 k = 1, kk
3395 DO 180 i = 1, m
3396 ct( i ) = ct( i ) + a( k, i )*b( j, k )
3397 g( i ) = g( i ) + abs1( a( k, i ) )*
3398 $ abs1( b( j, k ) )
3399 180 CONTINUE
3400 190 CONTINUE
3401 END IF
3402 END IF
3403 END IF
3404 DO 200 i = 1, m
3405 ct( i ) = alpha*ct( i ) + beta*c( i, j )
3406 g( i ) = abs1( alpha )*g( i ) +
3407 $ abs1( beta )*abs1( c( i, j ) )
3408 200 CONTINUE
3409
3410
3411
3412 err = zero
3413 DO 210 i = 1, m
3414 erri = abs1( ct( i ) - cc( i, j ) )/eps
3415 IF( g( i ).NE.rzero )
3416 $ erri = erri/g( i )
3417 err = max( err, erri )
3418 IF( err*sqrt( eps ).GE.rone )
3419 $ GO TO 230
3420 210 CONTINUE
3421
3422 220 CONTINUE
3423
3424
3425 GO TO 250
3426
3427
3428
3429 230 fatal = .true.
3430 WRITE( nout, fmt = 9999 )
3431 DO 240 i = 1, m
3432 IF( mv )THEN
3433 WRITE( nout, fmt = 9998 )i, ct( i ), cc( i, j )
3434 ELSE
3435 WRITE( nout, fmt = 9998 )i, cc( i, j ), ct( i )
3436 END IF
3437 240 CONTINUE
3438 IF( n.GT.1 )
3439 $ WRITE( nout, fmt = 9997 )j
3440
3441 250 CONTINUE
3442 RETURN
3443
3444 9999 FORMAT( ' ******* FATAL ERROR - COMPUTED RESULT IS LESS THAN HAL',
3445 $ 'F ACCURATE *******', /' EXPECTED RE',
3446 $ 'SULT COMPUTED RESULT' )
3447 9998 FORMAT( 1x, i7, 2( ' (', g15.6, ',', g15.6, ')' ) )
3448 9997 FORMAT( ' THESE ARE THE RESULTS FOR COLUMN ', i3 )
3449
3450
3451