• 0 Vote(s) - 0 Average
• 1
• 2
• 3
• 4
• 5
 Precision check on [pentation.gp] SOLVED Cherrina_Pixie Junior Fellow Posts: 6 Threads: 3 Joined: Jun 2011 06/28/2011, 07:17 AM (This post was last modified: 07/02/2011, 01:43 AM by Cherrina_Pixie.) Code:? init(exp(Pi()/2));loop    base          4.81047738096535165547304    fixed point   0.E-68 + 1.00000000000000000000000*I    Pseudo Period 3.69464335841375533580710 + 1.06216001044294092389502*I 4 strm(s) out of 12 sexp(z) generates 36 Riemann samples, scnt= 43 8 rtrm(s) out of 18 riemaprx(z) generates 12 sexp samples sexp(-0.5)= 0.44150846390775332819992735716755 6.315800093 Riemann/sexp binary precision bits I=0.1200000000*I 1=loopcount -0.005257090997 recenter/renorm 0.005021271902 18 strm(s) out of 18 sexp(z) generates 26 Riemann samples, scnt= 36 12 rtrm(s) out of 13 riemaprx(z) generates 19 sexp samples sexp(-0.5)= 0.44146829736572151824153868728668 15.39341456 Riemann/sexp binary precision bits I=0.1200000000*I 2=loopcount -0.000006612199794 recenter/renorm 0.00001187700815 22 strm(s) out of 28 sexp(z) generates 38 Riemann samples, scnt= 50 19 rtrm(s) out of 19 riemaprx(z) generates 22 sexp samples sexp(-0.5)= 0.44146826345089157095194288166600 23.36284166 Riemann/sexp binary precision bits I=0.1200000000*I 3=loopcount -0.00000002312286847 recenter/renorm 0.00000007284067824 33 strm(s) out of 33 sexp(z) generates 50 Riemann samples, scnt= 62 25 rtrm(s) out of 25 riemaprx(z) generates 29 sexp samples sexp(-0.5)= 0.44146826265404109196819015909251 31.45371224 Riemann/sexp binary precision bits I=0.1200000000*I 4=loopcount -1.197253818 E-10 recenter/renorm 1.162183490 E-10 40 strm(s) out of 43 sexp(z) generates 64 Riemann samples, scnt= 75 32 rtrm(s) out of 32 riemaprx(z) generates 30 sexp samples sexp(-0.5)= 0.44146826265346684519672624211813 39.80570360 Riemann/sexp binary precision bits I=0.1200000000*I 5=loopcount -3.150872206 E-13 recenter/renorm 7.227756425 E-13 45 strm(s) out of 45 sexp(z) generates 80 Riemann samples, scnt= 88 40 rtrm(s) out of 40 riemaprx(z) generates 35 sexp samples sexp(-0.5)= 0.44146826265345961129724099420262 48.14861781 Riemann/sexp binary precision bits I=0.1200000000*I 6=loopcount -1.151385006 E-15 recenter/renorm 1.259558322 E-15 52 strm(s) out of 52 sexp(z) generates 94 Riemann samples, scnt= 100 47 rtrm(s) out of 47 riemaprx(z) generates 38 sexp samples sexp(-0.5)= 0.44146826265345960793409283844758 56.61697584 Riemann/sexp binary precision bits I=0.1200000000*I 7=loopcount -2.738837756 E-18 recenter/renorm 6.799954147 E-18 57 strm(s) out of 57 sexp(z) generates 110 Riemann samples, scnt= 113 55 rtrm(s) out of 55 riemaprx(z) generates 42 sexp samples sexp(-0.5)= 0.44146826265345960787456889962878 65.05123526 Riemann/sexp binary precision bits I=0.1200000000*I 8=loopcount -7.703853099 E-21 recenter/renorm 1.516819744 E-20 63 strm(s) out of 63 sexp(z) generates 124 Riemann samples, scnt= 126 62 rtrm(s) out of 62 riemaprx(z) generates 46 sexp samples sexp(-0.5)= 0.44146826265345960787455655053558 73.11307204 Riemann/sexp binary precision bits I=0.1200000000*I 9=loopcount -2.263257063 E-23 recenter/renorm 8.389515433 E-23 69 strm(s) out of 69 sexp(z) generates 138 Riemann samples, scnt= 139 69 rtrm(s) out of 69 riemaprx(z) generates 49 sexp samples sexp(-0.5)= 0.44146826265345960787455568662140 81.22151804 Riemann/sexp binary precision bits I=0.1200000000*I 10=loopcount -1.148962259 E-25 recenter/renorm 1.542482215 E-25 73 strm(s) out of 73 sexp(z) generates 154 Riemann samples, scnt= 151 77 rtrm(s) out of 77 riemaprx(z) generates 52 sexp samples sexp(-0.5)= 0.44146826265345960787455568616458 89.35146101 Riemann/sexp binary precision bits I=0.1200000000*I 11=loopcount -2.948891538 E-28 recenter/renorm 1.148593582 E-27 78 strm(s) out of 78 sexp(z) generates 168 Riemann samples, scnt= 164 84 rtrm(s) out of 84 riemaprx(z) generates 56 sexp samples sexp(-0.5)= 0.44146826265345960787455568615250 92.04321478 Riemann/sexp binary precision bits I=0.1200000000*I 12=loopcount 2.894810870 E-29 recenter/renorm 2.344784572 E-28 84 strm(s) out of 84 sexp(z) generates 172 Riemann samples, scnt= 168 86 rtrm(s) out of 86 riemaprx(z) generates 56 sexp samples sexp(-0.5)= 0.44146826265345960787455568615249 92.01407805 Riemann/sexp binary precision bits I=0.1200000000*I 13=loopcount 3.046137094 E-29 recenter/renorm 2.376219052 E-28 UNEXPECTED LOSS: curprecision of base ~4.81 and somehow, the indeterminacy of the 13th loop exceeded that of the 12th loop. Is there a way to solve this issue without changing the base? I'm not exactly certain about the possible consequences of this, other than the likely limit of precision... is a result like this 'normal' for larger bases? sheldonison Long Time Fellow Posts: 684 Threads: 24 Joined: Oct 2008 06/28/2011, 02:33 PM (This post was last modified: 06/28/2011, 07:01 PM by sheldonison.) (06/28/2011, 07:17 AM)Cherrina_Pixie Wrote: .... I ran through the loops for of base ~4.81 and somehow, the indeterminacy of the 13th loop exceeded that of the 12th loop. Is there a way to solve this issue without changing the base? I'm not exactly certain about the possible consequences of this, other than the likely limit of precision... is a result like this 'normal' for larger bases?The problem is with the tetration/sexp code (which is used by pentation.gp). I haven't seen this bug before in kneser.gp, despite having tried a large number of different bases from 1.45 to 100,000 or so. But this particular base, exp(pi/2), seems to be a problem. I don't yet know why. Although in earlier versions of kneser.gp I saw it when pari-gp wasn't carrying out calculations with enough precision -- although I've fixed all of those bugs by manually updating the precision for the samples used to generate theta(z). I tried a few things, but this particular base seems to stop converging around 80 or 90 binary bits of precision (depending on which algorithm I use). I don't understand it yet, and I'm somewhat concerned that this means many other bases don't converge either. I know the error term is approximately a linear dot product from one iteration to the next -- I'm going to try some experiments and I'll report back results later. Small update: The fixed point, L, for $b=\exp(\pi/2)$, has L=i, with $\Re(L)=0$. I don't know what is special about the fixed point of L=i. I tried literally hundred of random bases between 1.5 and 10, (mostly using a faster version of kneser.gp, which uses a Taylor series for the Schroeder function and its inverse which is 3x faster than the version here). All of the random bases converged, so the problem is isolated. Then I tried bases nearby $b=\exp(\pi/2)+\delta$. Using the default kneser.gp precision (\p 67), I get failures $b=\exp(\pi/2)+10^{-25}$, but not for larger values of delta. Using a higher precision, \p 134, I get failures for $\delta=10^{-28}$, but where the precision plateaus depends on how small delta is. I'll report some graphs later; hopefully I'll start to make some sense of what's going, and why this one particular base may not converge .... - Sheldon JmsNxn Ultimate Fellow Posts: 889 Threads: 110 Joined: Dec 2010 06/29/2011, 02:42 AM (This post was last modified: 06/29/2011, 02:44 AM by JmsNxn.) That's funny, I was just thinking about $exp^{\pi/2}$. Naively I wanna say it's because $\lim_{n\to\infty} ^n(e^{\pi/2}) = \lim_{n\to\infty} ^n(i^{1/i})= i$ which is what usually occurs in the self root function. However this would imply real to imaginary which isn't very pretty. Funny that the code fails around it tho. I think there is something special about $exp^{\pi/2}$, I'm just not sure what. sheldonison Long Time Fellow Posts: 684 Threads: 24 Joined: Oct 2008 06/29/2011, 05:05 AM (This post was last modified: 06/29/2011, 10:24 PM by sheldonison.) (06/29/2011, 02:42 AM)JmsNxn Wrote: That's funny, I was just thinking about $exp^{\pi/2}$. Naively I wanna say it's because $\lim_{n\to\infty} ^n(e^{\pi/2}) = \lim_{n\to\infty} ^n(i^{1/i})= i$ which is what usually occurs in the self root function. However this would imply real to imaginary which isn't very pretty. Funny that the code fails around it tho. I think there is something special about $exp^{\pi/2}$, I'm just not sure what.Hey James, You are correct and I think you stated the problem more clearly then I can. The Schroder function limit equations don't seem to quite match the inverse Schroder function limit equations for this base, and using more iterations doesn't seem to help, and using a smaller radius doesn't seem to help. So that in turn means my superfunction and its inverse are only approximately inverses of each other. And that in turn messes up the theta calculation so it stops converging. And all of this is somehow because in the limit, L=i, and in the limit, the function is rotating around approximately pi/2 per iteration. But I really don't understand it and I'm not sure how to fix it. When I graph $h^{-1}(h(i+z))-i-z$, around the unit circle, with a small radius for z, I get this ugly noisy error term. Here, h is the Schroeder function, which is used to generate the inverse superfunction. Alternatively, graphing $h(h^{-1}(z))$ gives similar ugly results. This plot was calculated with \p 134, for 134 significant digits of precision, which should have left 66 or 67 digits of accuracy...     I think a closed form for the coefficients of the Schroeder and inverse Schroeder series might work better.... At least this problem seems to effect only this one base, exp(Pi/2). Here is a code fragment, to generate the Schroeder and inverse Schroeder functions, and then generate the plot above. In the real code, of course these numbers would be slightly more approximate, with very small error terms, but that doesn't seem to matter. - Sheldon Code:\p 134 L=I;          /* fixed point */ B=exp(Pi/2);  /* base        */ lnB=Pi/2;     /* log(base)   */ rlnB=1/lnB;   /* 1/log(base) */ logLxlnB = log(I*lnB); /* =log(log(L)) */ scnt = 336;   /* 336 iteration count, good for 67 digits of precision in the superfunction */ superm = real(exp(-scnt*logLxlnB)); /* multiplier, for superfunction; real valued number */ isuperm = real(exp(scnt*logLxlnB)); /* multiplier, inverse superfunction, also real number */ superf(z) = { /* complex valued superfunction for base B, generated from the fixed point I */ /* This is the inverse Schroeder function */   local (i,y);   y=L+z*superm;   for (i=1,scnt, y=exp(y*lnB));   return(y); } isuperf(z) = { /* complex valued inverse superfunction for base B, generated from the fixed point L */ /* This is the Schroeder function */   local (y,i);   y=z;   for (i=1,scnt, y=log(y)*rlnB); /* rlnB=1/log(B) */   y=y-L;   y=y*isuperm;   return(y); } ploth(t=0,2*Pi,x=0.1*exp(t*I);z=superf(isuperf(L+x))-L-x;[real(z),imag(z)]); ploth(t=0.001,1.999*Pi,x=0.1*exp(t*I);z=isuperf(superf(x))-x;[real(z),imag(z)]); Also, there's a big error spike in the imaginary of: isuperf(superf(0.1))-0.1 2.123527622415261486036 E-29 + 2.929305537712787288382 E-24*I sheldonison Long Time Fellow Posts: 684 Threads: 24 Joined: Oct 2008 06/29/2011, 10:36 PM (This post was last modified: 06/29/2011, 10:42 PM by sheldonison.) (06/29/2011, 05:05 AM)sheldonison Wrote: .... my superfunction and its inverse are only approximately inverses of each other. And that in turn messes up the theta calculation so it stops converging..... But I really don't understand it and I'm not sure how to fix it. So, its the inverse Schroeder equation that is not working for $b=\exp(\pi/2)$. The Schroeder equation, which has the iterated logarithms, and corresponds to the inverse superfunction, is converging very nicely. But the inverse Schroeder equation which I use for the superfunction seems to conerge -- only up to a point -- and then stops converging. The defining characteristic of the inverse Schroeder equation, which I will call g(z), is $b=\exp(\pi/2)$ $L=i\;\;\; b^L=L$ $\lambda = \log(L)=\frac{\pi i}{2}$ $b^{g(z)} = g(\lambda z)$, for all z. The g(z) function should be entire. $g(z) = \lim_{n \to \infty} \exp_b^{on} (L+\lambda^{-n})$ The trouble I'm having, is that increasing the number of iterations and the precision doesn't seem to improve the results. I can sample g(z) around a unit circle with a reasonable radius so the function is well behaved, and yet there seems to be a lot of "noise" in the resulting Taylor series for the function. I was using \p 134, for 134 decimal digits of precision, and 336 iterations for n, and still only getting 25 digits of precision for g(z). And I can get the same or even better results, with 67 decimal digits precision and 165 iteration. The code for g(z) is in the previous post, where I graphed $g(g^{-1}(L+z))-L-z$. I verified that the taylor series approximation for g(z) isn't following the defining equation definition too well, $b^{g(z)} = g(\lambda z)$. I guess this could be some kind of pari-gp artifact for this particular base, but I still don't understand it, or it could be mathematical, but I still don't understand it. I also verified that the Schroeder equation (inverse superfunction) is very well behaved for this base, and has nothing to do with the problem. Also, for all other bases I've tried, both the superfunction and the inverse superfunction are well behaved. Only bases very close to this base are affected, within a 10^-15 radius. - Sheldon sheldonison Long Time Fellow Posts: 684 Threads: 24 Joined: Oct 2008 06/30/2011, 06:00 AM (This post was last modified: 06/30/2011, 06:16 AM by sheldonison.) (06/29/2011, 10:36 PM)sheldonison Wrote: (06/29/2011, 05:05 AM)sheldonison Wrote: .... my superfunction and its inverse are only approximately inverses of each other. And that in turn messes up the theta calculation so it stops converging..... But I really don't understand it and I'm not sure how to fix it. So, its the inverse Schroeder equation that is not working for $b=\exp(\pi/2)$.... the inverse Schroeder equation which I use for the superfunction seems to conerge -- only up to a point -- and then stops converging. .... I guess this could be some kind of pari-gp artifact for this particular base, but I still don't understand it.... Also, for all other bases I've tried, both the superfunction and the inverse superfunction are well behaved. Only bases very close to this base are affected, within a 10^-15 radius.It is a pari-gp precision artifact in the superfunction code that affects this one and only base. I ran into similar problems, so that in general, I needed to update the precision of the function manually once per iteration, before taking the inverse superfunction. But this problem is unique to this base, $\exp(\pi/2)$. Pari-gp gets confused on how much precision the real and imaginary parts of the superfunction have, because the real part of the fixed point is zero, and the imaginary part of the fixed point is non-zero. For the test code case in my previous post, initially, the imaginary part has 134 decimal digits of precision, and pari-gp quickly decides the real part has only 67 digits of precision. We need 67 decimal digits of precision at the end, for the superfunction to work. Pari-gp gradually loses more and more precision in the calculations, while iterating through the loop, and 2/3rds of the way through the iterations, the real and imag both have only 38 decimal digits of precision remaining, which apparently leaves about 29 decimal digits of precision in the result. For other bases, pari-gp keeps the same precision throughout the entire loop. Here is a code patch, which I think slows down the code significantly enough that I'm not ready to make the default. Maybe I need special code for this one case when $\Re(L)=0$? - Sheldon Code:/* patch for superf(z) loss of precision for base exp(pi/2) */ superf(z) = { /* complex valued superfunction for base B, generated from the fixed point L */   local (y);   y=z-scnt;   y=L+((LxlnB)^y); /* LxlnB=L*log(B) */   for (i=1,scnt,     y=exp(y*lnB);     y=precision(y,precis); /* this is the patch */   );   return(y); } /* for comparison, current superf function, without the patch */ superf(z) = { /* complex valued superfunction for base B, generated from the fixed point L */   local (y);   y=z-scnt;   y=L+((LxlnB)^y); /* LxlnB=L*log(B) */   for (i=1,scnt, y=exp(y*lnB));   return(y); } sheldonison Long Time Fellow Posts: 684 Threads: 24 Joined: Oct 2008 07/01/2011, 10:56 PM The new version of pentation.gp includes the fix for $b=\exp(pi/2)$, and also runs faster. - Sheldon Cherrina_Pixie Junior Fellow Posts: 6 Threads: 3 Joined: Jun 2011 07/02/2011, 01:39 AM Code:? init(exp(Pi()/2));loop;genpent    base          4.81047738096535165547304    fixed point   0.E-145 + 1.00000000000000000000000*I    Pseudo Period 3.69464335841375533580710 + 1.06216001044294092389502*I generating superf taylor series; inverse Schroder equation, scnt 361 generating isuperf taylor series; Schroder equation, scnt 361 sexp(-0.5)= 0.44149388556590800271258410632775 1=loopcnt  6.734438825 Riemann/sexp binary precision bits sexp(-0.5)= 0.4414683320304119199073161516641704293457800006586290888270406592 2=loopcnt  16.14879308 Riemann/sexp binary precision bits sexp(-0.5)= 0.4414682630880826737650093490296176463777662860412226388888800474 3=loopcnt  24.03560113 Riemann/sexp binary precision bits sexp(-0.5)= 0.4414682626543550892577935620363027361724384104152456844251132995 4=loopcnt  32.04840534 Riemann/sexp binary precision bits sexp(-0.5)= 0.4414682626534634768199249898188360294331037280968748020127519588 5=loopcnt  40.38278448 Riemann/sexp binary precision bits sexp(-0.5)= 0.4414682626534596146811726186051271337321983516231997071615402916 6=loopcnt  48.80427930 Riemann/sexp binary precision bits sexp(-0.5)= 0.4414682626534596079045772784636752654658548653394290591963824627 7=loopcnt  57.30135824 Riemann/sexp binary precision bits sexp(-0.5)= 0.4414682626534596078746214147072484385025198343669362617159556147 8=loopcnt  65.83172331 Riemann/sexp binary precision bits sexp(-0.5)= 0.4414682626534596078745561240878470035479673141387946649519338734 9=loopcnt  73.82815204 Riemann/sexp binary precision bits sexp(-0.5)= 0.4414682626534596078745556870620231393941232039016898103310735241 10=loopcnt  81.91862008 Riemann/sexp binary precision bits sexp(-0.5)= 0.4414682626534596078745556861590530560306188724303721881960943663 11=loopcnt  90.09474083 Riemann/sexp binary precision bits sexp(-0.5)= 0.4414682626534596078745556861541674585781220374176541575975256630 12=loopcnt  98.32748089 Riemann/sexp binary precision bits sexp(-0.5)= 0.4414682626534596078745556861541578659514474143180100763842212618 13=loopcnt  106.6292769 Riemann/sexp binary precision bits sexp(-0.5)= 0.4414682626534596078745556861541578179603766836663909037659999642 14=loopcnt  114.9625455 Riemann/sexp binary precision bits sexp(-0.5)= 0.4414682626534596078745556861541578178760347939093471360666938970 15=loopcnt  123.3730561 Riemann/sexp binary precision bits sexp(-0.5)= 0.4414682626534596078745556861541578178754814000844681999624069557 16=loopcnt  131.4033281 Riemann/sexp binary precision bits sexp(-0.5)= 0.4414682626534596078745556861541578178754801786292553376291713396 17=loopcnt  139.5050918 Riemann/sexp binary precision bits sexp(-0.5)= 0.4414682626534596078745556861541578178754801717399690722642964464 18=loopcnt  147.6454299 Riemann/sexp binary precision bits sexp(-0.5)= 0.4414682626534596078745556861541578178754801717245455070461469587 19=loopcnt  155.8374822 Riemann/sexp binary precision bits sexp(-0.5)= 0.4414682626534596078745556861541578178754801717244680794029640458 20=loopcnt  164.0688306 Riemann/sexp binary precision bits sexp(-0.5)= 0.4414682626534596078745556861541578178754801717244679136515673714 21=loopcnt  172.3506366 Riemann/sexp binary precision bits sexp(-0.5)= 0.4414682626534596078745556861541578178754801717244679134861786060 22=loopcnt  182.0804715 Riemann/sexp binary precision bits sexp(-0.5)= 0.4414682626534596078745556861541578178754801717244679134851850198 23=loopcnt  190.4876162 Riemann/sexp binary precision bits sexp(-0.5)= 0.4414682626534596078745556861541578178754801717244679134851843918 24=loopcnt  198.7628907 Riemann/sexp binary precision bits sexp(-0.5)= 0.4414682626534596078745556861541578178754801717244679134851843818 25=loopcnt  207.2180495 Riemann/sexp binary precision bits sexp(-0.5)= 0.4414682626534596078745556861541578178754801717244679134851843818 26=loopcnt  215.6331171 Riemann/sexp binary precision bits sexp(-0.5)= 0.4414682626534596078745556861541578178754801717244679134851843818 27=loopcnt  223.9743246 Riemann/sexp binary precision bits sexp(-0.5)= 0.4414682626534596078745556861541578178754801717244679134851843818 28=loopcnt  232.3153383 Riemann/sexp binary precision bits sexp(-0.5)= 0.4414682626534596078745556861541578178754801717244679134851843818 29=loopcnt  239.6076665 Riemann/sexp binary precision bits pentation base        4.81047738096535165547304 pentation(-0.5)       0.460975347014689190545301 sexp fixed point      -1.94648466297646768934135 sexp slope at fixed   11.8300204422121443717564 pentation period      2.54314035015181068719998*I pentation singularity -2.58922659360979627202546 + 1.27157017507590534359999*I pentation precision, via sexp(pent(-0.5))-pent(0.5)                       -6.42327456356856623983420 E-45` Excellent!! =) I used \p 144 and that took less than three minutes! « Next Oldest | Next Newest »

 Possibly Related Threads… Thread Author Replies Views Last Post Attempt to find a limit point but each step needs doubling the precision... Gottfried 15 36,648 11/09/2014, 10:25 PM Last Post: tommy1729 Expansion of base-e pentation andydude 13 42,306 07/02/2011, 01:40 AM Last Post: Cherrina_Pixie How to force precision in SAGE? jaydfox 2 8,603 08/18/2007, 11:08 PM Last Post: jaydfox

Users browsing this thread: 1 Guest(s)