Fres,

thanks for your explanation. However, I am confused with the way I should call the bindings then. Why is it necessary to allocate memory for U, if it is not going to be referenced (as stated in the comments of the fortran code).

Currently, I've 'circumvented' that requirement by calling gesvd as such:

gesvd('O', 'O', 'S', A, s, A, VT);

However, I am not sure if this is safe to do.

--nico

On Sat, Mar 15, 2008 at 4:03 AM, Kresimir Fresl <fresl@grad.hr> wrote:
Nico Galoppo wrote:

> Upon calling the lapack gesvd in numeric/bindings/lapack/gesvd.cpp,
> with JobU = 'O', I'm getting an assertion error on line 418:
>
> 417         assert ((jobu == 'N' && traits::leading_dimension (u) >= 1)
> 418                 || (jobu == 'O' && traits::leading_dimension (u) >= 1)
> 419                 || (jobu == 'A' && traits::leading_dimension (u) >= m)
> 420                 || (jobu == 'S' && traits::leading_dimension (u) >= m));
>
> I don't understand why this assertion is correct. jobu == 'O' means
> that the left-hand vectors should be written to input matrix A, not to
> u, hence, the u input to gesvd should essentially be able to be a null
> pointer (or zero-size matrix).
>
> Am I missing something?

This is the requirement in the "original" LAPACK xgesvd.

DGESVD subroutine is in http://www.netlib.org/lapack/double/dgesvd.f
declared as:

      SUBROUTINE DGESVD( JOBU, JOBVT, M, N, A, LDA, S, U, LDU, VT, LDVT,
     $                   WORK, LWORK, INFO )
*     ..
*     .. Scalar Arguments ..
      CHARACTER          JOBU, JOBVT
      INTEGER            INFO, LDA, LDU, LDVT, LWORK, M, N
*     ..
*     .. Array Arguments ..
      DOUBLE PRECISION   A( LDA, * ), S( * ), U( LDU, * ),
     $                   VT( LDVT, * ), WORK( * )

Description of the arguments says:

*  U       (output) DOUBLE PRECISION array, dimension (LDU,UCOL)
*          (LDU,M) if JOBU = 'A' or (LDU,min(M,N)) if JOBU = 'S'.
*          If JOBU = 'A', U contains the M-by-M orthogonal matrix U;
*          if JOBU = 'S', U contains the first min(m,n) columns of U
*          (the left singular vectors, stored columnwise);
*          if JOBU = 'N' or 'O', U is not referenced.
*
*  LDU     (input) INTEGER
*          The leading dimension of the array U.  LDU >= 1; if
*          JOBU = 'S' or 'A', LDU >= M.

And there is a check later in the code:

      ELSE IF( LDU.LT.1 .OR. ( WNTUAS .AND. LDU.LT.M ) ) THEN
         INFO = -9

Regards,

fres

_______________________________________________
ublas mailing list
ublas@lists.boost.org
http://lists.boost.org/mailman/listinfo.cgi/ublas



--
Nico Galoppo :: http://www.ngaloppo.org