2436
2437
2438
2439
2440
2441
2442
2443
2444
2445
2446
2447
2448 COMPLEX ZERO
2449 parameter( zero = ( 0.0, 0.0 ) )
2450 REAL RZERO, RONE
2451 parameter( rzero = 0.0, rone = 1.0 )
2452
2453 COMPLEX ALPHA, BETA
2454 REAL EPS, ERR
2455 INTEGER KK, LDA, LDB, LDC, LDCC, M, N, NOUT
2456 LOGICAL FATAL, MV
2457 CHARACTER*1 TRANSA, TRANSB
2458
2459 COMPLEX A( LDA, * ), B( LDB, * ), C( LDC, * ),
2460 $ CC( LDCC, * ), CT( * )
2461 REAL G( * )
2462
2463 COMPLEX CL
2464 REAL ERRI
2465 INTEGER I, J, K
2466 LOGICAL CTRANA, CTRANB, TRANA, TRANB
2467
2468 INTRINSIC abs, aimag, conjg, max, real, sqrt
2469
2470 REAL ABS1
2471
2472 abs1( cl ) = abs( real( cl ) ) + abs( aimag( cl ) )
2473
2474 trana = transa.EQ.'T'.OR.transa.EQ.'C'
2475 tranb = transb.EQ.'T'.OR.transb.EQ.'C'
2476 ctrana = transa.EQ.'C'
2477 ctranb = transb.EQ.'C'
2478
2479
2480
2481
2482
2483 DO 220 j = 1, n
2484
2485 DO 10 i = 1, m
2486 ct( i ) = zero
2487 g( i ) = rzero
2488 10 CONTINUE
2489 IF( .NOT.trana.AND..NOT.tranb )THEN
2490 DO 30 k = 1, kk
2491 DO 20 i = 1, m
2492 ct( i ) = ct( i ) + a( i, k )*b( k, j )
2493 g( i ) = g( i ) + abs1( a( i, k ) )*abs1( b( k, j ) )
2494 20 CONTINUE
2495 30 CONTINUE
2496 ELSE IF( trana.AND..NOT.tranb )THEN
2497 IF( ctrana )THEN
2498 DO 50 k = 1, kk
2499 DO 40 i = 1, m
2500 ct( i ) = ct( i ) + conjg( a( k, i ) )*b( k, j )
2501 g( i ) = g( i ) + abs1( a( k, i ) )*
2502 $ abs1( b( k, j ) )
2503 40 CONTINUE
2504 50 CONTINUE
2505 ELSE
2506 DO 70 k = 1, kk
2507 DO 60 i = 1, m
2508 ct( i ) = ct( i ) + a( k, i )*b( k, j )
2509 g( i ) = g( i ) + abs1( a( k, i ) )*
2510 $ abs1( b( k, j ) )
2511 60 CONTINUE
2512 70 CONTINUE
2513 END IF
2514 ELSE IF( .NOT.trana.AND.tranb )THEN
2515 IF( ctranb )THEN
2516 DO 90 k = 1, kk
2517 DO 80 i = 1, m
2518 ct( i ) = ct( i ) + a( i, k )*conjg( b( j, k ) )
2519 g( i ) = g( i ) + abs1( a( i, k ) )*
2520 $ abs1( b( j, k ) )
2521 80 CONTINUE
2522 90 CONTINUE
2523 ELSE
2524 DO 110 k = 1, kk
2525 DO 100 i = 1, m
2526 ct( i ) = ct( i ) + a( i, k )*b( j, k )
2527 g( i ) = g( i ) + abs1( a( i, k ) )*
2528 $ abs1( b( j, k ) )
2529 100 CONTINUE
2530 110 CONTINUE
2531 END IF
2532 ELSE IF( trana.AND.tranb )THEN
2533 IF( ctrana )THEN
2534 IF( ctranb )THEN
2535 DO 130 k = 1, kk
2536 DO 120 i = 1, m
2537 ct( i ) = ct( i ) + conjg( a( k, i ) )*
2538 $ conjg( b( j, k ) )
2539 g( i ) = g( i ) + abs1( a( k, i ) )*
2540 $ abs1( b( j, k ) )
2541 120 CONTINUE
2542 130 CONTINUE
2543 ELSE
2544 DO 150 k = 1, kk
2545 DO 140 i = 1, m
2546 ct( i ) = ct( i ) + conjg( a( k, i ) )*b( j, k )
2547 g( i ) = g( i ) + abs1( a( k, i ) )*
2548 $ abs1( b( j, k ) )
2549 140 CONTINUE
2550 150 CONTINUE
2551 END IF
2552 ELSE
2553 IF( ctranb )THEN
2554 DO 170 k = 1, kk
2555 DO 160 i = 1, m
2556 ct( i ) = ct( i ) + a( k, i )*conjg( b( j, k ) )
2557 g( i ) = g( i ) + abs1( a( k, i ) )*
2558 $ abs1( b( j, k ) )
2559 160 CONTINUE
2560 170 CONTINUE
2561 ELSE
2562 DO 190 k = 1, kk
2563 DO 180 i = 1, m
2564 ct( i ) = ct( i ) + a( k, i )*b( j, k )
2565 g( i ) = g( i ) + abs1( a( k, i ) )*
2566 $ abs1( b( j, k ) )
2567 180 CONTINUE
2568 190 CONTINUE
2569 END IF
2570 END IF
2571 END IF
2572 DO 200 i = 1, m
2573 ct( i ) = alpha*ct( i ) + beta*c( i, j )
2574 g( i ) = abs1( alpha )*g( i ) +
2575 $ abs1( beta )*abs1( c( i, j ) )
2576 200 CONTINUE
2577
2578
2579
2580 err = zero
2581 DO 210 i = 1, m
2582 erri = abs1( ct( i ) - cc( i, j ) )/eps
2583 IF( g( i ).NE.rzero )
2584 $ erri = erri/g( i )
2585 err = max( err, erri )
2586 IF( err*sqrt( eps ).GE.rone )
2587 $ GO TO 230
2588 210 CONTINUE
2589
2590 220 CONTINUE
2591
2592
2593 GO TO 250
2594
2595
2596
2597 230 fatal = .true.
2598 WRITE( nout, fmt = 9999 )
2599 DO 240 i = 1, m
2600 IF( mv )THEN
2601 WRITE( nout, fmt = 9998 )i, ct( i ), cc( i, j )
2602 ELSE
2603 WRITE( nout, fmt = 9998 )i, cc( i, j ), ct( i )
2604 END IF
2605 240 CONTINUE
2606 IF( n.GT.1 )
2607 $ WRITE( nout, fmt = 9997 )j
2608
2609 250 CONTINUE
2610 RETURN
2611
2612 9999 FORMAT(' ******* FATAL ERROR - COMPUTED RESULT IS LESS THAN HAL',
2613 $ 'F ACCURATE *******', /' EXPECTED RE',
2614 $ 'SULT COMPUTED RESULT' )
2615 9998 FORMAT( 1x, i7, 2( ' (', g15.6, ',', g15.6, ')' ) )
2616 9997 FORMAT( ' THESE ARE THE RESULTS FOR COLUMN ', i3 )
2617
2618
2619