中国剩余定理的算法实现（基于c语言miracl库）

某电信安数基实验

要求：三个（老师要求的）方程组成的一次同余方程，文件中大数的顺序依次为a1,a2,a3,m1,m2,m3

函数：

multiply
函数原型: void multiply(big x, big y, big z);
功能说明: 两个大数相乘，z=x*y。

fdiv
函数原型：void fdiv(x,y,z)；
功能说明：将两个大数相除，z=x/y。

xgcd
函数原型: int xgcd(bigx, big y, big xd, big yd, big z);
功能说明: 计算两个大数的扩展最大公约数，也可以用来计算模逆。这个实验我们是用它来计算模逆。
例子: xgcd(x,p,x,x,x); //计算x^-1 mod p
/* x = 1/x mod p (p is prime) */

add
函数原型: void add(big x, big y, big z);
功能说明: 两个大数相加，z=x+y。
Example: add(x,x,x); // This doubles the value of x.

powmod
函数原型: void powmod(big x, big y,big z, big w);
功能说明: 模幂运算，w=xy(modz)。

cinstr
函数原型: int cinstr(big x, char*s);
功能说明: 将大数字符串转换成大数
返回值: 输入字符数的个数

egcd(x, y, z);//求最大公约数g=(x,y)

直接贴代码啦

#include<stdio.h>
#include<math.h>
#include<stdlib.h>
#include "miracl.h"
#include<Windows.h>

int main()
{
	FILE *fp;
	char fpname[100];//文件名
	int i, j, k;//备用
	char ch, datachar[7][500];//datachar为字符型的大数,使用7是为了防止最后一个换行符导致的数组越界
	miracl *mip = mirsys(500, 10);
	big a1 = mirvar(0);
	big a2 = mirvar(0);
	big a3 = mirvar(0);
	big m1 = mirvar(0);
	big m2 = mirvar(0);
	big m3 = mirvar(0);
	big m = mirvar(0);
	big Mi1 = mirvar(0);//为M1，M2，M3
	big Mi2 = mirvar(0);
	big Mi3 = mirvar(0);
	big Mj1 = mirvar(0);//为M1，M2，M3各自的逆
	big Mj2 = mirvar(0);
	big Mj3 = mirvar(0);
	big g1 = mirvar(0);//g1,g2,g3先用于判断m1,m2,m3是否互素,而后用于计算g = Mi*Mj*a 和 最后加法运算的中间变量
	big g2 = mirvar(0);
	big g3 = mirvar(0);
	big one = mirvar(1);//one=1

	printf("输入文件完整路径\n");
	scanf("%s", fpname);
	if ((fp = fopen(fpname, "r")) == NULL)
	{
		printf("fail to open the file\n");
		system("pause");
		exit(0);
	}
	else
	{
		i = 0; j = 0;
		while (!feof(fp))
		{
			fscanf(fp, "%c", &ch);//从文件中读入一个字符
			if (ch == '\n')
			{
				datachar[i][j] = '\0';//datacha[i][]输入完毕
				i++; j = 0; continue;//开始输入下一个数
			}
			datachar[i][j] = ch;
			j++;
		}
		datachar[i][j-1] = '\0';//datacha[i][]输入完毕,使用\0结尾
		//文件录入结束

		//把字符串转化为big型大数
		cinstr(a1, datachar[0]);
		cinstr(a2, datachar[1]);
		cinstr(a3, datachar[2]);
		cinstr(m1, datachar[3]);
		cinstr(m2, datachar[4]);
		cinstr(m3, datachar[5]);

		//判断mm1,m2,m3是否互素
		egcd(m1, m2, g1);//求最大公约数
		egcd(m2, m3, g2);
		egcd(m3, m1, g3);
		if (!((compare(g1, one) == 0) && (compare(g2, one) == 0) && (compare(g3, one)==0)))//当最大公约数都为1的时候，两两互素，则不进入if语句
		{
			printf("三个m并不两两互素，不能直接使用中国剩余定理\n");
			//结束
		}
		else {

			multiply(m1, m2, m);
			multiply(m, m3, m);//计算m=m1*m2*m3

			fdiv(m, m1, Mi1);//计算Mi1
			fdiv(m, m2, Mi2);//计算Mi2
			fdiv(m, m3, Mi3);//计算Mi3

			xgcd(Mi1, m1, Mj1, Mj1, Mj1);//Mj1为Mi1的逆
			xgcd(Mi2, m2, Mj2, Mj2, Mj2);//Mj2为Mi2的逆
			xgcd(Mi3, m3, Mj3, Mj3, Mj3);//Mj3为Mi3的逆

			multiply(Mi1, Mj1, g1);
			multiply(g1, a1, g1);//g1=Mi1*Mj1*a1

			multiply(Mi2, Mj2, g2);
			multiply(g2, a2, g2);//g2=Mi2*Mj2*a2

			multiply(Mi3, Mj3, g3);
			multiply(g3, a3, g3);//g3=Mi3*Mj3*a3

			/*cotnum(g1, stdout);
			cotnum(g2, stdout);
			cotnum(g3, stdout);*/

			add(g1, g2, g1);//g1 = Mi1*Mj1*a1+Mi2*Mj2*a2
			add(g1, g3, g1);//g1 = Mi1*Mj1*a1+Mi2*Mj2*a2+Mi3*Mj3*a3
			powmod(g1, one, m, g1);// g1 = g1的一次方 mod m

			printf("\n");
			printf("X is : \n");
			cotnum(g1, stdout);
			printf("m is : \n");
			cotnum(m, stdout);
			printf("\n");
			//即 X ≡ g1 ( mod m )
		}
	}
	fclose(fp);
	mirexit();
	system("pause");
}