Far more interesting is the question of computing
sin y (and cos y) to arbitrary precision.

a. Since sin y is periodic the argument that needs
to be handled by the core routine is in the range 
0 to 2p1 .

b. By repeated application of the sum formula the 
argument y can be further restricted to the range 
0 to 0.25p1 .

c. The Taylor series for sin and cos are:

(sin y) = -/ (y ^ i) % ! i=. 1 + 2 * i. n
(cos y) = -/ (y ^ i) % ! i=.     2 * i. n

Since 0.25p1 is less than 1, the series converge 
rapidly.  The error is bounded by (y^1+2*n)%!1+2*n
for sin and (y^2*n)%!2*n for cos.

For example, what is sin 0.8 to 40 digits?  sin 0.8
is equivalent to cos 0.5p1-0.8 .  

   y=: 0.5p1 - 0.8
   y
0.770796
   1e_40 (> i. 1:) | (y^2*n)%!2*n=: i.30
17
   (y^2*n)%!2*n=: 17
4.85046e_43

So 17 terms suffice.  First some sanity checks:

   sin 0.8
0.717356
   cos 0.5p1 - 0.8
0.717356
   -/ (y^i) % ! i=: 2 * i. 17
0.717356

Now the actual calculation:

   y1=: 4r5 -~ (2 * 10^40x) %~ <[EMAIL PROTECTED] 10^40x
   y1
15415926535897932384626433832795028841971r20000000000000000000000000000000000000000
   0j40 ": y1
0.7707963267948966192313216916397514420986
   0j40 ": -/ (y1^i) % ! i=: 2 * i. 17x
0.7173560908995227616271746105813853661928

In comparison, Abramowitz & Stegun says (Table 4.6,
page 158) that sin 0.8 is

0.71735 60908 99522 76162 718

I believe the last digit in the A&S figure is 
incorrect.


----------------------------------------------------------------------
For information about J forums see http://www.jsoftware.com/forums.htm

Reply via email to