数论基础模板

一部分是直接copy别人的，主要是为了之后方便用。

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;

///素数筛选
const int MAXN=1000010;
int prime[MAXN+5];
void getPrime()
{
    memset(prime,0,sizeof(prime));
    for(int i=2;i<=MAXN;i++)
    {
        if(!prime[i])prime[++prime[0]]=i;
        for(int j=1;j<=prime[0]&&prime[j]<=MAXN/i;j++)
        {
            prime[prime[j]*i]=1;
            if(i%prime[j]==0) break;
        }
    }
}

///大区间素数筛选（POJ 2689）//来自kuangbin
/*
* POJ 2689 Prime Distance
* 给出一个区间 [L,U]，找出区间内容、相邻的距离最近的两个素数和距离最远的两个素数。
* 1<=L<U<=2,147,483,647 区间长度不超过 1,000,000
* 就是要筛选出 [L,U] 之间的素数
*/
#include<bits/stdc++.h>
const int MAXN=100010;
int prime[MAXN+1];
void getPrime(){
    memset(prime,0,sizeof(prime));
    for(int i=2;i<=MAXN;i++){
        if(!prime[i])prime[++prime[0]]=i;
        for(int j=1;j<=prime[0]&&prime[j]<=MAXN/i;j++){
            prime[prime[j]*i]=1;
            if(i%prime[j]==0)break;
        }
    }
}
bool notprime[1000010];
int prime2[1000010];
void getPrime2(int L,int R){
    memset(notprime,false,sizeof(notprime));
    if(L<2)L=2;
    for(int i=1;i<=prime[0]&&(long long)prime[i]*prime[i]<=R;i++){
        int s=L/prime[i]+(L%prime[i]>0);
        if(s==1)s=2;
        for(int j=s;(long long)j*prime[i]<=R;j++)
        if((long long)j*prime[i]>=L)
            notprime[j*prime[i]-L]=true;
    }
    prime2[0]=0;
    for(int i=0;i<=R-L;i++)
    if(!notprime[i])
    prime2[++prime2[0]]=i+L;
}
int main(){
    getPrime();
    int L,U;
    while(scanf("%d%d",&L,&U)==2)
    {
        getPrime2(L,U);///区间内的素数都放在prime2中
        if(prime2[0]<2)
        printf("There are no adjacent primes.\n");
        else
        {
            int x1=0,x2=100000000,y1=0,y2=0;
            for(int i=1;i<prime2[0];i++)
            {
                if(prime2[i+1]-prime2[i]<x2-x1)
                {
                    x1=prime2[i];
                    x2=prime2[i+1];
                }
                if(prime2[i+1]-prime2[i]>y2-y1)
                {
                    y1=prime2[i];
                    y2=prime2[i+1];
                }
            }
            printf("%d,%dare closest,%d,%dare most distant.\n",x1,x2,y1,y2);
        }
    }
}

///合数分解
ll factor[100][2];
int fatCnt;
int getFactors(ll x)
{
    fatCnt=0;
    ll tmp=x;
    for(int i=1;prime[i]<=tmp/prime[i];i++)
    {
        factor[fatCnt][1]=0;
        if(tmp%prime[i]==0)
        {
            factor[fatCnt][0]=prime[i];
            while(tmp%prime[i]==0)
            {
                factor[fatCnt][1]++;
                tmp/=prime[i];
            }
            fatCnt++;
        }
    }
    if(tmp!=1)
    {
        factor[fatCnt][0]=tmp;
        factor[fatCnt++][1]=1;
    }
    return fatCnt;
}

///Miller_Rabin判断大数是否为素数
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
#define MT 5
int prime[] = {2, 3, 7, 61, 24251};
ll mulmod(ll a, ll b, ll k)///a*b相乘爆long long的话就要用快速乘了
{
    return (a*b)%k;
}
ll powermod(ll a, ll b, ll k)
{
    ll re = 1, temp = a;
    while(b)
    {
        if (b & 1) re = mulmod(re, temp, k);
        temp = mulmod(temp, temp, k);
        b >>= 1;
    }
    return re;
}
int TwiceDetect(ll a, ll b, ll k) {
    int t = 0;
    ll x, y;
    while ((b & 1) == 0)
    {
        b >>= 1;
        t++;
    }
    y = x = powermod(a, b, k);
    while (t--)
    {
        y = mulmod(x, x, k);
        if (y == 1 && x != 1 && x != k - 1)
            return 0;
        x = y;
    }
    return y;
}
bool Miller_Rabin(ll n)
{
    int i;
    ll tmp;
    for (i = 0; i < MT; i++) {
        tmp = prime[i];
        if (n == prime[i]) return true;
        if (TwiceDetect(tmp, n - 1, n) != 1)
            break;
    }
    return (i == MT);
}
int main()
{
    ll n;
    while (scanf("%lld", &n) == 1)
    {
        if ((n > 1) && Miller_Rabin(n)) {
            printf("YES\n");
        }else {
            printf("NO\n");
        }
    }
    return 0;
}

/// pollard_rho 算法进行质因素分解
///kuangbin
//*********************************************
ll factor[100];//质因素分解结果（刚返回时时无序的）
int tol;//质因素的个数，编号 0--tol-1

ll gcd(ll a,ll b){
    ll t;
    while(b){
        t = a;
        a = b;
        b = t%b;
    }
    if(a >= 0)return a;
    else return -a;
}
//找出一个因子
ll pollard_rho(ll x,ll c){
    ll i = 1, k = 2;
    srand(time(NULL));
    ll x0 = rand()%(x-1) + 1;
    ll y = x0;
    while(1){
        i ++;
        x0 = (mult_mod(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;}
    }
}
//对 n 进行素因子分解，存入 factor. k 设置为 107 左右即可
void findfac(ll n,int k){
    if(n == 1)return;
    if(Miller_Rabin(n))
    {
        factor[tol++] = n;
        return;
    }
    ll p = n;
    int c = k;
    while( p >= n)p = pollard_rho(p,c--);//值变化，防止死循环 k
    findfac(p,k);
    findfac(n/p,k);
}
//POJ 1811
//给出一个N(2 <= N < 2^54),如果是素数, 输出"Prime", 否则输出最小的素因子
int main(){
    int T;
    ll n;
    scanf("%d",&T);
    while(T--){
        scanf("%I64d",&n);
        if(Miller_Rabin(n))printf("Prime\n");
        else{
            tol = 0;
            findfac(n,107);
            ll ans = factor[0];
            for(int i = 1;i < tol;i++)
            ans = min(ans,factor[i]);
            printf("%I64d\n",ans);
        }
    }
    return 0;
}

///筛法求欧拉函数
const int N=3e6;
int euler[N+10];
void getEuler()
{
    memset(euler,0,sizeof euler);
    euler[1]=1;
    for(int i=2;i<=N;i++)
    if(!euler[i])
    {
        for(int j=i;j<=N;j+=i)
        {
            if(!euler[j])
                euler[j]=j;
            euler[j]=euler[j]/i*(i-1);
        }
    }
}
///欧拉线性筛
const int N = 2e7+5;
bool flag[N];
int phi[N];
int p[N];
int cnt = 0;
void Get_phi()
{
    cnt = 0;
    memset(flag, true, sizeof(flag));
    phi[1] = 1;
    for(int i=2; i< N; i++)
    {
        if(flag[i])
        {
            p[cnt++] = i;
            phi[i] = i-1;
        }
        for(int j=0; j<cnt; j++)
        {
            if(i*p[j] > N)
                break;
            flag[i*p[j]] = false;
            if(i%p[j] == 0)
            {
                phi[i*p[j]] = p[j] * phi[i];
                break;
            }
            else
                phi[i*p[j]] = (p[j]-1) * phi[i];
        }
    }
}

///求单个欧拉函数
ll eular(ll n)
{
    ll ans=n;
    for(ll i=2;i*i<=n;i++)
    if(n%i==0)
    {
        ans-=ans/i;
        while(n%i==0)
            n/=i;
    }
    if(n>1) ans-=ans/n;
    return ans;
}

int eular(int x)//效率更高点的求单个欧拉函数，需要先筛出质数
{
    int r=x;
    for(int i=1;prime[i]*prime[i]<=x;i++)
    {
        if(x%prime[i]==0)
        {
            r=r-r/prime[i];
            while(x%prime[i]==0)
                x/=prime[i];
        }
    }
    if(x>1)r=r-r/x;
    return r;
}

//小于等于n,与n互质的数的和是  euler(n)*n/2;
//欧拉函数是积性函数——若m,n互质，euler(m*n)=euler(n)*euler(m) 
//当n为奇数时,euler(2*n)=euler(n) 
//若n为质数则,euler(n)=n-1
//若a,m互质，a^phi(m) %m = 1 (欧拉定理)。当m为质数p时，a^(p-1) % p = 1 (费马小定理)。通常利用结论降幂
//p是素数，且n = p^k，那么φ(n) = (p-1)*p^(k-1)

#define MOD 9973
///快速乘
//计算 (a*b)%c， a,b都是long long的数，直接相乘可能溢出的，a,b,c <2^63
ll mult_mod(ll a,ll b,ll c)
{
    a%=c;
    b%=c;
    ll ret=0;
    while(b)
    {
        if(b&1){ret+=a;ret%=c;}
        a<<=1;
        if(a>=c)a%=c;
        b>>=1;
    }
    return ret;
}
///快速幂
ll qmod(ll x,ll p)
{
    ll ans=1;
    while(p)
    {
        if(p&1) ans=ans*x%MOD;
        x=x*x%MOD;
        p/=2;
    }
    return ans;
}

gcd(x^a−1,x^b−1)=x^gcd(a,b)−1
gcd(fib[a],fib[b])=fib[gcd(a,b)]  ///斐波那契数列

///欧几里得x=__gcd(a,b);
///扩展欧几里得
//返回d=gcd(a,b);和对应于等式ax+by=d中的x,y
ll extend_gcd(ll a,ll b,ll &x,ll &y)
{
    if(a==0&&b==0) return -1;//无最大公约数
    if(b==0){x=1;y=0;return a;}
    ll d=extend_gcd(b,a%b,y,x);
    y-=a/b*x;
    return d;
}

//扩展欧几里德求ax+by=n 非负解的个数
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;

ll gcd(ll a, ll b) {
    return b ? gcd(b, a%b) : a;
}

ll lcm(ll a, ll b) {
    return a / gcd(a,b) * b;
}

ll extend_gcd(ll a,ll b,ll&x,ll&y) {
    if(!b) {
        x = 1;
        y = 0;
        return a;
    }
    ll xt = 0, yt = 0;
    ll d = extend_gcd(b, a % b, xt, yt);
    x = yt;
    y = xt - yt * (a/b);
    return d;
}

ll cal(ll a,ll b,ll n) {
    ll x = 0,y = 0,d;
    d = extend_gcd(a,b,x,y);
    if(n % d != 0) {
        return 0;
    }
    x *= n / d, y *= n / d;
    ll LCM = lcm(a,b);
    ll t1 = LCM / a, t2 = LCM / b;
    if(x<1) {
        ll num = (1-x) / t1;
        x += num * t1;
        y -= num * t2;
        if(x<1) {
            y -= t2;
            x += t1;
        }
    }
    if(y<1) {
        ll num = (1-y) / t2;
        y += num * t2;
        x -= num * t1;
        if(y<1) {
            y += t2;
            x -= t1;
        }
    }
    ll ans = x > 0 && y > 0;
    if(ans) {
        ans += min((x-1) / t1, ((n-1) / b - y) / t2);
        ans += min((y-1) / t2, ((n-1) / a - x) / t1);
    }
    return ans;
}
int main()
{
    int t;
    scanf("%d",&t);
    while(t--)
    {
        ll n,a,b;
        scanf("%lld%lld%lld",&n,&a,&b);
        printf("%lld\n",cal(a,b,n));
    }
    return 0;
}

//求ax+by=c的最小正整数解(x,y)
int main()
{
    LL a, b, c, gcd;

    scanf("%lld%lld%lld", &a, &b, &c);
    gcd = Gcd(a, b);

    if(c % gcd != 0)// 判断是否有解
        printf("Impossible\n");
    else
    {
        ex_gcd(a, b);
        LL x1, y1, b1;
        b1 = b/gcd;
        x1 = (x + b1) * (c/gcd);
        x1 = (x1 % b1 + b1) % b1;// 求解出 x 的 最小正整数解
        y1 = (c - a*x1) / b;

        printf("x = %lld , y = %lld\n", x1, y1);
    }

    return 0;
}

///求逆元
//ax = 1(mod n) 求x
ll mod_reverse(ll a,ll n)
{
    ll x,y;
    ll d=extend_gcd(a,n,x,y);
    if(d==1) return (x%n+n)%n;
    else return -1;
}
int main()
{
    ll B;
    cin>>B;
    cout<<mod_reverse(B,MOD)<<endl;///扩展欧几里得求B的逆元
    cout<<qmod(B,MOD-2)<<endl;///快速幂求B的逆元
    return 0;
}

///中国剩余定理
//一般的要求m[i]两两互质
//x%a[i]=m[i] x解
int CRT(int a[],int m[],int n)
{
    int M=1,ans=0;
    for(int i=1;i<=n;i++)
        M*=m[i];
    for(int i=1;i<=n;i++)
    {
        int x,y;
        int Mi=M/m[i];
        extend_gcd(Mi,m[i],x,y);
        ans=(ans+Mi*x*a[i])%M;
    }
    if(ans<0) ans+=M;
    return ans;
}

//m[i]不满足两两互质的时候
ll Mod;
ll gcd(ll a, ll b)
{
    if(b==0)
        return a;
    return gcd(b,a%b);
}
//a在模n乘法下的逆元，没有则返回-1
ll inv(ll a, ll n)
{
    ll x,y;
    ll t = extend_gcd(a,n,x,y);
    if(t != 1)
        return -1;
    return (x%n+n)%n;
}
//将两个方程合并为一个
bool merge(ll a1, ll n1, ll a2, ll n2, ll& a3, ll& n3)
{
    ll d = gcd(n1,n2);
    ll c = a2-a1;
    if(c%d)
        return false;
    c = (c%n2+n2)%n2;
    c /= d;
    n1 /= d;
    n2 /= d;
    c *= inv(n1,n2);
    c %= n2;
    c *= n1*d;
    c += a1;
    n3 = n1*n2*d;
    a3 = (c%n3+n3)%n3;
    return true;
}
//求模线性方程组x=ai(mod ni),ni可以不互质
ll CRT(int len, ll* a, ll* n)
{
    ll a1=a[0],n1=n[0];
    ll a2,n2;
    for(int i = 1; i < len; i++)
    {
        ll aa,nn;
        a2 = a[i],n2=n[i];
        if(!merge(a1,n1,a2,n2,aa,nn))
            return -1;
        a1 = aa;
        n1 = nn;
    }
    Mod = n1;
    return (a1%n1+n1)%n1;
}

///组合数
const int N=1e5;
const ll mod=100000007;
ll qmod(ll x,ll p)
{
    ll ans=1;
    while(p)
    {
        if(p&1) ans=ans*x%mod;
        x=x*x%mod;
        p>>=1;
    }
    return ans;
}

ll fac[N+10],inv[N+10];
void init()
{
    fac[0]=1;
    for(int i=1;i<=N;i++)
    {
        fac[i]=fac[i-1]*i%mod;
        if(fac[i]==0)
        fac[i]=fac[i]+mod;
    }
    inv[N]=qmod(fac[N],mod-2);//N!的逆元
    for(int i=N-1;i;i--)
        inv[i]=(i+1)*inv[i+1]%mod;
}
ll C(ll n,ll m)
{
    if(m>n) return 0;
    return fac[n]*inv[n-m]%mod*inv[m]%mod;
}

inv[0]=inv[1]=1;
for(int i = 2; i < N; ++i)
    inv[i] = (mod - mod/i) * inv[mod%i] % mod;//线性求i的逆元

//大数组合数，卢卡斯定理
typedef long long ll;
ll  n, m, p;
ll extend_gcd(ll a,ll b,ll &x,ll &y)
{
   if(b==0) { x=1, y=0; return a; }
   ll ret= extend_gcd(b,a%b,y,x);
   y-= a/b*x;
   return ret;
}
ll Inv(ll a,int m)
{
   ll d,x,y,t= (ll)m;
   d= extend_gcd(a,t,x,y);
   if(d==1) return (x%t+t)%t;
   return -1;
}
ll C(ll n, ll m, ll p)
{
    ll a=1, b=1;
    if(m>n) return 0;
    while(m)
    {
        a=(a*n)%p;
        b=(b*m)%p;
        m--;
        n--;
    }
    return (ll)a*Inv(b,p)%p;
}
int Lucas(ll n, ll m, ll p)
{
    if(m==0) return 1;
    return (ll)C(n%p,m%p,p)*(ll)Lucas(n/p,m/p,p)%p;
}

int main()
{
    int  T;
    cin >> T;
    while(T--)
    {
        scanf("%lld%lld%lld",&n,&m,&p);
        printf("%d\n",Lucas(n,m,p));
    }
    return 0;
}