Parallelization of element-wise matrix multiplication

Question

I use Armadillo as an interface to OpenBLAS. In my current program, I have a loop, in which I do multiplications of the form

for(long t = t0; t < t1; t+=tStep)
{
    stateMatrix %= elementWiseEvolutionMatrix;
}

The operator % is an element-wise multiplication operator. The problem here is that for matrices of side length of 500+ (the ones I have at hand), I can see that there is no parallelization whatsoever. Now I would like to point out that normal matrix multiplication is parallelized. But this kind of element-wise multiplication is not parallelized.

How do I know that? Because I go to htop in my linux system, and I see that only one core is busy, while if I do the same with normal matrix multiplication, I see that all cores go busy.

Now I tried to manually parallelize this with OpenMP, but no luck. I tried:

for(long t = t0; t < t1; t+=tStep)
{
    #pragma omp parallel for
    for(long i = 0; i < static_cast<long>(stateMatrix.n_rows); i++)
    {
        stat1eMatrix.row(i) %= elementWiseEvolutionMatrix.row(i);
    }
}

But this got all the cores busy, but the program became about a factor of 10 slower.

My question: How can I get the element-wise multiplication to be as fast as possible with parallelization?

Thanks.

EDIT: I would like to point out that I'm more than happy to use another library for the element-wise multiplication, if necessary.

Answer 1

There's no element wise multiplication operation in the BLAS library. Your best approach is probably to just implement the operation yourself using (e.g.) OpenMP threading.

Before you do this, you should consider Amdahl's law and whether speeding up this bit of your code is really going to help- chances are that these elementwise multiplications are not where your code is spending most of its time, and as a result you probably won't see much speedup from parallelizing this part of your code.

Answer 2

First off, I agree with Brian Borchers comments about profiling to make sure these element-wise multiplications are where your performance issue lies. However, since you are convinced that is your problem, here is another suggestion.

Before trying to exploit multiple CPU, I would make sure that you have an implementation that exploits vectorization. The SSE2 instruction set (available in most modern processors) has an operation to multiply vectors of double precision floating point numbers. Your code may or may not allow your compiler to exploit this instruction.

As far as I know, Armadillo does not have any direct support for SSE2. But since you indicated a willingness to switch libraries, the Eigen library (http://eigen.tuxfamily.org/index.php?title=Main_Page) definitely does generate code using the SSE2 instructions. This could give you a 4x improvement for single precision multiplies and 2x improvement for double precision. If you are fortunate enough to have a CPU that supports the AVX instruction set, the development version of Eigen supports this to provide additional speedup.