# [R] A slightly unorthodox matrix product.

Berry, Charles ccberry @end|ng |rom uc@d@edu
Sat Aug 4 21:31:31 CEST 2018

```
> On Aug 4, 2018, at 11:43 AM, Jeff Newmiller <jdnewmil using dcn.davis.ca.us> wrote:
>
> Sometimes a good old for loop performs best, even if it doesn't look sexy:
>
>

Fair enough, but a vectorized solution beats them all (see below).

Also,

[SNIP]

> # Charles
> ans1b <- function( a, b )
> {

The lapply you put here was from Eric's solution:

>  xxx <- lapply( seq.int( nrow( A ) )
>               , function( i ) {
>                    A[ i, ] %o% B[ i, ]
>                 }

This is what I had in mind:

ans1b.corrected <- function( a, b ) {
yyy <- sapply( seq.int( nrow( a ) )
, function( i ) a[ i, ] %o% b[ i, ]
, simplify = "array"
)
zzz <- aperm( yyy, c( 3, 1, 2 ) )
zzz
}

On my system it is slower than a for loop but a lot faster than your benchmark showed with the superfluous code from Eric's solution.

For speed, a vectorized solution is faster than a for loop by a factor of 3 on my laptop:

ans0 <- function(A,B){
nca <- ncol(A)
ncb <- ncol(B)
j.index <- rep(1:nca, times = ncb)
k.index <- rep(1:nca, each = ncb)
res <- array(A[, j.index] * B[, k.index], c(nrow(A), nca, ncb))
res
}

> microbenchmark(
+   res0 <- ans0(A, B),
+   res1 <- ans1(A, B),
+   res1b <- ans1b.corrected(A, B),
+   res2 <- ans2(A, B),
+   res3 <- ans3(A, B)
+ )
Unit: microseconds
expr     min       lq      mean   median       uq     max neval   cld
res0 <- ans0(A, B)  13.281  18.4960  21.52723  19.9905  23.4750  61.556   100 a
res1 <- ans1(A, B) 353.121 369.8635 409.77788 381.5840 444.3290 701.256   100     e
res1b <- ans1b.corrected(A, B)  82.816  89.4185 101.52321  95.4275 107.1700 217.357   100   c
res2 <- ans2(A, B)  49.674  54.4825  61.78278  58.7540  65.5265 172.625   100  b
res3 <- ans3(A, B) 317.772 342.4220 392.25065 360.4675 436.2125 602.346   100    d
>

FWIW, if there are dimnames on A and B, sapply( row names(A), ..., simplify="array") preserves them without further ado.

Chuck

```