Проблемы при подстановке матрицы в полином

Question

Пример: пусть

M = Matrix([[1,2],[3,4]]) # and 
p = Poly(x**3 + x + 1)    # then
p.subs(x,M).expand()

выдает ошибку:

TypeError: невозможно добавить <class'sympy.matrices.immutable.ImmutableDenseMatrix '> и <class' sympy.core.numbers.One '>

что очень правдоподобно, поскольку два первых члена становятся матрицами, но последний член (постоянный член) является не матрицей, а скаляром. Чтобы исправить эту ситуацию, я изменил многочлен на

p = Poly(x**3 + x + x**0)    # then

та же ошибка сохраняется. Обязан ли я набирать выражение вручную, заменяя x на M? В этом примере многочлен имеет только три члена, но на самом деле я сталкиваюсь (с многомерными многочленами) с десятками членов.

Answer 1

Поэтому я думаю, что вопрос в основном вращается вокруг концепции полинома матрицы :

(где P - многочлен, а A - матрица)

Я думаю, это говорит о том, что свободный член - это число, и его нельзя сложить с остальным, которое является матрицей, фактически операция сложения не определена между этими двумя типами.

TypeError: невозможно добавить <class'sympy.matrices.immutable.ImmutableDenseMatrix '> и <class' sympy.core.numbers.One '>

Однако это можно обойти, определив функцию, которая оценивает матричный полином для конкретной матрицы. Разница здесь в том, что мы используем возведение в степень , поэтому мы правильно вычисляем свободный член матричного полинома, a_0 * Iгде I=A^0- единичная матрица требуемой формы:

from sympy import *
x = symbols('x')
M = Matrix([[1,2],[3,4]])
p = Poly(x**3 + x + 1)

def eval_poly_matrix(P,A):
    res = zeros(*A.shape)
    for t in enumerate(P.all_coeffs()[::-1]):
        i, a_i = t
        res += a_i * (A**i)
    return res

eval_poly_matrix(p,M)

Выход:

В этом примере многочлен имеет только три члена, но на самом деле я сталкиваюсь (с многомерными многочленами) с десятками членов.

eval_poly_matrixВышеупомянутая функция может быть расширена для работы с многомерными многочленами, используя .monoms()метод извлечения одночленов с ненулевыми коэффициентами , например:

from sympy import *
x,y = symbols('x y')
M = Matrix([[1,2],[3,4]])
p = poly( x**3 * y + x * y**2 + y )

def eval_poly_matrix(P,*M):
    res = zeros(*M[0].shape)
    for m in P.monoms():
        term = eye(*M[0].shape)
        for j in enumerate(m):
            i,e = j
            term *= M[i]**e
        res += term
    return res

eval_poly_matrix(p,M,eye(M.rows))

Примечание. Возможны некоторые проверки работоспособности, обработка крайних случаев и оптимизации:

1. Количество переменных, присутствующих в полиноме, относится к количеству матриц, переданных в качестве параметров (первое никогда не должно быть больше второго, и если оно ниже, чем требуется какая-то логика для обработки этого, я рассмотрел только случай когда двое равны)
2. Все матрицы должны быть квадратными в соответствии с определением матричного полинома.
3. Обсуждение многомерной версии правила Хорнера содержится в комментариях к этому вопросу . Это может быть полезно для минимизации количества умножений матриц.
4. Учитывайте тот факт, что в матрице многочлен x*yотличается от того, y*xчто умножение матриц некоммутативно . По-видимому, поли-функции в sympy не поддерживают некоммутативные переменные, но вы можете определять символы с помощью, commutative=Falseи, похоже, там есть путь вперед

Что касается четвертого пункта выше, в SymPy есть поддержка выражений Matrix, и это может помочь здесь:

from sympy import *
from sympy.matrices import MatrixSymbol

A = Matrix([[1,2],[3,4]])
B = Matrix([[2,3],[3,4]])

X = MatrixSymbol('X',2,2)
Y = MatrixSymbol('Y',2,2)
I = eye(X.rows)

p = X**2 * Y + Y * X ** 2 + X ** 3 - I
display(p)

p = p.subs({X: A, Y: B}).doit()
display(p)

Выход:

Для получения дополнительной информации об этой функции следуйте # 18555