Re: [Rd] Numerical accuracy of matrix multiplication
> peter dalgaard > on Fri, 16 Sep 2016 13:33:11 +0200 writes: > On 16 Sep 2016, at 12:41 , Alexis Sarda wrote: >> Hello, >> >> while testing the crossprod() function under Linux, I noticed the following: >> >> set.seed(883) >> x <- rnorm(100) >> x %*% x - sum(x^2) # equal to 1.421085e-14 >> >> Is this difference normal? It seems to be rather large for double precision. >> > It's less than .Machine$double.eps, relative (!) to x %*% x ~= 100. indeed! Still, it gives exactly 0 on my platform(s), where I'm using R's own version of BLAS / Lapack. Are you perhaps using an "optimized" BLAS / LAPACK , i.e, one that is fast but slightly less so accurate ? Martin Maechler, ETH Zurich > -pd >> Regards, >> Alexis. >> >> [[alternative HTML version deleted]] >> >> __ >> R-devel@r-project.org mailing list >> https://stat.ethz.ch/mailman/listinfo/r-devel > -- > Peter Dalgaard, Professor, > Center for Statistics, Copenhagen Business School > Solbjerg Plads 3, 2000 Frederiksberg, Denmark > Phone: (+45)38153501 > Office: A 4.23 > Email: pd@cbs.dk Priv: pda...@gmail.com > __ > R-devel@r-project.org mailing list > https://stat.ethz.ch/mailman/listinfo/r-devel __ R-devel@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-devel
Re: [Rd] Numerical accuracy of matrix multiplication
> Alexis Sarda > on Tue, 20 Sep 2016 17:33:49 +0200 writes: > I just realized that I was actually using a different random number > generator, could that be a valid reason for the discrepancy? > The code should be: > RNGkind("L'Ecuyer") > set.seed(883) > x <- rnorm(100) > x %*% x - sum(x^2) # equal to 1.421085e-14 Yes, now I get the same result so my story on "BLAS / LAPACK" is not relevant here. But do note the main point from Peter Dalgaard that this is well within Machine epsilon precision. More precisely, here it is really one bit difference in the least significant bit : > print(rbind( x%*%x, crossprod(x), sum(x^2)), digits= 19) [,1] [1,] 103.5096830356289814 [2,] 103.5096830356289814 [3,] 103.5096830356289672 > cbind(sprintf("%a", c(x%*%x, crossprod(x), sum(x^2 [,1] [1,] "0x1.9e09ea598568fp+6" [2,] "0x1.9e09ea598568fp+6" [3,] "0x1.9e09ea598568ep+6" > > Regards, > Alexis Sarda. > On Tue, Sep 20, 2016 at 5:27 PM, Martin Maechler > wrote: >> > peter dalgaard >> > on Fri, 16 Sep 2016 13:33:11 +0200 writes: >> >> > On 16 Sep 2016, at 12:41 , Alexis Sarda >> wrote: >> >> >> Hello, >> >> >> >> while testing the crossprod() function under Linux, I noticed the >> following: >> >> >> >> set.seed(883) >> >> x <- rnorm(100) >> >> x %*% x - sum(x^2) # equal to 1.421085e-14 >> >> >> >> Is this difference normal? It seems to be rather large for double >> precision. >> >> >> >> > It's less than .Machine$double.eps, relative (!) to x %*% x ~= 100. >> >> indeed! >> >> Still, it gives exactly 0 on my platform(s), where I'm using R's >> own version of BLAS / Lapack. >> >> Are you perhaps using an "optimized" BLAS / LAPACK , i.e, one >> that is fast but slightly less so accurate ? >> >> Martin Maechler, >> ETH Zurich >> >> >> > -pd >> >> >> Regards, >> >> Alexis. >> >> >> >> [[alternative HTML version deleted]] >> >> >> >> __ >> >> R-devel@r-project.org mailing list >> >> https://stat.ethz.ch/mailman/listinfo/r-devel >> >> > -- >> > Peter Dalgaard, Professor, >> > Center for Statistics, Copenhagen Business School >> > Solbjerg Plads 3, 2000 Frederiksberg, Denmark >> > Phone: (+45)38153501 >> > Office: A 4.23 >> > Email: pd@cbs.dk Priv: pda...@gmail.com >> >> > __ >> > R-devel@r-project.org mailing list >> > https://stat.ethz.ch/mailman/listinfo/r-devel >> > [[alternative HTML version deleted]] __ R-devel@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-devel
Re: [Rd] Numerical accuracy of matrix multiplication
I just realized that I was actually using a different random number generator, could that be a valid reason for the discrepancy? The code should be: RNGkind("L'Ecuyer") set.seed(883) x <- rnorm(100) x %*% x - sum(x^2) # equal to 1.421085e-14 Regards, Alexis Sarda. On Tue, Sep 20, 2016 at 5:27 PM, Martin Maechler wrote: > > peter dalgaard > > on Fri, 16 Sep 2016 13:33:11 +0200 writes: > > > On 16 Sep 2016, at 12:41 , Alexis Sarda > wrote: > > >> Hello, > >> > >> while testing the crossprod() function under Linux, I noticed the > following: > >> > >> set.seed(883) > >> x <- rnorm(100) > >> x %*% x - sum(x^2) # equal to 1.421085e-14 > >> > >> Is this difference normal? It seems to be rather large for double > precision. > >> > > > It's less than .Machine$double.eps, relative (!) to x %*% x ~= 100. > > indeed! > > Still, it gives exactly 0 on my platform(s), where I'm using R's > own version of BLAS / Lapack. > > Are you perhaps using an "optimized" BLAS / LAPACK , i.e, one > that is fast but slightly less so accurate ? > > Martin Maechler, > ETH Zurich > > > > -pd > > >> Regards, > >> Alexis. > >> > >> [[alternative HTML version deleted]] > >> > >> __ > >> R-devel@r-project.org mailing list > >> https://stat.ethz.ch/mailman/listinfo/r-devel > > > -- > > Peter Dalgaard, Professor, > > Center for Statistics, Copenhagen Business School > > Solbjerg Plads 3, 2000 Frederiksberg, Denmark > > Phone: (+45)38153501 > > Office: A 4.23 > > Email: pd@cbs.dk Priv: pda...@gmail.com > > > __ > > R-devel@r-project.org mailing list > > https://stat.ethz.ch/mailman/listinfo/r-devel > [[alternative HTML version deleted]] __ R-devel@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-devel
Re: [Rd] Handlers in setGraphicsEventHandlers() can recursively call getGraphicsEvent(). Intended behavior?
Hi Is the correct patch to remove the setting of the gettingEvent flag or would it be better to flip the TRUE/FALSE setting (set to TRUE before handling then reset to FALSE after handling) ? Also, for this patch and for the other two you sent, one difficulty will be with testing the patches. I have no testing code for this, so would need at least a test or two from you (ideally someone would also have some regression tests, beyond ?getGraphicsEvent, to ensure continuity of previous behaviour). Paul On 18/09/2016 3:29 a.m., Richard Bodewits wrote: Hey all. As in general it's a bad idea to allow an event handler to generate an event, and as comments in the code seem to suggest this isn't the intention either, I was wondering about recursion in getGraphicsEvent(). In main/gevents.c, both doMouseEvent() and doKeybd() have the following line; dd->gettingEvent = FALSE; /* avoid recursive calls */ And they reset it to TRUE before returning. The effective result of this is that the event handlers on the R side are allowed to call getGraphicsEvent(), and recurse until they eventually would run out of stack space. Though R does catch and handle that cleanly with the error; Error: evaluation nested too deeply: infinite recursion / options(expressions=)? A quick scan of the SVN logs suggests this code has been untouched since its introduction in 2004, so I'm left to wonder if this is intended behavior. It stands out as a bit of a sore thumb due to the generic check for recursion in do_getGraphicsEvent() in the same file, which will error out with error(_("recursive use of 'getGraphicsEvent' not supported")) if dd->gettingEvent is already set to TRUE. Which would suggest recursively calling it is very much not intended to be possible. To me, setting gettingEvent to FALSE seems like an easy mistake to make if you temporarily interpret gettingEvent to mean that event(s) are allowed to still come in. But the actual interpretation in do_getGraphicsEvents() is the opposite, as it's interpreted as an indicator of whether or not an event is currently being processed. I've removed the gettingEvent altering lines from both doMouseEvent() and doKeybd() to no ill effect, and doing so disabled the ability to call getGraphicsEvent() from inside one of the assigned handlers as expected. But after 12 (!) years, it's conceivable that people have come to depend on this behavior in existing scripts. Is this something that should be left alone to minimize disruption? Or should this be fixed (if it is indeed unintended) for the sake of protecting people from infinite recursions? I've included a small patch as attachment that removes the offending lines. - Richard Bodewits __ R-devel@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-devel -- Dr Paul Murrell Department of Statistics The University of Auckland Private Bag 92019 Auckland New Zealand 64 9 3737599 x85392 p...@stat.auckland.ac.nz http://www.stat.auckland.ac.nz/~paul/ __ R-devel@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-devel
[Rd] Undocumented 'use.names' argument to c()
'c' has an undocumented 'use.names' argument. I'm not sure if this is a documentation or implementation bug. > c(a = 1) a 1 > c(a = 1, use.names = F) [1] 1 Karl __ R-devel@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-devel
Re: [Rd] Undocumented 'use.names' argument to c()
> On Sep 20, 2016, at 7:18 PM, Karl Millar via R-devel > wrote: > > 'c' has an undocumented 'use.names' argument. I'm not sure if this is > a documentation or implementation bug. It came up on stackoverflow a couple of years ago: http://stackoverflow.com/questions/24815572/why-does-function-c-accept-an-undocumented-argument/24815653#24815653 At the time it appeared to me to be a documentation lag. > >> c(a = 1) > a > 1 >> c(a = 1, use.names = F) > [1] 1 > > Karl > > __ > R-devel@r-project.org mailing list > https://stat.ethz.ch/mailman/listinfo/r-devel David Winsemius Alameda, CA, USA __ R-devel@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-devel