Re: complex number policy for beast

[Stefan Westerfeld <stefan space twc de>]

   Hi!

On Sat, Nov 04, 2006 at 06:41:58PM +0100, Tim Janik wrote:
> On Sat, 4 Nov 2006, Tim Janik wrote:
> 
> >- C++ source files can use BseComplex but may opt to include <complex>
> >  and use that. (they have to convert back-and-forth between std::complex
> >  and BseComplex for API though).
> 
> after reading some more C++ headers, i have to add that using C++ complex 
> numbers can not be recommended from an accuracy point of view, from
> /usr/include/c++/3.4/complex:
> 
>   // 26.2.5/13
>   // XXX: This is a grammar school implementation.
>   template<typename _Tp>
>     template<typename _Up>
>     complex<_Tp>&
>     complex<_Tp>::operator*=(const complex<_Up>& __z)
>     {
>       const _Tp __r = _M_real * __z.real() - _M_imag * __z.imag();
>       _M_imag = _M_real * __z.imag() + _M_imag * __z.real();
>       _M_real = __r;
>       return *this;
>     }
> 
> this has many cancellation/etc problems, described e.g. in Numerical Recipes
> and is usually worked around by in math packages, or e.g. the _Cdiv() from 
> the
> C9x specification.
> also our bse_complex_div() from bsemath.[hc]* has none of these problems
> (similar things hold for our other complex number functions), so if accuracy
> matters at all, we better use BseComplex over std::complex<double>.

This problem has been brought up on the gcc mailing lists, here is a
link:

   http://gcc.gnu.org/ml/gcc/2003-04/msg00814.html

The relevant information is that while this is the generic
implementation, there are template specializations for the
std::complex<float>, std::complex<double> and std::complex<long double>
case, which use compiler builtin __complex__ types / operations. So
as I understand it there should be not more numeric problems than in
gccs C9x implementation.

   Cu... Stefan
-- 
Stefan Westerfeld, Hamburg/Germany, http://space.twc.de/~stefan