On 10/09/2010 6:46 AM, Renaud Gaujoux wrote:
Hi,
suppose you have two versions of the same algorithm: one in pure R, the
other one in C/C++ called via .Call().
Assuming there is no bug in the implementations (i.e. they both do the
same thing), is there any well known reason why the C/C++ implementation
could return numerical results non identical to the one obtained from
the pure R code? (e.g. could it be rounding errors? please explain.)
Has anybody had a similar experience?
R often uses extended reals (80 bit floating point values on Intel
chips) for intermediate values. C compilers may or may not do that.
By not identical, I mean very small differences (< 2.4 e-14), but enough
to have identical() returning FALSE. Maybe I should not bother, but I
want to be sure where the differences come from, at least by mere curiosity.
Briefly the R code perform multiple matrix product; the C code is an
optimization of those specific products via custom for loops, where
entries are not computed in the same order, etc... which improves both
memory usage and speed. The result is theoretically the same.
Changing the order of operations will often affect rounding. For
example, suppose epsilon is the smallest number such that 1 + epsilon is
not equal to 1. Then 1 + (epsilon/2) + (epsilon/2) will evaluate to
either 1 or 1 + epsilon, depending on the order of computing the additions.
Duncan Murdoch
Thank you,
Renaud
______________________________________________
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
