[Insight-users] Code integration and code design - Question

[Steve M. Robbins]

Hi Alex,

Apologies for the big delay in responding -- real life got in the way (birth 
of second daughter).

So my memory of Shewchuk's code -- and this goes back to 2001, so it may be 
faulty -- is that the code requires the floating point calculations to strictly 
adhere to IEEE arithmetic in round-to-even mode.  I ran into two problems: the 
x86 machines extended precision, and the rounding mode.

The x86 machines used extended precision by default, which I had to disable in 
my code prior to using Shechuk's predicates.  Disabling extended precision 
used platform-dependent system calls; I don't recall all the details and 
eventually I used the GNU Scientific Library (GSL) method gsl_ieee_set_mode() 
to hide the system variations.  The other issue with this is that if you 
change the IEEE mode, other system functions -- e.g. math library functions 
like sin(), cos(), exp() -- may fail to work because they were coded assuming 
extended precision.  Or worse: calling sin() will turn on extended precision, 
do the calculation and leave the system in extended precision.   You might 
consider wrapping your predicates in code that saves the floating point mode, 
sets it for the calculation, then restores to the initial value.  However, I 
seem to recall that setting the IEEE mode took a non-trivial amount of time.

The second problem is similar: some machines I was using did not use "round-
to-even" as the default.  I was using Pentium III running linux and SGI MIPS 
machines (R10k, R12k) running IRIX.  I think it may have been the latter.  In 
fact, I think the MIPS compiler by default did not obey all IEEE rules and I 
had to use special compiler options to force strict IEEE compliance.  A 
portable library would have to ensure the correct rounding mode and possibly 
restore it after the computation was done.

A possible third issue involves the IEEE fault handling.  My initialization 
code disables faults for denormalized numbers and for underflow.  
Unfortunately, I can't recall for sure whether this was an issue for the 
Shewchuk code or for something else.

And finally, if necessary, one has to be sure the code is built with options to 
ensure strict IEEE compliance.

As noted, my experience was in 2001 and I have no idea whether any of this is 
an issue with today's hardware.



On November 13, 2012 10:42:47 PM Alexandre GOUAILLARD wrote:
> Dear steve,
> 
> That's a very good point.
> 
> Yes, we saw the warning. We have implemented a very simple test for the
> robustness, which actually test both normal ITK, and the new code with 4
> points lying on circle, and adding an arbitrary small number to the 4th
> point coordinate, and we could not find a platform where this was a
> problem. We tried linux with gcc 4.2, mac with both gcc and clang, and
> windows (xp, 7) with MSVC 2005, 2008 and 2010.
> 
> Now, that does not mean there is no problem. 

Agreed: I would not find this convincing enough to bet my exact computation 
results on it :-)

Regards,
-Steve