bzoj 3738: [Ontak2013]Kapitał_给出三个数字 n,m,k,求 c(n+m,n) 去掉所有末尾的 0 后对 10^k 取模的结果。k<-优快云博客

本文链接：https://blog.youkuaiyun.com/herobrine_tkj/article/details/78693441

本文介绍了一种利用中国剩余定理求解大组合数的方法，通过将问题分解为模2^9和模5^9两部分，分别计算后再合成最终答案。文章详细解释了如何通过阶乘模质数幂的方法来解决模数非质数的情况。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

题意已经够简短就不概括了吧。
题解：
k最大是9，那就直接当所有都 $\%10^9$ 吧。答案再保留k位就好。
把 $10^9$ 拆成 $2^9 \times 5^9$ ，假设我们知道 $C_{N+M}^N\%2^9$ 和 $C_{N+M}^N\%5^9$ 的值，我们就可以用中国剩余定理或扩展欧几里得求出答案。具体来说，就是设 $x=C_{N+M}^N,\begin{cases}x \equiv a\pmod{2^9}\\\ x \equiv b\pmod{5^9}\end{cases}$ ，那么 $x=k\cdot 2^9\cdot 5^9+a\cdot 5^9\cdot inv(5^9,2^9)+b\cdot 2^9\cdot inv(2^9,5^9)$ ，模掉 $10^9$ 就只剩后面两项了。

那么怎么求 $C_{N+M}^N\%2^9$ 和 $C_{N+M}^N\%5^9$ 呢？Lucas定理是行不通了，因为模数不是质数。分子分母分别求出来再除（乘逆元）的话，可能会出现这种情况：结果是 $\frac{a\cdot p^\alpha}{b\cdot p^\beta}$ ，当 $\alpha$ 和 $\beta$ 都 $\ge9$ 时上下都=0，而答案可能不等于0。所以我们要保留 $\alpha$ 和 $\beta$ ，而 $a$ 和 $b$ 就可以随便模。求 $C$ 就可以直接用阶乘模质数的幂来搞了。阶乘模质数的幂这样来做：

n! = 1 \times 2 \times 3 \times \dots \times n = 1 \times \dots \times (p - 1) \times (p + 1) \times (p + 2) \times \dots \times (2 p - 1) \times (2 p + 1) \times \dots \times n \times p ⌊ n / p ⌋ \times (1 \times 2 \times \dots \times ⌊ n / p ⌋) = 1 \times \dots \times (p - 1) \times (p + 1) \times (p + 2) \times \dots \times (2 p - 1) \times (2 p + 1) \times \dots \times n \times p ⌊ n / p ⌋ \times ⌊ n / p ⌋!

$n!\\\ =1\times2\times3\times\cdots\times n\\\ =1\times\cdots\times(p-1)\times(p+1)\times(p+2)\times\cdots\times(2p-1)\times(2p+1)\times\cdots\times n\times p^{\lfloor n/p\rfloor}\times (1\times2\times\cdots\times \lfloor n/p\rfloor)\\\ =1\times\cdots\times(p-1)\times(p+1)\times(p+2)\times\cdots\times(2p-1)\times(2p+1)\times\cdots\times n\times p^{\lfloor n/p\rfloor}\times \lfloor n/p\rfloor!$
前面的部分预处理，最后那部分可以递归下去搞，返回一个

a⋅px $a\cdot p^x$ 。除的时候把x相减，乘上a的逆元就好了。
好像已经做完了。
代码：

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

long long n,m,k,fac1[513],fac2[1953126];
struct num
{
    long long x,y;//x*p^y
};

num times(num x,num y,long long p)
{
    return (num){x.x*y.x%p,x.y+y.y};
}
long long pow(long long x,long long y,long long p)
{
    if(y==0)
    return 1;
    long long ans=pow(x,y>>1,p);
    ans=ans*ans%p;
    if(y&1)
    ans=ans*x%p;
    return ans;
}
long long exgcd(long long a,long long b,long long&x,long long&y)
{
    if(b==0)
    {
        x=1;
        y=0;
        return a;
    }
    long long tx,ty,d=exgcd(b,a%b,tx,ty);
    x=ty;
    y=tx-(a/b)*ty;
    return d;
}
long long inv(long long a,long long p)//ax%p=1
{
    long long x,y;
    exgcd(a,p,x,y);
    x=(x%p+p)%p;
    return x;
}
num inv(num x,long long p)
{
    x.x=inv(x.x,p);
    x.y*=-1;
    return x;
}
num div(num x,num y,long long p)
{
    return times(x,inv(y,p),p);
}
num wk(long long n,long long p,long long pk,long long* fac)//n!%p^k  pk=p^k
{
    if(n==0)
    return (num){1,0};
    long long hh=pow(fac[pk],n/pk,pk)*fac[n%pk]%pk;
    return times((num){hh,n/p},wk(n/p,p,pk,fac),pk);
}
int main()
{
    scanf("%lld%lld%lld",&n,&m,&k);
    fac1[0]=1;
    for(int i=1;i<=512;i++)
    fac1[i]=fac1[i-1]*((i&1)?i:1)%512;
    fac2[0]=1;
    for(int i=1;i<=1953125;i++)
    fac2[i]=fac2[i-1]*((i%5)?i:1)%1953125;
    num a=wk(n+m,2,512,fac1);
    a=div(a,wk(n,2,512,fac1),512);
    a=div(a,wk(m,2,512,fac1),512);
    num b=wk(n+m,5,1953125,fac2);
    b=div(b,wk(n,5,1953125,fac2),1953125);
    b=div(b,wk(m,5,1953125,fac2),1953125);
    long long hh=min(a.y,b.y);
    long long ans2=a.x*pow(2,a.y-hh,512)%512*pow(inv(5,512),hh,512)%512;
    long long ans5=b.x*pow(5,b.y-hh,1953125)%1953125*pow(inv(2,1953125),hh,1953125)%1953125;
    long long ans=(ans2*1953125%1000000000*inv(1953125,512)%1000000000+ans5*512*inv(512,1953125)%1000000000)%1000000000;
    char s[10];
    sprintf(s,"%%0%lldlld",k);
    printf(s,ans%pow(10,k,1000000001));
}