【模板】数论

公约数
gcd

LL gcd(LL a,LL b){
    return b==0? a:gcd(b,a%b);
}

ex_gcd

//返回值为最大公约数
LL ex_gcd(LL a,LL b,LL &x,LL &y){
    LL d = a;
    if(!b){x = 1,y = 0;}
    else{
      d = ex_gcd(b,a%b,y,x);
      y-=a/b*x;
    }
    return d;
}

素数

素因子分解

线性筛

void xxs(){
    int tot = 0;
    memset(check, 0, sizeof(check));
    for (int i = 2; i < MAXL; ++i){
        if (!check[i]){
            prime[tot++] = i;
        }
        for (int j = 0; j < tot; ++j){
            if (i * prime[j] > MAXL) break;
    	    check[i*prime[j]] = 1;
     	    if (i % prime[j] == 0) break;
    }
}

朴素算法\(O(\sqrt n)\)

void prime_factor(int n,map<int,int> &pf){
    for(int i =2 ; i*i<=n ; ++i){//n为素数时!
        while(n%i==0){
            ++pf[i];
            n/=i;
        }
    }
    if(n!=1)pf[n] = 1;
}

Eratosthenes筛法

void Eratosthenes(int n){
    memset(is_prime,true,sizeof(is_prime));

    for(int i = 2 ; i*i<=n; ++i)
        if(is_prime[i])
        for(int j=i*i ; j<=n ; j+=i)is_prime[j] = false;
}

区间筛法

void segment_sieve(LL a,LL b){//[a,b]
    memset(is_prime_ab,true,sizeof(is_prime_ab[0])*(b-a+1));
    memset(is_prime_sqrtb,true,sizeof(is_prime_sqrtb[0])*(sqrt(b)+2));
    for(LL i = 2 ; i*i<=b ; ++i)
    if(is_prime_sqrtb[i]){
        for( LL j = i*i ; j*j<=b ; j+=i)is_prime_sqrtb[j] = false;
        for(LL j = max(i*i,(a-1)/i+1)*i ; j<=b ; j+=i)is_prime_ab[j-a] = false;
    }
}

大素数分解与大素数测试

miller_rabin

已知最快的素数分解算法.\(O(lgV)\)

bool witness(LL a,LL n,LL u,LL t){
	LL x0 = power_mod(a,u,n),x1;
	for(int i=1 ;i<=t ; ++i){
		x1 = mulmod(x0,x0,n);
		if(x1==1 && x0!=1 && x0!=n-1)
			return false;
		x0 = x1;
	}
	if(x1 !=1)return false;
	return true;
}

bool miller_rabin(LL n, int times = 20){
	if(n<2)return false;
	if(n==2)return true;
	if(!(n&1))return false;
	LL u = n-1,t =0;
	while (u%2==0) {
		t++;u>>=1;
	}
	while (times--) {
		LL a = random(1,n-1);
		//if(a == 0)std::cout << a << " "<<n<< " "<<u<<" " << t<<'\n';
		if(!witness(a,n,u,t))return false;
	}
	return true;
}

pollard_rho 分解一个合数\(V\)的运行时间\(O(V^{1/4 })\)

/*
*pollard_rho分解n,
*c : 随机迭代器，每次运行设置为随机种子往往更快.
*/
LL pollard_rho(LL n,LL c = 1){
	LL x  = random(1,n);
	LL i =1,k =2,y = x;
	while (1) {
		i++;
		x = (mulmod(x,x,n)+c)%n;
		LL d = gcd(y-x>=0?y-x:x-y,n);
		if(d!=1 && d!=n)return d;//非平凡因子.
		if(y==x)return n;//重复.
		if(i==k){ y = x ; k<<=1;}//将x_1,2,4,8,16,..赋值为y.
	}
}

找出因子分解 $O(V^{1/4}lgV)$

void find_factor(LL n,std::map<LL, int>  & m){
	if(n<=1)return ;
	if(miller_rabin(n)){
		++m[n];
		return ;
	}
	LL p = n;

	while (p==n)p = pollard_rho(p,random(1,n));
	find_factor(p,m);
	find_factor(n/p,m);
}

euler phi函数

\[\phi(n) = n\prod_{i = 1}^{k}(1-\frac{1}{p_i}) \]

证明详见《初等数论及其应用》

int  euler_phi(int  n){
    int  ans = n;
    for(int i=2 ; i*i<=n ; ++i)
        if(ans%i ==0)
        {
            ans = ans/i*(i-1);
            while(n%i==0)n/=i;
        }
        if(n>1)ans = ans/n*(n-1);
    return ans;
}

phi_table,类似于线性筛的做法

void phi_table(int n){
    memset(phi,0,sizeof(phi[0])*(n+5));
    phi[1] = 1;
    for(int i = 2 ; i<=n ; ++i){
        if(!phi[i]){//素数标记
            for(int j = i ; j<=n ; j+=i){
                if(!phi[j])phi[j] = j;
                phi[j] = phi[j]/i*(i-1);
            }
        }
    }
}

模运算
power_mod

LL power_mod(LL x,LL n,LL mod){
    LL res = 1;
    while(n){
        if(n&1)res = (res*x)%mod;
        x = (x*x)%mod;
        n>>=1;
    }
    return res;
}

大数乘法取模

LL mulmod(LL a,LL b,LL mod){
	LL res =0,y = a%mod;
	while (b) {
		if(b&1)res = (res+y)%mod;
		b>>=1;
		y = (y<<1)%mod;
	}
	return res;
}

大整数取模

LL big_mod(string val,LL mod){
    LL res = 0;
    for(int i=0 ; i<val.length() ; ++i){
        res = ((res)*10+val[i]-'0')%mod;
    }
    return res;
}

模方程

LL MLE(LL a,LL b,LL n){
    LL d,x,y;
    d = ex_gcd(a,n,x,y);
    if(b%d !=0){
        return -1;
    }else{
        LL x0 = x*b/d%n+n;
        return x0%(n/d);//模（n/d)
    }
}

乘法逆元 a在模n意义下的逆

LL inv(LL a,LL n){
    LL x,y;
    LL d = ex_gcd(a,n,x,y);
    return d==1? (x+n)%n:-1;//非负性保证.
}

中国剩余定理

//x % m[i] = a[i]
LL china(int n,int *a,int *m){
    LL M = 1,x = 0,y,z;
    for(int i=0 ; i<n ; ++i)M*=m[i];
    for(int i=0 ; i<n ; ++i){
        LL M_i = M/m[i];
        ex_gcd(M_i,m[i],y,z);//M_i*y = 1(mod m[i])
        x = (x+M_i*a[i]*y)%M;
    }
    return (x+M)%M;
}

朴素模方程(\(m_i不两两互素的时候\))

LL MLE(int *r,int *mod,int n){
	LL lm = 0, lb = 1;
	for (int i = 0; i < n; i++)
	{
		LL k1,k2;
		LL d= exgcd(lb, mod[i],k1,k2);	// x=c1(mod r1)
		if ((lm - r[i]) % d) { return -1; }	// 联立x=r2(mod m2)，(r1-r2)=0(mod gcd)才有解
		lb = lb / d * mod[i];							// lcm
		LL z = k2 * ((lm - r[i]) / d);						// 求出k2
		lm = z * mod[i] + r[i];											// 得到方程组的一个最小解
		lm = ((lm % lb) + lb) % lb;										// 保证最小解大于0
	}
	return lm;
}

莫比乌斯反演 线性筛法\(O(n)\)

int prime[maxn],cnt;
int mu[maxn];

void init_mobius(){
    memset(prime,0,sizeof(prime));
    mu[1] =1;
    cnt =0;
    for(int i=2 ;i<maxn ; ++i){
        if(!prime[i]){
            prime[cnt++] = i;
            mu[i] = -1;
        }
        for(int j=0 ; j<cnt && i*prime[j]<maxn ; ++j){
            prime[i*prime[j]] = 1;
            if(i%prime[j])mu[i*prime[j]] = -mu[i];
            else {
                mu[i*prime[j]] = 0;
                break;
            }
        }
    }
}