数论
1、分數的四則運算：用結構體表示
2、高精度計算
3、數論
3.1最小公倍數 GCD
int gcd(int a,int b){
//輾轉相除
if(b==0) return a;
else return gcd(b,a%b);
}
int gcd(int a,int b){
return !b ? a:gcd(b,a%b);
}
最大公約數 LCM
lcm(a,b)=a*b/gcd(a,b)
3.2 快速冪
快速冪及擴展的矩陣快速冪，快速計算a^n的值
兩種方法
//分治
int fastpow(int a,int n){
if(n==1) return a;
int temp=fastpow(a,n/2);
if(n%2==1) return temp*temp*a;
else return temp*temp;
}
//位運算,加上取模
int fastpow(int a,int n,int mod){
int base=a;
int res=1;
while(n){
if(n&1) res=res*base%mod;
base=base*base%mod;
n>>=1;
}
return res;
}
矩陣快速冪，但是難點是把遞推關系轉化為矩陣的關系，例如：斐波那契數列
struct matrix{
int m[maxn][maxn];
matrix(){
memset(m,0,sizeof(m));
}
};
matrix multi(matrix a,matrix b){
matrix res;
for(int i=0;i<maxn;i++){
for(int j=0;j<maxn;j++){
for(int k=0;k<maxn;k++){
res.m[i][j]=(res.m[i][j]+a.m[i][k]*b.m[k][j])%mod;
}
}
}
return res;
}
matrix fastm(martix a,int n){
matrix res;
for(int i=0;i<maxn;i++)
res.m[i][i]=1; //初始化為單位矩陣 就像res=1
while(n){
if(n&1) res=multi(res,a);
a=multi(a,a);
n>>=1;
}
return res;
}
3.3 求素數
普通的做法(試除法），即從2~sqrt(n)之內判斷的----n<=10^12 以內可以接受 O(√N）
判斷：[2,sqrt(n)]內的所有數-----[2,sqrt(n)]內的所有素數
bool is_prime(int n){
if(n==1) return false;
for(int i=2;i*i<=n;i++){
if(n%i==0) return false;
}
return true;
}
埃氏篩法：O(NlongNlongN)
因為空間用到了is[maxn]，所以當maxn=1e7時是可以接受的，不然空間就太大了
優化（1）用來做篩的數到sqrt(n)就可以了（2）for(int j=i+i 可以改為j=i*i，因為前面已經篩過了
int prime[maxn];
bool is[maxn];
int ans=0;
void findd(){
for(int i=2;i*i<=maxn;i++){ //是<
if(is[i]==0){
prime[ans++]=i;
for(int j=i*i;j<maxn;j+=i){
p[j]=1;
}
}
}
}
線性篩模板
void prim(){
for(int i=2;i<=n;i++){
if(!vis[i]){
prime[++ans]=i;
}
for(int j=1;j<=ans&&prime[j]*i<=n;j++){
vis[prime[j]*i]=1;
if(i%prime[j]==0) break;
}
}
}
大區間素數
埃氏篩法能解決n<=1e7的問題，但是如果n更大，那么就可以擴展位大區間素數，例如求區間[a,b]內的素數，a<b<=1e12, b-a<=1e6
先用埃氏篩法求出[2,sqrt(b)]范圍內的素數，再用這些素數去篩空間[a,b]內的素數
題目：1619： 【例 1】Prime Distance
3.4 擴展歐幾里得算法與二元一次方程的整數解
給出整數a,b,n，問方程ax+by=n什么時候有整數解
！！！有解的充要條件是gcd(a,b)可以整除n
如果確定有解，那么解題方法是先找到一個可行解(x0,y0),然后用通解公式
x=x0+bt
y=y0-at求解，那么問題就是如何求(x0,y0)
擴展歐幾里得方法：當方程符合ax+by=gcd(a,b)時，可以用擴展歐幾里得算法求(x0,y0)
void extend_gcd(int a,int b,int &x,int &y){
if(b==0){
x=1;y=0;
return;
}
extend_gcd(b,a%b,x,y);
int temp=x;
x=y;
y=temp-(a/b)*y;
}
求出來的x,y就是所得到x0,y0,那么通解公式就是：
x=x0+b/gcd*t gcd=gcd(a,b)
y=y0-a/gcd *t
x對應的最小非負整數解為(x%(b/gcd))+b/gcd)%(b/gcd)
那么問題ax+by=n的整數解
利用上面的結果，步驟如下：
（1）判斷ax+by=n是否有整數解，也就是gcd(a,b)能不能整除n
（2）擴展歐幾里得求得x0,y0
（3）在ax0+by0=gcd(a,b)兩邊同時乘以n/gcd(a,b)
（4）對照ax+by=n，得到解為：
x0'=x0*n/gcd(a,b) y0'=y0*n/gcd(a,b)
那么通解公式為：
x=x0'+b/gcd*k
y=y0'-a/gcd *k
運用場合：（1）求解不定方程（2）求解模的逆元（3）求解同余方程
3.5 同余與逆元
a%m=b%m，那么a和b對m同余，m稱為同余的模，記為m|(a-b)，即a-b是m的整數倍，同余的符號記為a=b(mod m) (等號三橫）
逆元：給出a,m，求方程ax=1(mod m)，即ax除以m的余數為1
有解的條件是gcd(a,m)=1，即a，m互素，等價于求ax+my=1
方程ax=1(mod m)的一個解x，稱x為a模m的逆元。這樣的x很多，都稱為逆
int mod_inverse(int a,int m){
int x,y;
extend_gcd(a,m,x,y);
return (m+x%m)%m; //可能是負數，所以需要處理
}
逆元的應用：除法的模，求(a/b)mod m，a,b很大，但是求除法之后再取模會損失精度
把除法的模運算轉換成了乘法模運算
(a/b)modm=((a/b)modm)((bk)modm)=((a/b)*(bk))modm=akmodm
一元線性同余方程
ax=b(mod m)，即ax-b是m的整數倍，得到ax+my=b，當且僅當gcd(a,m)能整除b是有解
當gcd(a,m)=b時，能夠運用擴展歐幾里得方法直接求解ax+my=b
當gcd(a,m)!=b 時，就需要運用逆元
如果a'是a模m的逆，則a'a=1mod(m),那么在ax=b(mod m)兩邊同時乘以a',那么就是x=a'b(mod m)
步驟:求解ax+my=b 同余方程為ax=b(mod m)
（1）有解的條件：gcd(a,m)能夠整除b
（2）求ax=1(mod m)得逆元a',等價于擴展歐幾里得算法求出ax+my=1
（3）一個特解是x=a'b
（4）代入方程ax+my=b，求解y
中國剩余定理（孫子定理）
https://www.cnblogs.com/wkfvawl/p/9633188.html
但是這樣的中國剩余定理有一個要求就是m1,m2,m3...mn是互素的
#include<iostream>
#include<cstring>
#include<cmath>
#include<algorithm>
#include<stack>
#include<cstdio>
#include<queue>
#include<map>
#include<vector>
#include<set>
using namespace std;
const int maxn=1010;
const int INF=0x3fffffff;
typedef long long LL;
typedef unsigned long long ull;
LL n;
//中國剩余定理（孫子定理）
//https://www.cnblogs.com/wkfvawl/p/9633188.html
LL aa[12],bb[12];
void extend_gcd(LL a,LL b,LL &x,LL &y){
if(b==0){
x=1;
y=0;return;
}
extend_gcd(b,a%b,x,y);
LL tmp=x;
x=y;
y=tmp-(a/b)*y;
}
LL gcd(LL x,LL y){
if(y==0) return x;
else return gcd(y,x%y);
}
int main(){
scanf("%lld",&n);
LL mod=1;
for(int i=1;i<=n;i++){
scanf("%lld %lld",&aa[i],&bb[i]);
mod*=aa[i];
}
LL ans=0;
LL x,y;
for(int i=1;i<=n;i++){
LL p=mod/aa[i];
extend_gcd(p,aa[i],x,y);
//cout<<x<<endl;
ans=(ans+bb[i]*x*p)%mod;
}
cout<<(ans+mod)%mod;
/* 超時到懷疑人生哈哈哈哈
for(LL i=mm;;i++){
bool flag=0;
for(int j=1;j<=n;j++){
if(!judge(i-bb[j],aa[j])) {
flag=1;break;
}
}
if(!flag){
cout<<i<<endl;break;
}
}
*/
return 0;
}
View Code
擴展中國剩余定理
擴展中國剩余定理跟中國剩余定理沒半毛錢關系，一個是用擴展歐幾里得，一個是用構造
https://www.cnblogs.com/zwfymqz/p/8425731.html
判定有無解：（兩個一次處理）
#include<iostream>
#include<cstring>
#include<cmath>
#include<algorithm>
#include<stack>
#include<cstdio>
#include<queue>
#include<map>
#include<vector>
#include<set>
using namespace std;
const int maxn=1e5+10;
const int INF=0x3fffffff;
typedef long long LL;
typedef unsigned long long ull;
//變形是需要判斷有沒有解
//中國剩余定理的要求就是：m1,m2....mn兩兩互素
//但是如果不是兩兩互素的話，那么就需要擴展中國剩余定理
//https://www.cnblogs.com/zwfymqz/p/8425731.html
LL n;
//x=c1(mod m1)
LL c[maxn],m[maxn],x,y;
LL gcd(LL a,LL b){
if(b==0) return a;
return gcd(b,a%b);
}
LL extend_gcd(LL a,LL b,LL &x,LL &y){ //不僅算了一元二次方程的解，還算出來了公約數
if(b==0){
x=1;y=0;return a;
}
LL r=extend_gcd(b,a%b,x,y);
LL tmp=x;
x=y;
y=tmp-(a/b)*y;
return r;
}
LL inv(LL a,LL b){ //a在模b下的逆
LL r=extend_gcd(a,b,x,y); //最大公約數
while(x<0) x+=b; //取最小非負的逆
return x;
}
int main(){
while(~scanf("%lld",&n)){
for(int i=1;i<=n;i++){
scanf("%lld %lld",&m[i],&c[i]);
}
bool flag=0;
for(int i=2;i<=n;i++){ //兩兩合并位1個，并繼續擴展
LL m1=m[i-1],m2=m[i];
LL c1=c[i-1],c2=c[i];
LL t=gcd(m1,m2);
if((c2-c1)%t!=0){
flag=1;break; //必須余數為0 因為式子中要出現相除
}
m[i]=(m1*m2)/t;
c[i]=(inv(m1/t,m2/t)*(c2-c1)/t)%(m2/t)*m1+c1;
c[i]=(c[i]%m[i]+m[i])%m[i]; //還有這一步，別忘了
}
printf("%lld
",flag? -1:c[n]);
}
return 0;
}
View Code
3.6 質因子分解
#include<iostream>
#include<cstring>
#include<cmath>
#include<algorithm>
using namespace std;
const int maxn=100010;
bool is_prime(int x){
if(x==1) return 0;
for(int i=2;i*i<=x;i++){
if(x%i==0) return 0;
}
return 1;
}
//2147483646
int prime[maxn],ans=0;
struct node{
int x,num; //用結構體存
}fac[10]; //10個就夠了
void find_prime(){ //求素數表
for(int i=1;i<maxn;i++){
if(is_prime(i)) prime[ans++]=i;
}
}
int main(){
find_prime();
int n,num=0;
cin>>n;
if(n==1) {
cout<<"1"<<endl;
return 0; //特判
}
cout<<n<<"=";
int sqr=(int)(sqrt(1.0*n));
for(int i=0;prime[i]<=sqr&&i<ans;i++){
if(n%prime[i]==0){
fac[num].x=prime[i];
fac[num].num=0;
while(n%prime[i]==0){
fac[num].num++;
n/=prime[i];
}
num++;
}
if(n==1) break; //提前退出
}
if(n!=1){ //最后不等于1 ，說明有一個大于sqrt(n）的因子
fac[num].x=n;
fac[num++].num=1;
}
for(int i=0;i<num;i++){
if(i>0) cout<<"*";
cout<<fac[i].x;
if(fac[i].num!=1) cout<<"^"<<fac[i].num;
}
return 0;
}
ps.如果進行因子分解后，得到各因子的個數位e1,e2...ek,那么因子個數就是(e1+1)*(e2+1)*...*(ek+1)
因子之和就是（1+p1+p1^2+...+p1^e1)*(1+p2+p2^2+...p2^e2)*...*(1+pk+pk^2+...+pk^ek)
3.7 關于n！的問題
1、求n！中有多少個質因子p
可以枚舉，但是很慢
n！中有(n/p+n/p^2+.....)個質因子p.
eg.10！中有5個2^1,2個2^2,1個2^3，所以質因子2的個數位5+2+1=8個，O(logN)
int cal(int n,int p){
int ans=0;
while(n){
ans+=n/p;
n/=p;
}
return ans;
}
應用：計算n!的末尾有幾個0，末尾0的個數等于n！中因子10的個數，等于質因子5的個數，計算cal(n,5)就可以了
遞歸版本：
int cal(int n,int p){
if(n<p) return 0;
return n/p+cal(n/p,p);
}
3.8 組合數的計算
一、計算Cnm
方法一：通過定義計算，即使用Long long也只能承受n<=20;
方法二：遞推或者遞歸，單次計算不會超過O(N^2)，記憶化
res[n][m]=js(n-1,m)+js(n-1,m-1)
res[i][j]=res[i][i-j]
long long res[67][67]; //n=67,m=33時溢出
long long js(long long n,long long m){
if(m==0||m==n) return 1;
if(res[n][m]!=0) return res[n][m];
return res[n][m]=js(n-1,m)+js(n-1,m-1);
}
//計算整張表
void js(){
int n=60;
for(int i=0;i<=n;i++){
res[i][0]=res[i][i]=1;
}
for(int i=2;i<=n;i++){
for(int j=1;j<=i/2;j++){
res[i][j]=res[i-1][j]+res[i-1][j-1];
res[i][i-j]=res[i][j];
}
}
}
方法三：定義式變形計算
ans=ans*(n-m+i)/i i=1~m、
//根據定義式，劃分為計算m次
//n=62,m=31時溢出 long long js(long long n,long long m){ long long ans=0; for(long long i=1;i<=m;i++){ ans=ans*(n-m+i)/i; } return ans; }
二、計算Cnmmodp
第一種：遞推或者遞歸
int res[1010][1010]; //m<=n<=1000 ，對p的素性和大小沒有限制
int js(int n,int m,int p){
if(m==0||m==n) return 1;
if(res[n][m]!=0) return res[n][m];
return res[n][m]=(js(n-1,m)+js(n-1,m-1))%p;
}
第二種：定義式計算
將組合數進行質因子分解，那么Cn,m=p1^c1*p2^c2...*pn^cn%p，通過快速冪計算每一組pi^ci%p，然后相乘取模就可以了
對Cn,m進行質因子分解：遍歷不超過n的所有質數，計算n!,m!,(n-m)!中分別含質因子pi的個數x,y,z，那么Cn,m中含質因子pi的個數為x-y-z
這里就用到了前面的求n!里面的質因子p的個數cal(int n,int p)
m<=n<=10^6，對p的素性沒有限制
int prime[maxn]; //素數表（已經得到的，埃氏篩法）
int js(int n,int m,int p){
int ans=1;
for(int i=0;prime[i]<=n;i++){
int c=cal(n,prime[i])-cal(m,prime[i])-cal(n-m,prime[i]);
ans=ans*binpow(prime[i],c,p)%p; //快速冪計算prime[i]^c%p
}
}
高次同余方程 BSGS
#include<iostream>
#include<cstring>
#include<cmath>
#include<algorithm>
#include<stack>
#include<cstdio>
#include<queue>
#include<map>
#include<vector>
#include<set>
using namespace std;
const int maxn=1010;
const int INF=0x3fffffff;
typedef long long LL;
typedef unsigned long long ull;
//1.y^zmod p
//2. x*y=z(mod p)
//3. y^x=z(mod p)
LL js1(LL y,LL z,LL p){
LL tmp=1;
while(z){
if(z&1) tmp=tmp*y%p;
y=y*y%p;
z>>=1;
}
return tmp%p;
}
void extend_gcd(LL a,LL b,LL &x,LL &y){
if(b==0){
x=1;
y=0;return;
}
extend_gcd(b,a%b,x,y);
LL tmp=x;
x=y;
y=tmp-(a/b)*y;
}
LL gcd(LL a,LL b){
if(b==0 ) return a;
else return gcd(b,a%b);
}
//x*y=z(mod p) ax=b(mod m)
//ax+my=b
//y*x+p*k=z 已知y,p,z
//求出x0之后
//最小非負解：((x=x0*z/gcd)+p/gcd)%(p/gcd)
void js2(LL y,LL z,LL p){
LL A=y,B=z,M=p;
LL xx,yy;
extend_gcd(A,M,xx,yy);
//gcd(a,m)%b==0
LL d=gcd(y,p);
if(z%d){
printf("Orz, I cannot find x!
");return;
}
xx=((xx*z/d)%(p/d)+p/d)%(p/d);
printf("%lld
",xx);
//return (xx%M+M)%M;
}
//y^x=z(mod p)
//BSGS可以解決高次同余方程：
//https://www.cnblogs.com/kamimxr/p/11555986.html
void js3(LL a,LL ans,LL p){
map<LL,LL> myhash;
ans%=p;
int tmp=sqrt(p)+1;
for(int i=0;i<tmp;i++){
//js1其實就是快速冪
myhash[(ans)*js1(a,i,p)%p]=i;
}
a=js1(a,tmp,p)%p;
if(a==0&&ans==0){
printf("1
");return;
}
if(a==0&&ans!=0){
printf("Orz, I cannot find x!
");
return;
}
for(int i=0;i<=tmp;i++){
if(myhash.find(js1(a,i,p))!=myhash.end()&&(i*tmp-myhash[js1(a,i,p)]>=0)){
printf("%d
",i*tmp-myhash[js1(a,i,p)]);
return;
}
}
printf("Orz, I cannot find x!
");
}
int main(){
int t,k;
LL y,z,p;
while(~scanf("%d %d",&t,&k)){
if(k==1){
for(int i=1;i<=t;i++){
scanf("%lld %lld %lld",&y,&z,&p);
printf("%lld
",js1(y,z,p));
}
}
else if(k==2){
for(int i=1;i<=t;i++){
scanf("%lld %lld %lld",&y,&z,&p);
js2(y,z,p);
}
}
else if(k==3){
////y^x=z(mod p)
for(int i=1;i<=t;i++){
scanf("%lld %lld %lld",&y,&z,&p);
js3(y,z,p);
}
}
}
return 0;
}
View Code
總結
- 上一篇： SQL进阶随笔--case用法（一）
- 下一篇： 尼斯湖水怪流传 1000 多年了，现在依