Re: tensor() function and sets

About this list Date view Thread view Subject view Author view Other groups

Subject: Re: tensor() function and sets
From: Martin Maechler (
Date: Thu 22 Jul 1999 - 02:44:18 EST

Message-Id: <>

>>>>> On Tue, 20 Jul 1999 15:26:35 +0100 (BST), Jonathan Rougier <> said:

JonR> Hi Everyone,

JonR> To complete the outer() and kronecker() functions in the base, may I
JonR> suggest the following tensor() function, which allows the multiplication
JonR> of arrays through sets of conformable dimensions. I am happy to write a
JonR> help page if required.

JonR> The code also needs a setdiff() function which prompts me to ask:
JonR> what about simple set functions? I expect many of us have written our own
JonR> (Brian has a setdiff() in drop1.lm(), for example), which seems like a
JonR> good reason for putting versions in the base. I would be happy to provide
JonR> mine for general scrutiny.

Actually your version of "setdiff" is in

    "%w/o%" <- function(x,y) x[!x %in% y] #-- x without y

Maybe we should make this part of core R, since it allows the nicer syntax
 (1:10) %w/o% c(3,7,12).

One problem with all these might be that you forget to think about efficiency
{ %w/o% -> %in% -> match(.) }
but that might be ok.

are two other set functions, the first of which is also contained in
example(match) { == example("%in%") }:
    intersect <- function(x, y) y[match(x, y, nomatch = 0)]

What do others think? Should (some of) these be put into "core R" ?

JonR> Jonathan Rougier Science Laboratories
JonR> Department of Mathematical Sciences South Road
JonR> University of Durham Durham DH1 3LE

JonR> "setdiff" <-
JonR> function (x, y) x[!(x %in% y)]

JonR> "tensor" <-
JonR> function (A, B, da, db)
JonR> {
JonR> # tensor product of A and B through da and db
JonR> <- is.null(na <- dimnames(A <- as.array(A)))
JonR> dima <- dim(A)
JonR> no.nb <- is.null(nb <- dimnames(B <- as.array(B)))
JonR> dimb <- dim(B)
JonR> if (any(dima[da] != dimb[db]))
JonR> stop("Mismatched dimensions")
JonR> kpa <- setdiff(seq(along = dima), da)
JonR> kpb <- setdiff(seq(along = dimb), db)
JonR> # fix up the dimnames (see outer)
JonR> if ( && no.nb)
JonR> nms <- NULL
JonR> else {
JonR> if (
JonR> na <- vector("list", length(dima))
JonR> else if (no.nb)
JonR> nb <- vector("list", length(dimb))
JonR> nms <- c(na[kpa], nb[kpb])
JonR> }
JonR> A <- matrix(aperm(A, c(kpa, da)), ncol = prod(dima[da]))
JonR> B <- matrix(aperm(B, c(db, kpb)), nrow = prod(dimb[db]))
JonR> array(A %*% B, c(dima[kpa], dimb[kpb]), dimnames = nms)
JonR> }

JonR> # examples

JonR> A <- array(1:6, 2:3, list(LETTERS[1:2], LETTERS[3:5]))
JonR> B <- 1:3
JonR> tensor(A, B, 2, 1) # same as drop(A %*% B)
JonR> A <- outer(A, B) # A now 2 by 3 by 3
JonR> tensor(A, B, 2, 1); tensor(A, B, 3, 1) # both 2 by 3 (nb dimnames)
JonR> B <- outer(1:2, B) # B now 2 by 3
JonR> tensor(A, B, c(1, 3), 1:2) # must be 3 vector
JonR> B <- outer(1:3, B) # B now 3 by 2 by 3
JonR> tensor(A, B, 1:3, c(2, 1, 3)) # must be a scalar
JonR> sum(A * aperm(B, c(2, 1, 3))) # same scalar value

Looks interesting, however is it really useful
{for doing statistics, graphics, numerical computing, I mean} ?

Yes, if you'd provide a *.Rd help file (use ``prompt(tensor)''),
we'd be even more inclined to include the tensor function....
Also, I'd like to see a "real" use for it.

        tensor(A, B, 1:3, c(2, 1, 3))
really easier than
        sum(A * aperm(B, c(2, 1, 3)))

[probably yes, you have to think less, but still:
 couldn't "da" and "db" (3rd and 4th arg) get reasonable default values ?]

r-devel mailing list -- Read
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !)  To:

About this list Date view Thread view Subject view Author view Other groups

This archive was generated by hypermail 2b25 : Tue 04 Jan 2000 - 14:16:06 EST