[Rd] Non identical numerical results from R code vs C/C++ code?
murdoch.duncan at gmail.com
Fri Sep 10 13:00:25 CEST 2010
On 10/09/2010 6:46 AM, Renaud Gaujoux wrote:
> 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.
> Thank you,
More information about the R-devel