xdoj 1227 Godv的数列(lucas,扩展lucas,中国剩余定理模版)

来源:互联网 发布:php商城订单 编辑:程序博客网 时间:2024/05/29 03:18

首先lucas定理是指对于一个

Lucas(n,m,p)=C(n%p,m%p)*Lucas(n/p,m/p,p)(p是质数)
解决组合的取模问题(解决mod的数字过大问题,排列取模可不能随便取)

需要注意的是对于lucas,需要计算出0-n,0-m的可以计算的全部组合,包括C(0,0)

然后如果取模对象p不是质数,那么可以用中国剩余定理来合,也就是分别如1001,对每个C(n,m)都分解成对7,11,13取模,然后用中国剩余定理合成出最后一个数再*

ll Extended_Euclid(ll a,ll b,ll &x,ll &y)//扩展欧几里得算法{ll d;if(b==0){x=1;y=0;return a;}d=Extended_Euclid(b,a%b,y,x);y-=a/b*x;return d;}ll Chinese_Remainder(ll a[],ll w[],ll len)//中国剩余定理  a[]存放余数  w[]存放两两互质的数{ll i,d,x,y,m,n,ret;ret=0;n=1;for (i=0;i<len;i++)n*=w[i];for (i=0;i<len;i++){m=n/w[i];d=Extended_Euclid(w[i],m,x,y);ret=(ret+y*m*a[i])%n;}return (n+ret%n)%n;}


#include<cstdio>#include<cmath>#include<stdio.h>#include<string.h>#include<algorithm>#include<vector>#include<set>#include<queue>#include<iostream>using namespace std;typedef long long ll;const int maxn=1e5+7;ll N[maxn][3];ll niv[maxn][3];ll pow1(ll a,ll b,ll mod){    if(b==0)return 1;    ll res1=pow1(a,b/2,mod)%mod;    res1=res1*res1%mod;    if(b%2==1)res1=res1*a%mod;    return res1;}ll MM[3];ll C[100][100][3];ll q1[maxn];ll ni(ll x,ll mod){    return pow1(x,mod-2,mod);}//Lucas要注意考虑C(i,j)得预处理全部的1-i和1-jll Lucas(ll n, ll m,ll k) //预处理了c[i][j][k]表示C(i,j)%MM[k]的结果,直接调用即可{    ll mod=MM[k];    if(m ==0)  return 1;    else return  (C[n%mod][m%mod][k]*Lucas(n/mod,m/mod,k))%mod;}ll Extended_Euclid(ll a,ll b,ll &x,ll &y)//扩展欧几里得算法{ll d;if(b==0){x=1;y=0;return a;}d=Extended_Euclid(b,a%b,y,x);y-=a/b*x;return d;}ll Chinese_Remainder(ll a[],ll w[],ll len)//中国剩余定理  a[]存放余数  w[]存放两两互质的数{ll i,d,x,y,m,n,ret;ret=0;n=1;for (i=0;i<len;i++)n*=w[i];for (i=0;i<len;i++){m=n/w[i];d=Extended_Euclid(w[i],m,x,y);ret=(ret+y*m*a[i])%n;}return (n+ret%n)%n;}//对于lucas,C[0][0]=1还是要的,因为lucas是会把其分解的//int main(){    //freopen("in.txt","r",stdin);    ll i,j,k,f1,f2,f3,f4,t1,t2,t3,t4,n,m;    int T;    MM[0]=7;MM[1]=11;MM[2]=13;    N[0][0]=N[0][1]=N[0][2]=1ll;    for(i=0;i<=2;i++)        for(j=1;j<MM[i];j++){            N[j][i]=(N[j-1][i]*j)%MM[i];        }    for(i=0;i<3;i++){        niv[MM[i]-1][i]=ni(N[MM[i]-1][i],MM[i]);    }    for(j=0;j<3;j++)    for(i=MM[j]-2;i>=1;i--)        niv[i][j]=(niv[i+1][j]*(i+1))%MM[j];   ll a[5],b[5];   b[0]=7;b[1]=11;b[2]=13;   for(k=0;k<3;k++)   for(i=0;i<MM[k];i++)        for(j=0;j<=i;j++)        if(i!=j&&j!=0)        C[i][j][k]=((N[i][k]*niv[i-j][k])%MM[k])*niv[j][k]%MM[k];        else        C[i][j][k]=1;   while(scanf("%d",&T)==1){   while(T--){    cin >>n;    n--;    for(i=0;i<=n;i++){        cin >>q1[i];        q1[i]%=1001;    }    ll Sum=0;    for(i=0;i<=n;i++){    for(j=0;j<3;j++)    a[j]=Lucas(n,i,j);    Sum=(Sum+q1[i]*Chinese_Remainder(a,b,3))%1001;    }    cout <<Sum << endl;   }   }    return 0;}


include <cstdio>include <iostream>using namespace std;typedef long long ll;ll p=1000000007;ll n,m,ans;ll fast_pow(ll x,ll y,ll p){ll ans=1ll;while(y){if (y&1)ans=(ans*x)%p;x=(x*x)%p;y>>=1;}return ans;}void exgcd(ll a,ll b,ll &x,ll &y){if (!b) x=1ll,y=0ll;else exgcd(b,a%b,y,x),y-=(a/b)*x;}ll inv(ll q,ll p){if (!q) return 0;ll a=q,b=p,x=0ll,y=0ll;exgcd(a,b,x,y);x=((x%b)+b)%b;if (!x) x+=b;return x;}ll fact(ll n,ll pi,ll pk){if (!n) return 1ll;ll ans=1ll;if (n/pk){for (ll i=2;i<=pk;i++)if (i%pi) ans=ans*i%pk;ans=fast_pow(ans,n/pk,pk);}for (ll i=2;i<=n%pk;i++)if (i%pi)ans=ans*i%pk;return ans*fact(n/pi,pi,pk)%pk;}ll C(ll n,ll m,ll p,ll pi,ll pk){if (m>n) return 0ll;ll a=fact(n,pi,pk),b=fact(m,pi,pk),c=fact(n-m,pi,pk);ll k=0ll,ans;for (ll i=n;i;i/=pi) k+=i/pi;for (ll i=m;i;i/=pi) k-=i/pi;for (ll i=n-m;i;i/=pi) k-=i/pi;ans=a*inv(b,pk)%pk*inv(c,pk)%pk*fast_pow(pi,k,pk)%pk;return ans*(p/pk)%p*inv(p/pk,pk)%p; }main(){cin>>n>>m>>p;for (ll x=p,i=2;i<=p;i++)if (x%i==0){ll pk=1ll;while(x%i==0) pk*=i,x/=i;ans=(ans+C(n,m,p,i,pk))%p;}cout<<ans;}





阅读全文
0 0
原创粉丝点击