home *** CD-ROM | disk | FTP | other *** search
/ ProfitPress Mega CDROM2 …eeware (MSDOS)(1992)(Eng) / ProfitPress-MegaCDROM2.B6I / PROG / PASCAL / MISCTI10.ZIP / TI227.ASC < prev    next >
Encoding:
Text File  |  1988-04-18  |  6.5 KB  |  302 lines
  1. PRODUCT : TURBO PASCAL WITH BCD SUPPORT     NUMBER : 227
  2. VERSION : 3.01
  3.      OS : MS-DOS, PC-DOS, CP/M-86
  4.    DATE : August 4, 1986
  5.  
  6.   TITLE : TRANSCENDENTAL FUNCTIONS 
  7.  
  8. The  following example routines are public domain  programs  that 
  9. have  been uploaded to our Forum on CompuServe.  As a courtesy to 
  10. our  users  that  do not have  immediate  access  to  CompuServe, 
  11. Technical Support distributes these routines free of charge.
  12.  
  13. However,  because these routines are public domain programs,  not 
  14. developed by Borland International,  we are unable to provide any 
  15. technical support or assistance using these routines. If you need 
  16. assistance   using   these   routines,    or   are   experiencing 
  17. difficulties,  we  recommend  that you log  onto  CompuServe  and 
  18. request  assistance  from the Forum members that developed  these 
  19. routines.
  20.  
  21. Written by Randall A. Gacek
  22.  
  23. This  is  a  first approximation of a set of routines to  do  the 
  24. transcendental functions LOG,  LN, SQRT, ARCTAN, SIN, COS and EXP 
  25. in the BCD version of Turbo Pascal.   
  26. WARNING:  The following code is specific to the implementation of 
  27. Turbo  Pascal with BCD support.   These functions should only  be 
  28. used with this implementation. 
  29.  
  30. program checkfuncs;
  31.  
  32. function sqrt(x:real):real;
  33.  
  34. var
  35.   n,i,m :integer;
  36.   f,y   :real;
  37.   v     :record case boolean of
  38.             true:(y:real);
  39.            false:(z:array[1..10] of byte)
  40.           end;
  41. begin
  42.   if x = 0.0 then
  43.     sqrt:=0.0
  44.   else if x < 0.0 then
  45.     halt
  46.   else begin
  47.     v.y:=x;
  48.     n:=v.z[1]-63;
  49.     v.z[1]:=63;
  50.     f:=v.y;
  51.     y:=0.580661+f/2.0-0.086462/(f+0.175241);
  52.     for i:=1 to 2 do
  53.       y:=0.5*(y+f/y);
  54.     y:=y+0.5*(f/y-y);
  55.     if odd(n) then
  56.     begin
  57.       y:=y*0.316227766016837933;
  58.       n:=n+1;
  59.     end;
  60.  
  61. checkfuncs Con't.
  62.  
  63.     m:=n div 2;
  64.     v.y:=y;
  65.     v.z[1]:=v.z[1]+m;
  66.     sqrt:=v.y;
  67.   end;
  68. end; { sqrt }
  69.  
  70. function log(x:real):real;
  71.   const
  72.     c0= 0.316227766016837933;
  73.     a0=-0.260447002405557636E+2;
  74.     a1= 0.554085912041205931E+2;
  75.     a2=-0.392737410203156250E+2;
  76.     a3= 0.103338571514793865E+2;
  77.     a4=-0.741010784161919239E+0;
  78.     b0=-0.899552077881033117E+2;
  79.     b1= 0.245347618868489348E+3;
  80.     b2=-0.244303035341829542E+3;
  81.     b3= 0.107109789115668009E+3;
  82.     b4=-0.193732345832854786E+2;
  83.     c=  0.868588963806503655;
  84.  
  85.   var
  86.     n:integer;
  87.     xn,f,s,w,aw,bw,rs2,rs:real;
  88.     v:record case boolean of
  89.       true:(y:real);
  90.       false:(z:array[1..10] of byte)
  91.     end;
  92.  
  93. begin
  94.   if x <= 0.0 then
  95.     halt;
  96.   v.y:=x;
  97.   n:=v.z[1]-63;
  98.   v.z[1]:=63;
  99.  
  100.   f:=v.y;
  101.   if f <= c0 then
  102.   begin
  103.     n:=n-1;
  104.     f:=f*10.0;
  105.   end;
  106.  
  107. checkfuncs Con't.
  108.  
  109.   s:=((f-0.5)-0.5)/(f+1.0);
  110.   w:=sqr(s);
  111.   aw:= (((a4*w+a3)*w+a2)*w+a1)*w+a0;
  112.   bw:=((((w+b4)*w+b3)*w+b2)*w+b1)*w+b0;
  113.   rs2:=w*aw/bw;
  114.   rs:=s*(c+rs2);
  115.   xn:=n;
  116.   log:=xn+rs;
  117. end; { log }
  118.  
  119. function ln(x:real):real;
  120.   const
  121.     c3=2.30258509299404568;
  122.  
  123. begin
  124.     ln:=log(x)*c3;
  125. end;
  126.  
  127. function exp(x:real):real;
  128.   const
  129.     bigx=147.365445951618923;
  130.     smallx=-145.062860858624878;
  131.     eps=5.0e-19;
  132.     onelnsqrt10=0.868588963806503655;
  133.     c1=1.151;
  134.     c2=2.92546497022842009e-4;
  135.     p0=0.333267029226801611e+6;
  136.     p1=0.100974148724273918E+5;
  137.     p2=0.420414268137450315E+2;
  138.     q0=0.666534058453603223E+6;
  139.     q1=0.757393346159883444E+5;
  140.     q2=0.841243584514154545E+3;
  141.     sqrt10=3.16227766016837933;
  142.   var
  143.     n:integer;
  144.     xn,g,z,gpz,qz,rg:real;
  145.     v:record case boolean of
  146.       true:(y:real);
  147.  
  148.       false:(z:array[1..10] of byte)
  149.       end;
  150.  
  151. checkfuncs Con't.
  152.  
  153. begin
  154.   if x > bigx then
  155.     halt;
  156.   if x < smallx then
  157.     halt;
  158.   if abs(x) < eps then
  159.     exp:=1.0
  160.   else begin
  161.     n:=round(x*onelnsqrt10);
  162.     xn:=n;
  163.     g:=(x-xn*c1)-xn*c2;
  164.     z:=sqr(g);
  165.     gpz:=((p2*z+p1)*z+p0)*g;
  166.     qz:= ((z+q2)*z+q1)*z+q0;
  167.     rg:=(0.5+gpz/(qz-gpz))*2.0;
  168.     if odd(n) then
  169.       if n >= 0 then
  170.         rg:=sqrt10*rg
  171.       else
  172.         rg:=rg/sqrt10;
  173.     n:=n div 2;
  174.     v.y:=rg;
  175.     v.z[1]:=v.z[1]+n;
  176.     exp:=v.y;
  177.   end;
  178. end; { exp }
  179.  
  180. function sincos(x,y,sgn:real):real;
  181.   const
  182.     ymax=3141592654.0;
  183.     onepi=0.318309886183790672;
  184.     c1= 3.141;                  { pi to 22 digits }
  185.     c2= 0.000592653589793238463;
  186.     eps=1.0e-9;
  187.     r1=-0.166666666666666651e+0;
  188.     r2= 0.833333333333316503E-2;
  189.     r3=-0.198412698412018405E-3;
  190.     r4= 0.275573192101527561E-5;
  191.     r5=-0.250521067982745845E-7;
  192.     r6= 0.160589364903715891E-9;
  193.     r7=-0.764291780689104677E-12;
  194.     r8= 0.272047909578888462E-14;
  195.  
  196. checkfuncs Con't.
  197.  
  198.   var
  199.     xn,f,t,g,rg:real;
  200. begin
  201.   if y >= ymax then
  202.     halt;
  203.   xn:=y*onepi;
  204.   xn:=int(xn+0.5);
  205.   if frac(xn / 2.0) <> 0.0 then
  206.     sgn:=-sgn;
  207.   if abs(x) <> y then { cos wanted }
  208.     xn:=xn-0.5;
  209.   f:=(abs(x)-xn*c1)-xn*c2;
  210.   if abs(f) < eps then
  211.     t:=f
  212.   else begin
  213.     g:=sqr(f);
  214.     rg:=(((((((r8*g+r7)*g+r6)*g+r5)*g+r4)*g+r3)*g+r2)*g+r1)*g;
  215.     t:=f+f*rg;
  216.   end;
  217.   sincos:=sgn*t;
  218. end; { sincos }
  219.  
  220. function sin(x:real):real;
  221. begin
  222.   if x < 0.0 then
  223.     sin:=sincos(x,-x,-1.0)
  224.   else
  225.     sin:=sincos(x,x,1.0);
  226. end; {sin}
  227.  
  228. function cos(x:real):real;
  229. begin
  230.   cos:=sincos(x,abs(x)+1.57079632679489662,1.0);
  231. end; {cos}
  232.  
  233. checkfuncs Con't.
  234.  
  235. function arctan(x:real):real;
  236.   const
  237.     twomsqrt3=0.267949192431122706;
  238.     sqrt3=1.73205080756887729;
  239.     a=0.732050807568877294;
  240.     eps=1e-9;
  241.     p0=-0.136887688941919269e+2;
  242.     p1=-0.205058551958616520e+2;
  243.     p2=-0.849462403513206835e+1;
  244.     p3=-0.837582993681500593e+0;
  245.     q0= 0.410663066825757813e+2;
  246.     q1= 0.861573495971302425e+2;
  247.     q2= 0.595784361425973445e+2;
  248.     q3= 0.150240011600285761e+1;
  249.  
  250.   var
  251.     n:integer;
  252.     f,result,g,gpg,qg,r:real;
  253. begin
  254.   f:=abs(x);
  255.   if f > 1.0 then
  256.   begin
  257.     f:=1.0/f;
  258.     n:=2;
  259.   end
  260.   else
  261.     n:=0;
  262.   if f > twomsqrt3 then
  263.   begin
  264.     f:=(((a*f-0.5)-0.5)+f)/(sqrt3+f);
  265.     n:=n+1;
  266.   end;
  267.   if abs(f) < eps then
  268.     result:=f
  269.   else begin
  270.     g:=sqr(f);
  271.     gpg:=(((p3*g+p2)*g+p1)*g+p0)*g;
  272.     qg:=(((g+q3)*g+q2)*g+q1)*g+q0;
  273.     r:=gpg/qg;
  274.     result:=f+f*r;
  275.   end;
  276.   if n > 1 then
  277.     result:=-result;
  278.  
  279. checkfuncs Con't.
  280.  
  281.   case n of
  282.     0:;
  283.     1:result:=0.523598775598298873+result;
  284.     2:result:=1.57079632679489662+result;
  285.     3:result:=1.04719755119659775+result;
  286.   end;
  287.   if x < 0.0 then
  288.     result:=-result;
  289.   arctan:=result;
  290. end; { arctan }
  291.  
  292. begin
  293.   writeln('sqrt= ',sqrt(25));
  294.   writeln('ln = ',ln(25));
  295.   writeln('exp = ',exp(25));
  296.   writeln('cos = ',cos(25));
  297.   writeln('sin = ',sin(25));
  298.   writeln('log = ',log(25));
  299.   writeln('arctan = ',arctan(25));
  300. end.
  301.  
  302.