【信息安全数学基础】中国剩余定理
算法原理:
设正整数 𝒎𝟏,𝒎𝟐,⋯,𝒎𝒌 两两互素,对任意整数 𝒂𝟏,𝒂𝟐,⋯,𝒂𝒌,一次同余方程组
𝒙 ≡ 𝒂𝟏 (𝒎𝒐𝒅 𝒎𝟏)
𝒙 ≡ 𝒂𝟐 (𝒎𝒐𝒅 𝒎𝟐)
⋮
𝒙 ≡ 𝒂𝒌 (𝒎𝒐𝒅 𝒎𝒌)
在模𝒎意义下有唯一解,该解可表示为 𝒙 ≡ 𝑴𝟏(𝑴𝟏^-1)𝒂𝟏 +𝑴𝟐(𝑴𝟐^−𝟏)𝒂𝟐 +⋯+𝑴𝒌(𝑴𝒌^−𝟏)𝒂𝒌( 𝒎𝒐𝒅 𝒎) 其中𝒎 = 𝒎𝟏𝒎𝟐 ⋯𝒎𝒌,𝑴𝒋 = 𝒎/ 𝒎𝒋,𝑴𝒋(𝑴𝒋−𝟏) ≡ 𝟏( 𝒎𝒐𝒅 𝒎𝒋) ,𝒋 = 𝟏,𝟐,⋯𝒌。
算法步骤:
(1) 判断正整数 𝒎𝟏,𝒎𝟐,⋯,𝒎𝒌 是否两两互素;是,则继续, 否则跳出,输出“不能直接利用中国剩余定理”
(2) 计算 (𝑴𝒋^−𝟏) (𝒎𝒐𝒅 𝒎𝒋)
(3) 计算 𝒙𝒋≡ 𝑴𝒋(𝑴𝒋^−𝟏)𝒂𝒋 (𝒎𝒐𝒅 𝒎)
(4) 计算 𝒙 ≡ ∑𝒋=𝟏 𝒌 𝒙𝒋 (𝒎𝒐𝒅 𝒎)。
代码设计及注释:
1 extern "C" 2 { 3 #include<miracl.h> 4 #include<mirdef.h> 5 } 6 #include<stdio.h> 7 #include<stdlib.h> 8 9 #define NUM 3 //从文本读入方程组个数,数据以分行存储 10 11 int main() 12 { 13 //中国剩余定理:x≡M1M1^1a1+M2M2^-1a2+……+MkMk^-1ak(mod m)其中m=m1m2m……mk,Mj=m/mj,MjMj^-1≡1(mod mj),j=1,2,……k 14 int i, j;//设定循环变量 15 FILE *fp; 16 int flag = 0; 17 big a[NUM], m[NUM], x[NUM], Mj[NUM], Mj1[NUM]; //设定大数数组,a存放同余号右边,m存放模值,x存放同余号左边,Mj与字面同义,Mj1存放Mj^-1 18 big t0, t1, M, m1, X, y, W; 19 miracl *mip = mirsys(4000, 10); //开辟空间 20 mip->IOBASE = 10; 21 for (i = 0; i < NUM; i++)//初始化数组 22 { 23 a[i] = mirvar(0); 24 m[i] = mirvar(0); 25 x[i] = mirvar(0); 26 Mj[i] = mirvar(0); 27 Mj1[i] = mirvar(0); 28 } 29 t0 = mirvar(0);//初始化各变量 30 t1 = mirvar(1); 31 M = mirvar(1); 32 m1 = mirvar(0); 33 X = mirvar(0); 34 y = mirvar(1); 35 W = mirvar(0); 36 37 fp = fopen("2.txt", "r"); 38 char tmp[4000]; 39 while (1) 40 { 41 fscanf(fp, "%s", tmp); 42 if (flag<NUM)//从文件中读入a列并存入数组 43 { 44 cinstr(a[flag], tmp); 45 cotnum(a[flag], stdout); 46 flag++; 47 } 48 else//从文件中读入m列存入数组 49 { 50 cinstr(m[flag - NUM], tmp); 51 cotnum(m[flag - NUM], stdout); 52 flag++; 53 } 54 if (flag == 2 * NUM) break; 55 } 56 printf("结果为:"); 57 fclose(fp); 58 59 for (i = 0; i<NUM; i++) //判断m列是否两两互素 60 { 61 for (j = i + 1; j<NUM; j++) 62 { 63 egcd(m[i], m[j], t0);//t0=存放公约数,若gcd(m[i], m[j])=1即m两两互素进行下一次判断,否则跳出 64 if (mr_compare(t0, t1)==0) continue; 65 else 66 { 67 printf("不能直接应用中国剩余定理\n"); 68 exit(0); 69 } 70 } 71 } 72 73 for (i = 0; i<NUM; i++) 74 { 75 multiply(m[i], M, M); //计算M=m1m2……mk 76 } 77 78 copy(M, W);//W=M 79 for (i = 0; i<3; i++) 80 { 81 divide(M, m[i], Mj[i]); //计算Mj 82 xgcd(Mj[i], m[i], Mj1[i], m1, t1); //计算Mj^1(mod mj) 83 copy(W, M);//M=W,数值还原 84 } 85 86 for (i = 0; i<NUM; i++) //计算xj及其累加结果X 87 { 88 multiply(Mj[i], Mj1[i], t0); 89 multiply(t0, a[i], x[i]); 90 add(x[i], X, X); 91 } 92 93 94 powmod(X, y, M, X);//X mod M 95 96 cotnum(X, stdout); 97 mirexit(); 98 99 return 0; 100 101 }