Re: [eigen] Eigen is slow with complex matrices

[Rohit Garg <rpg.314@xxxxxxxxx>]

AFAIR, some work was done to add complex number support, compatible
with std::complex. It is sitting in the FFT fork. What happened to it?
Also, when can we see the FFT work being merged back into eigen
mainline?

On Mon, Oct 19, 2009 at 7:11 PM, Benoit Jacob <jacob.benoit.1@xxxxxxxxx> wrote:
> Hi,
>
> 2009/10/19 Robert Lupton the Good <rhl@xxxxxxxxxxxxxxxxxxx>:
>> We're using eigen on a large astronomical project;  it's great -- no need to
>> find the proper combination du jour of optimised blas and lapack.
>>
>> Unfortunately, in the context of complex matrices we are having to
>> reconsider the use of BLAS. The reason appears to be in the complex class
>> itself (see below)  Is there any way that Eigen can work around this?
>
> There are several ways in which Eigen is currently not optimally fast
> with complex numbers, some are Eigen's fault, some are the GNU STL's
> fault as you mention. Perhaps in the particular case of your project,
> it reduces to one precise reason;
>
>>> A better implementation would be to have op* directly calculate both
>>> components of the result, and return it as a new complex object:
>>> {
>>>   return complex<_Tp>(
>>>       _x.real() * _y.real() - _x.imag() * _y.imag(),
>>>       _x.real() * _y.imag() + _x.imag() * _y.real());
>>> }
>
> Have you tried it? Does it fix your performance issues? If yes, you
> have a pretty clear path from here: just copy this standard complex
> implementation as a new class, say "ourcomplex" and make that change
> there; Then add support for "ourcomplex" as explained here,
> http://eigen.tuxfamily.org/dox/CustomizingEigen.html#CustomScalarType
>
> Otherwise, here are some thoughts on the performance issues with complex.
>
> A. Not our fault
> 0) What you mention --- I didn't know about it.
> 1) The basic operations are not vectorized by GCC as they should. For
> example, an addition of two complex<double> should be trivial to
> vectorize.
> 2) Because of that, we haven't added vectorization in Eigen (though
> for complex<float> we should be able to do something)
> 3) So we've been considering adding our own complex class, see the
> TODO. Needs someone with time to do it.
>
> B. Our fault
> 4) Eigen should be modified to avoid taking abs() of complex numbers
> too often, because that means a sqrt(). The idea would be to introduce
> ei_abs_or_abs2 that does abs for real and abs-squared for complex, so
> that no sqrt() is needed.
>
> To check if 4) is affecting you, edit Core/MathFunctions.h and put
> EIGEN_DONT_INLINE in front of the ei_abs() specializations for
> complex. Then profile your program, is ei_abs() accounting for a large
> part of the execution time?
>
> Cheers,
> Benoit
>
>
>



-- 
Rohit Garg

http://rpg-314.blogspot.com/

Senior Undergraduate
Department of Physics
Indian Institute of Technology
Bombay