Python Parallel Matrix Vector Multiplication

Question

I am trying to multiply each column of a matrix by a vector element-wise.

I have a serial solution that works correctly.

    for i in range(0,temp.shape[0]):
        for j in range(0,temp.shape[1]):
             temp[i,j] = temp[i,j] * h[i,0]

Below is the parallel solution that works for what i am trying to do, but does not return the same matrix as the above code.

def mv_mult(start_col,end_col,M,V):
for i in range(0,M.shape[0]):
    for j in range(start_col,end_col):
        M[i,j] = M[i,j] * V[i,0]

num_threads = multiprocessing.cpu_count()
threads = []
extra = temp.shape[1] % num_threads
start_col = 0
jump = temp.shape[1] / num_threads
for i in range(0,num_threads):
    print 'starting thread ', i
    if (i <  extra)
        args = start_col, start_col+jump+1,temp,h
        p = multiprocessing.Process(target=mv_mult,args=args)
        p.start()
        threads.append(p)
        start_col += jump+1
    else:
        args = start_col, start_col+jump,temp,h
        p = multiprocessing.Process(target=mv_mult,args=args)
        p.start()
        threads.append(p)
        start_col += jump

for i in threads:
    i.join()

I am relatively new to Python, but as far as I know everything is pass by reference, so the temp matrix passed to each new process is the same object, so it should work the same as the serial solution, minus the fact that it is split up by column.

Any ideas as to why it is not working?

Best by far to do this with numpy, and make sure your numpy is linked against a threaded BLAS (Atlas, OpenBLAS). You don't want to be writing your own linear algebra routines. — Jonathan Dursi

Jonathan Dursi Jonathan Dursi · Accepted Answer · 2014-07-30T21:18:38

This is just the numpy multiply operation:

import numpy

a = numpy.matrix(numpy.ones((5,5)))
b = numpy.arange(5).reshape(5,1)
print a
print b    
#[[ 1.  1.  1.  1.  1.]
# [ 1.  1.  1.  1.  1.]
# [ 1.  1.  1.  1.  1.]
# [ 1.  1.  1.  1.  1.]
# [ 1.  1.  1.  1.  1.]]
#[[0]
# [1]
# [2]
# [3]
# [4]]

c = numpy.multiply(a,b)
print c
#[[ 0.  0.  0.  0.  0.]
# [ 1.  1.  1.  1.  1.]
# [ 2.  2.  2.  2.  2.]
# [ 3.  3.  3.  3.  3.]
# [ 4.  4.  4.  4.  4.]]

(Note that in your serial version you're multiplying rows, not columns).

Comparing the time of this to your serial version:

#!/usr/bin/env python

import numpy
import numpy.random
import numpy.linalg
import sys
import time

def init(n):
    a = numpy.matrix(numpy.random.rand(n,n))
    b = numpy.random.rand(n,1)
    b = b.reshape((n,1))
    return a,b

def np_mult(a, b):
    c = numpy.multiply(a,b)
    return c

def manual_mult(a,b):
    c = numpy.matrix(numpy.zeros_like(a))
    for i in range(0,a.shape[0]):
        for j in range(0,a.shape[1]):
             c[i,j] = a[i,j] * b[i]
    return c

def main(argv=None):
    if argv is None:
        argv = sys.argv

    n = int(argv[1])
    niters = int(argv[2])
    a,b = init(n)

    start = time.time()
    for i in xrange(niters):
        cnp = np_mult(a,b)
    end = time.time()

    nptime = (end-start)/niters
    print 'Numpy: ', nptime

    start = time.time()
    for i in xrange(niters):
        cm = manual_mult(a,b)
    end = time.time()

    manualtime = (end-start)/niters
    print 'Manual: ', manualtime
    print 'Speed difference = ', manualtime/nptime
    print 'Diff = ', numpy.linalg.norm(cm - cnp)/numpy.linalg.norm(cnp)

    return 0

if __name__ == "__main__":
    sys.exit(main())

Running with a few small sizes gives:

$ ./mult.py 50 5
Numpy:  3.48091125488e-05
mManual:  0.021101808548
Speed difference =  606.215068493
Diff =  0.0

$ ./mult.py 100 5
Numpy:  6.59942626953e-05
Manual:  0.0830503940582
Speed difference =  1258.44869942
Diff =  0.0

$ ./mult.py 500 5
Numpy:  0.000931406021118
Manual:  2.08903641701
Speed difference =  2242.88481032
Diff =  0.0

with speed differences completely dwarfing anything you could make up for with parallelism.

Python Parallel Matrix Vector Multiplication

3 Answers