数论总结
2013-09-11 20:55
232 查看
1.欧几里得 求最大公约数,最小公倍数
(1)递归的写法:
(2)最小公倍数:
2.扩展欧几里得 求ax=b (mod m)=> ax+my=b
如果d=gcd(a,m)且b%d==0,
则同余方程有解,其最小解为x=x0*(b/d);
ax+by=c 如d=gcd(a,b),则存在x,y,使xa+yb=d;
当x+=b,y-=a后仍然成立
因为xa+yb+ab-ab=d;==>(x+b)a+(y-a)b=d
3.素数判定
(1)试除法:
(2)miller-rabin 算法
bool witness(i64 a,i64 n)
{
i64 x,d=1,i=ceil(log(n-1.0)/log(2.0))-1;
for(;i>=0;i--)
{
x=d; d=(d*d)%n;
if(d==1&&x!=1&&x!=n-1) return 1;
if(((n-1)&(1<<i))>0) d=(d*a)%n;
}
return d==1?0:1;
}
bool miller_rabin(i64 n)
{
if(n==2) return 1;
if(n==1||(n&1)==0) return 0;
i64 j,a;
for(j=0;j<50;j++)
{
a=rand()*(n-2)/RAND_MAX+1;
if(witness(a,n)) return 0;
}
return 1;
}
另一种写法,更好理解
bool witness(i64 a,i64 n)
{
int i,j=0;
i64 m=n-1,x,y;
while(m%2==0)
{
m>>=1;
j++;
}
x=pow(a,m,n);///快速幂取模
for(i=1;i<=j;i++)
{
y=pow(x,2,n);
if(y==1&&x!=1&&x!=n-1) return true;
x=y;
}
return y==1?false:true;
}
bool miller_rabin(i64 n)
{
if(n==2) return true;
if(n==1||n%2==0) return false;
for(int i=1;i<=10;i++)
{
i64 a=rand()%(n-1)+1;
if(witness(a,n)) return false;
}
return true;
}
4.素数筛法
(1) O(n)
(2) O(nlog(n))
(3)
素数的二次筛选
5.整数分解
(1)
(2)Pollard-pho大数分解
i64 Pollard(i64 n,int c)
{
i64 i=1,k=2,x=rand()%n,y=x,d;
srand(time(NULL));
while(true)
{
i++;
x=(mod_mult(x,x,n)+c)%n;
d=gcd(y-x,n);
if(d>1&&d<n) return d;
if(y==x) return n;
if(i==k) { y=x; k<<=1; }
}
}
6.
任意正整数都有且只有一种方式写出其素因子的乘积表达式:
A=(p1^k1)*(p2^k2)*(p3^k3)*....*(pn^kn) ,其中pi均为素数.
有A的所有因子之和为
S = (1+p1+p1^2+p1^3+...p1^k1) *(1+p2+p2^2+p2^3+….p2^k2) * (1+p3+ p3^3+…+p3^k3) * .... *(1+pn+pn^2+pn^3+...pn^kn)
(1).用递归二分求等比数列1+pi+pi^2+pi^3+...+pi^n:
1)若n为奇数,一共有偶数项,则:
1 + p + p^2 + p^3 +...+ p^n
= (1+p^(n/2+1))+ p * (1+p^(n/2+1)) +...+ p^(n/2) * (1+p^(n/2+1))
= (1+ p + p^2 +...+ p^(n/2)) * (1 +p^(n/2+1))
上式中颜色加粗的前半部分恰好就是原式的一半,那么只需要不断递归二分求和就可以了,后半部分为幂次式,将在下面第4点讲述计算方法。
2)若n为偶数,一共有奇数项,则:
1 + p + p^2 + p^3 +...+ p^n
= (1+p^(n/2+1))+ p * (1+p^(n/2+1)) +...+ p^(n/2-1) * (1+p^(n/2+1)) + p^(n/2)
= (1+ p + p^2 +...+ p^(n/2-1)) *(1+p^(n/2+1)) + p^(n/2);
上式颜色加粗的前半部分恰好就是原式的一半,依然递归求解
Code:
2) 求出1~n中所有数的欧拉函数值.
欧拉函数的递推关系:
假设素数p能整除n,那么
如果(n%p==0 && (n/p)%p==0) (p还能整除n / p) , PHI(n) = PHI(n / p) * p;
如果(n%p==0 && (n/p)%p==1) (p不能整除n / p), PHI(n) = PHI(n / p) * (p - 1);
可以用与筛选求素数非常类似的方法,在O(nloglogn)时间内计算完毕。
8.求逆元ax=1 (mod m) , x是a的逆元
=>ax+my=1=> gcd(a,m)=1
/*
int mod(int x,int n)
{
x%=n; //如x是负数,进行处理
if(x<0) x+=n;
return x; //return (x%n+n)%n;会越界
}
*/
用扩展欧几里得求:
9.
(1)幂取模( a^n mod m)
法一:
法二:用快速幂取模(a^n mod m )
(2)快速幂( a^k)
10.快速模乘 a*b%p (位运算的思路)
11.求解模线性方程组(中国剩余定理)
x=a1 mod m1
x=a2 mod m2
......
x=an mod mn 其中,a[],m[]已知,m[i]>0且m[i]与m[j]互质,求x.
同余模公式:
(a+b)%m=(a%m+b%m)%m
(a*b)%m=(a%m*b%m)%m
设m1,m2,...,mn是两两互素的正数,则对任意的整数a1,a2,...,an,同余方程组
其解为:X=((M_1*M1*a1)+(M_2*M2*a2)+...+(M_n*Mn*an)) mod m;
其中M=m1*m2*...*mn; Mi=M/mi; M_i是Mi的逆元(mod mi).
12.枚举出一个集合的子集。设原集合为mask,则下面的代码就可以列出它的所有子集:
for (i = mask ; i ; i = (i - 1) & mask) ;
(1)递归的写法:
int gcd(int a,int b) {return b?gcd(b,a%b):a;}
(2)最小公倍数:
int lcm(int a,int b) {return a/gcd(a,b)*b;}
2.扩展欧几里得 求ax=b (mod m)=> ax+my=b
如果d=gcd(a,m)且b%d==0,
则同余方程有解,其最小解为x=x0*(b/d);
ax+by=c 如d=gcd(a,b),则存在x,y,使xa+yb=d;
当x+=b,y-=a后仍然成立
因为xa+yb+ab-ab=d;==>(x+b)a+(y-a)b=d
long long ext_gcd(long long a, long long b, long long & x, long long &y) { if(!b) { x=1; y=0; return a; //d=a,x=1,y=0,此时等式d=ax+by成立 } longlong d=ext_gcd(b,a%b,x,y); long long xt=x; x=y; y=xt-a/b*y; return d; //对于拓展的欧几里德算法有:x=y0;y=x0-a/b*y0; }
3.素数判定
(1)试除法:
bool isprime(int n) { int temp=(int)(sqrt(n*1.0)+1); for(int i=2;i<=temp;i++) if(n%i==0) return false; return true; }
(2)miller-rabin 算法
bool witness(i64 a,i64 n)
{
i64 x,d=1,i=ceil(log(n-1.0)/log(2.0))-1;
for(;i>=0;i--)
{
x=d; d=(d*d)%n;
if(d==1&&x!=1&&x!=n-1) return 1;
if(((n-1)&(1<<i))>0) d=(d*a)%n;
}
return d==1?0:1;
}
bool miller_rabin(i64 n)
{
if(n==2) return 1;
if(n==1||(n&1)==0) return 0;
i64 j,a;
for(j=0;j<50;j++)
{
a=rand()*(n-2)/RAND_MAX+1;
if(witness(a,n)) return 0;
}
return 1;
}
另一种写法,更好理解
bool witness(i64 a,i64 n)
{
int i,j=0;
i64 m=n-1,x,y;
while(m%2==0)
{
m>>=1;
j++;
}
x=pow(a,m,n);///快速幂取模
for(i=1;i<=j;i++)
{
y=pow(x,2,n);
if(y==1&&x!=1&&x!=n-1) return true;
x=y;
}
return y==1?false:true;
}
bool miller_rabin(i64 n)
{
if(n==2) return true;
if(n==1||n%2==0) return false;
for(int i=1;i<=10;i++)
{
i64 a=rand()%(n-1)+1;
if(witness(a,n)) return false;
}
return true;
}
4.素数筛法
(1) O(n)
void get_prime() { memset(tag,0,sizeof(tag)); for(int i=2;i<N;i++) { if(!tag[i]) p[cnt++]=i; for(intj=0;j<cnt&&p[j]*i<N;j++) { tag[i*p[j]]=1;//tag[x]=1表示x不是素数 if(i%p[j]==0) break; } } }
(2) O(nlog(n))
void get_prime() //打表法求出所有的素数 { int m=sqrt(SIZE)+1; memset(yes,0,sizeof(yes)); for(int i=2;i<=m;i++) if(!yes[i]) { for(int j=i*i;j<SIZE;j+=i) yes[j]=1; } }
(3)
素数的二次筛选
/*********************************************** *intyes[SIZE]; *long longprim[SIZE]; *int k=0; //prim[] 中存放着k个素数 *intadp[1001000] //二次需要筛选的(U-L)< 1001000 *************************************************/ //获得基准素数 voidgetPrim() { memset(yes,0,sizeof(yes)); for(long long i=2;i<SIZE;i++) if(!yes[i]) { prim[k++]=i; for(long longj=i*i;j<SIZE;j+=i) yes[j]=1; } } //二次筛选范围为:L~U voidtwice_choose() { memset(adp,0,sizeof(adp)); for(int i=0;i<k&&prim[i]*prim[i]<=U;i++) { long longbeg=L/prim[i]+(L%prim[i]>0); if(beg==1) beg++; for(long longj=beg*prim[i];j<=U;j+=prim[i]) //调整j,使得L<=j<=U adp[j-L]=1; } }
5.整数分解
(1)
/******************************************** * prime[]数组保存筛选出的素数,长度为N; * p[]数组保存n的素数因子; * t[]保存的素数因子的个数,共tot个 *********************************************/ int tot=0; void split(intn,int *p,int *t) { for(int i=1;i<=N;i++) { int s=0; while(n%prime[i]==0) { s++; n/=prime[i]; } if(s) {p[tot]=prime[i]; t[tot]=s;tot++;} if(n==1) break; if(n<prime[i]*prime[i]) { p[tot]=n; t[tot]=1; n=1; break; } } }
(2)Pollard-pho大数分解
i64 Pollard(i64 n,int c)
{
i64 i=1,k=2,x=rand()%n,y=x,d;
srand(time(NULL));
while(true)
{
i++;
x=(mod_mult(x,x,n)+c)%n;
d=gcd(y-x,n);
if(d>1&&d<n) return d;
if(y==x) return n;
if(i==k) { y=x; k<<=1; }
}
}
6.
任意正整数都有且只有一种方式写出其素因子的乘积表达式:
A=(p1^k1)*(p2^k2)*(p3^k3)*....*(pn^kn) ,其中pi均为素数.
有A的所有因子之和为
S = (1+p1+p1^2+p1^3+...p1^k1) *(1+p2+p2^2+p2^3+….p2^k2) * (1+p3+ p3^3+…+p3^k3) * .... *(1+pn+pn^2+pn^3+...pn^kn)
(1).用递归二分求等比数列1+pi+pi^2+pi^3+...+pi^n:
1)若n为奇数,一共有偶数项,则:
1 + p + p^2 + p^3 +...+ p^n
= (1+p^(n/2+1))+ p * (1+p^(n/2+1)) +...+ p^(n/2) * (1+p^(n/2+1))
= (1+ p + p^2 +...+ p^(n/2)) * (1 +p^(n/2+1))
上式中颜色加粗的前半部分恰好就是原式的一半,那么只需要不断递归二分求和就可以了,后半部分为幂次式,将在下面第4点讲述计算方法。
2)若n为偶数,一共有奇数项,则:
1 + p + p^2 + p^3 +...+ p^n
= (1+p^(n/2+1))+ p * (1+p^(n/2+1)) +...+ p^(n/2-1) * (1+p^(n/2+1)) + p^(n/2)
= (1+ p + p^2 +...+ p^(n/2-1)) *(1+p^(n/2+1)) + p^(n/2);
上式颜色加粗的前半部分恰好就是原式的一半,依然递归求解
long longMOD=9901; long longp_sum(int x, long long y) { long long k1, k2, t; if(y==0) return 1; if(y==1) return (1+x); if(y%2==1) //求x的最高次数为y的等比数列。当y为奇数,可直接用公式。。。 { k1=(p_mi(x,y/2+1)+1)%MOD; k2=p_sum(x,y/2)%MOD; return (k1*k2)%MOD; } else //y为偶数时,则y-1为奇数,求出x^y+p_sum(x,y-1) { k1=p_mi(x,y)%MOD; k2=p_sum(x,y-1)%MOD; return (k1+k2)%MOD; } }
//递归求x^y次幂 long longp_mi(int x, long long y) { long long a; if(y==0) return 1; if(y==1) return x; if(y%2==1) //y为奇数时有,x^y=((x^((y-1)/2))^2)*x { a=p_mi(x,(y-1)/2)%MOD; return ((a*a%9901)*x)%MOD; } else { a=p_mi(x,y/2)%MOD; return (a*a)%MOD; } }
Code:
int euler_phi(int n) { int m=(int) sqrt(n+0.5); int ans=n; for(int i=2;i<=m;i++) if(n%i==0) { ans=ans/i*(i-1); while(n%i==0) n/=i; } if(n>1) ans=ans/n*(n-1); }
2) 求出1~n中所有数的欧拉函数值.
欧拉函数的递推关系:
假设素数p能整除n,那么
如果(n%p==0 && (n/p)%p==0) (p还能整除n / p) , PHI(n) = PHI(n / p) * p;
如果(n%p==0 && (n/p)%p==1) (p不能整除n / p), PHI(n) = PHI(n / p) * (p - 1);
可以用与筛选求素数非常类似的方法,在O(nloglogn)时间内计算完毕。
void phi_table(int n, int * phi) { for(int i=2;i<=n;i++)phi[i]=0; phi[1]=1; for(int i=2;i<=n;i++) if(!phi[i]) { for(intj=i;j<=n;j+=i) { if(!phi[j])phi[j]=j; phi[j]=phi[j]/i*(i-1); } } }
8.求逆元ax=1 (mod m) , x是a的逆元
=>ax+my=1=> gcd(a,m)=1
/*
int mod(int x,int n)
{
x%=n; //如x是负数,进行处理
if(x<0) x+=n;
return x; //return (x%n+n)%n;会越界
}
*/
用扩展欧几里得求:
int Inv(int a,int m) { int r,x,y; r=exgcd(a,m,x,y); if(r==1) return (x%m+m)%m; return -1; }
9.
(1)幂取模( a^n mod m)
法一:
int pow_mod(inta, int n, int m) { int x=pow_mod(a,n/2,m); long long ans=(long long)x*x%m; if(n%2==1) ans=ans*a%m; return (int)ans; }
法二:用快速幂取模(a^n mod m )
int ModPow(inta,int n,int m) { int rec=1; while(n) { if (n & 1) rec = (rec * a) % m; a = (a * a) % m; n >>= 1; } return rec % m; }
(2)快速幂( a^k)
int pow(int a,int k) { int rec = 1; while( k ) { if (k & 1) rec *= a; a *= a; k >>= 1; } return rec; }
10.快速模乘 a*b%p (位运算的思路)
int mul(int a, int b, int p) { intr=0; while(b) { if(b&1) r=(r+a)%p; a=(a<<1)%p; b>>=1; } return r; }
11.求解模线性方程组(中国剩余定理)
x=a1 mod m1
x=a2 mod m2
......
x=an mod mn 其中,a[],m[]已知,m[i]>0且m[i]与m[j]互质,求x.
同余模公式:
(a+b)%m=(a%m+b%m)%m
(a*b)%m=(a%m*b%m)%m
设m1,m2,...,mn是两两互素的正数,则对任意的整数a1,a2,...,an,同余方程组
其解为:X=((M_1*M1*a1)+(M_2*M2*a2)+...+(M_n*Mn*an)) mod m;
其中M=m1*m2*...*mn; Mi=M/mi; M_i是Mi的逆元(mod mi).
int china(int*a,int *m,int n) { int M=1,ans=0,Mi,i,x,y; for(i=0;i<n;i++) M*=m[i]; for(i=0;i<n;i++) { Mi=M/m[i]; exgcd(m[i],Mi,x,y); ans=(ans+Mi*y*a[i])%M; } return (ans%M+M)%M; }
12.枚举出一个集合的子集。设原集合为mask,则下面的代码就可以列出它的所有子集:
for (i = mask ; i ; i = (i - 1) & mask) ;