Ravi Varadhan
2015-Jul-03 19:01 UTC
[R] : Ramanujan and the accuracy of floating point computations - using Rmpfr in R
Thank you all. I did think about declaring `pi' as a special constant, but for some reason didn't actually try it. Would it be easy to have the mpfr() written such that its argument is automatically of extended precision? In other words, if I just called: mpfr(exp(sqrt(163)*pi, 120), then all the constants, 163, pi, are automatically of 120 bits precision. Is this easy to do? Best, Ravi ________________________________________ From: David Winsemius <dwinsemius at comcast.net> Sent: Friday, July 3, 2015 2:06 PM To: John Nash Cc: r-help; Ravi Varadhan Subject: Re: [R] : Ramanujan and the accuracy of floating point computations - using Rmpfr in R> On Jul 3, 2015, at 8:08 AM, John Nash <John.Nash at uottawa.ca> wrote: > > > > > Third try -- I unsubscribed and re-subscribed. Sorry to Ravi for extra > traffic. > > In case anyone gets a duplicate, R-help bounced my msg "from non-member" (I changed server, but not email yesterday, > so possibly something changed enough). Trying again. > > JN > > I got the Wolfram answer as follows: > > library(Rmpfr) > n163 <- mpfr(163, 500) > n163 > pi500 <- mpfr(pi, 500) > pi500 > pitan <- mpfr(4, 500)*atan(mpfr(1,500)) > pitan > pitan-pi500 > r500 <- exp(sqrt(n163)*pitan) > r500 > check <- "262537412640768743.99999999999925007259719818568887935385..." > savehistory("jnramanujan.R") > > Note that I used 4*atan(1) to get pi.RK got it right by following the example in the help page for mpfr: Const("pi", 120) The R `pi` constant is not recognized by mpfr as being anything other than another double . There are four special values that mpfr recognizes. ? Best; David> It seems that may be important, > rather than converting. > > JN > > On 15-07-03 06:00 AM, r-help-request at r-project.org wrote: > >> Message: 40 >> Date: Thu, 2 Jul 2015 22:38:45 +0000 >> From: Ravi Varadhan <ravi.varadhan at jhu.edu> >> To: "'Richard M. Heiberger'" <rmh at temple.edu>, Aditya Singh >> <aps6dl at yahoo.com> >> Cc: r-help <r-help at r-project.org> >> Subject: Re: [R] : Ramanujan and the accuracy of floating point >> computations - using Rmpfr in R >> Message-ID: >> <14ad39aaf6a542849bbf3f62a0c2f38f at DOM-EB1-2013.win.ad.jhu.edu> >> Content-Type: text/plain; charset="utf-8" >> >> 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. >> ------------------------------ > > > ______________________________________________ > 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.David Winsemius, MD Alameda, CA, USA
David Winsemius
2015-Jul-04 00:08 UTC
[R] : Ramanujan and the accuracy of floating point computations - using Rmpfr in R
It doesn?t appear to me that mpfr was ever designed to handle expressions as the first argument. ? David> On Jul 3, 2015, at 12:01 PM, Ravi Varadhan <ravi.varadhan at jhu.edu> wrote: > > Thank you all. I did think about declaring `pi' as a special constant, but for some reason didn't actually try it. > Would it be easy to have the mpfr() written such that its argument is automatically of extended precision? In other words, if I just called: mpfr(exp(sqrt(163)*pi, 120), then all the constants, 163, pi, are automatically of 120 bits precision. > > Is this easy to do? > > Best, > Ravi > ________________________________________ > From: David Winsemius <dwinsemius at comcast.net> > Sent: Friday, July 3, 2015 2:06 PM > To: John Nash > Cc: r-help; Ravi Varadhan > Subject: Re: [R] : Ramanujan and the accuracy of floating point computations - using Rmpfr in R > >> On Jul 3, 2015, at 8:08 AM, John Nash <John.Nash at uottawa.ca> wrote: >> >> >> >> >> Third try -- I unsubscribed and re-subscribed. Sorry to Ravi for extra >> traffic. >> >> In case anyone gets a duplicate, R-help bounced my msg "from non-member" (I changed server, but not email yesterday, >> so possibly something changed enough). Trying again. >> >> JN >> >> I got the Wolfram answer as follows: >> >> library(Rmpfr) >> n163 <- mpfr(163, 500) >> n163 >> pi500 <- mpfr(pi, 500) >> pi500 >> pitan <- mpfr(4, 500)*atan(mpfr(1,500)) >> pitan >> pitan-pi500 >> r500 <- exp(sqrt(n163)*pitan) >> r500 >> check <- "262537412640768743.99999999999925007259719818568887935385..." >> savehistory("jnramanujan.R") >> >> Note that I used 4*atan(1) to get pi. > > RK got it right by following the example in the help page for mpfr: > > Const("pi", 120) > > The R `pi` constant is not recognized by mpfr as being anything other than another double . > > > There are four special values that mpfr recognizes. > > ? > Best; > David > > >> It seems that may be important, >> rather than converting. >> >> JN >> >> On 15-07-03 06:00 AM, r-help-request at r-project.org wrote: >> >>> Message: 40 >>> Date: Thu, 2 Jul 2015 22:38:45 +0000 >>> From: Ravi Varadhan <ravi.varadhan at jhu.edu> >>> To: "'Richard M. Heiberger'" <rmh at temple.edu>, Aditya Singh >>> <aps6dl at yahoo.com> >>> Cc: r-help <r-help at r-project.org> >>> Subject: Re: [R] : Ramanujan and the accuracy of floating point >>> computations - using Rmpfr in R >>> Message-ID: >>> <14ad39aaf6a542849bbf3f62a0c2f38f at DOM-EB1-2013.win.ad.jhu.edu> >>> Content-Type: text/plain; charset="utf-8" >>> >>> 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. >>> ------------------------------ >> >> >> ______________________________________________ >> 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. > > David Winsemius, MD > Alameda, CA, USA >David Winsemius, MD Alameda, CA, USA
David Winsemius
2015-Jul-04 01:45 UTC
[R] : Ramanujan and the accuracy of floating point computations - using Rmpfr in R
> On Jul 3, 2015, at 5:08 PM, David Winsemius <dwinsemius at comcast.net> wrote: > > > > It doesn?t appear to me that mpfr was ever designed to handle expressions as the first argument.This could be a start. Obviously one would wnat to include code to do other substitutions probably using the all.vars function to pull out the other ?constants? and ?numeric? values to make them of equivalent precision. I?m guessing you want to follow the parse-tree and then testing the numbers for integer-ness and then replacing by paste0( ?mpfr(?, val, ?L, ?, prec,?)? ) Pre <- function(expr, prec){ sexpr <- deparse(substitute(expr) ) spi <- paste0("Const('pi',",prec,?)?) sexpr <- gsub("pi", spi, sexpr) print( eval(parse(text=sexpr))) } # Very minimal testing> Pre(pi,120)1 'mpfr' number of precision 120 bits [1] 3.1415926535897932384626433832795028847> Pre(exp(pi),120)1 'mpfr' number of precision 120 bits [1] 23.140692632779269005729086367948547394> Pre(log(exp(pi)),120)1 'mpfr' number of precision 120 bits [1] 3.1415926535897932384626433832795028847 At the moment it still needs to ahve other numbers mpfr-ified to succeed:> Pre(pi-4*atan(1),120)1 'mpfr' number of precision 120 bits [1] 1.2246467991473531772315719253130892016e-16> Pre(pi-4*atan(mpfr(1,120)),120)1 'mpfr' number of precision 120 bits [1] 0 Best; David.>> ? > David > >> On Jul 3, 2015, at 12:01 PM, Ravi Varadhan <ravi.varadhan at jhu.edu> wrote: >> >> Thank you all. I did think about declaring `pi' as a special constant, but for some reason didn't actually try it. >> Would it be easy to have the mpfr() written such that its argument is automatically of extended precision? In other words, if I just called: mpfr(exp(sqrt(163)*pi, 120), then all the constants, 163, pi, are automatically of 120 bits precision. >> >> Is this easy to do? >> >> Best, >> Ravi >> ________________________________________ >> From: David Winsemius <dwinsemius at comcast.net> >> Sent: Friday, July 3, 2015 2:06 PM >> To: John Nash >> Cc: r-help; Ravi Varadhan >> Subject: Re: [R] : Ramanujan and the accuracy of floating point computations - using Rmpfr in R >> >>> On Jul 3, 2015, at 8:08 AM, John Nash <John.Nash at uottawa.ca> wrote: >>> >>> >>> >>> >>> Third try -- I unsubscribed and re-subscribed. Sorry to Ravi for extra >>> traffic. >>> >>> In case anyone gets a duplicate, R-help bounced my msg "from non-member" (I changed server, but not email yesterday, >>> so possibly something changed enough). Trying again. >>> >>> JN >>> >>> I got the Wolfram answer as follows: >>> >>> library(Rmpfr) >>> n163 <- mpfr(163, 500) >>> n163 >>> pi500 <- mpfr(pi, 500) >>> pi500 >>> pitan <- mpfr(4, 500)*atan(mpfr(1,500)) >>> pitan >>> pitan-pi500 >>> r500 <- exp(sqrt(n163)*pitan) >>> r500 >>> check <- "262537412640768743.99999999999925007259719818568887935385..." >>> savehistory("jnramanujan.R") >>> >>> Note that I used 4*atan(1) to get pi. >> >> RK got it right by following the example in the help page for mpfr: >> >> Const("pi", 120) >> >> The R `pi` constant is not recognized by mpfr as being anything other than another double . >> >> >> There are four special values that mpfr recognizes. >> >> ? >> Best; >> David >> >> >>> It seems that may be important, >>> rather than converting. >>> >>> JN >>> >>> On 15-07-03 06:00 AM, r-help-request at r-project.org wrote: >>> >>>> Message: 40 >>>> Date: Thu, 2 Jul 2015 22:38:45 +0000 >>>> From: Ravi Varadhan <ravi.varadhan at jhu.edu> >>>> To: "'Richard M. Heiberger'" <rmh at temple.edu>, Aditya Singh >>>> <aps6dl at yahoo.com> >>>> Cc: r-help <r-help at r-project.org> >>>> Subject: Re: [R] : Ramanujan and the accuracy of floating point >>>> computations - using Rmpfr in R >>>> Message-ID: >>>> <14ad39aaf6a542849bbf3f62a0c2f38f at DOM-EB1-2013.win.ad.jhu.edu> >>>> Content-Type: text/plain; charset="utf-8" >>>> >>>> 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. >>>> ------------------------------ >>> >>> >>> ______________________________________________ >>> 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. >> >> David Winsemius, MD >> Alameda, CA, USA >> > > David Winsemius, MD > Alameda, CA, USA > > ______________________________________________ > 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.David Winsemius, MD Alameda, CA, USA
Gabor Grothendieck
2015-Jul-04 14:42 UTC
[R] : Ramanujan and the accuracy of floating point computations - using Rmpfr in R
You can do that with bc if you pass the entire expression to bc(...) in quotes but in that case you will have to use bc notation, not R notation so, for example, exp is e and atan is a.> library(bc) > bc("e(sqrt(163)*4*a(1))")[1] "262537412640768743.9999999999992500725971981856888793538563373369908627075374103782106479101186073116295306145602054347" On Fri, Jul 3, 2015 at 3:01 PM, Ravi Varadhan <ravi.varadhan at jhu.edu> wrote:> Thank you all. I did think about declaring `pi' as a special constant, > but for some reason didn't actually try it. > Would it be easy to have the mpfr() written such that its argument is > automatically of extended precision? In other words, if I just called: > mpfr(exp(sqrt(163)*pi, 120), then all the constants, 163, pi, are > automatically of 120 bits precision. > > Is this easy to do? > > Best, > Ravi > ________________________________________ > From: David Winsemius <dwinsemius at comcast.net> > Sent: Friday, July 3, 2015 2:06 PM > To: John Nash > Cc: r-help; Ravi Varadhan > Subject: Re: [R] : Ramanujan and the accuracy of floating point > computations - using Rmpfr in R > > > On Jul 3, 2015, at 8:08 AM, John Nash <John.Nash at uottawa.ca> wrote: > > > > > > > > > > Third try -- I unsubscribed and re-subscribed. Sorry to Ravi for extra > > traffic. > > > > In case anyone gets a duplicate, R-help bounced my msg "from non-member" > (I changed server, but not email yesterday, > > so possibly something changed enough). Trying again. > > > > JN > > > > I got the Wolfram answer as follows: > > > > library(Rmpfr) > > n163 <- mpfr(163, 500) > > n163 > > pi500 <- mpfr(pi, 500) > > pi500 > > pitan <- mpfr(4, 500)*atan(mpfr(1,500)) > > pitan > > pitan-pi500 > > r500 <- exp(sqrt(n163)*pitan) > > r500 > > check <- "262537412640768743.99999999999925007259719818568887935385..." > > savehistory("jnramanujan.R") > > > > Note that I used 4*atan(1) to get pi. > > RK got it right by following the example in the help page for mpfr: > > Const("pi", 120) > > The R `pi` constant is not recognized by mpfr as being anything other than > another double . > > > There are four special values that mpfr recognizes. > > ? > Best; > David > > > > It seems that may be important, > > rather than converting. > > > > JN > > > > On 15-07-03 06:00 AM, r-help-request at r-project.org wrote: > > > >> Message: 40 > >> Date: Thu, 2 Jul 2015 22:38:45 +0000 > >> From: Ravi Varadhan <ravi.varadhan at jhu.edu> > >> To: "'Richard M. Heiberger'" <rmh at temple.edu>, Aditya Singh > >> <aps6dl at yahoo.com> > >> Cc: r-help <r-help at r-project.org> > >> Subject: Re: [R] : Ramanujan and the accuracy of floating point > >> computations - using Rmpfr in R > >> Message-ID: > >> <14ad39aaf6a542849bbf3f62a0c2f38f at DOM-EB1-2013.win.ad.jhu.edu> > >> Content-Type: text/plain; charset="utf-8" > >> > >> 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. > >> ------------------------------ > > > > > > ______________________________________________ > > 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. > > David Winsemius, MD > Alameda, CA, USA > > > ______________________________________________ > 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. >-- Statistics & Software Consulting GKX Group, GKX Associates Inc. tel: 1-877-GKX-GROUP email: ggrothendieck at gmail.com [[alternative HTML version deleted]]