home *** CD-ROM | disk | FTP | other *** search
/ C!T ROM 2 / ctrom_ii_b.zip / ctrom_ii_b / PROGRAM / PASCAL / TPL60N19 / TESTPRGS / DATAN.PAS < prev    next >
Pascal/Delphi Source File  |  1993-02-14  |  9KB  |  312 lines
  1. PROGRAM DATan;     { ported from Fortran original 05-01-92 Norbert Juffa }
  2.  
  3. {$A+,B-,D-,E+,F-,G-,I-,L-,N-,O-,R-,S-,V-,X-}
  4.  
  5. USES MachArit, Power;
  6.  
  7. {
  8. C     PROGRAM TO TEST DATAN, DATAN2
  9. C
  10. C     DATA REQUIRED
  11. C
  12. C        NONE
  13. C
  14. C     SUBPROGRAMS REQUIRED FROM THIS PACKAGE
  15. C
  16. C        MACHAR - AN ENVIRONMENTAL INQUIRY PROGRAM PROVIDING
  17. C                 INFORMATION ON THE FLOATING-POINT ARITHMETIC
  18. C                 SYSTEM.  NOTE THAT THE CALL TO MACHAR CAN
  19. C                 BE DELETED PROVIDED THE FOLLOWING SIX
  20. C                 PARAMETERS ARE ASSIGNED THE VALUES INDICATED
  21. C
  22. C                 IBETA  - THE RADIX OF THE FLOATING-POINT SYSTEM
  23. C                 IT     - THE NUMBER OF BASE-IBETA DIGITS IN THE
  24. C                          SIGNIFICAND OF A FLOATING-POINT NUMBER
  25. C                 IRND   - 0 IF FLOATING-POINT ADDITION CHOPS,
  26. C                          1 IF FLOATING-POINT ADDITION ROUNDS
  27. C                 MINEXP - THE LARGEST IN MAGNITUDE NEGATIVE
  28. C                          INTEGER SUCH THAT DFLOAT(IBETA)**MINEXP
  29. C                          IS A POSITIVE FLOATING-POINT NUMBER
  30. C                 XMIN   - THE SMALLEST NON-VANISHING FLOATING-POINT
  31. C                          POWER OF THE RADIX
  32. C                 XMAX   - THE LARGEST FINITE FLOATING-POINT NO.
  33. C
  34. C        REN(K) - A FUNCTION SUBPROGRAM RETURNING RANDOM REAL
  35. C                 NUMBERS UNIFORMLY DISTRIBUTED OVER (0,1)
  36. C
  37. C     STANDARD FORTRAN SUBPROGRAMS REQUIRED
  38. C
  39. C         DABS, DLOG, DMAX1, DATAN, DATAN2, DFLOAT, DSQRT
  40. C
  41. C
  42. C     LATEST REVISION - DECEMBER 6, 1979
  43. C
  44. C     AUTHOR - W. J. CODY
  45. C              ARGONNE NATIONAL LABORATORY
  46. C
  47. C
  48. }
  49.  
  50.  
  51.  
  52. FUNCTION REN (K: LONGINT): REAL;
  53.  
  54. {
  55.       DOUBLE PRECISION FUNCTION REN(K)
  56. C
  57. C     RANDOM NUMBER GENERATOR - BASED ON ALGORITHM 266 BY PIKE AND
  58. C      HILL (MODIFIED BY HANSSON), COMMUNICATIONS OF THE ACM,
  59. C      VOL. 8, NO. 10, OCTOBER 1965.
  60. C
  61. C     THIS SUBPROGRAM IS INTENDED FOR USE ON COMPUTERS WITH
  62. C      FIXED POINT WORDLENGTH OF AT LEAST 29 BITS.  IT IS
  63. C      BEST IF THE FLOATING POINT SIGNIFICAND HAS AT MOST
  64. C      29 BITS.
  65. C
  66. }
  67.  
  68. VAR   J:  LONGINT;
  69. CONST IY: LONGINT = 100001;
  70.  
  71. BEGIN
  72.    J  := K;
  73.    IY := IY * 125;
  74.    IY := IY - (IY DIV 2796203) * 2796203;
  75.    REN:= 1.0 * (IY) / 2796203.0e0 * (1.0e0 + 1.0e-6 + 1.0e-12);
  76. END;
  77.  
  78.  
  79.  
  80. FUNCTION MAX1 (A, B:REAL): REAL;
  81. BEGIN
  82.    IF A > B THEN
  83.       MAX1 := A
  84.    ELSE
  85.       MAX1 := B;
  86. END;
  87.  
  88.  
  89.  
  90. VAR   I,IBETA,IEXP,IOUT,IRND,II,IT,I1,J,K1,
  91.       K2,K3,MACHEP,MAXEXP,MINEXP,N,NEGEP,NGRD: LONGINT;
  92.  
  93.       A,AIT,ALBETA,B,BETA,BETAP,DEL,EM,EPS,
  94.       EPSNEG,EXPON,HALF, OB32,ONE,R6,R7,SUM,
  95.       TWO,W,X,XL,XMAX,XMIN,XN,XSQ,X1,Y,T,Z,
  96.       THREE,FIVE,SIXTEEN,SEVENTEEN,ZERO,ZZ,
  97.       THREEFORTH,SIXTEENTH: REAL;
  98.  
  99. LABEL 100,105,110,120,200;
  100.  
  101. BEGIN
  102.  
  103.    N := 1000000;   { number of trials }
  104.  
  105.    MACHAR (IBETA,IT,IRND,NGRD,MACHEP,NEGEP,IEXP,MINEXP,MAXEXP,
  106.            EPS,EPSNEG,XMIN,XMAX);
  107.    PRINTPARAM (IBETA,IT,IRND,NGRD,MACHEP,NEGEP,IEXP,MINEXP,MAXEXP,
  108.                EPS,EPSNEG,XMIN,XMAX);
  109.    BETA      := IBETA;
  110.    ALBETA    := LN (BETA);
  111.    AIT       := IT;
  112.    ZERO      := 0;
  113.    ONE       := 1;
  114.    TWO       := 2;
  115.    THREE     := 3;
  116.    FIVE      := 5;
  117.    SIXTEEN   := 16;
  118.    SEVENTEEN := 17;
  119.    SIXTEENTH := 0.0625;
  120.    HALF      := 0.5;
  121.    THREEFORTH:= 0.75;
  122.    A         := -SIXTEENTH;
  123.    B         := -A;
  124.    OB32      := B * HALF;
  125.    XN        := N;
  126.    I1        := 0;
  127.  
  128. {----------------------------------------------------------------- }
  129. {     RANDOM ARGUMENT ACCURACY TESTS                               }
  130. {----------------------------------------------------------------- }
  131.  
  132.    FOR J := 1 TO 4 DO BEGIN
  133.       K1 := 0;
  134.       K3 := 0;
  135.       X1 := ZERO;
  136.       R6 := ZERO;
  137.       R7 := ZERO;
  138.       DEL:= (B - A) / XN;
  139.       XL := A;
  140.  
  141.       FOR I := 1 TO N DO BEGIN
  142.          X := DEL * REN(I1) + XL;
  143.          IF J = 2 THEN
  144.             X := ((ONE+X*A)-ONE)*SIXTEEN;
  145.          Z := ARCTAN (X);
  146.          IF J <> 1 THEN
  147.             GOTO 100;
  148.          XSQ := X * X;
  149.          EM  := SEVENTEEN;
  150.          SUM := XSQ / EM;
  151.          FOR II := 1 TO 7 DO BEGIN
  152.              EM := EM - TWO;
  153.              SUM:= (ONE/EM - SUM) * XSQ;
  154.          END;
  155.          SUM := -X * SUM;
  156.          ZZ  :=  X + SUM;
  157.          SUM := (X - ZZ) + SUM;
  158.          IF IRND = 0 THEN
  159.             ZZ := ZZ + (SUM + SUM);
  160.          GOTO 110;
  161. 100:     IF J <> 2 THEN
  162.             GOTO 105;
  163.          Y := (X - SIXTEENTH) / (ONE + X * A);
  164.          ZZ:= ARCTAN (Y) - 8.1190004042651526021e-5;  {  arctan (1/16) - 1/16 }
  165.          ZZ:= ZZ + OB32;
  166.          ZZ:= ZZ + OB32;
  167.          GOTO 110;
  168. 105:     Z := Z + Z;
  169.          Y := X / ((HALF + X * HALF)*((HALF - X) + HALF));
  170.          ZZ:= ARCTAN (Y);
  171. 110:     IF Z <> ZERO THEN
  172.             W := (Z - ZZ) / Z
  173.          ELSE IF ZZ <> ZERO THEN
  174.             W := ONE;
  175.          IF W > ZERO THEN
  176.             K1 := K1 + 1;
  177.          IF W < ZERO THEN
  178.             K3 := K3 + 1;
  179.          W := ABS(W);
  180.          IF W <= R6 THEN
  181.             GOTO 120;
  182.          R6 := W;
  183.          X1 := X;
  184. 120:     R7 := R7 + W * W;
  185.          XL := XL + DEL;
  186. 200:  END;
  187.  
  188.       K2 := N - K3 - K1;
  189.       R7 := SQRT(R7/XN);
  190.  
  191.       IF J = 1 THEN BEGIN
  192.          WRITELN;
  193.          WRITELN ;
  194.          WRITELN ('TEST OF ARCTAN(X) VS TRUNCATED TAYLOR SERIES');
  195.          WRITELN;
  196.          END;
  197.       IF J = 2 THEN BEGIN
  198.          WRITELN;
  199.          WRITELN;
  200.          WRITELN ('TEST OF ARCTAN(X) VS ARCTAN(1/16) + ARCTAN((X-1/16)/(1+X/16))');
  201.          WRITELN;
  202.          END;
  203.       IF J > 2 THEN BEGIN
  204.          WRITELN;
  205.          WRITELN;
  206.          WRITELN ('TEST OF 2*ARCTAN(X) VS ARCTAN(2X((1-X*X))');
  207.          WRITELN;
  208.          END;
  209.  
  210.       WRITELN (N, ' RANDOM ARGUMENTS WERE TESTED FROM THE INTERVAL');
  211.       WRITELN ('(', A, ',', B, ')');
  212.       WRITELN;
  213.  
  214.       WRITELN ('ARCTAN (X) WAS LARGER', K1:6, ' TIMES');
  215.       WRITELN ('               AGREED', K2:6, ' TIMES');
  216.       WRITELN ('      AND WAS SMALLER', K3:6, ' TIMES');
  217.  
  218.       WRITELN;
  219.       WRITELN ('THERE ARE ', IT, ' BASE ', IBETA,
  220.                ' SIGNIFICANT DIGITS IN A FLOATING-POINT NUMBER');
  221.       WRITELN;
  222.  
  223.       W := -999;
  224.       IF R6 <> ZERO THEN
  225.          W := LN (ABS(R6))/ALBETA;
  226.       WRITELN ('THE MAXIMUM RELATIVE ERROR OF          ', R6:12,
  227.                ' = ', IBETA, ' **', W:7:2);
  228.       WRITELN ('OCCURED FOR X = ', X1);
  229.       W := MAX1 (AIT+W,ZERO);
  230.       WRITELN;
  231.       WRITELN ('THE ESTIMATED LOSS OF BASE ', IBETA,
  232.                ' SIGNIFICANT DIGITS IS        ', W:7:2);
  233.       W := -999.0;
  234.       IF R7 <> ZERO THEN
  235.          W := LN (ABS(R7))/ALBETA;
  236.       WRITELN;
  237.       WRITELN ('THE ROOT MEAN SQUARE RELATIVE ERROR WAS', R7:12,
  238.                ' = ', IBETA, ' **' , W:7:2);
  239.       W := MAX1 (AIT+W,ZERO);
  240.       WRITELN;
  241.       WRITELN ('THE ESTIMATED LOSS OF BASE ', IBETA,
  242.                ' SIGNIFICANT DIGITS IS        ', W:7:2);
  243.       A := B;
  244.       IF J = 1 THEN
  245.          B := TWO - SQRT (THREE);
  246.       IF J = 2 THEN
  247.          B := SQRT (TWO) - ONE;
  248.       IF J = 3 THEN
  249.          B := ONE;
  250.    END;
  251.  
  252. {-----------------------------------------------------------------}
  253. {     SPECIAL TESTS                                               }
  254. {-----------------------------------------------------------------}
  255.  
  256.    WRITELN;
  257.    WRITELN;
  258.    WRITELN ('SPECIAL TESTS');
  259.    WRITELN;
  260.    WRITELN ('THE IDENTITY  ARCTAN(-X) = -ARCTAN(X)  WILL BE TESTED');
  261.    WRITELN;
  262.    WRITELN ('          X           F(X) + F(-X)');
  263.    WRITELN;
  264.  
  265.    A := FIVE;
  266.  
  267.    FOR I := 1 TO 5 DO BEGIN
  268.       X := REN(I1) * A;
  269.       Z := ARCTAN(X) + ARCTAN(-X);
  270.       WRITELN (X:18, Z:18);
  271.    END;
  272.  
  273.    WRITELN;
  274.    WRITELN;
  275.    WRITELN ('THE IDENTITY ARCTAN(X) = X , X SMALL, WILL BE TESTED.');
  276.    WRITELN;
  277.    WRITELN ('          X                 X-F(X)');
  278.    WRITELN;
  279.    BETAP := POW (BETA,IT);
  280.    X := REN(I1) / BETAP;
  281.  
  282.    FOR I := 1 TO 5 DO BEGIN
  283.       Z := X - ARCTAN (X);
  284.       WRITELN (X:18, Z:18);
  285.       X := X / BETA;
  286.    END;
  287.  
  288.    WRITELN;
  289.    WRITELN;
  290.    WRITELN ('TEST OF UNDERFLOW FOR VERY SMALL ARGUMENT');
  291.    WRITELN;
  292.    EXPON := MINEXP * THREEFORTH;
  293.    X := POW (BETA,EXPON);
  294.    Y := ARCTAN(X);
  295.    WRITELN ('ARCTAN (', X:18, ') = ', Y:18);
  296.  
  297. {-----------------------------------------------------------------}
  298. {     TEST OF ERROR RETURNS                                       }
  299. {-----------------------------------------------------------------}
  300.  
  301.    WRITELN;
  302.    WRITELN;
  303.    WRITELN ('TEST OF ERROR RETURNS');
  304.    WRITELN;
  305.    WRITELN ('ARCTAN WILL BE CALLED WITH THE ARGUMENT ', XMAX:18);
  306.    WRITELN ('THIS SHOULD NOT TRIGGER AN ERROR MESSAGE');
  307.    Z := ARCTAN(XMAX);
  308.    WRITELN ('ARCTAN (', XMAX:18, ') = ', Z:18);
  309.    WRITELN;
  310.    WRITELN ('THIS CONCLUDES THE TESTS');
  311. END. { DATan }
  312.