CODEWITHSUNDEEP
pythonnumpylinear-algebra
user2654096
user2654096

Reputation: 71

python-Interpolate polynomial where coefficients are matrices

I have a polynomial of the form:

p(y) = A + By + Cy^2 ... + Dy^n

Here, each of the coefficients A,B,..,D are matrices (and therefore p(y) is also a matrix). Say I interpolate the polynomial at n+1 points. I should now be able to solve this system. I'm trying to do this in Numpy. I have the following code right now:

a = np.vander([0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1.0,1.1,1.2]) #polynomial degree is 12, a -> (12x12)
b = np.random.rand(12,60,60) #p(x) is a 60x60 matrix that I have evaluated at 12 points
x = np.linalg.solve(a,b)

I get the following error:

ValueError: solve: Input operand 1 has a mismatch in its core dimension 0, with gufunc signature (m,m),(m,n)->(m,n) (size 60 is different from 12)

How can I solve this system in Numpy to get x? Is there a general mathematical trick to this?

Upvotes: 1

Views: 390

Answers (1)

Daniel F
Daniel F

Reputation: 14399

Essentially, you're just doing 3600 12d polynomial regressions and composing the coefficients into matrices. For instance, the component p(y)[0,0] is just:

p(y)[0, 0] = A[0, 0] + B[0, 0] * y + C[0, 0] * y**2 ... + D[0, 0] * y**n

The problem is that np.linalg.solve can only take one dimension of coefficients. But since your matrix elements are all independent (y is scalar), you can ravel them and you can do the calulation with the form (m,m),(m,n**2) -> (m,n**2) and reshape back to a matrix. So try:

a = np.vander([0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1.0,1.1,1.2]) #polynomial degree is 12, a -> (12x12)
b = np.random.rand(12,60,60) #p(x) is a 60x60 matrix that I have evaluated at 12 points
s = b.shape
x = np.linalg.solve(a, b.reshape(s[0], -1))
x = x.reshape(s)

Upvotes: 1

Related Questions