信息学竞赛中数论常见问题.doc

Prime
概述:
>1的数,除了1和本身没有其他因子.
1既不是素数也不是合数,0和所有的负整数同样如此.
100以内的素数
{2,3,5,7,11,13,17,19,23,29,31,37,41,43,47,53,59,61,67,71,73,79,83,89,97}
1.打表(1)
#define n 1000003
int prime[n];

//要判断的范围是多大,就要设置多大的数组,对应于某个数,看它下标对应的数是1，则是素数,否则,不是素数
void getPrime(){
int i,j;
int bound=sqrt((double)n);
for(i=2;i for (j=0;j0)return false;
double bound=sqrt((double)x);
for(int i=3;i<=bound;i+=2){
if(x%i==0){
return false;
}
}
return true;

}

GCD(最大公约数)---欧几里得算法
a,b为正整数，设集合A = {x*a+y*b|x,y是整数}，则A中最小正元素是gcd(a,b)

long kgcd(long a,long b)
{
if (a==0)return b;
if (b==0)return a;
if (!(a&1)&&!(b&1))return kgcd(a>>1,b>>1)<<1;
else if (!(b&1))return kgcd(a,b>>1);
else if (!(a&1))return kgcd(a>>1,b);
else return kgcd(abs(a-b),min(a,b));

}
LCM (最小公倍数)
LCM ( a, b ) = a * b / GCD ( a, b )
实际上最好写成a/GCD(a,b)*b
long lcm(long a,long b)
{
long c,dsw;
c=(a>=b)?a:b;
d=(a<=b)?a:b;
while (c%d!=0)
{
sw=c%d;
c=d;
d=sw;
}
return (a/d)*b;
}
求多个数的lcm,需要将res初始化为1
相关题目:
1. Least Common Multiple
2. wolf and rabbit

扩展欧几里得(可以求解模线性方程)
a*x+b*y=gcd(a,b)一系列解是x=x+b,y=y-a
int extgcd(int a,int b,int &x,int &y){
if(b==0){
x=1;
y=0;
return a;
}
int d=extgcd(b,a%b,x,y);
int t=x;
x=y;
y=t-a/b*y;
return d;
}
a*x+n*y=b==ax=b(mod n)模线性方程
定理：方程ax=b(mod n)对于未知量x有解，当且仅当gcd(a,n)|b.
令d=gcd(a,n).
设为x' ,y'为所求 满足ax'+ny'=gcd(a,n).
则原方程有一解x0=x'*(b/d)modn.
定理:假设ax=b(mod n) 有解,x0是该方程的任意一个解，则该方程模n恰有d个不同的解，分别为：xi=(x0+i*(n/d))mod n (其中i=1,2,...d-1).
// a*x=b(%n)
void modeq(__int64 a,__int64 b,__int64 n)
{
__int64 e,d,x,y,t;
d=extgcd(a,n,x,y);
if(b%d)
printf("Impossible\n");
else{
e=(x*(b/d))%n+n;
t=n/d;
t=t>0?t:-t;
e%=t;
if(e<0)
e+=t;
printf("%I64d\n",e);
}
}题目:
1.

若n=p1e1p2e2…prer，则
n的因数个数为
(1+e1)*(1+e2)*……(1+er)
n所有因数的和为(1+p1+p12+…+p1e1)*(1+p2+p22+…+p2e2) *…*(1+pr+pr2+…+prer)

数的各位之和
int sum( int number )
{
int sum = 0;

while( number != 0 )
{
sum += number % 10;
number /= 10;
}
return sum;
}//不知道数有几位,但是可以每次都取个位

质因数分解
对于i,从2到sqrt((double)m),如果m%i==0,则i是m的一个质因数,
然后m/i=m,i=2,重新判断,直至最后,剩下的那个数同样是一个质因数,要加上的.

int zhi(int n,int *x){
int count=0;
int i=2;
while(i<=(int)sqrt((double)n)){
ifn%i==0){
x[count++]=i;
n=n/i;
}else{
i++;
}

}
x[count++]=n;
return count;
}//返回质因数的个数,数组x是所有的因子

//n是要分解的数,pn是100范围内质数的个数,prim是100范围内的所有质数,q存储分解出来的质数相应的个数
int zhi_num(int n,int pn,int *prim,int *q)
{
for (int j=0;j if (n%prim[j]==0)
{
__int64 temp=0;
while (n%prim[j]==0)
{
temp++;//相应的prim[j]的个数
n/=prim[j];
}
q[j]+=temp;
}
}
}
调用:
prim[25]={2,3,5,7,11,13,17,19,23,29,31,37,41,43,4753,59,61,67,71,73,79,83,89,97};
int q[25]={0};
int pn=25;
zhi_num(20,pn,prim,q);
欧拉函数
欧拉函数值是小于或等于n的数中与n互质的数的个数 Phi(N)同样如此
公式:
φ(x)=x(1-1/p1)(1-1/p2)(1-1/p3)(1-1/p4)…..(1-1/pn),
其中p1, p2……pn为x的所有质因数，x是不为0的整数。
性质有:
1. φ(1)=1（唯一和1互质的数就是1本身）。
2. 若n是素数,则 φ(n)=n-1;
3. 若n是合数,则 φ(n)n为奇数时，φ(2n)=φ(n),
5. 若a,b互质,则, φ(a*b)=φ(a)*φ(b);
代码:
int euler(int x){
int i,res=x;
for(i=2;i<(int)sqrt(x*1.0)+1;i++){
if(x%i==0){
res=res/i*(i-1);
while(x%i==0)x/=i;
}
}
if(x>1)res=res/x*(x-1);
return res;
}

中国剩余定理
不管互质或者是不互质,只要在输入时进行互质处理,即可.而且,是_int64显然通用性更强
下面程序的输入输出格式:
3//每组有几组数据
3 2//除数 余数
5 3
7 2
#include
using namespace std;
__int64 x,y,t;
__int64 extend_gcd(__int64 a,__int64 b) //return gcd(a,b) 得到x,y
{
if(b == 0)
{
x = 1;
y = 0;
return a;
}
else
{
__int64 temp = extend_gcd(b,a%b);
t = x; x = y;
y = t - a/b*y;
return temp;
}
}
int main()
{
__int64 a1,a2,r1,r2,c,gcd,temp;
bool yes;
int n;
freopen("G://in.txt","r",stdin);
while(scanf("%d",&n) !=EOF)//每组有一组数据
{
yes = false;
scanf("%I64d %I64d",&a1 ,&r1);//a是除数,r是余数
for(int i=0 ;i yes = true;
continue;
}
temp = a2/gcd;
x = (c/gcd*x % temp + temp)%temp;
r1 = a1*x + r1;
a1 = a1*a2/gcd;
}
if(yes)
printf("-1\n");
else
printf("%I64d\n",r1);

}
return 0;
}

闰年是否并返回一年天数
/* 判断是否闰年 */
bool isleap(int& year)
{
if (year % 4 == 0 && year % 100 != 0 || year % 400 == 0) return true;
else return false;
}

/* 返回一年的最大天数 */
int maxday(int& year)
{
if (isleap(year)) return 366;
else return 365;
}
Days是指距离某个日期是多少天,应该均可以的,只是最终结果可能有所变化的.
string getweek(int& days)
{
return week[days % 7];
}
一个关于日历的题目
#include
#include

using namespace std;

string week[] = {"Saturday", "Sunday", "Monday", "Tuesday", "Wednesday", "Thursday", "Friday"};
int day[] = {0, 31, 28, 31, 30, 31, 30, 31, 31, 30, 31, 30, 31};

/* 判断是否闰年 */
bool isleap(int& year)
{
if (year % 4 == 0 && year % 100 != 0 || year % 400 == 0) return true;
else return false;
}

/* 返回一年的最大天数 */
int maxday(int& year)
{
if (isleap(year)) return 366;
else return 365;
}

/* 得到年份 */
int getyear(int& days)
{
int year = 2000;
while (days > maxday(year))
{
days -= maxday(year);
year++;
}
return year;
}

/* 得到月份 */
int getmonth(int year, int& days)
{
int month = 1;
if (isleap(year)) day[2] = 29;
else day[2] = 28;
while (days > day[month])
{
days -= day[month];
month++;
}
return month;
}

/* 得到天数 */
int getday(int& days)
{
return days;
}

/* 得到星期 */
string getweek(int& days)
{
return week[days % 7];
}

/* 主函数 */
int main()
{
int days;
int y, m, d;
string wk;
while(cin >> days && days != -1)
{
wk = getweek(days);
days++;
y = getyear(days);
m = getmonth(y, days);
d = getday(days);

cout << y << "-";
if (m < 10) cout << 0;
cout << m << "-";
if (d < 10) cout << 0;
cout << d << " " << wk << endl;
}
//system("pause");
return 0;
}
生成两个大的素数P,Q,乘起来得N=P*Q.然后算出N的欧拉函数Phi(N)=(P-1)*(Q-1).然后我们取一个范围在[1,phi(N)]中且与phi(N)互质的正整数E.它就是所谓的公钥。得到公钥之后，我们再算出E关于phi(N)的逆元D,即E*D mod phi(N)=1.这个D就是私钥。在得到这些数据以后，P,Q被丢弃，E,N做为公钥被公开，D做为私钥被解密人私人保存。
假设有一个明文M,那么它所对应的密文就是C=M^E mod N.
如果我们现在得到一个密文C,那么它所对应的明文就是M=C^D mod N
也就是说，任何人都可以用公钥对数据进行加密，但是只有拥有私钥的人才可以对数据进行解密。

将N分解成P*Q的乘积。那么就可以直接利用公式phi(N)=(P-1)*(Q-1)绕开暴力求解欧拉函数的过程，从而实现RSA的破解。

这道题就是模拟这个破解过程，下面来说说具体的做法：
1.首先用miller-rabin,pollard_rho做大整数的质因数分解，得到两个素数P,Q，pollard_rho的复杂度在N^0.25次方，那么一个64位的整数 要计算的次数为 2^64^0.25=2^16 =65536次，可以瞬间出解。
2.求出phi(N)=(P-1)*(Q-1)
3.然后用ext_gcd求出E关于phi(N)的逆元。
4.用得到的私钥对数据C进行解密即可。
PS:对这题而言，仅仅完成上述步骤还是不够的。因为N达到2^62次方，即使是使用无符号long long ,也很容易因为出乘法操作而溢出。这也是网上说要避免使用扩展欧几里德的原因。其实实现的时候，我们可以自己写一个特殊的乘法（内部用加法实现），由于使用的无符号的long long ，第64位刚好可以用来保存两个数加过之后的进位位，再模除又可以保证其在2^62范围内，即可避免溢出

#include
#include
#include
#include
#include
using namespace std;
#define bigint __int64//这道题目不能用无符号数，不然wa
//欧几里德，求最大公约数
bigint gcd(bigint a,bigint b)
{
while (b)
{
bigint c = a % b;
a = b;
b = c;
}
return a;
}
//a*b%n
bigint product_mod(bigint a,bigint b,bigint n)
{
bigint tmp = 0;
while (b)
{
if (b&1)
{
tmp += a;
if (tmp>=n) tmp-=n;
}
a<<=1;
if (a>=n) a-=n;
b>>=1;
}
return tmp;
}
//a^m%n
bigint power_mod(bigint a,bigint m,bigint n)
{
bigint tmp = 1;
a%=n;
while (m)
{
if (m&1) tmp = product_mod(tmp,a,n);
a = product_mod(a,a,n);
m>>=1;
}
return tmp;
}
int pri[] = {2,3,5,7,11,13,17,19,23,29};
//MillerRabin大素数判断
bool Miller_Rabin(bigint n)
{//n,s.取s个随机数值a,进行a是n为和数的证明判断
if (n<2) return false;
if (n==2) return true;
if (!(n&1)) return false;
bigint k = 0 , i ,j ,m ,a;
m = n - 1;
while (!(m&1)) m>>=1,k++;
for (i = 0 ; i < 10 ; i ++)
{
if (pri[i]>=n) return true;
a = power_mod(pri[i],m,n);//幂取模
if (a==1) continue;
for (j = 0 ; j < k ; j ++)
{
if (a==n-1) break;
a = product_mod(a,a,n);
}
if (j==k) return false;
}
return true;
}
//pollard_rho随机质因数分解
bigint pollard_rho(bigint C, bigint N) //返回一个平凡因子
{
bigint I, X, Y, K, D;
I = 1;
X = Y = rand() % N;
K = 2;
do
{
I++;
D = gcd(N + Y - X, N);//这里为了防止负数，先加上一个N
if (D > 1 && D < N) return D;//如果D不是非平凡因子
if (I == K) Y = X, K <<= 1;
X = (product_mod(X, X, N) + N - C) % N;//随机一个增量来求X,Y，使得gcd(Y - X,N)不是非平凡因子(不是1，跟N)
}
while (Y != X);
return N;
}
//二分，分解N的质因数,这里返回最小质因子，如果想要其他的话，可以用个数组存起来
bigint rho(bigint N)
{
if (Miller_Rabin(N))
return N;
bigint T = N;
while (T>=N)
T = pollard_rho(rand()%(N - 1) + 1,N);
bigint A = rho(T);
bigint B = rho(N / T);
return Ad,y,x);
y -= x * (a/b);
}
//求a关于n的逆
bigint inv(bigint a,bigint n)
{
bigint d,x,y;
Gcd(a,n,d,x,y);
if (d==1) return ( x % n + n ) % n;
else return -1;
}
bigint P,Q;// N = P * Q
bigint T;//T = (P - 1) * (Q - 1)
bigint C,E,N;

int main()
{
freopen("G://in","r",stdin);
srand((unsigned int)time(NULL));
//已知的是密文,公钥,N (= P * Q).
while (scanf("%I64d%I64d%I64d",&C,&E,&N)!=EOF)
{
P = rho(N);//找两个大素数,使得P*Q=N.现在是已知N,去求这两个大素数
Q = N / P;
T = (P - 1) * (Q - 1);//其实T就是phi[N],N的欧拉函数，因为P,Q都是素数，所以(P-1)*(Q-1)就是他的欧拉函数
bigint D = inv(E,T);//求E关于T的逆元,因为这里是gcd(E,T)==1,而又要E*D%T==1
bigint M = power_mod(C,D,N);//M = (C^D)%N,如果是求C,就是C=(M^E)%N
printf("%I64d\n",M);
}
}