[Rd] Suggest new outer for R-1.3

Jonathan Rougier J.C.Rougier@durham.ac.uk
Wed, 21 Mar 2001 11:16:10 +0000 (GMT)


Hi everyone,

Can I suggest the following modification of outer for R-1.3, in the
interests of speed and size of calculation:

*******************************************************************

"outer" <-
function (X, Y, FUN = "*", ...) 
{
    no.nx <- is.null(nx <- dimnames(X <- as.array(X)))
    dX <- dim(X)
    no.ny <- is.null(ny <- dimnames(Y <- as.array(Y)))
    dY <- dim(Y)
    if (is.character(FUN) && FUN=="*") {
      robj <- as.vector(X) %*% t(as.vector(Y))
      dim(robj) <- c(dX, dY)
    } else {
      FUN <- match.fun(FUN)
      Y <- rep(Y, rep(length(X), length(Y)))
      X <- rep(X, length.out = length(Y))
      robj <- array(FUN(X, Y, ...), c(dX, dY))
    }
    ## no dimnames if both don't have ..
    if (no.nx) 
        nx <- vector("list", length(dX))
    else if (no.ny) 
        ny <- vector("list", length(dY))
    if (!(no.nx && no.ny)) 
        dimnames(robj) <- c(nx, ny)
    robj
}

"%o%" <- .Alias(outer)

*******************************************************************

In other words, branch in the default case to a simple matrix outer
product.  For R-1.1 --vanilla on sparc-sun-solaris2.7 I have

> system.time(outer(1:1000, 1:300, "*"))
Error: heap memory (6144 Kb) exhausted [needed 1171 Kb more]
       See "help(Memory)" on how to increase the heap size.
Timing stopped at: 0.49 0 0.49 0 0 
> system.time(new.outer(1:1000, 1:300, "*"))
[1] 0.38 0.00 0.38 0.00 0.00
> system.time(outer(1:500, 1:300, "*"))
[1] 0.47 0.00 0.47 0.00 0.00
> system.time(new.outer(1:500, 1:300, "*"))
[1] 0.25 0.00 0.25 0.00 0.00

i.e. a saving of about 50% in time.  There will be a similar speed-up in
"kronecker".

Cheers, Jonathan.

Jonathan Rougier                       Science Laboratories
Department of Mathematical Sciences    South Road
University of Durham                   Durham DH1 3LE
tel: +44 (0)191 374 2361, fax: +44 (0)191 374 7388
http://www.maths.dur.ac.uk/stats/people/jcr/jcr.html

-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
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
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._