其他数论模板

miller-rabin素性判断+Pollard-rho算法分解质因数

下面这个代码不仅能判断素数还能给出质因数分解的结果（每一个质因子次数也能求出来），直接套用即可

#include <iostream> 
#include <cstdio> 
#include <algorithm>  
#include <cmath>  
#include <cstring>  
#include <map>  
#define LL long long
using namespace std;
 
const int times=20;
int number=0;
LL factor[100];
int tol;
map<LL,int> m;
                                                                
/*----------------------------------------------------------*/
/*----------------------------------------------------------*/
LL Random(LL n)			//生成[ 0 , n ]的随机数
{
	return ((double)rand() / RAND_MAX*n + 0.5);
}
 
LL q_mul(LL a, LL b, LL mod) //快速计算 (a*b) % mod
{
	LL ans=0;
	while (b) {
		if (b&1) {
			b--;
			ans=(ans+a)%mod;
		}
		b/=2;
		a=(a+a)%mod;
	}
	return ans;
}
 
LL q_pow(LL a,LL b,LL mod) //快速计算 (a^b) % mod
{
	LL ans=1;
	while (b) {
		if (b&1)  ans = q_mul(ans, a, mod);
		b/=2;
		a=q_mul(a, a, mod);
	}
	return ans;
}

LL gcd(LL a,LL b)
{
    if (a==0) return 1;
    if (a<0) return gcd(-a,b);
    while (b)
    {
        LL t=a%b;
        a=b;
        b=t;
    }
    return a;
}   //  允许负数 

/*----------------------------------------------------------*/
/*----------------------------------------------------------*/
bool witness(LL a, LL n)   //miller_rabin算法的精华 用检验算子a来检验n是不是素数
{
	LL tem=n-1;
	int j=0;
	while (tem%2 == 0) {
		tem/=2;
		j++;
	}
 
	LL x = q_pow(a,tem,n); 
	if (x==1 || x==n-1) return 1;
	while (j--) {
		x = q_mul(x, x, n);
		if (x==n-1) return 1;
	}
	return 0;
}
 
bool Miller_Rabin(LL n) 
{
	if (n==2) return 1;
	if (n<2 || n%2==0)  return 0;			
 
	for (int i=1; i<=times; i++) {
		LL a = Random(n-2)+1; 
		if (!witness(a,n))  return 0;
	}
	return 1;
}

/*----------------------------------------------------------*/
/*----------------------------------------------------------*/
LL Pollard_rho(LL x,LL c)
{
    LL i=1,k=2;
    LL x0=rand()%x;
    LL y=x0;
    while (1) {
        i++;
        x0 = (q_mul(x0,x0,x)+c)%x;
        LL d=gcd(y-x0,x);
        if (d!=1 && d!=x) return d;
        if (y==x0) return x;
        if (i==k) { y=x0; k+=k; }
    }
}

void findfac(LL n)
{
    if (Miller_Rabin(n)) {
        factor[++tol]=n;
        return;
    }
    LL p=n;
    while (p>=n) p = Pollard_rho(p,rand()%(n-1)+1);
    findfac(p);
    findfac(n/p);
}
 

int main()
{
	LL tar;
	while (cin >> tar)
	{
		if (Miller_Rabin(tar))	cout << "Yes, Prime!" << endl;
		else cout << "No, not prime.." << endl;
		
		tol=0;
		findfac(tar);
		for (int i=1; i<=tol; i++ ) cout<<factor[i]<<' '; cout<<endl;
	}
	return 0;
}

 

原根

下面这个代码能求出所有原根，直接套用即可

#include <vector>
#include <cstdio>
#include <algorithm>
using namespace std;
typedef long long LL;
const int N=1000001;
int m[N],phi[N],p[N],tot,prime[N];
//m[i]是i的最小素因数

void initPhi(){
    phi[1]=1;
    int k;
    for(int i=2;i<N;i++){
        if(!m[i]){
            p[tot++]=m[i]=i;
            prime[i]=1;
            phi[i]=i-1;
        }
        for(int j=0;j<tot&&(k=p[j]*i)<N;j++){
            m[k]=p[j];
            if(m[i]==p[j]){
                phi[k]=phi[i]*p[j];break;   
            }
            else phi[k]=phi[i]*(p[j]-1);
        }
    }
}

int root[N]={0};
void initRoot(){
    root[2]=root[4]=1;
    //从第二个质数3开始
    for(int i=1;i<tot;i++){
        for(LL j=p[i];j<N;j*=p[i]){
            root[j]=1;
        }
        for(LL j=2*p[i];j<N;j*=p[i]){
            root[j]=2;
        }
    }
}
 
vector<int> getfac(int n){
    vector<int>fac;
    LL tmp=n;
    for(int i=0;tmp>1&&i<tot;i++){
        if(tmp%p[i]==0){
            fac.push_back(p[i]);
            while(tmp%p[i]==0)tmp/=p[i];
        }
    }
    if(tmp>1)fac.push_back(tmp);
    return fac;
}

int gcd(int a,int b){return b==0?a:gcd(b,a%b);}
int ans[N],ia;
//已知n的一个原根x求n的所有phi(phi(n))个原根
void getRoot(int n,int x){
    ia=0;
    ans[ia++]=x;
    int y=x;
    for(int i=2;i<phi[n];i++){
        y=(y*x)%n;
        if(gcd(i,phi[n])==1)ans[ia++]=y;
    }
    sort(ans,ans+ia);
}

LL pow_mod(LL a,LL b,LL p){
    LL ret=1;
    while(b){
        if(b&1)ret=(ret*a)%p;
        a=(a*a)%p;
        b>>=1;
    }
    return ret;
}

//求n的一个原根
int oneRoot(int n){
    vector<int> fac=getfac(phi[n]);
    for(int i=1;i<n;i++){
        int flag=1;
        if(pow_mod(i,phi[n],n)!=1)flag=0;
        else
        for(int j=0;j<fac.size();j++){
            if(pow_mod(i,phi[n]/fac[j],n)==1){
                flag=0;
                break;
            }
        }
        if(flag)return i;
    }
    return 0;
}

void getRoot(int n){
    if(n==1||n==2||n==3)printf("%d\n",n-1);
    else if(root[n]){
        getRoot(n,oneRoot(n));
        for(int i=0;i<ia;i++){
            if(i)printf(" ");
            printf("%d",ans[i]);
        }
        printf("\n");
    }
    else printf("-1\n");
}

int main(){
    initPhi();
    initRoot();
    int n;
    while(~scanf("%d",&n)){
        getRoot(n);
    }
    return 0;
} 

离散对数

下面这个代码解决了A^x=B(mod C)的x，直接套用即可

#include<iostream>  
#include<cstdio>  
#include<cstring>  
#include<algorithm>  
#include<cmath>  
#define LL __int64  
#define N 1000000  
using namespace std;  
struct Node{  
    int idx;  
    int val;  
}baby[N];  

bool cmp(Node n1,Node n2){  
    return n1.val!=n2.val?n1.val<n2.val:n1.idx<n2.idx;  
}  

int gcd(int a,int b){  
    return b==0?a:gcd(b,a%b);  
}
  
int extend_gcd(int a,int b,int &x,int &y){  
    if(b==0){  
        x=1;  
        y=0;  
        return a;  
    }  
    int gcd=extend_gcd(b,a%b,x,y);  
    int t=x;  
    x=y;  
    y=t-a/b*y;  
    return gcd;  
}
  
int inval(int a,int b,int n){  
    int e,x,y;  
    extend_gcd(a,n,x,y);  
    e=((LL)x*b)%n;  
    return e<0?e+n:e;  
}
  
int PowMod(int a,int b,int MOD){  
    LL ret=1%MOD,t=a%MOD;  
    while(b){  
        if(b&1)  
            ret=((LL)ret*t)%MOD;  
        t=((LL)t*t)%MOD;  
        b>>=1;  
    }  
    return (int)ret;  
}
  
int BinSearch(int num,int m){  
    int low=0,high=m-1,mid;  
    while(low<=high){  
        mid=(low+high)>>1;  
        if(baby[mid].val==num)  
            return baby[mid].idx;  
        if(baby[mid].val<num)  
            low=mid+1;  
        else  
            high=mid-1;  
    }  
    return -1;  
}
  
int BabyStep(int A,int B,int C){  
    LL tmp,D=1%C;  
    int temp;  
    for(int i=0,tmp=1%C;i<100;i++,tmp=((LL)tmp*A)%C)  
        if(tmp==B)  
            return i;  
    int d=0;  
    while((temp=gcd(A,C))!=1){  
        if(B%temp) return -1;  
        d++;  
        C/=temp;  
        B/=temp;  
        D=((A/temp)*D)%C;  
    }  
    int m=(int)ceil(sqrt((double)C));  
    for(int i=0,tmp=1%C;i<=m;i++,tmp=((LL)tmp*A)%C){  
        baby[i].idx=i;  
        baby[i].val=tmp;  
    }  
    sort(baby,baby+m+1,cmp);  
    int cnt=1;  
    for(int i=1;i<=m;i++)  
        if(baby[i].val!=baby[cnt-1].val)  
            baby[cnt++]=baby[i];  
    int am=PowMod(A,m,C);  
    for(int i=0;i<=m;i++,D=((LL)(D*am))%C){  
        int tmp=inval(D,B,C);  
        if(tmp>=0){  
            int pos=BinSearch(tmp,cnt);  
            if(pos!=-1)  
                return i*m+pos+d;  
        }  
    }  
    return -1;  
}
  
int main(){  
    int A,B,C;  
    while(scanf("%d%d%d",&A,&C,&B)!=EOF){  
        if(B>=C){  
            puts("Orz,I can’t find D!");  
            continue;  
        }  
        int ans=BabyStep(A,B,C);  
        if(ans==-1)  
            puts("Orz,I can’t find D!");  
        else  
            printf("%d\n",ans);  
    }  
    return 0;  
}  

二次剩余

占坑，有空补