[R] : Ramanujan and the accuracy of floating point computations - using Rmpfr in R

[Richard M. Heiberger]

precedence does matter in this example.  the square
root was taken of a doubleprecision (53 bit) number.  my revision
takes the square root of a 120 bit number.

> sqrt(mpfr(pi, 120))
1 'mpfr' number of precision  120   bits
[1] 1.7724538509055159927515191031392484397
> mpfr(sqrt(pi), 120)
1 'mpfr' number of precision  120   bits
[1] 1.772453850905515881919427556567825377
> print(sqrt(pi), digits=20)
[1] 1.7724538509055158819
>

Sent from my iPhone

> On Jul 3, 2015, at 00:38, Ravi Varadhan <ravi.varadhan at jhu.edu> wrote:
>
> Hi Rich,
>
> The Wolfram answer is correct.
> http://mathworld.wolfram.com/RamanujanConstant.html
>
> There is no code for Wolfram alpha.  You just go to their web engine and plug in the expression and it will give you the answer.
> http://www.wolframalpha.com/
>
> I am not sure that the precedence matters in Rmpfr.  Even if it does, the answer you get is still wrong as you showed.
>
> Thanks,
> Ravi
>
> -----Original Message-----
> From: Richard M. Heiberger [mailto:rmh at temple.edu]
> Sent: Thursday, July 02, 2015 6:30 PM
> To: Aditya Singh
> Cc: Ravi Varadhan; r-help
> Subject: Re: [R] : Ramanujan and the accuracy of floating point computations - using Rmpfr in R
>
> There is a precedence error in your R attempt.  You need to convert
> 163 to 120 bits first, before taking
> its square root.
>
>> exp(sqrt(mpfr(163, 120)) * mpfr(pi, 120))
> 1 'mpfr' number of precision  120   bits
> [1] 262537412640768333.51635812597335712954
>
> ## just the last four characters to the left of the decimal point.
>> tmp <- c(baseR=8256, Wolfram=8744, Rmpfr=8333, wrongRmpfr=7837)
>> tmp-tmp[2]
>     baseR    Wolfram      Rmpfr wrongRmpfr
>      -488          0       -411       -907
>
> You didn't give the Wolfram alpha code you used.  There is no way of verifying the correct value from your email.
> Please check that you didn't have a similar precedence error in that code.
>
> Rich
>
>
>> On Thu, Jul 2, 2015 at 2:02 PM, Aditya Singh via R-help <r-help at r-project.org> wrote:
>>
>> Ravi
>>
>> I am a chemical engineer by training. Is there not something like law of corresponding states in numerical analysis?
>>
>> Aditya
>>
>>
>>
>> ------------------------------
>>> On Thu 2 Jul, 2015 7:28 AM PDT Ravi Varadhan wrote:
>>>
>>> Hi,
>>>
>>> Ramanujan supposedly discovered that the number, 163, has this interesting property that exp(sqrt(163)*pi), which is obviously a transcendental number, is real close to an integer (close to 10^(-12)).
>>>
>>> If I compute this using the Wolfram alpha engine, I get:
>>> 262537412640768743.99999999999925007259719818568887935385...
>>>
>>> When I do this in R 3.1.1 (64-bit windows), I get:
>>> 262537412640768256.0000
>>>
>>> The absolute error between the exact and R's value is 488, with a relative error of about 1.9x10^(-15).
>>>
>>> In order to replicate Wolfram Alpha, I tried doing this in "Rmfpr" but I am unable to get accurate results:
>>>
>>> library(Rmpfr)
>>>
>>>
>>>> exp(sqrt(163) * mpfr(pi, 120))
>>>
>>> 1 'mpfr' number of precision  120   bits
>>>
>>> [1] 262537412640767837.08771354274620169031
>>>
>>> The above answer is not only inaccurate, but it is actually worse than the answer using the usual double precision.  Any thoughts as to what I am doing wrong?
>>>
>>> Thank you,
>>> Ravi
>>>
>>>
>>>
>>>      [[alternative HTML version deleted]]
>>>
>>> ______________________________________________
>>> R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see
>>> https://stat.ethz.ch/mailman/listinfo/r-help
>>> PLEASE do read the posting guide
>>> http://www.R-project.org/posting-guide.html
>>> and provide commented, minimal, self-contained, reproducible code.
>>
>> ______________________________________________
>> R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide
>> http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.