[R] problem with the precision of numbers

[(Ted Harding)]

On 19-Jan-10 17:55:43, Ben Bolker wrote:
> kayj <kjaja27 <at> yahoo.com> writes:
>> Hi All,
>> 
>> I was wodering if it is possible to increase the precision using R.
>> I ran the script below in R and MAPLE and I got different results
>> when k is large.
>> Any idea how to fix this problem? thanks for your help
>> 
>> for (k in 0:2000){
>>  s=0
>>  for(i in 0:k){
>>  s=s+((-1)^i)*3456*(1+i*1/2000)^3000
>>  }
>> }
> 
> (1) see
> http://wiki.r-project.org/rwiki/doku.php?id=misc:r_accuracy:high_precisi
> on_arithmetic
> 
> (2) consider whether there is more accurate algorithm you
> could use. I don't recognize the series, but perhaps it
> has a closed form solution, maybe as a special function?
> How much accuracy do you really need in the solution?
> 
>   Ben Bolker

I suspect this is an invented computation -- the "3456" strikes
me as "unlikely" (it reminds me of my habitual illustrative use
of set.seed(54321)).

There is a definite problem with the development given by kayj.
When k=2000 and i=k, the formula requires evaluation of

  3456*(2^3000)

on a log10 scale this is

  log10(3456) + 3000*log10(2) = 906.6286

Since R "gives up" at 10^308.25471 = 1.79767e+308
(10^308.25472 => Inf), this algorithm is going to be tricky to
evaluate!

I don't know how well Rmpfr copes with very large numbers (the
available documentation seems cryptic). However, I can go along
with the recommendation in the URL the Ben gives, to use 'bc'
("Berkeley Calculator"), available on unix[oid] systems since
a long time ago. That is an old friend of mine, and works well
(it can cope with exponents up to X^2147483647 in the version
I have). It can eat for breakfast the task of checking whether
Kate Bush can accurately sing pi to 117 significant figures:

  http://www.absolutelyrics.com/lyrics/view/kate_bush/pi

(Try it in R).

Ted.

--------------------------------------------------------------------
E-Mail: (Ted Harding) <Ted.Harding at manchester.ac.uk>
Fax-to-email: +44 (0)870 094 0861
Date: 19-Jan-10                                       Time: 18:41:27
------------------------------ XFMail ------------------------------