home *** CD-ROM | disk | FTP | other *** search
/ Oakland CPM Archive / oakcpm.iso / sigm / sigmv082.ark / LN.PAS < prev    next >
Encoding:
Pascal/Delphi Source File  |  1984-04-29  |  1.3 KB  |  75 lines
  1. extern
  2. function ln ( x : real ): real;
  3.  
  4. const
  5. index = 16;
  6. ln4 = 0.13862943611196e+01;
  7. r3 = 0.33333333333333e+00;
  8. r5 = 0.20000000000000e+00;
  9. r7 = 0.14285714285714e+00;
  10. r9 = 0.11111111111111e+00;
  11. r11 = 0.90909090909091e-01;
  12. r13 = 0.76923076923077e-01;
  13. r15 = 0.66666666666667e-01;
  14. r17 = 0.58823529411765e-01;
  15. r19 = 0.52631578947368e-01;
  16. r21 = 0.47619047619048e-01;
  17. r23 = 0.43478260869565e-01;
  18. r25 = 0.40000000000000e-01;
  19. r27 = 0.37037037037037e-01;
  20. r29 = 0.34482758620690e-01;
  21. r31 = 0.32258064516129e-01;
  22.  
  23. var
  24. div_count,i : integer;
  25. result,term,term2 : real;
  26. p : array [1..index] of real;
  27.  
  28. begin (* ln - natural logarithm *)
  29. if x <= 0.0 then ln:= -0.99999999999999e+63
  30. else
  31. begin
  32. (* x must be in range 0.7 to 2.85 *)
  33. div_count:=0;
  34. while x < 0.7 do
  35.   begin
  36.   x:=x*4;
  37.   div_count:=div_count-1;
  38.   end;
  39. while x > 2.85 do
  40.   begin
  41.   x:=x/4;
  42.   div_count:=div_count+1;
  43.   end;
  44.  
  45. term:=(x-1.0)/(x+1.0);
  46. term2:=sqr(term);
  47. for i:=1 to index do
  48.   begin
  49.   p[i]:=term;
  50.   term:=term*term2;
  51.   end;
  52.  
  53. result:= 2 * (
  54.     p[1]
  55.     +(p[2] * r3)
  56.     +(p[3] * r5)
  57.     +(p[4] * r7)
  58.     +(p[5] * r9)
  59.     +(p[6] * r11)
  60.     +(p[7] * r13)
  61.     +(p[8] * r15)
  62.     +(p[9] * r17)
  63.     +(p[10] * r19)
  64.     +(p[11] * r21)
  65.     +(p[12] * r23)
  66.     +(p[13] * r25)
  67.     +(p[14] * r27)
  68.     +(p[15] * r29)
  69.     +(p[16] * r31)
  70.     );
  71. ln:=result + div_count * ln4;
  72. end; (* else *)
  73.  
  74. end;. (* ln *)
  75.