Fast sparse matrix dot multiplication in GF(256) with Scipy.Sparse

Question

I need to enhance the speed of my algorithm. The method takes the two matrices as argument and performs the dot multiplication. The only issue is that the elements have to be multiplied as octets in GF(256) and then added as XOR. Since I'm dealing with very large sparse matrices the performances are horrible.

def matrix_oct_mult(U, V, OCT_EXP, OCT_LOG):
temp_sum = 0
if shape(U)[1] == None and shape(V)[1] == None:
    for i in range(len(U)):
        temp_sum = oct_sum(temp_sum, oct_mult(U[i], V[i], OCT_EXP, OCT_LOG))
    return temp_sum
assert shape(U)[1] == shape(V)[0], "Wrong size requirements for matrix dot multiplication"
temp_sum = 0
W = sp.lil_matrix((shape(U)[0], shape(V)[1]))
for i in range (shape(U)[0]):
    for z in range(shape(V)[1]):
        for j in range (shape(U)[1]):
             temp_sum = oct_sum(temp_sum, oct_mult(U[i, j], V[j, z], OCT_EXP, OCT_LOG))

        W[i, z] = temp_sum
        temp_sum = 0
return W

As you may see I tried with the lil class, but the performance is still bad.

Is there any efficient why to solve this?

Python is not a compiled language and nested loops for matrix operations will always be exceptionally slow. You should restructure your problem to work with standard matrix operations or you should write this in a compiled language of your choice. — CJR, Nov 09 '20 at 19:11
I don't see where you are taking advantage of the sparsity of `U` or `V` or the result `W`. For iteration like this a normal `ndarray` (or even lists of lists) would be better. Even with `lil` element access of a sparse matrix is slow. — hpaulj, Nov 09 '20 at 19:59

score 0 · Answer 1 · edited Aug 11 '21 at 17:50

Since Python is interpreted, nested for loops notoriously have poor performance. Whereas the equivalent for loops in C would be fast. So best performance will come with compiled code.

For that reason, I wrote a Python library called galois that extends NumPy arrays to operate in Galois fields. I wrote the code in Python but JIT compiled it using Numba so the Galois field arithmetic is as fast, or nearly as fast, as native NumPy array arithmetic, see this performance comparison.

The library supports linear algebra on Galois field matrices using standard binary operators (+, @, etc) and normal NumPy linear algebra functions. Most of those routines are also JIT compiled for speed.

I believe you are trying to do matrix multiplication on two matrices ((M,N) x (N,K)). Here is an example doing that in galois.

In [1]: import galois                                                                                                                                                                          

# Create the GF(2^8) field using AES's irreducible polynomial -- your's may be different
In [2]: GF = galois.GF(2**8, irreducible_poly=0x11b)                                                                                                                                           

In [3]: print(GF.properties)                                                                                                                                                                   
GF(2^8):
  characteristic: 2
  degree: 8
  order: 256
  irreducible_poly: x^8 + x^4 + x^3 + x + 1
  is_primitive_poly: False
  primitive_element: x + 1

In [4]: A = GF.Random((3,4)); A                                                                                                                                                                
Out[4]: 
GF([[ 35, 130, 167, 111],
    [ 58, 161, 194, 200],
    [244,  65, 160,   8]], order=2^8)

In [5]: B = GF.Random((4,5)); B                                                                                                                                                                
Out[5]: 
GF([[ 50,  59,  23,  34,  38],
    [ 16,  59, 162,  80, 250],
    [205, 145, 114,   9,  40],
    [212, 250, 162, 210,  72]], order=2^8)

In [6]: A @ B                                                                                                                                                                                  
Out[6]: 
GF([[144, 236, 142,  90,  89],
    [ 83, 171,  34,   2, 117],
    [192,   1,  20, 208, 127]], order=2^8)

Fast sparse matrix dot multiplication in GF(256) with Scipy.Sparse

1 Answers1