ran 3 trials of 100K to check the recovery of a rational number from its decimal
American Citizen <[email protected]>
| Newsgroups | gmane.comp.mathematics.pari.user |
|---|---|
| Message-ID | <[email protected]> |
Found some interesting results. I used the random(10^30) function to
create 30 digits rational number then set their decimal value and sought
to recover the original
Listed are those programs and program results, using
contfracpnqn(contfrac(decimal)), then bestappr(decimal, N) and finally
algdep(decimal,1)
Using contfracpnqn(contfrac(decimal)) method
N=10^30;
default(realprecision,77);
default(seriesprecision,77);
for(s=1,100000,\
n=random(N);d=random(N);\
if(d!=0,a=n/d,a=1);\
b=1.0*a;\
c=contfracpnqn(contfrac(b));\
d=c[1,1]/c[2,1];\
write("trial1",a-d);\
);
##
*** last result: cpu time 16,722 ms, real time 16,722 ms.
owner@localhost:~> sort trial1 | uniq -c
100000 0 (all are correct)
CAVEAT: precision of at least twice digits in N must be used
2nd Method: bestappr(decimal, bounded)
N=10^30;
default(realprecision,77);
for(s=1,100000,\
n=random(N);d=random(N);\
if(d!=0,a=n/d,a=1);\
b=1.0*a;\
c=bestappr(b,N);\
write("trial2",a-c);\
);
##
owner@localhost:~> gp -q < 001
*** last result: cpu time 16,768 ms, real time 16,768 ms.
sort trial2 | uniq -c
100000 0
(all 100,000 okay)
CAVEAT: N has to be used to set the denominator size in the bestappr
command when not set, all 100,000 trials failed. This requires knowledge
ahead of time of the approximate digits in the rational to be recovered.
2nd CAVEAT: precision has to be a least twice that as N's digits (57
failed here)
3rd method - algdep(decimal,1)
N=10^30;
default(realprecision,77);
for(s=1,100000,\
n=random(N);d=random(N);\
if(d!=0,a=n/d,a=1);\
b=1.0*a;\
c=Vec(algdep(b,1));\
d=-c[2]/c[1];\
write("trial3",a-d);\
);
##
*** last result: cpu time 17,551 ms, real time 17,551 ms.
owner@localhost:~> sort trial3 | uniq -c
100000 0
CAVEAT: precision must be at least twice digits of N
So the contfracpnqn(contfrac(decimal)) method is slightly faster than
bestappr(decimal,N) and algdep takes about 1 second more, for 100K trials.
This confirms what I have observed, but now I have tested this for 30
digit rationals. (at least 100K of them)
And somewhat irritatingly, I got a zero from the random(N) function
call, which resulted in a zero in the denominator, so I had to watch out
for that.
Randall