08/20/2010, 01:05 AM
(08/19/2010, 02:35 AM)sheldonison Wrote: Later, I will probably post results showing how the Taylor and Fourier series terms converge, and how many terms it takes to get appropriate precision in the theta(z) fourier series so as to generate sufficiently accurate samples for the sexp(z) series......Here's an example, settings which I like as the next step up from the term term double precision results, which take less than a minute. Type "init(exp(1));loop(14);" and this next example runs for 7-8 minutes on my laptop, and reports that resulting sexp(z) Taylor series is accurate to 108 binary bits, or just a bit over 32 decimal digits. That precision would be within a radius of abs(z)<=0.5.
\p 67;
ctrm = 75;
cvec = 150;
rtrm = 1000;
rvec = 2000;
sfuncprec = 1E-32;
Here, gp's precision is set to 67 digits. The goal is 32 decimal digits of precision.
The sfuncprec (superfunction precision) is set to 1E-32 digits, to give superf/isuperf results accurate to 32 digits. That is, the isuperf/superf iteration count was set to 235 iterations. If you take the
67-32 leaves 35 digits of precision for gp to work with, or 3 spare bits for rounding errors, which surprisingly seems sufficient! We're asking for a Taylor series with 75 terms. By the Nyquist sampling criteria, we need at least 150 sample points. So cvec (the number of sample terms used to generate the sexp(z) Taylor series) =150. There will be a small sampling error term, due to higher frequencies in the superfunction being sampled. But, for base e sexp(z), the higher frequency terms roll off pretty fast, with each term being roughly half the magnitude of the term before. We want results accurate to 32 decimal digits, so we need 150 sample points, accurate to 32 decimal digits. The sample points will be in the upper half of the complex plane, centered around z=0. The points closest to the real axis are
rtrm = 1000 is driven by the requirement that the Riemann approximation, which is
We also need to show that our 75 term sexp(z) Taylor series supports 32 decimal digits of precision at a radius with abs(z)<=0.5. This is for the unit length we use to generate the RiemannCircle/Theta function from. The unit length will be from -0.5+idelta to +0.5+idelta. We want all terms over the unit length accurate to 32 decimal digits. The 75th term is -7.15 E-25. 7.15E-25*(0.5^75)=1.9E-47, which is much much less then 1E-32. The sexp(z) Taylor series has 32 decimal digits of precision out to a radius of about 0.78, which is greater then 1/2.
So, all of the parameters are sufficient to get around 32 bits of precision. Observationally, when iterating the "loop(14)", we get about 108 binary bits equal in the sexp(z) algorithm versus the riemaprx(z) algorithm. This is about 33 decimal digits. Graphing the error gainst a golden copy with higher precision shows similar results.
- Sheldon