Вычисление одного элемента сопряженной или обратной символьной двоичной матрицы

Я пытаюсь получить один элемент вспомогательного A_adj матрицы A, оба из которых должны быть символическими выражениями, где символы x_i двоичные, а матрица A симметричная и разреженная. Python sympy отлично подходит для небольших задач:

from sympy import zeros, symbols

size = 4
A = zeros(size,size)
x_i = [x for x in symbols(f'x0:{size}')]

for i in range(size-1):
    A[i,i] += 0.5*x_i[i]
    A[i+1,i+1] += 0.5*x_i[i]
    A[i,i+1] = A[i+1,i] = -0.3*(i+1)*x_i[i]

A_adj_0 = A[1:,1:].det()

A_adj_0

Это вычисляет первый элемент A_adj_0 матрицы кофакторов (который является соответствующим второстепенным) и правильно дает мне 0,125x_0x_1x_2 - 0,28x_2x_2 ^ 2 - 0,055x_1 ^ 2x_2 - 0,28x_1x_2 ^ 2, это выражение, которое мне нужно, но есть две проблемы:

  1. Это совершенно невозможно для больших матриц (мне это нужно для sizes из ~100).
  2. x_i являются двоичными переменными (то есть 0 или 1), и, похоже, sympy не может упростить выражения двоичных переменных, то есть упростить многочлены x_i^n = x_i.

Первую проблему можно частично решить, вместо этого решая систему линейных уравнений Ay = b, где b задается первым базисным вектором [1, 0, 0, 0], так что y является первым столбцом обратного A. Первая запись y является первым элементом, обратным A:

b = zeros(size,1)
b[0] = 1

y = A.LUsolve(b)

s = {x_i[i]: 1 for i in range(size)}

print(y[0].subs(s) * A.subs(s).det())
print(A_adj_0.subs(s))

Проблема здесь в том, что выражение для первого элемента y чрезвычайно сложно, даже после использования simplify() и так далее. Это будет очень простое выражение с упрощением двоичных выражений, как указано в пункте 2 выше. Это более быстрый метод, но все же неприемлемый для больших матриц.

Это сводится к моему актуальному вопросу:

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

Я также открыт для использования другого программного обеспечения.

Приложение 1:
Кажется, упрощение бинарных выражений в sympy возможно с помощью простой пользовательской замены, о которой я не знал:

A_subs = A_adj_0

for i in range(size):
    A_subs = A_subs.subs(x_i[i]*x_i[i], x_i[i])

A_subs

sympy binary matrix sparse-matrix inverse
person candlejack    schedule 25.01.2021    source источник
comment
Это отличные подсказки, спасибо! Это, безусловно, станет проблемой, которая возникнет в дальнейшем, но, к сожалению, в моей текущей задаче коэффициенты реальны. Следовательно, если я правильно понимаю, вычисления не могут быть выполнены над конечными полями.   -  person candlejack    schedule 25.01.2021

Ответы (1)


arrow_upward
2
arrow_downward

Вы должны убедиться, что используете Rational, а не числа с плавающей запятой в sympy, поэтому S(1)/2 или Rational(1, 2), а не 0.5.

Существует новая (недокументированная и на данный момент внутренняя) реализация матриц в sympy под названием DomainMatrix. Это, вероятно, будет намного быстрее для такой задачи и всегда дает полиномиальные результаты в полностью расширенной форме. Я ожидаю, что это будет намного быстрее для такого рода проблем, но все еще кажется довольно медленным для этого, потому что он не является разреженным внутри (пока - это, вероятно, изменится в следующем выпуске), и он не использует преимущества упрощения из символов, являющихся двоичными значениями. Его можно заставить работать над GF(2), но не с символами, которые, как предполагается, находятся в GF(2), что является чем-то другим.

В случае, если это полезно, вот как вы могли бы использовать его в sympy 1.7.1:

from sympy import zeros, symbols, Rational
from sympy.polys.domainmatrix import DomainMatrix


size = 10
A = zeros(size,size)
x_i = [x for x in symbols(f'x0:{size}')]

for i in range(size-1):
    A[i,i] += Rational(1, 2)*x_i[i]
    A[i+1,i+1] += Rational(1, 2)*x_i[i]
    A[i,i+1] = A[i+1,i] = -Rational(3, 10)*(i+1)*x_i[i]

# Convert to DomainMatrix:
dM = DomainMatrix.from_list_sympy(size-1, size-1, A[1:, 1:].tolist())

# Compute determinant and convert back to normal sympy expression:
# Could also use dM.det().as_expr() although it might be slower
A_adj_0 = dM.charpoly()[-1].as_expr()

# Reduce powers:
A_adj_0 = A_adj_0.replace(lambda e: e.is_Pow, lambda e: e.args[0])

print(A_adj_0)
person Oscar Benjamin    schedule 25.01.2021
comment
Я только что проверил это на более крупной проблеме с реальными значениями, и похоже, что это работает НАМНОГО быстрее, чем другие методы. Это не что иное, как магия. Не могу не выразить признательности за ваше время. Благодарю вас! - person candlejack; 25.01.2021