> Last night when I thought Eigen couldn't handle complex SVD (based on
> apparently outdated comments in SVD and JacobiSvd.h),
> I held my nose and cra**ed out some C++ code to use LAPACK's CGESDD function
> for complex<float> SVD.
> The lapack routines are all in Fortran with no C/C++ declarations and thus
> have some ugly calling conventions
>
> extern "C" {
> typedef void* fcomplex_ptr;
> void cgesdd_( const char* jobz, const int* m,const int* n, fcomplex_ptr
> a,
> const int* lda, float* s, fcomplex_ptr u,const int* ldu,
> fcomplex_ptr vt,const int* ldvt,
> fcomplex_ptr work, const int* lwork, float* rwork, int* iwork,
> int* info );
> #define cgesdd cgesdd_
> }
> ... lots of initialization ...
> cgesdd( "S", &m, &n, &_a[0], &lda, &_s[0], &_u[0], &ldu, &_vt[0], &ldvt,
> &_work[0], &lwork, &_rwork[0], &_iwork[0], &info );
> ...
>
>
> The Good News:
>
> This morning I saw Gael's note and was overjoyed that I could ditch the ugly
> code in favor of much more elegant Eigen code:
>
> _svd = JacobiSVD<typename
> MatType::PlainObject>(_rows,_cols,ComputeThinU|ComputeThinV);
> ...
> memcpy(_A.data(),inblock.data(),inblock.bytes() );
> _svd.compute(_A);
> _singVals = _svd.singularValues();
>
>
> The Bad News:
>
> The elegant code consistently runs much slower than LAPACK. :(
> E.g. ~5x longer to do a 128x1000 complex<float> SVD (see system details
> below)
>
> Both are computing full set of singular values and the thinned subset of the
> left&right singular vectors -- (This is the "S" argument for cgesdd;
> ComputeThinU|ComputeThinV for JacobiSvd) .
>
> Both are initialized outside the loop and then called several times.
> (calling cgesdd with lwork=-1 to find size, reusing JacobiSvd object )
>
> I double checked that Eigen is using SIMD ( i.e. lots of code in the
> assembly like 'mulps %xmm7,%xmm6')
>
> EIGEN_NO_DEBUG is defined
>
> Is there anything else I could be doing wrong that would cause such a big
> speed gap?
>
> -- Mark
>
> P.S.
> System details:
> ubuntu 11.10 x86_64
> Core i7 CPU Q 820 @ 1.73GHz
> g++ 4.6.1
> eigen 3.0.3
> liblapack-dev 3.3.1-1
>
>
>
>
> On 12/20/2011 02:45 AM, Gael Guennebaud wrote:
>>
>> JacobiSVD does work for complex matrices.
>>
>> gael
>>
>> On Tue, Dec 20, 2011 at 3:08 AM, Mark Borgerding<
mark@xxxxxxxxxxxxxx>
>> wrote:
>>>
>>> Is there an SVD decomposition for complex matrices?
>
>
>