bzoj1951: [Sdoi2010]古代猪文

数论综合。费马小定理，lucas定理，中国剩余定理，exgcd，快速幂，乘法逆元。

首先要计算出n的每个约数，简单的sqrt(n)枚举即可。

然后计算C(i,m)(m个中挑i个的组合数，ps:因为网上正反俩种都有，所以标注一下。。)

设s=sum(C(i,m))

题目要求g^(s)%mod,

由费马小定理得 g^(s)=g^(s%(mod-1))。所以这里需要特判一下g是否等于mod。

可是组合数计算因为用到除法,所以模数必须为质数。999911659-1=2×3×4679×35617。

这四个数设为d[i]，对这四个数分别进行计算。

由因为组合数很大，直接计算会超时。

需要用到lucas定理。 lucas(n,m)=lucas(n/mod,m/mod)*C(n%mod,m%mod)。这里mod不是999911659,而是那4个数。

所以就得到4个答案m[i]。

组合数的计算可以根据费马小定理转换成快速幂计算。

然后再用中国剩余定理合并。

这时得到的答案可以写成如下形式：

s%d[0]=m[0],s%d[1]=m[1];

改写为 s%a1=b1,s%a2=b2 (1)

则有 a1*x+b1=a2*y+b2 （x,y分别为俩个未知数，和上面的任何数没有任何关系）

设a=a1,b=a2,c=b2-b1

则上式变为 ax+by=c (正负号并没有关系，因为不会用到y)

计算t=extended_gcd(a,b,x,y)

如果c%t ！=0 的话，说明此方程是无解的。在本题中t=1。（一下都基于本题t=1的情况）

然后我们计算出的是 ax+by=t。 所以我们将每一项乘以c。

x=(x*c%b+b)%b。

当然这样的x有无数组为 x1=x+k*b。

带入(1)式有  s=a1(x+k*b)+b1。

s=a1*b*k+a1*x+b1.

所以有 s%(a1*b)=a1*x+b1。

所以 b1=b1+a1*x,a1=a1*b。

在本题中这样操作完以后a1>b1，实际上即使t！=1，b1也会小于1。

最后合并完返回b1即可。

EN TARO Tassadar

#include<cstdio>
#include<algorithm>
#include<cstring>
#include<cmath>
#define LL long long 
using namespace std;
const int mod = 999911659;
const int d[]={2,3,4679,35617};
int fac[4][40000],m[4];
int n,g;

int power(LL a,int e,int mod) {
    LL res=1;
    while(e) {
        if(e&1) res=res*a%mod;
        a=a*a%mod;
        e>>=1;
    }
    return res;
}

int C(int n,int m,int x) {
    if(n<m) return 0;    
    return fac[x][n]*power(fac[x][m]*fac[x][n-m]%d[x],d[x]-2,d[x])%d[x];
}

int lucas(int n,int m,int x) {
    if(m==0) return 1;
    return lucas(n/d[x],m/d[x],x)*C(n%d[x],m%d[x],x)%d[x];    
}

int exgcd(int a,int b,int &x,int &y) {
    if(b==0) {
        x=1; y=0; return a;
    }
    else {
        int t=exgcd(b,a%b,y,x);
        y-=(a/b)*x;
        return t;
    }
}

int solve() {
    int a1,b1,a2,b2,a,b,c,x,y;    
    a1=d[0]; b1=m[0];
    for(int i=1;i<4;i++) {
        a2=d[i];b2=m[i];
        a=a1,b=a2,c=b2-b1;
        exgcd(a,b,x,y);
        x=(c*x%b+b)%b;//if(c%t) printf("Test\n");
        b1=b1+a1*x;
        a1=a1*b;
    }
    return b1;
}

int main() {
    for(int i=0;i<4;i++) {
        fac[i][0]=1;
        for(int j=1;j<=40000;j++) fac[i][j]=fac[i][j-1]*j%d[i];
    }
    scanf("%d%d",&n,&g);
    if(g==mod) printf("0\n");
    else {
        g=g%mod;    
        int n2=(int)sqrt(n+0.5);
        for(int i=1,j;i<=n2;i++) 
        if(n%i==0) {
            j=n/i;
            for(int x=0;x<4;x++) {
                m[x]=(m[x]+lucas(n,i,x))%d[x];
                if(i!=j) m[x]=(m[x]+lucas(n,j,x))%d[x];    
            }
        }
        printf("%d\n",power(g,solve(),mod));
    }
    return 0;    
}