[Rd] DGESDD from Lapack for R-1.4.0?

Prof Brian D Ripley ripley@stats.ox.ac.uk
Sat, 17 Nov 2001 08:44:18 +0000 (GMT)


On Fri, 16 Nov 2001, Keith Richards-Dinger wrote:

>
> Hi,
>
> I'm just wondering if it is planned to include the Lapack
> routine DGESDD (and friends) in R-1.4.0?  This is faster
> (supposedly by a factor of ~6 for large matrices) than
> DGESVD which is currently (R-1.3.1) called by La.svd.

No plans.

> And if it is not in the plans yet, is there a chance it
> could be?  I've added it to my local version of R-1.3.1 and

If you want to contribute tested code to do this it could be incorporated,
but the timescale for 1.4.0 is short: if that is to be released before
Christmas feature-freeze is imminent (but not yet announced).

We don't see R as a linear algebra system, and the aim of introducing
La.svd was accuracy not speed (the LINPACK version of svd having some
algorithmic limitations, as well as some problems on machines with
extra-precision FPU registers).

> so far see a factor of 4 improvement over La.svd and a
> factor of 3 over svd (see output below).  I'd be glad to
> hand off what I've done so far to someone else or to work
> more on it (make it behave as svd and La.svd do if you don't
> want all the singular vectors; and add ZGESDD also).
> Though what I've done so far would probably be trivial for
> someone who knew what they were doing...
>
> Mac G4, 450 MHz, LinuxPPC, pre-compiled ATLAS-3.2.1
> libraries for G4 from www.netlib.org, R-1.3.1:
>
> > m <- matrix(rnorm(25e4), 5e2, 5e2)
> > system.time(svd(m))
> [1] 22.31  0.10 23.01  0.00  0.00
> > system.time(La.svd(m))
> [1] 29.60  0.10 29.78  0.00  0.00
> > system.time(La.sdd(m))
> [1] 7.28 0.12 8.00 0.00 0.00

I'd want to see a wider comparison, including on non-ATLAS systems.
In particular, without ATLAS my 1GHZ PIII gives

> m <- matrix(rnorm(25e4), 5e2, 5e2)
> system.time(svd(m))
[1] 22.76  0.06 22.80  0.00  0.00
> system.time(La.svd(m))
[1] 16.32  0.15 16.95  0.00  0.00

and on Doug Bates' franz (with ATLAS)

> system.time(svd(m))
[1] 14.70  0.03 17.96  0.00  0.00
> system.time(La.svd(m))
[1] 15.00  0.03 17.71  0.00  0.00

which suggests something peculiar on your system, the only one I have ever
seen that has La.svd appreciably slower than svd.  We've had a couple of
instances where changing to use BLAS slowed down the code on the basic
system, so we know that the precise BLAS is important.


-- 
Brian D. Ripley,                  ripley@stats.ox.ac.uk
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford,             Tel:  +44 1865 272861 (self)
1 South Parks Road,                     +44 1865 272860 (secr)
Oxford OX1 3TG, UK                Fax:  +44 1865 272595

-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-devel mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !)  To: r-devel-request@stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._