Andreas Klein <[EMAIL PROTECTED]> writes:
| Hello
|
| I have found a bug in the implementation of the complex library of
| g++ and the complex.h library of the gcc (c99 support).
|
| The simplest program that demonstrates the bug is:
|
| --------------------------------
|
| #include<iostream>
| #include<complex>
|
| using namespace std;
|
| int main() {
| complex<double> a = 1.0/0.0; // + infinity
| complex<double> b = 1.0;
| complex<double> c = b/a; // should be 0 but is (NaN,NaN)
|
| cout << a << endl;
| cout << b << endl;
| cout << c << endl;
| }
|
| --------------------------------
|
| Since all values are real one should expect the result 0
| (the IEEE floating-point value for 1/infinity). More serious is that
| if a contains a large but representable value for which |a|^2 is not
| representable the implementation will yields bad results. For example
| let a be the largest representable number and b be a/2 then b/a should
| be 0.5 and not NaN.
|
| The bug comes from the naive implementation
|
| a + b I ac + bd bc - ad
| ------- = ------- + ------- I
| c + d I c^2+d^2 c^2+d^2
|
| We get unnecessary overflow errors if c^2+d^2 is too large.
|
As a matter of fact, the implementation of <complex> is criticized,
once in a while, because it does NOT use the grammar school rule you
present above. However, for float, double, long double it specializes
to __complex__ T which is what the compiler uses to implement complex
numbers in C99 and Fortran. So, the problem is a compiler problem not
libstdc++ problem.
[...]
| So my suggestion is to add template specialization for the
| operator/= for the types complex<float> complex<double> and
| complex<long double> that use Smith formula. (For integer types I
| think we should stay with the old implementation.)
Did you look at the actual implementation?
-- Gaby
