Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

How to exploit symmetry in outer product in Numpy (or other Python solutions)?

Tags:

Suppose that we want to compute the outer product of a vector with itself:

import numpy as np

a = np.asarray([1, 1.5, 2, 2.5, 3])
A = np.outer(a, a)

print(A)

which results in:

[[ 1.    1.5   2.    2.5   3.  ]
 [ 1.5   2.25  3.    3.75  4.5 ]
 [ 2.    3.    4.    5.    6.  ]
 [ 2.5   3.75  5.    6.25  7.5 ]
 [ 3.    4.5   6.    7.5   9.  ]]

This results in a symmetric matrix. Prior knowledge of the fact that the two vectors in the outer product are the same could be exploited by only computing one triangle of the matrix, and filling in the remaining entries from the corresponding entries in the triangle.


Question: are there any easy ways to exploit such knowledge in numpy (or other solutions in Python)? Of course, it wouldn't be too difficult to write a custom solution in Python, but if that comes at the cost of not using a solid BLAS it would be unlikely to be worth it.

I am primarily concerned with computation time, not so much necessarily with RAM usage.

like image 518
Dennis Soemers Avatar asked Jul 23 '18 13:07

Dennis Soemers


1 Answers

If you can deal with the result only being valid in the lower triangle, this Numba solution is the best I've come up with:

import numba

@numba.njit                 
def outer(arr):
    n = len(arr)
    res = np.empty((n,n), arr.dtype)
    for ii in range(n):
        for jj in range(ii+1):
            res[ii,jj] = arr[ii] * arr[jj]
    return res

It is twice as fast as np.outer() for large vectors, and five times as fast for small ones. If you need the fully populated solution you can set res[jj,ii] = res[ii,jj] in the inner loop and it is still somewhat faster than np.outer().

For some reason, np.multiply.outer() is faster for small vectors than np.outer() (and no slower for large vectors).

like image 75
John Zwinck Avatar answered Sep 28 '22 17:09

John Zwinck