home *** CD-ROM | disk | FTP | other *** search
/ Programmer 7500 / MAX_PROGRAMMERS.iso / PASCAL / BPL70N12.ZIP / TESTPRGS.ZIP / DSIN.PAS < prev    next >
Encoding:
Pascal/Delphi Source File  |  1993-02-14  |  9.3 KB  |  329 lines
  1. PROGRAM DSin;   { 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 DSIN/DCOS
  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 FIVE
  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                 MINEXP - THE LARGEST IN MAGNITUDE NEGATIVE
  26. C                          INTEGER SUCH THAT  DFLOAT(IBETA)**MINEXP
  27. C                          IS A POSITIVE FLOATING-POINT NUMBER
  28. C                 EPS    - THE SMALLEST POSITIVE FLOATING-POINT
  29. C                          NUMBER SUCH THAT 1.0+EPS .NE. 1.0
  30. C                 EPSNEG - THE SMALLEST POSITIVE FLOATING-POINT
  31. C                          NUMBER SUCH THAT 1.0-EPSNEG .NE. 1.0
  32. C
  33. C        REN(K) - A FUNCTION SUBPROGRAM RETURNING RANDOM REAL
  34. C                 NUMBERS UNIFORMLY DISTRIBUTED OVER (0,1)
  35. C
  36. C
  37. C     STANDARD FORTRAN SUBPROGRAMS REQUIRED
  38. C
  39. C         DABS, DLOG, DMAX1, DCOS, DLOAT, DSIN, 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. FUNCTION REN (VAR K: LONGINT): REAL;
  52.  
  53. {
  54.       DOUBLE PRECISION FUNCTION REN(K)
  55. C
  56. C     RANDOM NUMBER GENERATOR - BASED ON ALGORITHM 266 BY PIKE AND
  57. C      HILL (MODIFIED BY HANSSON), COMMUNICATIONS OF THE ACM,
  58. C      VOL. 8, NO. 10, OCTOBER 1965.
  59. C
  60. C     THIS SUBPROGRAM IS INTENDED FOR USE ON COMPUTERS WITH
  61. C      FIXED POINT WORDLENGTH OF AT LEAST 29 BITS.  IT IS
  62. C      BEST IF THE FLOATING POINT SIGNIFICAND HAS AT MOST
  63. C      29 BITS.
  64. C
  65. }
  66.  
  67. VAR   J:  LONGINT;
  68. CONST IY: LONGINT = 100001;
  69.  
  70. BEGIN
  71.    J  := K;
  72.    IY := IY * 125;
  73.    IY := IY - (IY DIV 2796203) * 2796203;
  74.    REN:= 1.0 * (IY) / 2796203.0e0 * (1.0e0 + 1.0e-6 + 1.0e-12);
  75. END;
  76.  
  77.  
  78.  
  79. FUNCTION MAX1 (A, B:REAL): REAL;
  80. BEGIN
  81.    IF A > B THEN
  82.       MAX1 := A
  83.    ELSE
  84.       MAX1 := B;
  85. END;
  86.  
  87.  
  88.  
  89. VAR   I,IBETA,IEXP,IOUT,IRND,IT,I1,J, K1,K2,
  90.       K3,MACHEP,MAXEXP,MINEXP,N,NEGEP,NGRD: LONGINT;
  91.  
  92.       A,AIT,ALBETA,B,BETA,C,DEL,EPS,EPSNEG,HALF,ONE,
  93.       T,THREE,BETAP,TEMP,EXPON,R6,R7,TENTH,W,X,XL,PI,
  94.       XMAX,XMIN,XN,X1,Y,Z,ZERO,ZZ,THREEFOURTH,FOUR,SIX: REAL;
  95.  
  96. LABEL 100, 110, 120, 150, 160, 210, 220, 230, 240, 300;
  97.  
  98.  
  99. BEGIN
  100.  
  101.    N := 1000000;     { number of trials }
  102.  
  103.    MACHAR (IBETA,IT,IRND,NGRD,MACHEP,NEGEP,IEXP,MINEXP,MAXEXP,
  104.            EPS,EPSNEG,XMIN,XMAX);
  105.    PRINTPARAM (IBETA,IT,IRND,NGRD,MACHEP,NEGEP,IEXP,MINEXP,MAXEXP,
  106.                EPS,EPSNEG,XMIN,XMAX);
  107.    BETA       := IBETA;
  108.    ALBETA     := LN (BETA);
  109.    AIT        := IT;
  110.    ZERO       := 0;
  111.    ONE        := 1;
  112.    THREE      := 3;
  113.    FOUR       := 4;
  114.    SIX        := 6;
  115.    HALF       := 0.5;
  116.    THREEFOURTH:= 0.75;
  117.    PI         := 3.14159265358979323846264338327950288;
  118.    A          := ZERO;
  119.    B          := HALF*PI;
  120.    C          := B;
  121.    XN         := N;
  122.    I1         := 0;
  123.  
  124. {-----------------------------------------------------------------}
  125. {     RANDOM ARGUMENT ACCURACY TESTS                              }
  126. {-----------------------------------------------------------------}
  127.  
  128.    FOR J := 1 TO 3 DO BEGIN
  129.       K1 := 0;
  130.       K3 := 0;
  131.       X1 := ZERO;
  132.       R6 := ZERO;
  133.       R7 := ZERO;
  134.       DEL:= (B - A) / XN;
  135.       XL := A;
  136.  
  137.       FOR I := 1 TO N DO BEGIN
  138.          X := DEL * REN (I1) + XL;
  139.          Y := X / THREE;
  140.          Y := (X + Y);
  141.          Y := Y - X;
  142.          X := THREE * Y;
  143.          IF J = 3 THEN
  144.             GOTO 100;
  145.          Z := SIN (X);
  146.          ZZ:= SIN (Y);
  147.          W := ONE;
  148.          IF Z <> ZERO THEN BEGIN
  149.             TEMP := FOUR * ZZ * ZZ;
  150.             TEMP := THREE - TEMP;
  151.             TEMP := ZZ * TEMP;
  152.             TEMP := Z - TEMP;
  153.             W := TEMP / Z;
  154.             END;
  155.          GOTO 110;
  156. 100:     Z  := COS (X);
  157.          ZZ := COS (Y);
  158.          W  := ONE;
  159.          IF Z <> ZERO THEN BEGIN
  160.             TEMP := FOUR * ZZ * ZZ;
  161.             TEMP := THREE - TEMP;
  162.             TEMP := ZZ * TEMP;
  163.             TEMP := Z + TEMP;
  164.             W := TEMP / Z;
  165.             END;
  166. 110:     IF W > ZERO THEN
  167.             K1 := K1 + 1;
  168.          IF W < ZERO THEN
  169.             K3 := K3 + 1;
  170.          W := ABS(W);
  171.          IF W <= R6 THEN
  172.             GOTO 120;
  173.          R6 := W;
  174.          X1 := X;
  175. 120:     R7 := R7 + W * W;
  176.          XL := XL + DEL;
  177.       END;
  178.  
  179.       K2 := N - K3 - K1;
  180.       R7 := SQRT (R7/XN);
  181.  
  182.       IF J = 3 THEN
  183.          GOTO 210;
  184.       WRITELN;
  185.       WRITELN;
  186.       WRITELN ('TEST OF SIN(X) VS 3*SIN(X/3)-4*SIN(X/3)**3');
  187.       WRITELN;
  188.       WRITELN (N, ' RANDOM ARGUMENTS WERE TESTED FROM THE INTERVAL');
  189.       WRITELN ('(', A, ',', B, ')');
  190.       WRITELN;
  191.       inline ($fa/$fb);
  192.       WRITELN ('SIN (X) WAS LARGER', K1:6, ' TIMES');
  193.       WRITELN ('            AGREED', K2:6, ' TIMES');
  194.       WRITELN ('   AND WAS SMALLER', K3:6, ' TIMES');
  195.       GOTO 220;
  196. 210:  WRITELN;
  197.       WRITELN;
  198.       WRITELN ('TEST OF COS(X) VS 4*COS(X/3)**3-3*COS(X/3)');
  199.       WRITELN;
  200.       WRITELN (N, ' RANDOM ARGUMENTS WERE TESTED FROM THE INTERVAL');
  201.       WRITELN ('(', A, ',', B, ')');
  202.       WRITELN;
  203.       WRITELN ('COS (X) WAS LARGER', K1:6, ' TIMES');
  204.       WRITELN ('            AGREED', K2:6, ' TIMES');
  205.       WRITELN ('   AND WAS SMALLER', K3:6, ' TIMES');
  206. 220:  WRITELN;
  207.       WRITELN ('THERE ARE ', IT, ' BASE ', IBETA,
  208.                ' SIGNIFICANT DIGITS IN A FLOATING-POINT NUMBER');
  209.       WRITELN;
  210.       W := -999.0;
  211.       IF R6 <> ZERO THEN
  212.          W := LN (ABS(R6))/ALBETA;
  213.       WRITELN ('THE MAXIMUM RELATIVE ERROR OF          ', R6:12,
  214.                ' = ', IBETA, ' **', W:7:2);
  215.       WRITELN ('OCCURED FOR X = ', X1);
  216.       W := MAX1 (AIT+W, ZERO);
  217.       WRITELN;
  218.       WRITELN ('THE ESTIMATED LOSS OF BASE ', IBETA,
  219.                ' SIGNIFICANT DIGITS IS        ', W:7:2);
  220.       W := -999.0;
  221.       IF R7 <> ZERO THEN
  222.          W := LN (ABS(R7))/ALBETA;
  223.       WRITELN;
  224.       WRITELN ('THE ROOT MEAN SQUARE RELATIVE ERROR WAS', R7:12,
  225.                ' = ', IBETA, ' **' , W:7:2);
  226.       W := MAX1 (AIT+W,ZERO);
  227.       WRITELN;
  228.       WRITELN ('THE ESTIMATED LOSS OF BASE ', IBETA,
  229.                ' SIGNIFICANT DIGITS IS        ', W:7:2);
  230.       A := SIX * PI;
  231.       IF J = 2 THEN
  232.          A := B + C;
  233.       B := A + C;
  234.    END;
  235.  
  236. {-----------------------------------------------------------------}
  237. {     SPECIAL TESTS                                               }
  238. {-----------------------------------------------------------------}
  239.  
  240.    WRITELN;
  241.    WRITELN;
  242.    WRITELN ('SPECIAL TESTS');
  243.    WRITELN;
  244.    C := ONE / POW (BETA, (IT/2));
  245.    Z := SIN (A+C);
  246.    Z := Z - SIN (A-C);
  247.    TEMP := C + C;
  248.    Z := Z / TEMP;
  249.    WRITELN ('IF ',Z:18,' IS NOT ALMOST 1.0, SIN HAS THE WRONG PERIOD.');
  250.  
  251.    WRITELN;
  252.    WRITELN;
  253.    WRITELN ('THE IDENTITY   SIN(-X) = -SIN(X)   WILL BE TESTED.');
  254.    WRITELN;
  255.    WRITELN ('        F(X)           F(X)+F(-X)');
  256.    WRITELN;
  257.  
  258.    FOR I := 1 TO 5 DO BEGIN
  259.       X := REN(I1) * A;
  260.       Z := SIN(X) + SIN(-X);
  261.       WRITELN (X:18, Z:18);
  262.    END;
  263.  
  264.    WRITELN;
  265.    WRITELN;
  266.    WRITELN ('THE IDENTITY SIN(X) = X , X SMALL, WILL BE TESTED.');
  267.    WRITELN;
  268.    WRITELN ('         X                 X-F(X)');
  269.    WRITELN;
  270.    BETAP := POW (BETA, IT);
  271.    X := REN(I1) / BETAP;
  272.  
  273.    FOR I := 1 TO 5 DO BEGIN
  274.       Z := X - SIN(X);
  275.       WRITELN (X:18, Z:18);
  276.       X := X / BETA;
  277.    END;
  278.  
  279.    WRITELN;
  280.    WRITELN;
  281.    WRITELN ('THE IDENTITY   COS(-X) = COS(X)   WILL BE TESTED.');
  282.    WRITELN;
  283.    WRITELN ('         X              F(X)-F(-X)');
  284.    WRITELN;
  285.    FOR I := 1 TO 5 DO BEGIN
  286.       X := REN (I1) * A;
  287.       Z := COS (X) - COS (-X);
  288.       WRITELN (X:18, Z:18);
  289.    END;
  290.  
  291.    WRITELN;
  292.    WRITELN;
  293.    WRITELN ('TEST OF UNDERFLOW FOR VERY SMALL ARGUMENT.');
  294.    WRITELN;
  295.    EXPON := MINEXP * THREEFOURTH;
  296.    X := POW (BETA, EXPON);
  297.    Y := SIN (X);
  298.    WRITELN ('SIN (', X:15, ') =                   ', Y:15);
  299.    WRITELN;
  300.    WRITELN ('THE FOLLOWING THREE LINES ILLUSTRATE THE LOSS IN SIGNIFICANCE');
  301.    WRITELN ('FOR LARGE ARGUMENTS.  THE ARGUMENTS ARE CONSECUTIVE.');
  302.    WRITELN;
  303.    Z := SQRT (BETAP);
  304.    X := Z * (ONE - EPSNEG);
  305.    Y := SIN(X);
  306.    WRITELN ('SIN (', X:15, ') =                   ', Y:15);
  307.    Y := SIN(Z);
  308.    WRITELN ('SIN (', X:15, ') =                   ', Y:15);
  309.    X := Z * (ONE + EPS);
  310.    Y := SIN(X);
  311.    WRITELN ('SIN (', X:15, ') =                   ', Y:15);
  312.  
  313. {-----------------------------------------------------------------}
  314. {     TEST OF ERROR RETURNS                                       }
  315. {-----------------------------------------------------------------}
  316.  
  317.    WRITELN;
  318.    WRITELN;
  319.    WRITELN ('TEST OF ERROR RETURNS');
  320.    WRITELN;
  321.    X := BETAP;
  322.    WRITELN ('SIN WILL BE CALLED WITH THE ARGUMENT ',  X:15);
  323.    WRITELN ('THIS SHOULD TRIGGER AN ERROR MESSAGE');
  324.    Y := SIN(X);
  325.    WRITELN ('SIN RETURNED THE VALUE ', Y:15);
  326.    WRITELN;
  327.    WRITELN ('THIS CONCLUDES THE TESTS');
  328. END. { DSin }
  329.