[Rd] "LAPACK routine DGESDD gave error code -12" with Debian (PR#2822)

[camm at enhanced.com]

Greetings!

Douglas Bates <bates at stat.wisc.edu> writes:

> Peter Dalgaard BSA <p.dalgaard at biostat.ku.dk> writes:
> 
> > Camm Maguire <camm at enhanced.com> writes:
> > 
> > > 2) Given the error code, and the scaling behavior with matrix size,
> > >    I'd say the lwork parameter (size of the work array) passed to
> > >    dgesdd is not always large enough, i.e. is not scaling properly
> > >    with n,m.  Please see 'man dgesdd' for interpretations of the error
> > >    code.  It is the responsibility of the calling routine to allocate
> > >    and pass the work array to dgesdd.  With most lapack routines, one
> > >    can make a 'workspace query' call first by setting lwork to -1, or
> > >    some such.  check the man page for details.  This of course would
> > >    have to be done with each change in n,m.  Alternatively, you could
> > >    take the minimum workspace requirements from the manpage.  
> > 
> > Right, but that's actually what we do, use the workspace query. It's
> > all very weird, because the -12 value indicates that the lwork
> > parameter is wrong, but it is computed from an exactly identical call,
> > except lwork=-l:
> > 
> >         lwork = -1;
> > 
> >         F77_CALL(dgesdd)(CHAR(STRING_ELT(jobu, 0)),
> >                          &n, &p, xvals, &n, REAL(s),
> >                          REAL(u), &ldu,
> >                          REAL(v), &ldvt,
> >                          &tmp, &lwork, iwork, &info);
> >         lwork = (int) tmp;
> 
> It looks like a problem in ilaenv, the Lapack routine that returns the
> tuning parameters, like the optimal temporary storage size, for
> various Lapack routines.  The value returned in tmp will depend upon
> the results of several calls to ilaenv.  These results can vary between
> different implementations of the blas (or atlas).
> 

Just a note that you should check the returned info value on the
workspace call to make sure tmp has been filled in.  (the man page
says this is done only if the other parameters are valid.)

Does this problem vary with blas, and if so, how?  You can run under
whatever blas you want, including reference, via the LD_LIBRARY_PATH
variable.  Asuming you've installed all the i386 atlas versions and the blas
package:

LD_LIBRARY_PATH     used blas:

/usr/lib                reference
/usr/lib/atlas          base vanilla i386 atlas
/usr/lib/{sse,sse2,3dnow}/atlas    atlas with ISA extensions


ilaenv might possibly be an issue, but only realistically if the
problem blas is coming from the 3dnow atlas.  When I put together the
lapack package, I read in the lapack notes how many of the ilaenv
constants can be hardwired, saving a certain amount of time on the
first call.  I chose not to do this only because of the existence of
the 3dnow, non-ieee compliant blas option at runtime, as one of the
parameters pertains directly to ieee.  So ilaenv is calculating its
parameters n first call at runtime on Debian.  I'm dubious as to the
relevance of this, though, as this should only kick in for single
precision if at all, and should not affect integer values in any case.

Take care,

> 
> 
> 

-- 
Camm Maguire			     			camm at enhanced.com
==========================================================================
"The earth is but one country, and mankind its citizens."  --  Baha'u'llah