Multiplication of Matrices composed of polynomials

Question

Would it be possible to use numpy/scipy to multiply matrices composed of polynomials?

Specifically I wish to multiply a 120 by 120 sparse matrix who's entries can look like a+7*b+c by itself.

Honestly, I haven't tried very hard to do this. I see that there is a polynomial module in numpy but I have no experience with it. I am just hoping that someone sees this and says "obviously it's possible, do this".

There is one relevant question asked before from what I've seen: Matrices whose entries are polynomials

Answer 1

I don't know about sparse, but numpy object arrays work fine.

In [1]: from numpy.polynomial import Polynomial as P

In [2]: a = np.array([[P([1,2]), P([3,4])]]*2)

In [3]: a
Out[3]: 
array([[Polynomial([ 1.,  2.], [-1,  1], [-1,  1]),
        Polynomial([ 3.,  4.], [-1,  1], [-1,  1])],
       [Polynomial([ 1.,  2.], [-1,  1], [-1,  1]),
        Polynomial([ 3.,  4.], [-1,  1], [-1,  1])]], dtype=object)

In [4]: np.dot(a, a)
Out[4]: 
array([[Polynomial([  4.,  14.,  12.], [-1.,  1.], [-1.,  1.]),
        Polynomial([ 12.,  34.,  24.], [-1.,  1.], [-1.,  1.])],
       [Polynomial([  4.,  14.,  12.], [-1.,  1.], [-1.,  1.]),
        Polynomial([ 12.,  34.,  24.], [-1.,  1.], [-1.,  1.])]], dtype=object)