[R] logical inconsistency

[John C Nash]

This actually goes back a very long way. Peter is right to remind us 
that "optimizers"  (in the sense of compilers) can corrupt algorithms 
that are well-designed. Optimizing in tests is something some of us have 
fought for nearly 40 years, but compiler writers don't do much 
floating-point computation -- they sometimes save a few microseconds of 
computer time for many days and weeks of human time designing around 
such silliness. It is worse than a "bug" in that the code to work around 
the optimization can be very complicated and arcane -- and may not port 
across architectures. Clearly a program should execute as the code 
instructs, and not as a compiler designer decides to reinterpret it. 
Nevertheless, with a fairly large offset of 16.0 or 256.0 I have never 
seen such a test fail -- and it can port across different precision and 
different radix arithmetics.  However, I must say I tend to prefer to 
turn off compile optimization and try to keep my code clean i.e, 
manually optimized if possible. I also save the two sides of the test to 
separate variables, though I can guess that some compilers would corrupt 
that pretty easily.

Note that historically this issue didn't upset early codes where 
optimization was manual, then became an issue to watch out for on each 
implementation when we computed locally, and now could be a nasty but 
hidden trap with grid and cloud computing when algorithms may be running 
on a number of different systems, and possibly controlling critical 
systems.

JN


Peter Dalgaard wrote:
> nashjc at uottawa.ca wrote:
>> This comment is orthogonal to most of the others. It seems that folk 
>> often
>> want to test for equality of "real" numbers. One important one is for
>> convergence tests. When writing my Compact Numerical Methods book I 
>> had to
>> avoid lots of logical tests, but wanted to compare two REALs. I found 
>> that
>> the following approach, possibly considered a trick, is to use an offset
>> and compare
>>
>>     xnew + offset
>>
>> to
>>
>>     xold + offset
>>
>> This works on the examples given to motivate the current thread with an
>> offset of 10, for example.
>>
>> Motivation: Small xold, xnew  compare offset with itself. Large xold and
>> xnew are compared bitwise. Essentially we change from using a 
>> tolerance to
>> using 1/tolerance.
>>
>> Perfect? No. But usable? Yes. And I believe worth keeping in mind for
>> those annoying occasions where one needs to do a comparison but wants to
>> get round the issue of knowing the machine precision etc.
>
> Hmm. Echos of some early battles with R's qbeta() in this. I don't 
> think it can be recommended.
>
> The problem is that you can end up in a situation where xnew=xold-1ulp 
> and xnewnew is xnew+1ulp. I.e. in two iterations you're back at xold.
>
> Even in cases where this provably cannot happen, modern optimizers may 
> make it happen anyway...
>