home *** CD-ROM | disk | FTP | other *** search
/ norge.freeshell.org (192.94.73.8) / 192.94.73.8.tar / 192.94.73.8 / pub / computers / cpm / alphatronic / JRTPAS40.ZIP / LN.PAS < prev    next >
Pascal/Delphi Source File  |  1999-04-03  |  2KB  |  67 lines

  1.  
  2. EXTERN
  3.  
  4.  
  5. {=================================================================}
  6. FUNCTION ln (x : real) : real;
  7. CONST
  8. index = 16;
  9. ln4 = 0.13862943611196e+01;
  10. r3 = 0.33333333333333e+00;
  11. r5 = 0.20000000000000e+00;
  12. r7 = 0.14285714285714e+00;
  13. r9 = 0.11111111111111e+00;
  14. r11 = 0.90909090909091e-01;
  15. r13 = 0.76923076923077e-01;
  16. r15 = 0.66666666666667e-01;
  17. r17 = 0.58823529411765e-01;
  18. r19 = 0.52631578947368e-01;
  19. r21 = 0.47619047619048e-01;
  20. r23 = 0.43478260869565e-01;
  21. r25 = 0.40000000000000e-01;
  22. r27 = 0.37037037037037e-01;
  23. r29 = 0.34482758620690e-01;
  24. r31 = 0.32258064516129e-01;
  25. VAR
  26. div_count, i : integer;
  27. result, term, term2 : real;
  28. p :  ARRAY [1..index] OF real;
  29.  
  30. BEGIN  (* ln - natural logarithm *)
  31.  
  32. IF x <= 0.0 THEN
  33.        ln :=  - 0.99999999999999e+63
  34. ELSE
  35.        BEGIN
  36.          (* x must be in range 0.7 to 2.85 *)
  37.        
  38.        div_count := 0;
  39.        WHILE x < 0.7 DO
  40.               BEGIN
  41.               x := x * 4;
  42.               div_count := div_count - 1;
  43.               END;
  44.        WHILE x > 2.85 DO
  45.               BEGIN
  46.               x := x / 4;
  47.               div_count := div_count + 1;
  48.               END;
  49.        term := (x - 1.0) / (x + 1.0);
  50.        term2 := sqr(term);
  51.        FOR i := 1 TO index DO
  52.               BEGIN
  53.               p[i] := term;
  54.               IF (abs(term) <= 0.9e-51) THEN
  55.                      term := 0.0
  56.               ELSE
  57.                      term := term * term2;
  58.               END;
  59.        result := 2.0 * (p[1] + (p[2] * r3) + (p[3] * r5) + (p[4] * 
  60.               r7) + (p[5] * r9) + (p[6] * r11) + (p[7] * r13) + (
  61.               p[8] * r15) + (p[9] * r17) + (p[10] * r19) + (p[11] * 
  62.               r21) + (p[12] * r23) + (p[13] * r25) + (p[14] * r27) + (
  63.               p[15] * r29) + (p[16] * r31));
  64.        ln := result + div_count * ln4;
  65.        END;  (* else *)
  66. END;  (* ln *)  .
  67.