基于GMP大数库用 “ 快速幂取模算法 ” + “ 扩展欧几里得算法 ”实现RSA加密算法（C++实现）

RSA涉及到大数运算，所以我本能地想到了之前打CTF了解到的一个第三方库gmp，这次RSA算法索性用gmp来实现了，顺便入门一下gmp的API函数。

基于“VS2017+GMP”实现RSA算法

• 一、RSA算法简述
• 二、子算法
• • （一）扩展欧几里得算法
  • • 1.思想
    • 2.案例
    • 3.代码构思与实现
  • （二）快速幂取模算法
  • • 1.算法设计
    • 2.代码实现
• 三、环境
• • （一）配置GMP
  • （二）API
  • • 1.初始化与赋值函数
    • • void mpz_init (mpz_t x)
      • void mpz_set (mpz_t rop, const mpz_t op)
      • int mpz_set_str (mpz_t rop, const char *str, int base)
    • 2.运算函数
    • • void mpz_sub (mpz_t rop, const mpz_t op1, const mpz_t op2)
      • void mpz_add (mpz_t rop, const mpz_t op1, const mpz_t op2)
      • void mpz_mod (mpz_t r, const mpz_t n, const mpz_t d)
      • void mpz_mul(mpz_t rop, const mpz_t op1, const mpz_t op2)
      • void mpz_div (mpz_t rop, const mpz_t op1, const mpz_t op2)
      • void mpz_gcd (mpz_t rop, const mpz_t op1, const mpz_t op2)
    • 3.变量释放
    • • void mpz_clear (mpz_t x)
    • 4.判断函数
    • • int mpz_probab_prime_p (const mpz_t n, int reps)
      • int mpz_odd_p (const mpz_t op) 和 int mpz_even_p (const mpz_t op)
      • int mpz_cmp (const mpz_t op1, const mpz_t op2)
    • 4.随机数生成函数
    • • void gmp_randinit_default (gmp_randstate_t state)
      • void gmp_randseed_ui (gmp_randstate_t state, unsigned long int seed)
      • void mpz_urandomb (mpz_t rop, gmp_randstate_t state, mp_bitcnt_t n)
    • 5.其他
    • • void mpz_nextprime (mpz_t rop, const mpz_t op)
• 三、原创源码
• 四、实现效果的呈现与加密结果的验证
• • （一）随机生成密钥对并自动加密
  • （二）自定义密钥，对明文加密

一、RSA算法简述

关于rsa算法的思想，网上很多资料可以参考，这边不再赘述，用一图可概述如下。

二、子算法

在实现RSA算法之前我们需要先学会两个必要的算法

（一）扩展欧几里得算法

在上图RSA算法简述的例图中，算法第三条要求出公钥e关于n的欧拉函数的乘法逆元，即私钥d，而求乘法逆元就需要用到扩展欧几里得算法。

1.思想

已知整数a、b就可以在求得a、b的最大公约数的同时，能找到整数x、y，使它们满足贝祖等式 ax+by=gcd(a,b)

2.案例

以求二元一次不定方程 47x+30y=1 的整数解为例。

我们很容易得到下面式子

47=30*1+17 //被除数==47，除数==30，余数==17，商==1
30=17*1+13 //被除数==30，除数==17，余数==13，商==1
17=13*1+4  //被除数==17，除数==13，余数==4，商==1
13=4*3+1   //被除数==13，除数==4 ，余数==1，商==3

然后把它们余数都放到左边，右边放其他项。

17=47*1+30*(-1)
13=30*1+17*(-1)
4=17*1+13*(-1)
1=13+4*(-3)

因为第n个式子的除数是第n-1个式子的余数，所以可以从最后一个式子开始用上一个式子带入，也就是用第n-1个式子代替第n个式子的除数，别的不要动，也不要计算，带入一次后各个除数和余数的系数合并整理，进入下一次的替换，这就形成了递归，直到最后n == 1，就会出现出现形如 1 = 47x + 30y 。

1=13+4*(-3)
1=13+[17*1+13*(-1)]*(-3)     //4=17*1+13*(-1)替换除数4
1=17*(-3)+13*4               //整理
1= 17*(-3)+[30*1+17*(-1)]*4  //13=30*1+17*(-1)替换除数13
1=30*4+17*(-7)               //整理
1=30*4+[47*1+30*(-1)]*(-7)   //17=47*1+30*(-1)替换除数17
1=47*(-7)+30*11              //整理

回到求公钥e关于n的乘法逆元，$$(d\ast e)mod\phi (n)=1$$
上面的式子不就可以写成$$d\ast e-\phi (n)y=1$$
d和n的欧拉函数都已知，e是要求的未知数 x ，同样的，用上面案例中的方法。

3.代码构思与实现

还是先以求二元一次不定方程 47x+30y=1 的整数解为例。

原理知道了，但我们说到底要求的是系数，具体代码实现肯定需要两个变量存储每一次递归的系数，观察每一次替换整理后得到的式子如下
$$1=13\ast 1+4\ast (-3)$$
$$1=17\ast (-3)+13\ast 4$$
$$1=30\ast 4+17\ast (-7)$$
$$1=47\ast (-7)+30\ast 11$$
（以 前两个式子 为例）罗列一下表达式系数间的关系。

为了更容易理解，这里的代码实现就不用gmp的语法了。

//x和y都是系数
void Extended_Euclid(int a,int b,int &temp,int &x,int &y)
{
    if (b==0)                    //初始化
    {
        temp=a;x=1;y=0;
        return;
    }
    Extended_Euclid(b,a%b,d,y,x);//回代，x,y交叉代入
    y-=a/b*x;                    //y需要减去商
}
int main()
{
    scanf("%d%d",&m,&n);         //求m关于n的逆元
    Extended_Euclid(m,n,d,x,y);
    printf("%d\n",(x+n)%n);
    return 0;
}

（二）快速幂取模算法

我们在对明文m加密的时候，用到的公式是
$$c=(m^e)modn$$我们知道，公钥e一般是大数，开次方无疑会产生的数据很大，用long long型恐怕都不够存储，所以我们需要对这个加密算法进行优化，快速幂取模算法无疑是最好的选择。

1.算法设计

首先e其实可以写成这样的形式
$$e=\sum _{i=0}^n{x}_i\ast 2^i=x_n\ast 2^n+...+x_i\ast 2^i+...+x_2\ast 2^2+x_1\ast 2^1+x_0\ast 2^0$$
其中 x_i 都是 e 转化成二进制之后各位上的值，为0或1。

那么不难得出$$c=(m^e)modn=m^{x_n\ast 2^n}\ast ...\ast m^{x_i\ast 2^i}\ast ...\ast m^{x_2\ast 2^2}\ast m^{x_1\ast 2^1}\ast m^{x_0\ast 2^0}modn$$
而对于某一项

$$m^{x_i\ast 2^i}modn$$如果

• x_i=0：无效项
• x_i=1：有效项

所以得出结论，表达式$$c=(m^e)modn=m^{x_n\ast 2^n}\ast ...\ast m^{x_i\ast 2^i}\ast ...\ast m^{x_2\ast 2^2}\ast m^{x_1\ast 2^1}\ast m^{x_0\ast 2^0}modn$$是可以简化成所有有效项的乘积的，而某一项是否是有效项，就取决于对应的二进制位是否为1.

2.代码实现

int PowerMod(int a, int b, int c)
{
    int ans = 1;
    a = a % c;
    while(b>0) //取出的b==0时也就代表取二进制位取完了，循环也就结束
    {
        if(b % 2 = = 1) //根据二进制位是否是1决定要不要乘进去
            ans = (ans * a) % c;
        b = b/2;   //b就是x_i，存储每一个二进制位
        a = (a * a) % c;//更新a的值,a就代表每一项
    }
    return ans;
}

三、环境

（一）配置GMP

关于windows下如何配置第三方库gmp，强推一个博主的写的教程《windows环境下的gmp大数运算库的配置》，亲测有效。

（二）API

关于GMP内置函数的使用，csdn上目前没有博主做全面的整理，还是要去看英文版 GMP官方手册
下面我参考官方手册对本次RSA算法用到的gmp内置函数作简单的整理概述

1.初始化与赋值函数

void mpz_init (mpz_t x)

对mpz_t类型的变量x初始化值为0（每定义一个mpz_t都需要初始化变量）

void mpz_set (mpz_t rop, const mpz_t op)

用已经初始化的mpz_t变量op给另一个变量rop赋值

int mpz_set_str (mpz_t rop, const char *str, int base)

用一个base进制的字符串给rop变量赋值，base是进制数
例如：

mpz_set_str (rop, "1234567890", 10)
//相当于rop=1234567890

2.运算函数

void mpz_sub (mpz_t rop, const mpz_t op1, const mpz_t op2)

rop = op1 * op2

void mpz_add (mpz_t rop, const mpz_t op1, const mpz_t op2)

rop = op1 + op2

void mpz_mod (mpz_t r, const mpz_t n, const mpz_t d)

r = n % d

void mpz_mul(mpz_t rop, const mpz_t op1, const mpz_t op2)

rop = op1 * op2

void mpz_div (mpz_t rop, const mpz_t op1, const mpz_t op2)

rop = op1 / op2

void mpz_gcd (mpz_t rop, const mpz_t op1, const mpz_t op2)

rop = gcd(op1 , op2)

3.变量释放

void mpz_clear (mpz_t x)

释放之前定义的x变量

4.判断函数

int mpz_probab_prime_p (const mpz_t n, int reps)

判断n是不是素数，

• 如果确定是就返回2 ，
• 如果可能是素数就返回1，
• 如果不是素数就返回0，

判断误差取决于reps的值，理想的reps值在15-50之间

int mpz_odd_p (const mpz_t op) 和 int mpz_even_p (const mpz_t op)

判断是不是奇数（偶数），是就返回1，不是就返回0。

int mpz_cmp (const mpz_t op1, const mpz_t op2)

判断op1和op2的大小

• op1 > op2 : 返回正值
• op1 = op2 : 返回0
• op1 < op2 : 返回负值

4.随机数生成函数

void gmp_randinit_default (gmp_randstate_t state)

在随机生成一个数之前，需要定义一个gmp_randstate_t类型的状态变量，这个函数是对状态变量进行初始化的函数。

void gmp_randseed_ui (gmp_randstate_t state, unsigned long int seed)

初始化状态之后，要给状态变量state设置种子seed

void mpz_urandomb (mpz_t rop, gmp_randstate_t state, mp_bitcnt_t n)

在0~2ⁿ范围内随机生成一个数

5.其他

void mpz_nextprime (mpz_t rop, const mpz_t op)

生成一个比op大的素数rop。

三、原创源码

//基于GMP大数运算库的环境下运行
#define _CRT_SECURE_NO_WARNINGS
#include<stdio.h>
#include<gmp.h>
#include<string>
#include<stdlib.h>
#include<time.h>
#include<windows.h>
#define N 1000000
#define M 2000
using namespace std;

/*
  gmp语法参考：
  https://www.cnblogs.com/y3w3l/p/5947450.html
  https://www.cnblogs.com/ECJTUACM-873284962/p/8350320.html
  https://blog.csdn.net/iteye_17686/article/details/82310869?ops_request_misc=&request_id=&biz_id=102&utm_term=mpz_init_set_str&utm_medium=distribute.pc_search_result.none-task-blog-2~all~sobaiduweb~default-5-.nonecase&spm=1018.2226.3001.4187
*/

/*求欧拉函数*/
void GetEulerValue(mpz_t p, mpz_t q,mpz_t *euler)
{
	mpz_t i, j,temp;
	mpz_init(i);                        //初始化i,j
	mpz_init(j);
	mpz_init_set_str(temp, "1", 10);    //temp初始化赋值为1(把字符串初始化为gmp大整数)
	mpz_sub(i, p, temp);                //i = p - temp，即i=p-1
	mpz_sub(j, q, temp);                //j = q - temp，即j=q-1
	mpz_mul(*euler, j, i);              //*euler=j * i            
	mpz_clear(i);
	mpz_clear(j);
	mpz_clear(temp);
}

/*随机生成大素数，参考：https://www.cnblogs.com/y3w3l/p/5947450.html */
void CreatePrime(mpz_t *rand_p, mpz_t *rand_q, int judge, mpz_t euler,int range)
{
	unsigned long seed = time(NULL);
	mpz_t temp, i;

	mpz_init(temp);
	mpz_init_set_str(i, "1", 10);

	//对一个随机状态的定义即初始化
	gmp_randstate_t state;
	gmp_randinit_default(state);

	gmp_randseed_ui(state, seed);           //设置随机状态的种子

	int BOOL = 0;
	switch (judge)
	{
	case 0:
		while (BOOL == 0)
		{
			mpz_urandomb(*rand_p, state, range);   //生成随机数	
			BOOL = mpz_probab_prime_p(*rand_p, 15);//判断是否是素数，返回BOOL==1代表是素数
		}
		BOOL = 0;
		while (BOOL == 0)
		{
			mpz_nextprime(*rand_q, *rand_p);
			BOOL = mpz_probab_prime_p(*rand_q, 15);//判断是否是素数，返回BOOL==1代表是素数
		}
		break;
	case 1:
		while (BOOL == 0)
		{
			//生成随机数
			mpz_urandomm(*rand_p, state, euler);
			//求rand和euler的最大公约数
			mpz_gcd(temp, *rand_p, euler);
			//判断是否满足互素条件
			BOOL = (mpz_cmp(temp, i)==0);
		}
		break;
	}
	mpz_clear(temp);
	mpz_clear(i);
}

/*快速幂算法学习：https://blog.csdn.net/chen_zan_yu_/article/details/90522763?ops_request_misc=&request_id=&biz_id=102&utm_term=%E5%BF%AB%E9%80%9F%E5%B9%82%E5%8F%96%E6%A8%A1%E7%AE%97%E6%B3%95&utm_medium=distribute.pc_search_result.none-task-blog-2~all~sobaiduweb~default-3-.pc_search_result_before_js&spm=1018.2226.3001.4187
int PowerMod(int a, int b, int c)
{
	int ans = 1;
	a = a % c;
	while(b>0) {
		if(b % 2 = = 1)
		ans = (ans * a) % c;
		b = b/2;
		a = (a * a) % c;
	}
	return ans;
}*/

void Encrypt(mpz_t *result,mpz_t n, mpz_t m, mpz_t e)
{
	mpz_t i,j,k,l;                  //为了参与gmp库的大数运算，需要一些高精度变量存储一些数据

	mpz_init_set_str(i, "2", 10);     //用i代表常数2
	mpz_init_set_str(j, "1", 10);     //用j代表常数1
	mpz_init_set_str(l, "0", 10);     //用l代表常数0
	mpz_init(k);
	mpz_set(*result, j);

	mpz_mod(m, m, n);             //(a*b) mod c = ((a mod c)*(b mod c)) mod c
	while (mpz_cmp(e, l)>0) 
	{
		
		if (mpz_odd_p(e))             //判断是否是奇数 
		{
			/* result=(result*temp) mod euler */
			mpz_mul(*result, *result, m);
			mpz_mod(*result, *result, n);
		}
		mpz_div(e, e, i);
		mpz_mul(m, m, m);
		mpz_mod(m, m, n);
	}
	mpz_clear(i);
	mpz_clear(j);
	mpz_clear(k);
	mpz_clear(l);
}

/*
扩展欧几里得算法参考学习：
https://blog.csdn.net/qq_16964363/article/details/79251433?ops_request_misc=&request_id=&biz_id=102&utm_term=%E6%B1%82%E4%B9%98%E6%B3%95%E9%80%86%E5%85%83&utm_medium=distribute.pc_search_result.none-task-blog-2~all~sobaiduweb~default-6-.nonecase&spm=1018.2226.3001.4187
https://blog.csdn.net/my_acm/article/details/38423553?ops_request_misc=&request_id=&biz_id=102&utm_term=%E6%B1%82%E4%B9%98%E6%B3%95%E9%80%86%E5%85%83&utm_medium=distribute.pc_search_result.none-task-blog-2~all~sobaiduweb~default-3-.nonecase&spm=1018.2226.3001.4187
https://blog.csdn.net/qq_21120027/article/details/52921404?ops_request_misc=%257B%2522request%255Fid%2522%253A%2522162201937316780262575134%2522%252C%2522scm%2522%253A%252220140713.130102334..%2522%257D&request_id=162201937316780262575134&biz_id=0&utm_medium=distribute.pc_search_result.none-task-blog-2~all~sobaiduend~default-1-52921404.pc_search_result_before_js&utm_term=%E6%B1%82%E4%B9%98%E6%B3%95%E9%80%86%E5%85%83&spm=1018.2226.3001.4187
https://www.jianshu.com/p/9a5f60efeb0d
*/

void Get_d(mpz_t e, mpz_t euler,mpz_t *temp,mpz_t *a,mpz_t *b)
{
	mpz_t zero,one,i;
	mpz_init_set_str(zero, "0", 10);
	mpz_init_set_str(one, "1", 10);
	mpz_init(i);
	if (mpz_cmp(zero, euler) == 0)
	{
		mpz_set(*temp, e);
		mpz_set(*a, one);
		mpz_set(*b, zero);
		return;
	}
	mpz_mod(i, e, euler);
	Get_d(euler,i,temp,b,a);
	mpz_div(i, e, euler);
	mpz_mul(i, i, *a);
	mpz_sub(*b, *b, i);
	mpz_clear(zero);
	mpz_clear(one);
	mpz_clear(i);
}
int main()
{
	system("title MADE BY YQC");
	char *temp_p = (char*)malloc(sizeof(char*)*N),
		 *temp_q = (char*)malloc(sizeof(char*)*N),
		 *temp_m = (char*)malloc(sizeof(char*)*N), 
		 *temp_n = (char*)malloc(sizeof(char*)*N),
		 *temp_e = (char*)malloc(sizeof(char*)*N);
	mpz_t p, q, e, d, euler, c, m, x, y, temp,n;
	int base=10;                                   //默认10进制输入

	mpz_init(euler);
	mpz_init(c);
	mpz_init(d);
	mpz_init(x);
	mpz_init(y);
	mpz_init(temp);
	mpz_init(n);

	int operation,range,ch='\0';
	while (1)
	{
		system("CLS");
		printf("（注意：本程序读入明文密文密钥默认都是以数据读入，如果想对字符串加密，需要额外转化成ASCII值再输入，例如想对'0'加密，输入的明文应该是'0'对应的ASCII值48）\n1.随机生成大素数p,q,e\n2.手动输入大素数p,q,e\n3.退出程序\n请选择要进行的操作:");
		scanf("%d", &operation);
		switch (operation)
		{
		case 1:
			system("CLS");
			mpz_init_set_str(p,"0",10);
			mpz_init_set_str(q,"0",10);
			mpz_init_set_str(e,"0",10);
			printf("请输入明文m:");
			scanf("%s", temp_m);
			mpz_init_set_str(m, temp_m, base);

			CreatePrime(&p, &q, 0, euler,M);
			gmp_printf("\n┑(￣Д ￣)┍:已生成测试大素数p:\n%Zd\n", p);
			gmp_printf("┑(￣Д ￣)┍:已生成测试大素数q:\n%Zd\n", q);
			mpz_mul(n, p, q);
			gmp_printf("┑(￣Д ￣)┍:计算出n=\n%Zd\n", n);
			GetEulerValue(p, q, &euler);
			CreatePrime(&e, &q, 1, euler,M);
			gmp_printf("┑(￣Д ￣)┍:已生成测试公钥e:\n%Zd\n", e);
			Get_d(e, euler, &temp, &x, &y);
			mpz_add(d, euler, x);
			mpz_mod(d, d, euler);
			gmp_printf("┑(￣Д ￣)┍:计算出私钥d为:\n%Zd\n", d);
			Encrypt(&c, n, m, e);
			gmp_printf("┑(￣Д ￣)┍:计算出密文c为:\n%Zd\n", c);
			label_3:
			gmp_printf("┑(￣Д ￣)┍:是否继续(Y/N):\n");
			scanf_s("%c", &ch);
			if(ch=='y'||ch=='Y') break;
			else goto label_3;
		case 2:
			system("CLS");
			printf("请输入明文m:");
			scanf("%s", temp_m);
			printf("请输入全局数据采用的进制:");
			scanf("%d", &base);
			mpz_init_set_str(m, temp_m, base);

			printf("1.（私钥获取模式）输入p和q，这将帮你算出私钥和密文\n2.（加密机模式）输入n，这将直接算出密文\n请选择模式:");
			scanf("%d",&operation);
			switch (operation)
			{
			case 1:
				label_1:
				gmp_printf("请输入大素数p:");
				scanf("%s", temp_p);
				mpz_init_set_str(p, temp_p, base);
				if (!mpz_probab_prime_p(p, 10))  goto label_p;

				label_2:
				gmp_printf("请输入大素数q:");
				scanf("%s", temp_q);
				mpz_init_set_str(q, temp_q, base);
				if (!mpz_probab_prime_p(q, 10))  goto label_q;

				mpz_mul(n, p, q);
				gmp_printf("┑(￣Д ￣)┍:已计算出n=%Zd\n", n);

				GetEulerValue(p, q, &euler);
				gmp_printf("┑(￣Д ￣)┍:已计算出n的欧拉函数Euler=%Zd\n", euler);

				gmp_printf("请输入公钥e【确保gcd(Euler,e)==1】:");
				scanf("%s", temp_e);
				mpz_init_set_str(e, temp_e, base);

				Get_d(e, euler, &temp, &x, &y);
				mpz_add(d, euler, x);
				mpz_mod(d, d, euler);
				gmp_printf("┑(￣Д ￣)┍:已计算出私钥d为：%Zd\n", d);

				Encrypt(&c, n, m, e);
				gmp_printf("┑(￣Д ￣)┍:最终计算出密文c为:%Zd\n", c);
				break;
			case 2:
				gmp_printf("请输入公钥n:");
				scanf("%s", temp_n);
				mpz_init_set_str(n, temp_n, base);

				gmp_printf("请输入公钥e【确保gcd(Euler,e)==1】:");
				scanf("%s", temp_e);
				mpz_init_set_str(e, temp_e, base);

				Encrypt(&c, n, m, e);
				gmp_printf("┑(￣Д ￣)┍:最终计算出密文c为:%Zd\n", c);
				break;
			}
		label_4:
			gmp_printf("是否继续（Y\N）:");
			scanf("%c", &ch);
			if (ch == 'y' || ch == 'Y') break;
			else goto label_4;
		case 3:
			mpz_clear(p);
			mpz_clear(q);
			mpz_clear(e);
			mpz_clear(d);
			mpz_clear(euler);
			mpz_clear(c);
			mpz_clear(m);
			mpz_clear(x);
			mpz_clear(y);
			mpz_clear(temp);
			mpz_clear(n);
			free(temp_m);
			free(temp_p);
			free(temp_q);
			free(temp_e);
			free(temp_n);
			exit(0);
			break;
		default:
			printf("\n输入有误！！！");
			break;
		}
	}
label_p:
	printf("输入的p不是素数，请重新输入\n");
	goto label_1;
label_q:
	printf("输入的q不是素数，请重新输入\n");
	goto label_2;
	return 0;
}

四、实现效果的呈现与加密结果的验证

（一）随机生成密钥对并自动加密

（二）自定义密钥，对明文加密

测试明文m：48（'0'）
测试私钥p：176845155422953166325517120209549292739
测试私钥q：309450490878749684523977649682267771163
生成$n=p\cdot q=54724820155161639103996662503795892597711455603713962373761530076820549485457$
测试公钥e：13

自定义输入包含以下两种模式：

• 选择自行设计p和q，然后进行加密
• 选择直接输入n加密，也就是加密机
  
  得到密文数据7180192468708211294208。

用RSA-Tool验证如下，