Gerhard Niklasch on Tue, 15 Jan 2002 14:14:34 +0100 (MET)


[Date Prev] [Date Next] [Thread Prev] [Thread Next] [Date Index] [Thread Index]

Re: easy precision problem ( I hope ;)


In response to:
> Message-Id: <001c01c19da3$d6103120$f436e0d4@schnippie>
> From: "Judith Gatz" <gatz@in.tum.de>
> To: "pari-users list" <pari-users@list.cr.yp.to>
> References: <Pine.GSO.4.33.0201081824480.19895-100000@geo>
> Date: Tue, 15 Jan 2002 10:05:16 +0100

> I want to make a number of exact factorisations like:
> 
> ? for(x=65,75,print(factor(truncate(precision(exp(x),1000)),0)))

First, you don't need 1000 digits to get the exact values of
truncate(exp(x)) in this range - 50 digits are more than enough.

Second, precision(exp(x),1000) does the wrong thing:  it computes
exp(x) at the current default realprecision  (28 digits)  and
_then_ tweaks the resulting object to look like a number of
greater precision  (by padding the significand with 0 bits,
but without changing the already-rounded value).

> the results cant be exact, there are too many 2^x:
> [2, 1; 3, 2; 941605135783518730078768673, 1]
> [2, 2; 19, 2; 43, 1; 1217, 1; 119723, 1; 5092511140946699, 1]

These two are correct.

> [2, 3; 401, 1; 1213, 1; 20929, 1; 1537753234015292017, 1]

This one isn't;  28 digits (default realprecision) are insufficient
to represent the integer part of exp(67).  Normally, truncate()
would complain:
  ***   precision loss in truncation
but the precision(_,1000) has fooled it into believing it had been
passed a sufficiently precise value.

> before that, I tried to increase the precision with:  (I don´t know if this
> makes sence...)
> ? primelimit=1000000000000000000000000000000000000000000000000000000000000
> %1 = 1000000000000000000000000000000000000000000000000000000000000
> ? parisize=100000000000000
> %2 = 100000000000000
> ? realprecision= 500
> %3 = 500

These introduce variables of the indicated names, but do not affect
the PARI defaults in any way.  Use
  default(realprecision, 80)
to increase the realprecision default, and leave the primelimit well
alone  (you'd need a computer about as large as the observable universe
to store a table of primes up to that number!  but factorization works
reliably at the default primelimit).

And then run the for loop without the precision(_,1000) call.

(Increase the realprecision gradually for larger x, guided by the size
of exp(x) - using ln(10)=2.302585..., it is easy to estimate how many
significant digits you are going to need.  truncate() will warn you,
but you may lose a bit or two just before you get that warning, so it's
better to work out in advance what's needed.)

(parisize is also excessive - the above would amount to a few hundred
Terabytes of memory.  The present computation does not use a lot of
stack space, and will happily proceed in the default few MBy.)

Hope this helps,
Gerhard