P2183 [国家集训队]礼物 题解


以下题解仅供学习参考使用。

抄袭、复制题解,以达到刷AC率/AC数量或其他目的的行为,在洛谷是严格禁止的。

洛谷非常重视学术诚信。此类行为将会导致您成为作弊者。具体细则请查看洛谷社区规则

评论

  • cellur925
    %%%黑题dalao
  • nosta
    讲的钛好了!%%%
  • cellur925
    讲的钛好了!%%%
  • 可爱即是正义
    tql%%%
  • noi题库正在
    tql
  • Hola
    ~~抄的~~
  • cplusplus 
    您佩服自己干嘛?
  • baoyu
    我佩服我自己(
  • TLE自动机
    大佬可以用有个东西叫LATEX
作者: ___new2zy___ 更新时间: 2018-08-15 16:25  在Ta的博客查看 举报    8  

题解 P2183 [国家集训队]礼物

题目传送门:

https://www.luogu.org/problemnew/show/P2183

==========================================

不得不说这是我做过的最正经的一道黑题了= =

在这里先吐槽一下:数论真恶心

很是佩服这题的第一个题解,就是有点小蒙,我来给大家讲一下我自己的理解方式

希望同学们能理解本题的核心数论算法

扩展卢卡斯定理(ExLucas)

(这题其实可以当个模板题来做= =)

前置知识:

扩展欧几里得(Exgcd),乘法逆元,

中国剩余定理(CRT),Lucas定理(不必备)

(其实理解了做法也没什么难的对吧)(快逃)

==========================================

简化版题面:

你手中一共有n件礼品,你有m个好(基)

你打算送给每个人wi件礼品

请你求出送出礼品的方案数,并让答案对p取模

请注意,礼物之间两两不同

==========================================

可能本人废话比较多,欢迎来吐槽

看完题目emmm。。。让我冷静一下= =

我们不妨假设所有朋友收到的礼物数之和sum:

sum=∑(i:1~m)wi

如果记答案为ans的话,那么有:

ans={C(n,sum)*C(sum,w1)*C(sum-w1,w2)*...*C(wi,wi)}%p

在上式中,C(n,m)代表组合数,每一项组合数乘在一起就是ans

那么如果n<sum显然无解即"impossible"

(礼物不够当然没法送出去啦QAQ)

emmm....继续冷静思考,发现好像很简单?(逃

(只要暴力递推计算就好啦)

但是,如果我们直接递推暴力求解ans的话,那么会发现时间复杂度不能接受

题目中给定的n的范围是1e9,预处理时直接递推求阶乘的话估计就要跑个很长时间吧= =所以这种做法显然是不可取的

那么这样就很难受,继续思考= =

"n很大,还要求组合数取模,那么用Lucas定理吧"

灵机一动,想到在数论问题中,对于类似于

C(n,m)%p(必须保证p为prime)的问题,

我们有一个定理:Lucas定理

这个定理大概长这个样子:

*Lucas(n,m)=C(n%p,m%p)Lucas(n/p,m/p)**

其中Lucas(n,m)其实就是要计算的式子C(n,m)%p

对于这个定理,在此就不证明了(我比较菜所以不会证明)

(注意:Lucas定理适用于n<=1e5)

回到这题上来。。。

读者会发现,,好像我之前说的是废话= =

这个题n<=1e9好吧!!你在逗我?没法做。。

先淡定一下,到了这里其实已经离正解不远了,毕竟求解有了一定的方向

再次陷入沉思= =发现我们的Lucas定理不可以直接使用,因为题目中还有一个很致命的限制条件:

p<=1e9(此外什么也没说)

那么就抛来一个很重要的问题:模数p不保证是质数

那么导致我们如果直接用逆元取模的话不行,因为可能在这些阶乘中出现取模,这样答案直接就变为0了,因此我们要尝试把这些数提取出来。

emmm。。。恐怖如斯= =根本没法做了啊

但是我们发现,题目中给了提示:

设P=p1^c1*p2^c2*p3^c3*…*pt^ct,pi为质数。

其实这是数论中的唯一分解定理

任意大于1的正整数N,存在唯一分解式N=p1^c1*p2^c2*...*pi^ci

其中p1~pi是质数,^是次方,c1~ci是次数

那么又可以很自然的想到:

我们可以将分子分母都对于p的唯一分解式中的每一项(即pi^ci)取模(就保证是prime了),最后将每一项用CRT合并就得到了解

我们不妨再来考虑一下组合数取模的计算公式:

*C(n,m)%p={fac(n)/[fac(m)fac(n-m)]}%p**

(注意,这里只是形式,实际上除法不能直接取模)

其中fac(x)表示x的阶乘,即x!

如果拆分的话,我们可以对于n!,m!,(n-m)!分别%p再进行合并即可求出答案

又可以注意到,其实整个问题就剩下了一个式子:

求解(x!)%(pi^ci),保证pi为质数

对于这个式子展开可以发现,它的求解分为三部分:

  1. 在x!中,若存在pi的倍数那么就可以拆出来,与pi^ci抵消,即可以转化为子问题:求(x/pi)!%pi^ci,只需要递归下去求解就行了

  2. n!中可能有会包含多个完整的1~pi^ci-1(在%pi^ci下的剩余系),这部分可以先预处理一下(pi^ci-1)! (注意是不包含pi的倍数的阶乘,因为我们还要记录一下阶乘中pi的次数最后一起处理),然后pi的若干次幂用快速幂处理。

  3. 那么进行完1,2两步之后,对于剩下的部分直接求一发逆元计算组合数,最后合并即可(扩欧求逆元,中国剩余定理合并)

有些懵逼?我来举个栗子食用更佳哈:

假设pi=5,展开阶乘可以发现:

n!= (1234678911…)5^(n/5)*(n/5)!

(乘号*没了555,变成了斜体请自行脑补)

那么我们直接提取出所有pi=5的倍数,正好分为三部分,直接按步骤上面计算就行了

其实以上就是(最恶心的)数论算法:扩展卢卡斯定理(ExLucas)

讲到这里。。。突然发现这东西好像和卢卡斯定理一点关系没有= =(所以Lucas不是必备的知识)

大概思路讲完了,我们在放代码之前可以来总结一下这次的求解方法,即扩展卢卡斯定理(ExLucas)

扩展卢卡斯定理适用于计算任意的C(n,m)%p问题,其中n,m,p为任意正整数

求解过程主要应用了乘法逆元扩欧Exgcd求)、中国剩余定理(CRT也就只有一行)、带取模快速幂(码量很小)

总体来说思路还是很清楚的吧

那么好啦,既然有了解决方案,题目里的式子就很好求了吧,放个代码code~~~

PS:代码里也有解释哦

#include<iostream>
#include<cstdio>
#include<cmath>
#include<algorithm>
#include<queue>
#include<cstring>
using namespace std;
typedef long long ll;
const int inf=1e9+7;
inline ll poww(ll a,ll b,ll mod)//快速幂 
{
    ll base=a,ans=1;
    while(b)
        {
          if(b&1)
             ans=(1ll*ans*base)%mod;
          base=(1ll*base*base)%mod;
          b>>=1;
        }
    return 1ll*ans;
}
inline ll read()//快读
{
    ll p=0,f=1;char c=getchar();
    while(c<'0'||c>'9'){if(c=='-')f=-1;c=getchar();}
    while(c>='0'&&c<='9'){p=p*10+c-'0';c=getchar();}
    return 1ll*f*p;}
ll n,m,sum,mod,ans=1,A[19];
inline void Exgcd(ll a,ll b,ll &x,ll &y)
//扩欧用来求逆元 
{
    if(!b){x=1,y=0;return ;}
    Exgcd(b,a%b,y,x);y-=a/b*x;
}
inline ll rev(ll k,ll p)
//求k在mod p下的逆元(转换一下变成正整数)
{
    if(!k)return 0;
    ll x=0,y=0,a=k,b=p;
    Exgcd(a,b,x,y);
    x=(x%b+b)%b;
    if(!x)x+=b;
    return 1ll*x;
}
inline ll mul(ll n,ll p,ll pk)
//求n!%pi^ci,pk=pi^ci
{
    if(!n)return 1;
    ll ans=1;
    for(ll i=2;i<=pk;i++)
        if(i%p)ans=ans*i%pk;
    ans=poww(ans,n/pk,pk);
    for(ll i=2;i<=n%pk;i++)
        if(i%p)ans=ans*i%pk;
    return 1ll*ans*mul(n/p,p,pk)%pk;
    //递归下去求解(n/pi)!%pi^ci 
}
inline ll C(ll n,ll m,ll mod,ll p,ll pk)
//求C(n,m)%mod,其中唯一分解之后质因子为p,总乘积为pk(pi^ci)
{
    if(m>n)return 0;
    ll a=mul(n,p,pk),b=mul(m,p,pk),c=mul(n-m,p,pk),k=0;
    //先求一下n!%pi^ci,m!%pi^ci,(n-m)!%pi^ci 
    for(ll i=n;i;i/=p)k+=i/p;
    for(ll i=m;i;i/=p)k-=i/p;
    for(ll i=n-m;i;i/=p)k-=i/p;
    //先除掉n!,m!,(n-m)!在%mod下的质因子p
    ll ans=1ll*a*rev(b,pk)%pk*rev(c,pk)%pk*poww(p,k,pk)%pk;
    //除去质因子p之后直接逆元求组合数(剩余部分) 
    return 1ll*ans*(mod/pk)%mod*rev(mod/pk,pk)%mod;
    //找到逆元了再乘回去(CRT合并)
}
int main()
{
    /*这是我在做这个题的时候写的解释
    思路:不妨设sum=∑wi(i:1~m)
    题目很简单,就是要求式子: 
    C(n,sum)*C(sum,w1)*C(sum-w1,w2)*...*C(wi,wi)
    但组合数要取模 
    自然想到逆元和Lucas定理 
    但是模数p不保证是质数,所以要用扩展Lucas
    不妨先把p唯一分解,对于每一个pi^ci都两两互质 
    可以求出C(n,m)%pi^ci
    对于求这个,除掉质因子pi然后求逆元 
    之后CRT一下乱搞就好啦 
    */
    mod=read(),n=read(),m=read();
    for(int i=1;i<=m;i++)
        A[i]=read(),sum+=A[i];
    if(sum>n)//要送出的礼物多于有的礼物显然无解 
        {
            printf("Impossible");
            return 0;
        }
    for(ll k=1;k<=m;k++)//枚举每一个人要的礼物 
        {
            n-=A[k-1];
            ll now=0,x=mod;
            for(ll i=2;i<=sqrt(mod);i++)
                //找到mod的每一个质因数p 
                if(!(x%i))
                 {
                      ll pk=1;
                      while(!(x%i))pk*=i,x/=i;//除掉质因数p
                      now=(now+C(n,A[k],mod,i,pk))%mod;
                      //求出C(n,A[k])%pi^ci累加
                 }
            if(x>1)now=(now+C(n,A[k],mod,x,x))%mod;
            ans=ans*now%mod;//统计答案 
        }
    printf("%lld",ans);//愉快的输出答案
    return 0;
}

好了,到这里大概就讲完了吧= =

真香,又是一道毒瘤呢(逃

感谢阅读~~~

最后来推广一下我的blog:

https://www.luogu.org/blog/new2zy/

拜拜~~~ >=<

评论

  • 没名字可被用
    顺带一提,这个程序分解质因数所用的数组太大了,35就够了.
  • 突然颓废
    %%%大佬
  • steven张
    %%%%%%tql
  • 萌新qyz
    自古蓝名出dalao
  • xsl666
    合数的合写成了和
  • zhangchengkai
    %%%lty
作者: 没名字可被用 更新时间: 2018-01-25 23:05  在Ta的博客查看 举报    8  

没人在这里发这题题解啊,我来发一篇。

明显答案是$C_n^{w_1}\times C_{n-w_1}^{w_2} \times \cdots$,如果要求的礼物总数超过了买的礼物,则直接输出Impossible

求解一个组合数$C_n^m$对$p_i^{k_i}$的过程可以分别求分子和分母,然后再进行运算,即,求出$n!\mod p_i^{k_i}$,$m!\mod p_i^{k_i}$,$(n-m)!\mod p_i^{k_i}$后合并。注意分子可以含有大于$k_i$个$p_i$,这样的话求出来就直接为0了,而除掉分母的时候是可能会除掉一些$p_i$的;同时分母本身也可能有含有$p_i$的约数,于是统一把$p_i$移除来,不乘进答案,接下来用分子有的$p_i$的质数减掉分母的,然后乘上去即可。不可能出现减掉后为负数的情况,因为那样意味着组合数除不尽;有可能出现减掉后仍然大于$k_i$的情况,那么这个组合数在模$p_i^{k_i}$下为$0$。

具体地,对于求$n! \mod p_i^{k_i}$可以分为几个部分求解。将每个$p_i$的倍数提出一个$p_i$来,组成$p_i^{\lfloor \frac{n}{p_i} \rfloor}$;对于剩余的部分,容易发现除掉了$p_i$的数(原来是$p_i$的倍数)组成了一个新的阶乘,即$\lfloor \frac{n}{p_i} \rfloor !$,这个我们可以递归求解。剩下的数形如$1,2,\cdots,p-1,p+1,p+2,\cdots,2p-1 \cdots$(因为$p_i$的倍数已经通过前两步瓜分到阶乘和次幂里去了),由于$a \times b \mod p=(a \mod p) \times (b \mod p)$,易得$\prod _{i=1}^{p_i^{k_i}-1}i[i \mod p \neq 0]$与$\prod _{i=p_i^{k_i}+1}^{2p_i^{k_i}-1}i[i \mod p \neq 0] $在$\mod p_i^{k_i}$下是相等的,那么只需求一个$p^k$的长度的乘积在模$p_i^{k_i}$下的结果即可,而这个阶乘里包含了$ \lfloor \frac{n}{p_i^{k_i}} \rfloor$个这个余数,快速幂即可;还有剩下来的、没有达到$p_i^{k_i}$的一段数,由于长度不会超过$p_i^{k_i}$(即$10^5$),直接计算即可。

考虑求解答案,因为模的是一个和数,所以不能直接用普通的$lucas$求解。将和数$P$分解,得到$P=\prod_{i=1}^{k}p_i^{k_i}$,我们要求的答案模$p_i^{k_i}$的结果,对于单独求组合数在模每个$p_i^{k_i}$的剩余系下的结果要相同。于是考虑对于每个组合数,单独求出对$p_i^{k_i}$的结果,然后通过中国剩余定理求解即可。

代码

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const int maxn=6;
const int maxp=(int)1e5+5;
int P,n,m,w[maxn],pri[maxp],cnt[maxp],tot=0,c[maxp],R[maxp],M[maxp];
map<int ,int > re;
inline int maybs(int x) {return x<0?-x:x;}
void div(int x)
{
    for(int i=2;i*i<=x;i++) {
        if(x%i==0) {
            pri[++tot]=i;re[i]=tot;
            while(x%i==0) x/=i,cnt[tot]++;
        }
    }
    if(x!=1) {
        if(re.find(x)==re.end()) pri[++tot]=x,re[x]=tot;
        cnt[re[x]]++;
    }
}
int quickpow(int x,int k)
{
    int ret=1;
    while(k>0) {
        if(k&1) ret=ret*x;
        x=x*x; k>>=1;
    }
    return ret;
}
int quickpow(int x,int k,int mo)
{
    int ret=1;
    while(k>0) {
        if(k&1) ret=1ll*ret*x%mo;
        x=1ll*x*x%mo; k>>=1;
    }
    return ret;
}
int gcd(int a,int b) {return b==0?a:gcd(b,a%b);}
int extgcd(int a,int b,int &x,int &y)
{
    if(b==0) {x=1,y=0; return a;}
    else {int d=extgcd(b,a%b,y,x);y-=a/b*x;return d;}
}
int inv(int q,int p)
{
    assert(gcd(p,q)==1);
    int x,y;
    extgcd(q,p,x,y);
    x=(x%p+p)%p; if(!x) x+=p;
    return x;
}
int Cnt(int x,int p) {
    if(x==0) return 0;
    return Cnt(x/p,p)+x/p;
}
typedef pair<int ,int > pii;
#define fir first
#define sec second
#define MP make_pair
pii Cal(int x,int p,int k)
{
    if(x==1 || x==0) return MP(1,0);
    int ret=1,del=0,tmp=quickpow(p,k);
    int cou=Cnt(x,p);
    del=cou;

    pii nw=Cal(x/p,p,k);
    ret=1ll*ret*nw.fir%tmp;

    if(x>=tmp) {
        int lev=1;
        for(int i=1;i<tmp;i++) {
            if(i%p==0) continue;
            lev=1ll*lev*i%tmp;
        }
        ret=1ll*ret*quickpow(lev,x/tmp,tmp)%tmp;
    }

    for(int i=x;i>=1;i--) {
        if(i%tmp==0) break;
        if(i%p==0) continue;
        ret=1ll*ret*i%tmp;
    }
    return MP(ret,del);
}
int getAns()
{
    for(int i=1;i<=tot;i++) M[i]=1;
    for(int i=1;i<=tot;i++) {
        for(int j=1;j<=tot;j++)
            if(i==j) continue;
            else M[i]=1ll*M[i]*R[j]%P;
    }
    int ret=0;
    for(int i=1;i<=tot;i++) ret=(ret+((1ll*c[i]*M[i]%P)*1ll*inv(M[i],R[i]))%P)%P;
    return ret;
}
int main()
{
    scanf("%d",&P);
    div(P);
    scanf("%d%d",&n,&m);
    int sum=0;
    for(int i=1;i<=m;i++) {scanf("%d",&w[i]);sum+=w[i];}
    if(sum>n) {printf("Impossible\n");exit(0);}
    int ans=1,N=n,M,res=1;
    pii up,dw1,dw2;
    for(int i=1;i<=m;i++) {
        M=w[i];
        for(int j=1;j<=tot;j++) R[j]=quickpow(pri[j],cnt[j]);
        for(int j=1;j<=tot;j++) {
            ans=1;
            up=Cal(N,pri[j],cnt[j]);
            dw1=Cal(M,pri[j],cnt[j]),dw2=Cal(N-M,pri[j],cnt[j]);
            up.sec-=dw1.sec; up.sec-=dw2.sec;
            assert(up.sec>=0);

            if(up.sec>=cnt[j]) ans=0;
            else {
                ans=1ll*ans*up.fir%R[j];
                ans=1ll*ans*inv(dw1.fir,R[j])%R[j]; ans=1ll*ans*inv(dw2.fir,R[j])%R[j];
                ans=1ll*ans*quickpow(pri[j],up.sec,R[j])%R[j];
            }
            c[j]=ans;
        }
        res=1ll*res*getAns()%P;
        N-=w[i];
    }
    printf("%d\n",res);
    return 0;
}

评论

  • EternHope
    个人感觉,这么考虑是最简单的!
  • ιχγббб
    式子很好懂
作者: da_AA 更新时间: 2018-09-27 16:04  在Ta的博客查看 举报    4  

看上去我的式子和别人的不一样所以我就来水一波233

我是这样考虑的:

对于这几个礼物,我们可以枚举排列,然后前$w_1$个给第一个人,之后的$w_2$个给第二个人……但是这样会有重复的,即有可能某个人拿的礼物的顺序不一样但是本质是相同的,这样的话,就要再除以$w_i!$,不要忘了剩下没分的那部分也有重的,所以最后的式子就是$\dfrac{n!}{(n-tot)!w_1!w_2!w_3!w_4!w_5!}$。

这样运用扩展lucas的思想就可以算了。

#include <cstdio>
#include <cmath>

typedef long long LL;

LL n, m, P;
LL w[6], tot;

void exgcd(LL a, LL b, LL &x, LL &y, LL &d) {
    if (!b) {
        d = a;
        x = 1;
        y = 0;
    } else {
        exgcd(b, a%b, y, x, d);
        y -= a / b * x;
    }
}

LL pow_mod(LL x, LL p, LL mod) {
    LL ans = 1;
    for (; p; p>>=1, x=x*x%mod) if (p&1) ans = ans * x % mod;
    return ans;
}

LL fac(LL n, LL p, LL pk) {
    if (!n) return 1;
    LL ans = 1;
    for (int i = 2; i <= pk; ++i) if (i%p) ans = ans * i % pk;
    ans = pow_mod(ans, n/pk, pk);
    for (int i = 2; i <= n%pk; ++i) if (i%p) ans = ans * i % pk;
    return ans * fac(n/p, p, pk) % pk;
}

LL inv(LL n, LL mod) {
    LL x, y, d;
    exgcd(n, mod, x, y, d);
    return (x%=mod) < 0 ? x+mod : x;
}

LL CRT(LL b, LL mod) {
    return b*(P/mod)%P*inv(P/mod, mod)%P;
}

LL calc(LL x, LL p) {
    LL k = 0;
    for (; x; x /= p) k += x / p;
    return k;
}

LL C(LL p, LL pk) {
    LL u = fac(n, p, pk), d = inv(fac(n-tot, p, pk), pk);
    for (int i = 1; i <= m; ++i) {
        d = d * inv(fac(w[i], p, pk), pk) % pk;
    }
    LL k = 0;
    k += calc(n, p);
    k -= calc(n-tot, p);
    for (int i = 1; i <= m; ++i) {
        k -= calc(w[i], p);
    }
    return u * d % pk * pow_mod(p, k, pk) % pk;
}

void work() {
    if (tot > n) {
        puts("Impossible");
        return;
    }
    LL ans = 0, tmp = P, pk;
    LL lim = sqrt(P) + 5;
    for (int i = 2; i <= lim; ++i) if (tmp % i == 0) {
        pk = 1;
        while (tmp % i == 0) pk*=i, tmp/=i;
        ans = (ans + CRT(C(i, pk), pk)) % P;
    }
    if (tmp > 1) {
        ans = (ans + CRT(C(tmp, tmp), tmp)) % P;
    }
    printf("%lld\n", ans);
}

int main() {
    scanf("%lld%lld%lld", &P, &n, &m);
    for (int i = 1; i <= m; ++i) scanf("%lld", &w[i]), tot += w[i];
    work();
    return 0;
}

评论

  • zrzzrz
    前排资瓷神仙
作者: eee_hoho 更新时间: 2019-06-26 18:22  在Ta的博客查看 举报    1  

题不是很难,黑题应该虚高

看到题之后,应该就能写出式子

$$ans=\prod_{i=1}^{m}C_{n-\sum_{j=1}^{i-1}w_j}^{w_i}(mod\ P)$$

式子看起来很麻烦,其实很好想

当我们给第一个人送礼物的时候,我们有$n$个礼物,要选$w_1$个,方案就有$C_n^{w_{1}}$种

然后给第二个人送礼物的时候,我们剩下$n-w_1$个礼物,要选$w_2$个,方案有$C_{n-w_1}^{w_2}$种

就这样一直送到第$m$个人,而根据乘法原理,就得到这个式子了

这道题似乎就做完了,但在算$C_n^m(mod\ p_i^{c_i})$时又遇到了瓶颈:$p_i^{c_i}$不是质数

既然不是质数,我们就不能用$Lucas$定理来求,这就需要用到$ExLucas$了

虽然是扩展的,但和$Lucas$完全沾不上边

我们观察到题目给了$P=\prod _{i=1}^{t}p_{i}^{c_i}$

而$p_i$是质数,那么所有的$p_i^{c_i}$都是互质的

那么我们只要对于每个$p_i^{k=c_i}$求出其$C_n^m\ mod\ p_i^k$的值,然后用中国剩余定理就可以算出最小正整数解

问题转化成了如何快速的求$C_n^m\ mod\ p^k$

把组合数展开,我们得到

$$\frac{n!}{m!(n-m)!}\ mod\ p^k$$

那么只要快速求出模意义下的阶乘就好了

如果我们要求$x!\ mod\ p^k$,假设$x=17,p=2,k=3$,观察一下式子

$$17!=1\times2\times3\times4\times5\times6\times7\times8\times9\times10\times11\times12\times13\times14\times15\times16\times17$$

化式子$17!=$

$$2\times4\times6\times8\times10\times12\times14\times16\times1\times3\times5\times7\times9\times11\times13\times15\times17$$

$$(2\times1)\times(2\times2)\times(2\times3)\times…\times(2\times8)\times(1\times3\times5\times…\times17)$$

$$2^8\times(1\times2\times3\times…\times8)\times(1\times3\times5\times…\times17)$$

$$2^8\times8!\times(1\times3\times5\times7\times9\times11\times13\times15\times17)$$

而在模$p^k=2^3=8$意义下$1\times3\times5\times7=9\times11\times13\times15$

所以式子就变成了

$$2^8\times8!\times(1\times3\times5\times7)^2\times17$$

那么$x!$就变成了$p^n\times t!\times(a_1\times a_2\times…a_m)^r\times a_{m+1}\times…a_{q}$

可以看出$n=t=\left \lfloor \frac{x}{p} \right \rfloor,r=\left \lfloor \frac{x}{p^k} \right \rfloor$,其中$t!$可以一直递归求解,然后后面的暴力算

这样这个题就做完啦QAQ

Code

#include <iostream>
#include <cstdio>
#include <algorithm>
#include <cstring>
#define int long long
using namespace std;
int p[200000],n,m,P,cnt,z[200000],a[200000],ans;
void exgcd(int a,int b,int &x,int &y)     //扩欧
{
    if (!b)x=1,y=0;
    else
    {
        exgcd(b,a%b,x,y);
        int t=x;
        x=y;
        y=t-a/b*y;
    }
}
int mypow(int a,int x,int p)      //快速幂
{
    int s=1;
    while (x)
    {
        if (x&1)s=s*a%p;
        a=a*a%p;
        x>>=1;
    }
    return s;
}
int fac(int n,int a,int b)    //求阶乘
{
    if (!n)return 1;       
    int s=1;
    for (int i=1;i<=b;i++)    //处理a1*a2*…am
        if (i%a!=0)
            s=s*i%b;
    s=mypow(s,n/b,b);
    for (int i=1;i<=n%b;i++)  //处理am+1*…aq
        if (i%a!=0)
            s=s*i%b; 
    return s*fac(n/a,a,b)%b;  //递归求解
}
int inv(int a,int b)    //求逆元
{
    int x,y;
    exgcd(a,b,x,y);
    return (x%b+b)%b;
}
int C(int m,int n,int a,int b)   //处理组合数
{
    int nn=fac(n,a,b),mm=fac(m,a,b),nm=fac(n-m,a,b),po=0;  //求阶乘
    for (int i=n;i;i/=a)      //处理n^p中的p
        po+=i/a;
    for (int i=m;i;i/=a)
        po-=i/a;
    for (int i=n-m;i;i/=a)
        po-=i/a;
    return nn*mypow(a,po,b)%b*inv(mm,b)%b*inv(nm,b)%b; 
}
signed main()
{
    cin>>n>>m>>P;
    int pp=P;
    for (int i=2;i*i<=P;i++)     //分解
        if (pp%i==0)
        {
            int x=i;
            z[++cnt]=i;
            while (pp%x==0)x*=i;
            x/=i;
            p[cnt]=x;
            pp/=x;
        }
    if (pp!=1)p[++cnt]=pp,z[cnt]=pp;
    for (int i=1;i<=cnt;i++)
        a[i]=C(m,n,z[i],p[i]);
    for (int i=1;i<=cnt;i++)    //中国剩余定理
    {
        int pi=P/p[i],x,y;
        exgcd(pi,p[i],x,y);
        ans=((ans+pi*a[i]*x%P)+P)%P;
    }
    cout<<ans<<endl;
    return 0;
}

评论

  • 还没有评论
作者: M_sea 更新时间: 2018-12-30 11:10  在Ta的博客查看 举报    2  

这题最多是紫题吧,式子我都能想出来qwq

前置姿势:扩展Lucas

不会的左转模板区

首先可以推出答案为$\large C_{n}^{w-1}\times C_{n-w_1}^{w_2}\times... \mod P$

组合意义很明确,就是在剩下的礼物中选出$w_i$个的方案数,再根据乘法原理乘起来。

显然我们要每步都对算出的组合数取模,但是这里的$P$不一定是质数,所以直接上$\color{red}{\text{扩展Lucas}}$即可。

细节见代码

 
反馈
如果你认为某个题解有问题,欢迎向洛谷反馈,以帮助更多的同学。



请具体说明理由,以增加反馈的可信度。