2

我正在使用 Python 中的有限域。我有一个包含多项式的矩阵,每个多项式都表示为一个整数。例如,多项式x^3 + x + 1表示为 11,因为:

x^3 + x + 1 ==> (1,0,1,1) ==> (8,0,2,1) ==> 8 + 0 + 2 + 1 ==> 11

给定一个矩阵,例如:

6, 2, 1
7, 0, 4
1, 7, 3

如何在 Python 中计算矩阵的逆?我已经实现了以下函数(对于多项式,而不是多项式矩阵):add(), sub(), mul(), div(),inv()

4

2 回答 2

1

我写了一个 python 包,它在 Galois 字段上实现了 numpy 数组。除了普通的数组算术,它还支持线性代数。https://github.com/mhostetter/galois

这是一个执行您所要求的示例的示例。

In [1]: import numpy as np                                                              

In [2]: import galois                                                                   

In [3]: GF = galois.GF(2**3)                                                            

In [4]: print(GF.properties)                                                            
GF(2^3):
  characteristic: 2
  degree: 3
  order: 8
  irreducible_poly: Poly(x^3 + x + 1, GF(2))
  is_primitive_poly: True
  primitive_element: GF(2, order=2^3)

In [5]: A = GF([[6,2,1],[7,0,4],[1,7,3]]); A                                            
Out[5]: 
GF([[6, 2, 1],
    [7, 0, 4],
    [1, 7, 3]], order=2^3)

In [6]: np.linalg.inv(A)                                                                
Out[6]: 
GF([[5, 5, 4],
    [3, 0, 1],
    [4, 3, 7]], order=2^3)

In [7]: np.linalg.inv(A) @ A                                                            
Out[7]: 
GF([[1, 0, 0],
    [0, 1, 0],
    [0, 0, 1]], order=2^3)
于 2021-05-01T15:40:23.020 回答
0
  1. 假设你的有限域是 $GF(8)$,你所有的计算函数(例如add(), mul())必须在同一个域 (GF(8)) 上进行计算。

  2. 在有限域上计算矩阵逆的方式与在其他域上计算矩阵逆的方式完全相同,而 Gauss-Jordan 消元法就可以做到。

以下是部分代码(来自:https ://github.com/foool/nclib/blob/master/gmatrixu8.cc ),用于计算 GF(256) 上的矩阵逆。希望对您有所帮助。

void GMatrixU8::Inverse(){
    GMatrixU8 mat_inv;
    int w = this->ww;
    int dim;
    int i, nzero_row, j;
    int temp, jtimes;

    dim = this->rr;
    mat_inv.Make_identity(this->rr, this->cc, this->ww);

    for(i = 0; i < dim; i++){
        nzero_row = i;
        if(0 == Get(i,i)){
            do{
                ++nzero_row;
                if(nzero_row >= dim){ ERROR("Non-full rank matrix!"); }
                temp = Get(nzero_row, i);
            }while((0 == temp)&&(nzero_row < dim));
            Swap_rows(i, nzero_row);
            mat_inv.Swap_rows(i, nzero_row);
        }
        for(j = 0; j < dim; j++){
            if(0 == Get(j,i))
                continue;
            if(j != i){
                jtimes = galois_single_divide(Get(j,i), Get(i,i), w);
                Row_plus_irow(j, i, jtimes);
                mat_inv.Row_plus_irow(j, i, jtimes);
            }else{
                jtimes = galois_inverse(Get(i, i), w);
                Row_to_irow(i , jtimes);
                mat_inv.Row_to_irow(i , jtimes);
            }
        }
    }
    this->ele = mat_inv.ele;
}
于 2017-09-14T14:47:47.947 回答