[Top][All Lists]

[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
## Re: [Toon-members] NonSymEigen for TooN, and (SVD full....)

**From**: |
E. Rosten |

**Subject**: |
Re: [Toon-members] NonSymEigen for TooN, and (SVD full....) |

**Date**: |
Thu, 28 Jan 2010 20:50:56 +0000 (GMT) |

**User-agent**: |
Alpine 2.00 (LSU 1167 2008-08-23) |

On Thu, 28 Jan 2010, Gabriel Nützi wrote:

Hello guys

`I have adapted the SymEigen for Non-Symetric Matrix, because there is no
``Complex Support in TooN we have just Imaginary and real parts for the
``EigenValues and Eigenvectors
`

`In principle, TooN does support Vector<std::complex<double> >. I think
``some minor modification will be required, but not much.
`

`One particular problem is that LAPACK seems to want to return the real and
``imaginary parts of the eigenvalues in different arrays. I find this a
``little odd seeing as FORTRAN supports complex types.
`

`To make that work with TooN, we would need a new Allocator type which
``stores two pointers, and upon indexing constructs and returns a
``complex<double>. This would require that Vector also be parameterized on
``the data pointer. Perhaps the underlying allocator type could provide this
``type. Such parameterization would allow a Vector to be constructed off a
``pair<double*,double*> rather than a complex<double>*.
`

`An intermediate solution would be to simply copy the two Vectors in to a
``single Vector<complex>. I can't see that this would be a problem since
``copying is a fast O(N) operation and eigendecomposition is a rather slow
``O(N^3) operation.
`

`Is returning a Vector<complex<double> > better than returning two
``vectors?
`

`Making that change would allow one to be able to decouple the underlying
``datatype from the datatype returned by operator[]. This would allow one to
``take the conjugate transpose of a matrix cleanly.
`
Thoughts?

`--> NonSymEigen.h (it's not finished (I leave the documentation and add ons
``for this class to other guys because I am not the pro here) but would be
``nice to have that in the later version of TooN :-). It has been tested and it
``works.... :-) uses Lapack dgeev_ and sgeev_
`

`We can leave out all the helper functions until someone needs them. Better
``to have missing features than to put in untested code.
`

TestCode :
Matrix<4,4> M;
...
NonSymEigen<> myEig(M);
Matrix<>a= myEig.get_evector();
Vector<>b= myEig.get_evaluesIm();
Vector<>c= myEig.get_evaluesRe();
====================================

`--> The SVD.h file contains an advanced version for full or not full
``decomposition, so we dont have to Pad with zeros....
`

`This seems to be directly supported by LAPACK, so it is probably more
``efficient this way than padding the matrix. If it is not more efficient,
``then the easiest solution is to have a separate FullSVD class which wraps
``SVD and pads the matrix.
`

--> Dont know what you think about that but its handy :-)

`I think it would be better to specify it using a boolean template
``parameter. This would ensure that the returned matrices are statically
``sized. This could be made compatible with existing code by renaming the
``SVD class, and making two wrapper classes, SVD and FullSVD.
`
That aside, I have a couple of qureies:

`Even with full_compute=true, it appears that get_vertical() depends on the
``input matrix. Given that, it appeats that get_U and get_VT are missing
``some if(full_compute) code paths.
`
-Ed

Thanks for the support and maintaining the class :-)
Gabriel

--
(You can't go wrong with psycho-rats.)(http://mi.eng.cam.ac.uk/~er258)
/d{def}def/f{/Times s selectfont}d/s{11}d/r{roll}d f 2/m{moveto}d -1
r 230 350 m 0 1 179{ 1 index show 88 rotate 4 mul 0 rmoveto}for/s 12
d f pop 235 420 translate 0 0 moveto 1 2 scale show showpage