From simon.urbanek at r-project.org  Tue Aug 24 04:03:20 2010
From: simon.urbanek at r-project.org (Simon Urbanek)
Date: Mon, 23 Aug 2010 22:03:20 -0400
Subject: [Rd] Speeding up matrix multiplies
In-Reply-To: <20100823233927.GA21529@cs.toronto.edu>
References: <20100823233927.GA21529@cs.toronto.edu>
Message-ID: <95C122CC-D17A-4007-ABDF-C77841E4000C@r-project.org>


On Aug 23, 2010, at 7:39 PM, Radford Neal wrote:

> Regarding my previous message on speeding up matrix multiplies, I've
> realized that the size of the result matrix is not really the right
> criterion for deciding to do the computation without using the Fortran
> routine.  A better criterion would be based on the ratio of the time
> for the matrix multiply to the time for the ISNAN check.  So I've
> replaced my previous check by
> 
>        /* Faster to just do it here if the time to do the matrix multiply
>           is not much greater than the time for the NA/NaN check below. */
> 
>        do_it_here = nrx*ncy <= 5*(nrx+ncy);
> 
> In particular, all matrix x vector and vector x matrix products will 
> in this version be done in the matprod routine, not the Fortran routine.
> And this is the right thing to do, since the time for the ISNAN check
> before calling the Fortan routine will be comparable to the time for
> the matrix multiply.  So avoiding it will be desirable unless the Fortran
> routine is really, really faster than the C code in matprod.
> 

It is, many times in fact, if you use threaded BLAS on a multi-core machine and large matrices.


> Of course, the real solution is to change the Fortran routine (which
> seems to be in src/extra/blas/blas.f, and isn't all that complicated)

Nope - it can be external and BLAS standard doesn't handle NAs.

Cheers,
Simon



> to handle NAs and NaNs as desired.  Or to write a carefully optimized
> C matrix multiply routine that does the right thing.  For this,
> knowing what the problem is supposed to be would help.  My quick
> attempt at using the bug reporting system didn't manage to find the
> PR#4582 that is mentioned in the comment in the code, but I don't
> really know how the bug report system is supposed to work.
> 
>    Radford Neal
> 
> ______________________________________________
> R-devel at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-devel
> 
> 



[Rd] Speeding up matrix multiplies
Simon Urbanek