如何在 Python 中求解线性丢番图同余?

问题描述 投票:0回答:1

我面临一个挑战,解决方案涉及求解一系列 14 个变量的线性模方程。以下是这些方程的一部分:

3a + 3b + 3c + 3d + 3e + 3f + 3g + h + i + j + k + l + m + n = 15
7a + 9b + 17c + 11d + 6e + 5f + g = 3
13a + 2b + 9c + 8d + 12f + 13g = 17
5a + 2b + 16c + 12d + 5e + 7f + g = 11
6a + 4b + 9c + 6d + 4e + 9f + 6g + h + 7i + 11j + k + 7l + 11m + n = 8
10a + 15b + 13c + 10d + 15e + 13f + 10g + 12h + 18i + 8j + 12k + 18l + 8m + 12n = 18
9a + 12b + 14c + 4d + 9e + 16f + 3g + 7h + 17i + 11j + 14k + 3l + 18m + n = 15
9a + 12b + 16c + 15d + e + 14f + 6g + 11h + 2i + 9j + 12k + 16l + 15m + n = 14

我最终将它们复制到这个模块化方程求解器中以获得参数化解。但是,我希望能够在 Python 程序中自动执行此操作,而不依赖于该网站。

最好,我应该能够使用任意一组(线性)方程来完成此操作,而不仅仅是这些特定的方程。我的挑战解决方案的一部分要求我为不同的场景编写一些不同的方程,并根据需要交换它们。

我知道 SymPy 有 丢番图方程求解器 几乎可以完成我想要它做的事情,但在文档中我没有看到一种方法让它强制执行某个模数(mod 19、mod 23、等)。

python sympy number-theory modular-arithmetic diophantine
1个回答
0
投票

丢番图这个词在这里有误导性。你想要做的是 Z_p 上的线性代数,即以素数 p 为模的整数。如果您使用 DomainMatrix 而不是 Matrix,SymPy 可以做到这一点。这是在 SymPy 1.12 中执行此操作的方法:

In [1]: from sympy import *

In [2]: from sympy.polys.matrices import DomainMatrix

In [3]: a, b, c, d, e, f, g, h, i, j, k, l, m, n = syms = symbols('a:n')

In [4]: eqs = '''
   ...: 3a + 3b + 3c + 3d + 3e + 3f + 3g + h + i + j + k + l + m + n = 15
   ...: 7a + 9b + 17c + 11d + 6e + 5f + g = 3
   ...: 13a + 2b + 9c + 8d + 12f + 13g = 17
   ...: 5a + 2b + 16c + 12d + 5e + 7f + g = 11
   ...: 6a + 4b + 9c + 6d + 4e + 9f + 6g + h + 7i + 11j + k + 7l + 11m + n = 8
   ...: 10a + 15b + 13c + 10d + 15e + 13f + 10g + 12h + 18i + 8j + 12k + 18l + 8m + 12n = 18
   ...: 9a + 12b + 14c + 4d + 9e + 16f + 3g + 7h + 17i + 11j + 14k + 3l + 18m + n = 15
   ...: 9a + 12b + 16c + 15d + e + 14f + 6g + 11h + 2i + 9j + 12k + 16l + 15m + n = 14
   ...: '''

In [5]: eqs = [parse_expr(eq, transformations='all').subs(I, i) for eq in eqs.strip().splitlines()]

In [6]: M, b = linear_eq_to_matrix(eqs, syms)

In [7]: M
Out[7]: 
⎡3   3   3   3   3   3   3   1   1   1  1   1   1   1 ⎤
⎢                                                     ⎥
⎢7   9   17  11  6   5   1   0   0   0  0   0   0   0 ⎥
⎢                                                     ⎥
⎢13  2   9   8   0   12  13  0   0   0  0   0   0   0 ⎥
⎢                                                     ⎥
⎢5   2   16  12  5   7   1   0   0   0  0   0   0   0 ⎥
⎢                                                     ⎥
⎢6   4   9   6   4   9   6   1   18  0  1   7   11  1 ⎥
⎢                                                     ⎥
⎢10  15  13  10  15  13  10  12  26  0  12  18  8   12⎥
⎢                                                     ⎥
⎢9   12  14  4   9   16  3   7   28  0  14  3   18  1 ⎥
⎢                                                     ⎥
⎣9   12  16  15  1   14  6   11  11  0  12  16  15  1 ⎦

In [8]: b
Out[8]: 
⎡15⎤
⎢  ⎥
⎢3 ⎥
⎢  ⎥
⎢17⎥
⎢  ⎥
⎢11⎥
⎢  ⎥
⎢8 ⎥
⎢  ⎥
⎢18⎥
⎢  ⎥
⎢15⎥
⎢  ⎥
⎣14⎦

In [9]: Z_17 = GF(17, symmetric=False)

In [12]: dM = DomainMatrix.from_Matrix(M).convert_to(Z_17).to_dense()

In [13]: db = DomainMatrix.from_Matrix(b).convert_to(Z_17).to_dense()

In [14]: dM
Out[14]: DomainMatrix([[3 mod 17, 3 mod 17, 3 mod 17, 3 mod 17, 3 mod 17, 3 mod 17, 3 mod 17, 1 mod 17, 1 mod 17, 1 mod 17, 1 mod 17, 1 mod 17, 1 mod 17, 1 mod 17], [7 mod 17, 9 mod 17, 0 mod 17, 11 mod 17, 6 mod 17, 5 mod 17, 1 mod 17, 0 mod 17, 0 mod 17, 0 mod 17, 0 mod 17, 0 mod 17, 0 mod 17, 0 mod 17], [13 mod 17, 2 mod 17, 9 mod 17, 8 mod 17, 0 mod 17, 12 mod 17, 13 mod 17, 0 mod 17, 0 mod 17, 0 mod 17, 0 mod 17, 0 mod 17, 0 mod 17, 0 mod 17], [5 mod 17, 2 mod 17, 16 mod 17, 12 mod 17, 5 mod 17, 7 mod 17, 1 mod 17, 0 mod 17, 0 mod 17, 0 mod 17, 0 mod 17, 0 mod 17, 0 mod 17, 0 mod 17], [6 mod 17, 4 mod 17, 9 mod 17, 6 mod 17, 4 mod 17, 9 mod 17, 6 mod 17, 1 mod 17, 1 mod 17, 0 mod 17, 1 mod 17, 7 mod 17, 11 mod 17, 1 mod 17], [10 mod 17, 15 mod 17, 13 mod 17, 10 mod 17, 15 mod 17, 13 mod 17, 10 mod 17, 12 mod 17, 9 mod 17, 0 mod 17, 12 mod 17, 1 mod 17, 8 mod 17, 12 mod 17], [9 mod 17, 12 mod 17, 14 mod 17, 4 mod 17, 9 mod 17, 16 mod 17, 3 mod 17, 7 mod 17, 11 mod 17, 0 mod 17, 14 mod 17, 3 mod 17, 1 mod 17, 1 mod 17], [9 mod 17, 12 mod 17, 16 mod 17, 15 mod 17, 1 mod 17, 14 mod 17, 6 mod 17, 11 mod 17, 11 mod 17, 0 mod 17, 12 mod 17, 16 mod 17, 15 mod 17, 1 mod 17]], (8, 14), GF(17))

In [15]: db
Out[15]: DomainMatrix([[15 mod 17], [3 mod 17], [0 mod 17], [11 mod 17], [8 mod 17], [1 mod 17], [15 mod 17], [14 mod 17]], (8, 1), GF(17))

现在你有了 dM 和 db 作为 Z_17 上的矩阵,你可以使用任何你想要的矩阵运算并使用

to_Matrix
转换回矩阵:

In [16]: dM.hstack(db).rref()[0].to_Matrix() # RREF of augmented matrix
Out[16]: 
⎡1  0  0  0  0  0  0  0  11  11  5   2   2   6   0 ⎤
⎢                                                  ⎥
⎢0  1  0  0  0  0  0  0  14  2   11  12  14  13  6 ⎥
⎢                                                  ⎥
⎢0  0  1  0  0  0  0  0  3   7   10  11  15  0   5 ⎥
⎢                                                  ⎥
⎢0  0  0  1  0  0  0  0  2   2   1   16  16  13  5 ⎥
⎢                                                  ⎥
⎢0  0  0  0  1  0  0  0  13  13  16  7   14  0   15⎥
⎢                                                  ⎥
⎢0  0  0  0  0  1  0  0  16  10  5   11  15  11  7 ⎥
⎢                                                  ⎥
⎢0  0  0  0  0  0  1  0  8   10  6   13  1   0   7 ⎥
⎢                                                  ⎥
⎣0  0  0  0  0  0  0  1  4   6   9   6   8   8   16⎦

In [17]: dM.nullspace().to_Matrix() # nullspace of M
Out[17]: 
⎡15  16  1   12  10  11  14  7   11  0   0   0   0   0 ⎤
⎢                                                      ⎥
⎢15  12  8   12  10  9   9   2   0   11  0   0   0   0 ⎥
⎢                                                      ⎥
⎢13  15  9   6   11  13  2   3   0   0   11  0   0   0 ⎥
⎢                                                      ⎥
⎢12  4   15  11  8   15  10  2   0   0   0   11  0   0 ⎥
⎢                                                      ⎥
⎢12  16  5   11  16  5   6   14  0   0   0   0   11  0 ⎥
⎢                                                      ⎥
⎣2   10  0   10  0   15  0   14  0   0   0   0   0   11⎦

您可以使用这些来构建您的参数化解决方案:

In [30]: p = dM.hstack(db).rref()[0].to_Matrix()[:,-1]

In [31]: p = Matrix.vstack(p, p.zeros(6, 1))

In [32]: p.transpose()
Out[32]: [0  6  5  5  15  7  7  16  0  0  0  0  0  0]

In [33]: (M*p - b).applyfunc(lambda e: e % 17).is_zero_matrix
Out[33]: True

In [34]: N = dM.nullspace().to_Matrix()

In [35]: N = dM.nullspace().transpose().to_Matrix()

In [36]: (M*N).applyfunc(lambda e: e % 17).is_zero_matrix
Out[36]: True

参数解是从特定解和零空间生成的:

In [38]: sol = p + N*Matrix(symbols('alpha:6'))

In [39]: sol
Out[39]: 
⎡  15⋅α₀ + 15⋅α₁ + 13⋅α₂ + 12⋅α₃ + 12⋅α₄ + 2⋅α₅  ⎤
⎢                                                ⎥
⎢16⋅α₀ + 12⋅α₁ + 15⋅α₂ + 4⋅α₃ + 16⋅α₄ + 10⋅α₅ + 6⎥
⎢                                                ⎥
⎢      α₀ + 8⋅α₁ + 9⋅α₂ + 15⋅α₃ + 5⋅α₄ + 5       ⎥
⎢                                                ⎥
⎢12⋅α₀ + 12⋅α₁ + 6⋅α₂ + 11⋅α₃ + 11⋅α₄ + 10⋅α₅ + 5⎥
⎢                                                ⎥
⎢   10⋅α₀ + 10⋅α₁ + 11⋅α₂ + 8⋅α₃ + 16⋅α₄ + 15    ⎥
⎢                                                ⎥
⎢11⋅α₀ + 9⋅α₁ + 13⋅α₂ + 15⋅α₃ + 5⋅α₄ + 15⋅α₅ + 7 ⎥
⎢                                                ⎥
⎢     14⋅α₀ + 9⋅α₁ + 2⋅α₂ + 10⋅α₃ + 6⋅α₄ + 7     ⎥
⎢                                                ⎥
⎢ 7⋅α₀ + 2⋅α₁ + 3⋅α₂ + 2⋅α₃ + 14⋅α₄ + 14⋅α₅ + 16 ⎥
⎢                                                ⎥
⎢                     11⋅α₀                      ⎥
⎢                                                ⎥
⎢                     11⋅α₁                      ⎥
⎢                                                ⎥
⎢                     11⋅α₂                      ⎥
⎢                                                ⎥
⎢                     11⋅α₃                      ⎥
⎢                                                ⎥
⎢                     11⋅α₄                      ⎥
⎢                                                ⎥
⎣                     11⋅α₅                      ⎦

让我们检查一下:

In [47]: (M*sol)._rep.convert_to(Z_17[symbols('alpha:6')]).to_Matrix() == b.applyfunc(lambda e: e % 17)
Out[47]: True
© www.soinside.com 2019 - 2024. All rights reserved.