您的位置:首页 > 其它

BZOJ 1951: [Sdoi2010]古代猪文(Lucas定理 &&中国剩余定理&&费马小定理)

2013-01-02 18:23 447 查看

1951: [Sdoi2010]古代猪文

求G^Sigma{C(N, i),i | N} mod M的值,其中M = 999911659。

M是一个素数,故根据费马小定理G^(M - 1) ≡ 1 (mod M)。

则 G^Sigma{C(N, i),i | N} ≡ G^(Sigma{C(N, i),i | N} mod (M-1)) ≡ G^Sigma{C(N, i) mod (M-1),i | N} (mod M)。
于是我们只需求 C(N, K) mod (M - 1) 然后累加即可。

注意到 :M - 1 = 2 * 3 * 4679 * 35617。 

于是我们可以得到
P ≡ a1 (mod 2)
P ≡ a2 (mod 3)
P ≡ a3 (mod 4679)
P ≡ a4 (mod 35617)

求得a1,a2,a3,a4后通过中国剩余定理合并即可。

欲求a1,a2,a3,a4我们只需求 C(N, K) mod P (P为素数)(n,k<= 10^10) 

n,k很大,我们可以通过lucas定理解决这一问题

calc(n,m,p) = C(n%p,m%p) * calc(n/p,m/p,p);

之后C(n,m)直接通过预处理阶乘 + 乘法逆元解决;

另外 :当G == M 时 飞马小定理不成立 要输出 0; 

#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <iostream>
#include <assert.h>
using namespace std;
typedef long long lint;
lint P,n,g;
lint a[10] = {2,3,4679,35617};
lint fac[40000];
inline lint pow(lint x,lint y,lint p)
{
lint res(1);
x %= p;
while(y)
{
if(y&1) res = (res*x)%p;
x = (x*x)%p;
y >>= 1;
}
return res%p;
}
inline lint C(lint n,lint m,lint p)
{
if(m > n) return 0;
if(!m) return 1;
lint res(1);
res = res * fac
% p;
res = res * pow(fac[m]*fac[n-m]%p,p-2,p) % p;
return res;
}
inline lint calc(lint n,lint m,lint p)
{
if(m > n) return 0;
if(m == 0 || m == n) return 1;
return calc(n/p,m/p,p)*C(n%p,m%p,p)%p;
}
inline lint get_ans(lint n,lint m)
{
lint res(0);
for(lint i = 0; i < 4; i++)
res = (res + (P/a[i])*pow(P/a[i],a[i]-2,P)%P*calc(n,m,a[i])%P)%P;
return (res%P+P)%P;
}
inline lint solve()
{
if(g%(P+1) == 0) return 0;
lint pw = 0;
for(lint i = 1; i * i <= n; i++)
if(n%i == 0)
{
if(i * i == n)
pw = (pw + get_ans(n,i))%P;
else
pw = (pw + get_ans(n,i) + get_ans(n,n/i))%P;
}
return pow(g,pw,P+1);
}
int main()
{
P = 999911658;
fac[0] = 1;
for(lint i = 1; i <= a[3]; i++)
fac[i] = fac[i-1] * i % P;
cin>>n>>g;
cout<<solve()<<endl;
return 0;
}
内容来自用户分享和网络整理,不保证内容的准确性,如有侵权内容,可联系管理员处理 点击这里给我发消息
标签: 
相关文章推荐