How to compute the sine of huge numbers

For several days, I've been wondering how it would be possible to compute the sine of huge numbers like 100000! (radians). I obviously don't use double but cpp_rational from the boost multiprecision library. But I can't simply do 100000! mod 2pi and then use the builtin function sinl (I don't need more than 10 decimal digits..) as I'd need several million digits of pi to compute this accurately.

Is there any way to achieve this?


I believe you may be able to calculate this without obscene numbers of digits of $\pi$ if you take advantage of the fact that these are factorials. To simplify the algebra, we can calculate $a_n=e^{i(n!)}$ instead (you want the imaginary part). Then $$a_{n+1}=e^{i(n!)(n+1)}=a_n^{n+1},$$ and it's perfectly reasonable to calculate $a_{100000}$ recursively with a high-precision library.

The downside is that to start the recursion you need a very good approximation of $e^i$, and I don't know if the error dependence works out any differently than in the $\pmod{2\pi}$ approach.

But to answer your actual question, Mathematica doesn't even break a sweat with the mere million digits needed for this:

> Block[{$MaxExtraPrecision = 1000000}, N[Sin[100000!], 10]]

-0.9282669319

takes about 15 ms on my computer.


For calculating the sine or cosine of a large arbitrary precision real number $x$, the gains of this method (which are tuned for $\sin n$ for integer $n$) are mostly lost, so I would recommend your original idea of reducing the argument $\bmod 2\pi$. As has been noted, the main bottleneck is a high-precision estimation of $\pi$. Your answer will be useless unless you can at least calculate $\frac{x}{\pi}$ to within $1$ (otherwise you may as well answer "somewhere between $-1$ and $1$"), so you need at least $\log_2(x/\pi^2)$ bits of precision for $\pi$. With $x\approx100000!$, that's about $1516701$ bits or $456572$ digits. Add to this the number $a$ of bits of precision you want in the result, so about $1516734$ digits of $\pi$ to calculate $33$ bits ($\approx 10$ digits) of $\sin x$ in the range $x\approx 100000!$.

Once you have an integer $n$ such that $y=2\pi n$ is close to $x$ (ideally $|x-2\pi n|\le2\pi$, it doesn't have to be perfectly rounded), calculate $\pi$ to precision $a+\log_2(n)$, so that $y$ is known to precision $a$, and then $x-y$ is precision $a$ and $\sin x=\sin (x-y)$ can be calculated to precision $a$ as well.