学习笔记 Stein算法

介绍

Stein算法是一种计算两个数最大公约数的算法，是针对欧几里德算法在对大整数进行运算时，需要试商导致增加运算时间的缺陷而提出的改进算法。

欧几里得算法（辗转相除法）缺陷

欧几里德算法是计算两个数最大公约数的传统算法，无论从理论还是从实际效率上都是很好的。但是却有一个致命的缺陷，这个缺陷在素数比较小的时候一般是感觉不到的，只有在大素数时才会显现出来。

一般实际应用中的整数很少会超过64位（当然已经允许128位了），对于这样的整数，计算两个数之间的模是很简单的。对于字长为32位的平台，计算两个不超过32位的整数的模，只需要一个指令周期，而计算64位以下的整数模，也不过几个周期而已。但是对于更大的素数，这样的计算过程就不得不由用户来设计，为了计算两个超过64位的整数的模，用户也许不得不采用类似于多位数除法手算过程中的试商法，这个过程不但复杂，而且消耗了很多CPU时间。对于现代密码算法，要求计算128位以上的素数的情况比比皆是，设计这样的程序迫切希望能够抛弃除法和取模。

基本原理

1.gcd（a，a）=a：一个数和其自身的最大公约数仍是其本身

2.gcd(ka,kb)=k*gcd(a,b):最大公约数运算和倍乘运算可以交换，取特殊值k=2时，说明两个偶数的最大公约数必然能被2整除

3.gcd(ka,b)=gcd(a,b):当ka与b互为质数时，约掉两个整数中只有其中一个整数含有的因子时，不会影响这两个数的最大公约数。取特殊值k=2时，说明求一个偶数和一个奇数的最大公约数时，可以先将偶数除以2

（插入概念：互为质数 一般指互质数。 互质数为数学中的一种概念，即两个或多个整数的公因数只有1的非零自然数。公因数只有1的两个非零自然数，叫做互质数。）

算法步骤

1.设整数A，B

2.若A=B，最大公约是就是其本身

3.若A=0, B就是最大公约数；反之亦然

4.若A,B都是偶数，最大公约数：2*gcd(A/2,B/2)

5.若A为even、B为odd，最大公约数为：gcd(A/2,B)；反之亦然

6.若AB都是odd，使用更相损减法：gcd(abs(a-b),min(a,b))；由于AB为odd,A-B必然为一个偶数，于是跳转第五步，优化为：gcd(abs(a-b)/2,min(a.b))

代码实现

#include<stdio.h>
#include<stdlib.h>                                //包含abs函数
#define MIN(a,b) a < b ? a : b
int gcd(int n, int m)
{
	if (n == m)
	{
		return n;
	}
	if (n == 0)
	{
		return m;
	}
	if (m == 0)
	{
		return n;
	}
	if (n % 2 == 1)
	{
		if (m % 2 == 1)
		{
			return gcd(abs(n-m) / 2, MIN(n, m));
		}
		else
		{
			return gcd(n, m / 2);
		}
	}
	else
	{
		if (m % 2 == 1)
		{
			return gcd(n / 2, m);
		}
		else
		{
			return 2 * gcd(n / 2, m / 2);
		}
	}
}
int main()
{
	int n = 0, m = 0;
	scanf("%d %d", &n, &m);
	printf("%d", gcd(n, m));
	return 0;
}

优化代码

通过位移运算符达到提高运算速度的效果：因为位移运算符直接作用于二进制补码，所以速度会提高很多

#include<stdio.h>
#include<stdlib.h>
#define MIN(a,b) a < b ? a : b
int gcd(int n, int m)
{
	if (n == m)
	{
		return n;
	}
	if (n == 0)
	{
		return m;
	}
	if (m == 0)
	{
		return n;
	}
	if (n % 2 == 1)
	{
		if (m % 2 == 1)
		{
			return gcd(abs(n - m) >> 1, MIN(n, m));
		}
		else
		{
			return gcd(n, m >> 1);
		}
	}
	else
	{
		if (m % 2 == 1)
		{
			return gcd(n >> 1, m);
		}
		else
		{
			return  gcd(n >> 1, m >> 1) << 1;
		}
	}
}
int main()
{
	int n = 0, m = 0;
	scanf("%d %d", &n, &m);
	printf("%d", gcd(n, m));
	return 0;
}