[R] Basic question on function "identical"

[Thomas Lumley]

On Mon, 15 Dec 2003 Ted.Harding at nessie.mcc.ac.uk wrote:

> On 15-Dec-03 Thomas Lumley wrote:
> >
> > One reason that which.max() exists is that we cannot guarantee
> > which(x==max(x)) to work. It is possible, though rather unlikely, for
> > there to be no x such that x==max(x).  One reason is the unpredictable
> > use of 10-byte wide floating point registers on Intel chips.
>
> Hmmm ...
>
> I'd be interested to learn more of what you mean. For instance,
> in C-speak, to find the maximum of an array of double x[], of
> length n, something like the following code could be written:
>
>   xmax=x[1];
>   for(i=1;i<n;i++) if(x[i+1]>x[i]) xmax=x[i+1];
>
> Regardless of the accuracy of the comparison, each assignment
> xmax = ... should make xmax an internally exact copy of the
> thing on the righthand side.

No, this is not guaranteed (though it is very likely).  Suppose that you
have just computed x as, say, sin(theta), and the value of x[n] happens to
still be stored in a floating point register on the CPU, with 64 bits of
precision (and 16 exponent bits, making a total of 10 bytes).

The compiler can use this value and decide it is the largest, and could
end up with xmax being this value in the register.  This value is not
necessarily the same as would be obtained by loading x[n] back in from
memory, and so need not compare equal.

In this case I think it is very unlikely that there is a problem, and I
think it can be guaranteed that
	which(x==max(x))
will return the right answer or nothing.

 However, we have really had a problem with the math library at one point
that was due to something vaguely of this sort.

This is why compilers provide options like --ffloat-store, and one reason
to use `volatile'.  It's also probably why so many Windows drivers seem
downgrade the floating point precision of the chip to 54 bits, to match
the precision of a double, although I seem to recall Brian Ripley saying
that this doesn't actually buy you as much as you would think.

>				However, your reply suggests that this
> may not happen, as a result of "unpredictable use of 10-byte wide
> floating point registers on Intel chips".
> Is this really the case? If so, how would a discrepancy arise?
> (I know C programmers like to use "register variables" where
> possible, for speed, but the copying should be faithful, surely).
>

This doesn't really have anything to do with the C type qualifier
"register", which is apparently ignored by most modern compilers, as they
think they can do register allocation better than the user.

	-thomas