数论及其应用——同余式定理
这篇文章我们将介绍数论当中几个很重要的定理：威尔逊定理、费马小定理以及欧拉定理，并讨论一些基于这些定理的算法。
首先我们给出费马小定理：如果p是素数，并且gcd(a,p) = 1 , 那么有a^(p-1) = 1(mod p)。
而关于这个定理的证明，也是不难理解的。
在证明之前，我们先需要知道这样两个引理。
引理1：若a、b、c为任意3个整数，且gcd(m,c) = 1,则当ac = bc(mod m)时，有a = b (mod m)。
引理2：设m是一个整数,且m>1,b是一个整数且gcd(m,b)=1.如果a1,a2,a3,a4,…am是模m的一个完全剩余系,那么ba[1],ba[2],ba[3],ba[4],…ba[m]也构成模m的一个完全剩余系。
这里需要解释一下的是，所谓完全剩余系，即{a1,a2,a3,a4,…am}中任意一个元素ai % m != 0。针对引理1，通过同余式本身的特点，都十分易证，随后再基于引理1，引理2也不难得证。这里便不再累述。 有了这两个引理做铺垫，下面便可以开始费马小引理的证明了。
构造素数p的完全剩余系P = {1,2,3,……p-1}。
因为gcd(a,p) = 1,由引理2得素数p的新的完全剩余系：A = {a,2a,3a,……(p-1)a}。
我们选出A、P中第i个元素，由引理1得，Ai = Pi(mod p)，因此我们容易得到如下等式。
1 * 2 * 3 *……(p-1) = a * 2a * 3a * ……(p-1)a (mod p)
等式两边同除(p - 1)!，得到a^(p-1) = 1 (mod p)，定理得证。
既然知道了费马小定理的内容和证明了，这里我们就用它进行一下拓展利用。
首先我们知道，在费马小定理的等式a^(p - 1) = 1 (mod p)中，如果p是素数，那么等式是成立的。
即 isprime(p) = 1 => a^(p - 1) = 1 (mod p) 。 那么我们根据逆否命题与原命题的等价性，可以得到，a^(p - 1) != 1 (mod p) => isprime(p) = 0。
那么逆命题和否命题与原命题是否等价呢？如果等价，我们岂不是可以通过判断是否满足费马小定理的等式来判断p是否为素数了？
实践表明，满足费马小定理的a、p组合中，p是可以为合数的，不妨自己举几个例子。 因此对于满足费马小定理等式的a、p，有这样的定义。如果p是合数，则称p是以a为基的伪素数。
我们不妨通过一个题目来理解一下所谓伪素数的概念。(Problem source:pku 3641)
Description
Fermat's theorem states that for any prime number p and for any integer a > 1, ap = a (mod p). That is, if we raise a to the pth power and divide by p, the remainder is a. Some (but not very many) non-prime values of p, known as base-a pseudoprimes, have this property for some a. (And some, known as Carmichael Numbers, are base-a pseudoprimes for all a.)
Given 2 < p ≤ 1000000000 and 1 < a < p, determine whether or not p is a base-a pseudoprime.
Input
Input contains several test cases followed by a line containing "0 0". Each test case consists of a line containing p and a.
Output
For each test case, output "yes" if p is a base-a pseudoprime; otherwise output "no".
题目大意：基于伪素数的概念，给出a、p，请你判断p是否是以a为基的伪素数。
数理分析：有了上文对伪素数概念的较少，我们知道，伪素数首先必须是合数，然后需要满足费马小定理。
编程实现：在判断是否是伪素数的时候，由于是判断单个数字，所以使用朴素的判断法即可。而判断是否满足费马小定理的时候，则用到了我们在之前讨论过的乘法快速幂算法。 参考代码如下。
#include<iostream> #include<math.h> using namespace std; bool ifprime(long long n) { long long i; int iffind = 0; if(n == 2 || n == 1) return 1; else { for(i = 2;i <= sqrt(n + 0.0);i++) if(n%i == 0) { iffind = 1; break; } if(iffind) return 0; else return 1; } } long long quick_mod(long long a , long long b , long long m) //快速幂取模，计算a^b(mod m)的值 { long long ans = 1; while(b) { if(b&1) { ans = (ans * a)%m; b--; } b /= 2; a = a * a %m; } return ans; } int main() { long long n , a, p ,b; while(cin >> p >> a) { if(p == 0 && a == 0) break; if(ifprime(p)) cout<<"no"<<endl; else { long long ans = quick_mod(a,p,p); if(ans == a) cout<<"yes"<<endl; else cout<<"no"<<endl; } } return 0; }
通过上文对费马小定理初步的理解，我们这里介绍在其基础上建立的一种更加高效的素数检测法——Miller-Rabin素数测试法。
通过上面的介绍，我们已经知道，满足费马小定理仅仅是一个数是素数的必要条件，那么对于一个满足费马小定理的数字，我们再通过什么方法来判断它是素数还是合数呢？
可能有人会想到之前学过的筛法，但是还是不够高效，在判断较大的数是否是素数的时候，只要一用筛法，时间复杂度就会激增，而避免这一情况的发生的办法就是想出线性时间复杂度的算法。
这里我们给出数论中一个线性时间复杂度判断素数的定理。
二次探测定理：如果p是一个素数，且0<x<p，则方程x^2 = 1(mod p)的解为x = 1或x = p - 1。
通过这个名字我们也可以看到，好像这个定理就是为了弥补利用费马小定理检测素数的诞生的。那么我们可以简单的概括这个新的检测素数的方法了。
step1：对于整数p，利用费马小定理判断。如果不满足，输出结果，如果满足，进入step2。
step2：利用二次探测定理判断，输出结果。 这就是基于费马小定理的改进版的Miller-Rabin素数检测法。由于两个定理中都包含未知量，所以再编程实现的时候我们都需要随机生成符合要求的参量，因此,Miller-Rabin算法本质上是一种随机算法。
那么下面我们开始讨论具体实现该算法的一些细节的处理。
虽然我们在描述算法步骤的时候，似乎将两个定理的应用分开了，但是通过观察两个定理的形式，会发现其实有着相似点，我们尝试将两个定理结合在一起。
假设我们想要判断n是否为素数，并随机生成了整数a(为了满足应用二次探测定理的条件),基于费马小定理，我们想要判断a^(n-1) = 1 (mod n)是否成立。
为了将此表达式更近得像二次探测定理靠拢，我们要找到含有平方的形式。那么我们就将n - 1化成 d * 2^r的形式(此处d为奇数)，这样费马小定理的左式就可以化成平方的形式了。
现在判断a^(d*2^r) = 1 (mod n)是否成立，如果成立，n是素数。 即(((a^d)^2)^2)^2)... = 1(mod n) ①
此时左式有r + 1层括号，并且这个等式同时包含了两个定理的形式。我们从最内层的a^d分析。
1.如果a^d = 1(mod n)，则①显然成立，n是素数。
2.如果a^d ≠ 1 (mod n)，结合二次探测定理，会有如下两种情况。
1'.a^d = n - 1,则(n-1)^2 = 1(mod n)，这就引导出了情况1.，n为素数。
2'.a^d ≠n - 1，等式不成立， n不是素数。
以上便完成了二次探测定理和费马小定理的巧妙结合。 我们通过一个题目来具体实现以下改进版的Miller-Rabin算法。(Problem source : hdu 2138)
题目大意：给出你一些数，让你判断素数的个数。
数理分析：由于题目给出的素数的范围很大，因此筛法就显得太过耗时，因此这里就可以用Miller-Rabin算法了。
编程实现：结合上文给出的分析，利用右移位运算来实现n - 1 = d*(2^r)的转换，并且还需要我们曾经讨论过的乘法快速幂的算法。
参考代码如下。
#include<stdio.h> #include<stdlib.h> #include<cmath> long long quick_mod(long long a,long long b,long long m) { long long ans = 1; while(b) { if(b&1) { ans = (ans * a) % m; b--; } b /= 2; a = a * a % m; } return ans; } bool witness(long long a,long long n) { long long m = n - 1; int j = 0; while(!(m&1)) //进行n-1 = d*(2^r)的转化 { j++; m >>= 1; } long long x = quick_mod(a, m , n); if(x == 1 || x == n -1) return false; while(j--) { x = x * x % n; if(x == n-1) return false; } return true; } bool miller_rabin(__int64 n) { if(n==2) return true; if(n==1 || ((n&1)==0)) return false; for(int i=0;i<50;i++)//随机生成次数为50次，该次数怎么得来的涉及一些概率方面的实验，在这里不作展开论述，会在以后关于随机算法的深入讨论中分析。 { __int64 a=rand()*(n-2)/RAND_MAX +1; //生成随机数a保证a<n if(witness(a, n)) return false; } return true; } int main() { int n,cnt; __int64 a; while(scanf("%d",&n)!=EOF) { cnt=0; while(n--) { scanf("%I64d",&a); if(miller_rabin(a)) cnt++; } printf("%d ",cnt); } return 0; }
基于上文对费马小定理的引入及其相关周边概念的介绍，我们再来通过一个题目来看看费马小定理在解决数论问题中有着怎样的应用。(Problem source : hdu 4704)
题目大意：定义函数s(k)，表示一个整数n被分成k个数的方案数，现在给定n，求解∑s(i) mod (1e9 + 7)的值，其中i∈[1,n]。
数理分析：这是一道比较综合的数论题目，其中涉及了很多的数论知识点。下面我们依次来击破。
首先我们分析函数s(i)的值，其功能是给出整数n被分成i个数字的方案数，我们不妨将整数n写成1+1+……+1的形式，这n个1中间形成了n-1个空间，我们插入i-1个隔板，那么我们整数n就被分成了i份，显然s(i) = C(n-1,i-1),这里运用到了组合数学中的鸽巢定理。
所以，我们所需要求的变成了∑C(n-1,i-1) , i∈[1,n]，由组合数学中二项式的相关性质，我们容易知道上式等于2^(n-1)。
但是基于n的取值太大，我们要考虑简化方案。
记m = 1e9 + 7。
当前我们所求的是：2^(n-1) = ? (mod m)。
通过费马小定理我们看到，a^(p-1) = 1(mod p)，当gcd(a,p)=1的时候。容易看到，gcd(a , m) = 1，那么有2^(m - 1) = 1(mod m)。
我们假设有n-1 = x(m-1) + q,那么基于同余等式的基本性质，2^(n-1) = 2^(x*(m - 1)) * 2^q (mod m) = 2^(x*(m-1)) (mod m) * 2^q(mod m) = 2^q(mod m)。
这样，一个很大的指数n-1，以上的等价转换，幂次就大大降低了，随后用快速幂取模来实现2^p (mod m)来得到最终结果即可。
总结起来，这道题目综合了组合数学的二项式的性质、鸽巢原理、费马小定理、快速幂取幂等多个算法，是一个很综合的数论题目。
参考代码如下。
#include<stdio.h> using namespace std; const int MAX = 100000 + 10; const int mod = 1000000000 + 7; char s[MAX]; long long MOD(char *a , int Mod) { long long sum = 0; for(int i = 0;a[i]!='