|
|
Ublas : |
From: Peter Melchior (pmelchior_at_[hidden])
Date: 2006-02-08 08:47:43
Hello,
On Thu, 2006-01-26 at 15:36 +0100, Karl Meerbergen wrote:
> All I can say is that the bindings are written for column major matrices.
> Your matrix is row_major. I do not believe it is sufficient to transpose it,
> since for a row_major matrix the leading dimension corresponds to the number
> of columns. Lapack assumes it corresponds to the number of rows.
It should be sufficient. At least this is sufficient when you call fortran routine with
a ordinary C 2D-array. For sure, the information about the ublas matrix dimensions is updated
when transposing it.
The LAPACK routine DGESVD takes both dimensions of the matrix as
arguments. They should be passed correctly from the bindings.
So, if I transpose the matrix, and pass LAPACK the correct values for N
and M, where is the problem?
Peter
> You could pass trans(M) as an argument to lapack. trans(M) creates an object
> of type matrix_unary2<M,scalar_identity<typename M::value_type> >. There are
> no bindings for this kind of object. These should be written.
> This would be the clean way of doing it.
>
> Karl
>
>
> On Thursday 26 January 2006 15:17, Peter Melchior wrote:
> > Hello,
> >
> > I want to do a SVD for non-square real matrices.
> > When I have a Matrix with more columns than rows, I get this runtime
> > error:
> > ** On entry to DGESVD parameter number 6 had an illegal value
> >
> > When I have a matrix with more rows than columns I get a wrong result
> > (checked against Mathematica and IDL).
> > But when using a square matrix, everything works fine.
> > I'm definitely sure, that the U, S and V matrix have the correct sizes.
> >
> > This is the code I use:
> >
> > template <class T>
> > inline void svd(boost::numeric::ublas::matrix<T>& M,
> > boost::numeric::ublas::matrix<T>& U, boost::numeric::ublas::matrix<T>&
> > S, boost::numeric::ublas::matrix<T>& Vt) {
> > // dimensions of matrices
> > int n = matrix.size1();
> > int m = matrix.size2();
> > int length = min(n,m);
> >
> > // since LAPACK works on column-major matrices, transpose it.
> > boost::numeric::ublas::matrix<T> transpose =
> > boost::numeric::ublas::trans(M);
> > boost::numeric::ublas::vector<T> s(length);
> > boost::numeric::bindings::lapack::gesvd('A','A',transpose,s,U,Vt);
> >
> > // transpose U and Vt
> > U = boost::numeric::ublas::trans(U);
> > Vt = boost::numeric::ublas::trans(Vt);
> >
> > // addressing the diagonal elements of S to store singular values
> >
> > boost::numeric::ublas::matrix_vector_range<boost::numeric::ublas::matrix<T>
> > > (S, boost::numeric::ublas::range (0, length),
> > boost::numeric::ublas::range (0, length)) = s; }
> >
> >
> >
> > Does anybody have a clue what is going on here?
> >
> >
> > Best regards,
>