On Dec 3, 2005, at 8:09 PM, John Maddock wrote:
ld = hi + low as above is not among the bit patterns he lists as valid. If you apply the definition of epsilon to the "normal" long doubles only, then you get a value like the one Paul Bristow computed for NTL's quad_float.
I agree: that explanation is quite explicit when it says it follows Kahan's "double double" and the IEEE spec. So the low part must be normalised so that it's bit's "follow on" from those in the high part. As it stands, 1+numeric_limits<>::epsilon() should evaluate to 1, but we really need to check this out. Does anyone of a numerical inclination, want to run some tests on Darwin?
Since i asked the students in my class to write a program calculating the machine epsilon as a homework, I just ran our solution (posted at http://www.itp.phys.ethz.ch/lectures/PT/Exercises/Uebung1/Solution/ eps.cpp). Please overlook that the code is ugly, and does not use any templates - they just had one week of C++ training. Anyway, this code: long double epsilon_long_double() { long double eps = 1.; long double onepluseps; do { eps /= 2; onepluseps = 1. + eps; } while ( onepluseps != 1.); return 2*eps; } returns 0, hence there is really the problem that numeric_limits<long double>::epsilon() == numeric_limits<long double>::min() Matthias