C++实现伽罗华域生成及四则运算（三）

最新推荐文章于 2026-01-03 19:30:00 发布

原创最新推荐文章于 2026-01-03 19:30:00 发布 · 885 阅读

19 ·

CC 4.0 BY-SA版权

文章标签：

#c++ #开发语言

C++ 同时被 3 个专栏收录

7 篇文章

订阅专栏

密码

4 篇文章

订阅专栏

数学

3 篇文章

订阅专栏

原来的求行列式代码

/*方阵行列式,r是阶数
代数余子式法
*/
template<uint8_t n>
uint16_t GFM<n>::delta(uint16_t* A, uint16_t r) {
    if (r == 0) throw TypeException("矩阵阶数出错");
    if (r == 1) return A[0];
    if (r == 2) return this->multi(A[0],A[3]) ^ this->multi(A[1],A[2]); //对角线法
    // 按第0行展开代数余子式
    uint16_t* B = (uint16_t*)malloc(sizeof(uint16_t) * (r-1) * (r-1));
    uint16_t res = 0;
    for (uint16_t omit = 0; omit < r; omit++) { //指定列号，求指定列的代数余子式
        for (uint16_t i = 1, k = 0; i < r; i++) { //从第1行开始读取，写入B[k]
            for (uint16_t j = 0; j < r; j++) { //逐列读取
                if (j == omit) continue; //跳过指定的列
                B[k++] = A[i*r+j];
            }
        }
        res ^= this->multi(A[omit],delta(B, r-1)); //按A的第0行代数余子式展开 
    }
    free(B); //释放内存
    return res;
}

是利用嵌套展开代数余子式的方法，uint16_t* B = (uint16_t*)malloc(sizeof(uint16_t) * (r-1) * (r-1));申请内存，用于存放比方阵A低1阶的子方阵，而对于该子方阵，只用再次调用delta函数就可以了，理解上应该没有困难，但仔细一想，太费内存了。
假如A是16阶，则B就是15阶，B的子方阵就是14阶，依次类推，直到最后到达2阶if (r == 2) return this->multi(A[0],A[3]) ^ this->multi(A[1],A[2]);才结束，这前后一共申请内存 $(16^2+15^2+\dots+2^2)*sizeof(uint16\_t)=2982Bytes$ ，足足约3K。

改进后的求行列式代码

/*方阵行列式,r是阶数
上三角对角线法
*/
template<uint8_t n>
uint16_t GFM<n>::delta(uint16_t* A, uint16_t r) {
    if (r == 0) throw TypeException("矩阵阶数出错");
    if (r == 1) return A[0];
    //复制A
    uint16_t* B = (uint16_t*)malloc(sizeof(uint16_t) * r * r);
    memcpy(B, A, sizeof(uint16_t) * r * r); 
    uint16_t res = 1;
    //将方阵转化成上三角形
    for (uint16_t i = 0; i < r-1; i++) { //指定与第i行异或
        if (B[i*r+i] == 0) { //第i行i列若已经是0，则要与不为0的下一行交换一整行，行列式变号，但GF(2^n)域内特殊，不用变号
            for (uint16_t j = i+1; j < r; j++) {
                if (B[j*r+i] == 0) continue;
                uint16_t temp[r];
                memcpy(temp, B+i*r, sizeof(uint16_t)*r);
                memcpy(B+i*r, B+j*r, sizeof(uint16_t)*r);
                memcpy(B+j*r, temp, sizeof(uint16_t)*r);
            }
            if (B[i*r+i] == 0) return 0; //全部查询一遍，结果仍为0，则结果必然为0
        };
        for (uint16_t j = i+1; j < r; j++) { //指定当前参与变换的行
            if (B[j*r+i] == 0) continue; //已经是0，刚好不用乘
            uint16_t row_times = this->multi(this->adver(B[i*r+i]),B[j*r+i]); //计算商=B[j][i]/B[i][i]
            for (uint16_t k = i; k < r; k++) { //当前行加上第i行乘row_times，行列式不变，省略前i个0
                B[j*r+k] ^= this->multi(B[i*r+k],row_times);
            }
        }
    }
    //计算对角线乘积
    for (uint16_t i = 0; i < r; i++) {
        res = this->multi(B[i*r+i], res);
    }
cout << res << endl;
    free(B); //释放内存
    return res;
}

看着有些长，但它不嵌套，采用的是变换上三角方法，最后直接把对角线作乘积就得出了其行列式值。这里用到2个定理：

一行的倍数加到另一行，行列式不变
$\begin{vmatrix} a & b & c\\ d & e & f\\ g & h & i \end{vmatrix} =\begin{vmatrix} a & b & c\\ d+a*k & e+b*k & f+c*k\\ g & h & i \end{vmatrix}$
交换两行，行列式变号
$\begin{vmatrix} a & b & c\\ d & e & f\\ g & h & i \end{vmatrix} =-\begin{vmatrix} d & e & f\\ a & b & c\\ g & h & i \end{vmatrix}$

GF(256)域幻方4阶矩阵求逆

幻方4阶矩阵形如下：
$\begin{bmatrix} a & b & c & d\\ d & a & b & c\\ c & d & a & b\\ b & c & d & a \end{bmatrix} 或 \begin{bmatrix} a & b & c & d\\ b & c & d & a\\ c & d & a & b\\ d & a & b & c \end{bmatrix}$
其特点是每行依次左移或右移1位，在很多加密算法里都比较常见。
那么如何求它的逆矩阵呢？用前2点的方法肯定可以，现在我介绍第3种方法，先看推导：
$\begin{bmatrix} a & b & c & d\\ d & a & b & c\\ c & d & a & b\\ b & c & d & a \end{bmatrix} ^4= \left( \begin{bmatrix} a & b & c & d\\ d & a & b & c\\ c & d & a & b\\ b & c & d & a \end{bmatrix} \times \begin{bmatrix} a & b & c & d\\ d & a & b & c\\ c & d & a & b\\ b & c & d & a \end{bmatrix} \right) ^2$

$=\begin{bmatrix} a*a+b*d+c*c+d*b & a*b+b*a+c*d+d*c & a*c+b*b+c*a+d*d & a*d+b*c+c*b+d*a\\ d*a+a*d+b*c+c*b & d*b+a*a+b*d+c*c & d*c+a*b+b*a+c*d & d*d+a*c+b*b+c*a\\ c*a+d*d+a*c+b*b & c*b+d*a+a*d+b*c & c*c+d*b+a*a+b*d & c*d+d*c+a*b+b*a\\ b*a+c*d+d*c+a*b & b*b+c*a+d*d+a*c & b*c+c*b+d*a+a*d & b*d+c*c+d*b+a*a \end{bmatrix} ^2$

$=\begin{bmatrix} a*a+c*c & 0 & b*b+d*d & 0\\ 0 & a*a+c*c & 0 & b*b+d*d\\ b*b+d*d & 0 & a*a+c*c & 0\\ 0 & b*b+d*d & 0 & a*a+c*c \end{bmatrix} ^2$

$=\begin{bmatrix} (a*a+c*c)^2+(b*b+d*d)^2 & 0 & (a*a+c*c)(b*b+d*d)+(b*b+d*d)(a*a+c*c) & 0\\ 0 & (a*a+c*c)^2+(b*b+d*d)^2 & 0 & (a*a+c*c)(b*b+d*d)+(b*b+d*d)(a*a+c*c)\\ (b*b+d*d)(a*a+c*c)+(a*a+c*c)(b*b+d*d) & 0 & (b*b+d*d)^2+(a*a+c*c)^2 & 0\\ 0 & (b*b+d*d)(a*a+c*c)+(a*a+c*c)(b*b+d*d) & 0 & (b*b+d*d)^2+(a*a+c*c)^2 \end{bmatrix}$

$=\begin{bmatrix} (a*a+c*c)^2+(b*b+d*d)^2 & 0 & 0 & 0\\ 0 & (a*a+c*c)^2+(b*b+d*d)^2 & 0 & 0\\ 0 & 0 & (a*a+c*c)^2+(b*b+d*d)^2 & 0\\ 0 & 0 & 0 & (a*a+c*c)^2+(b*b+d*d)^2 \end{bmatrix}$

$=\begin{bmatrix} a^4+c^4+b^4+d^4 & 0 & 0 & 0\\ 0 & a^4+c^4+b^4+d^4 & 0 & 0\\ 0 & 0 & a^4+c^4+b^4+d^4 & 0\\ 0 & 0 & 0 & a^4+c^4+b^4+d^4 \end{bmatrix}=A$

令 $x=a^4+b^4+c^4+d^4$ 因此
$\div x = \begin{bmatrix} 1 & 0 & 0 & 0\\ 0 & 1 & 0 & 0\\ 0 & 0 & 1 & 0\\ 0 & 0 & 0 & 1 \end{bmatrix}$

由此可以推断：
$\begin{bmatrix} a & b & c & d\\ d & a & b & c\\ c & d & a & b\\ b & c & d & a \end{bmatrix} ^3 \div x 与 \begin{bmatrix} a & b & c & d\\ d & a & b & c\\ c & d & a & b\\ b & c & d & a \end{bmatrix} 互为逆矩阵$

证：
$\begin{bmatrix} a & b & c & d\\ d & a & b & c\\ c & d & a & b\\ b & c & d & a \end{bmatrix} ^3 \div x = \begin{bmatrix} a*a+c*c & 0 & b*b+d*d & 0\\ 0 & a*a+c*c & 0 & b*b+d*d\\ b*b+d*d & 0 & a*a+c*c & 0\\ 0 & b*b+d*d & 0 & a*a+c*c \end{bmatrix} \times \begin{bmatrix} a & b & c & d\\ d & a & b & c\\ c & d & a & b\\ b & c & d & a \end{bmatrix} \div x$
$\begin{bmatrix} a^3+ac^2+b^2c+cd^2 & a^2b+bc^2+b^2d+d^3 & a^2c+c^3+ab^2+ad^2 & a^2d+c^2d+b^3+bd^2\\ a^2d+c^2d+b^3+bd^2 & a^3+ac^2+b^2c+cd^2 & a^2b+bc^2+b^2d+d^3 & a^2c+c^3+ab^2+ad^2\\ a^2c+c^3+ab^2+ad^2 & a^2d+c^2d+b^3+bd^2 & a^3+ac^2+b^2c+cd^2 & a^2b+bc^2+b^2d+d^3\\ a^2b+bc^2+b^2d+d^3 & a^2c+c^3+ab^2+ad^2 & a^2d+c^2d+b^3+bd^2 & a^3+ac^2+b^2c+cd^2 \end{bmatrix} \div x$
$\begin{bmatrix} \frac{a^3+ac^2+b^2c+cd^2}{x} & \frac{a^2b+bc^2+b^2d+d^3}{x} & \frac{a^2c+c^3+ab^2+ad^2}{x} & \frac{a^2d+c^2d+b^3+bd^2}{x}\\ \frac{a^2d+c^2d+b^3+bd^2}{x} & \frac{a^3+ac^2+b^2c+cd^2}{x} & \frac{a^2b+bc^2+b^2d+d^3}{x} & \frac{a^2c+c^3+ab^2+ad^2}{x}\\ \frac{a^2c+c^3+ab^2+ad^2}{x} & \frac{a^2d+c^2d+b^3+bd^2}{x} & \frac{a^3+ac^2+b^2c+cd^2}{x} & \frac{a^2b+bc^2+b^2d+d^3}{x}\\ \frac{a^2b+bc^2+b^2d+d^3}{x} & \frac{a^2c+c^3+ab^2+ad^2}{x} & \frac{a^2d+c^2d+b^3+bd^2}{x} & \frac{a^3+ac^2+b^2c+cd^2}{x} \end{bmatrix}$

$\begin{bmatrix} a & b & c & d\\ d & a & b & c\\ c & d & a & b\\ b & c & d & a \end{bmatrix} ^3 \div x \times \begin{bmatrix} a & b & c & d\\ d & a & b & c\\ c & d & a & b\\ b & c & d & a \end{bmatrix}$

在这里插入图片描述

$=\begin{bmatrix} \frac{a^4+d^4+c^4+b^4}{x} & \frac{0}{x} & \frac{0}{x} & \frac{0}{x}\\ \frac{0}{x} & \frac{a^4+d^4+c^4+b^4}{x} & \frac{0}{x} & \frac{0}{x}\\ \frac{0}{x} & \frac{0}{x} & \frac{a^4+d^4+c^4+b^4}{x} & \frac{0}{x}\\ \frac{0}{x} & \frac{0}{x} & \frac{0}{x} & \frac{a^4+d^4+c^4+b^4}{x} \end{bmatrix}=E$
证毕。

测试代码如下：

//测试4阶幻方矩阵
    GFM<8> gf8;
    // uint16_t a[16] = {0xe,0x9,0x8,0x7, 0x7,0xe,0x9,0x8, 0x8,0x7,0xe,0x9, 0x9,0x8,0x7,0xe};
    uint16_t a[16] = {0xe,0x9,0x8,0x7, 0x9,0x8,0x7,0xe, 0x8,0x7,0xe,0x9, 0x7,0xe,0x9,0x8};
    uint16_t res[16];
    uint16_t q[16];
    uint16_t x = gf8.adver(gf8.calc("0xe^4+0x9^4+0x8^4+0x7^4"));
    gf8.matrixPow(q, a, 4, 3);
    gf8.matrixNumMulti(q, q, 4, 4, x);
    gf8.matrixMulti(res, a, 4, 4, q, 4, 4);
    cout << "验证:" << endl;
    for (int i = 0; i < 4; i++) {
        for (int j = 0; j < 4; j++){
            cout << hex << setw(2) << res[i*4+j] << " ";
        }
        cout << endl;
    }