卢卡斯定理 $Lucas{~}Theory {~} O(log_{p}Nplogp)$

更好的阅读体验

$$ 先放结论 C_{a}^{b} (lucas)≡ C_{\frac{a}{p}}^{\frac{b}{p}}(lucas) C_{a{~}mod{~}p}^{b{~}mod{~}p}(mod{~}p)\\\\ $$

$$ a=a_{k}p^{k}+a_{k-1}p^{k-1}+…+a_{0}p^{0}\\\\ b=b_{k}p^{k}+b_{k-1}p^{k-1}+…+b_{0}p^{0}–①\\\\ $$
$$ \begin{align} (1+x)^{p^{k}}&=C_{p^k}^{0} \times 1+C_{p^k}^{1} \times x+…+C_{p^k}^{p^k} \times x^{p^k}\\\\ (C_{p^k}^{1} \sim C_{p^k}^{p^k-1} mod p = 0)&=1+x^{p^{k}} (mod{~}p)-② \end{align} $$
$$ \begin{align} (1+x)^{a}&=(1+x)^{a_{0}}((1+x)^{p})^{a_{1}}((1+x)^{p^{2}})^{a_{2}}…((1+x)^{p^{k}})^{a_{k}}\\\\ ((1+x)^{p^{k}}=1+x^{p^{k}} (mod{~}p)-式②)&=(1+x)^{a_{0}}(1+x^{p})^{a_{1}}(1+x^{p^{2}})^{a_{2}}…(1+x^{p^{k}})^{a_{k}} (mod{~}p)\\\\ &=C_{a_{0}}^{b_{0}} x^{b_{0}p^{0}}…C_{a_{k-1}}^{b_{k-1}}x^{b_{k-1} p^{k-1}}C_{a_{k}}^{b_{k}}x^{b_{k}p^{k}}+其他项\\\\ &=C_{a_{0}}^{b_{0}}… C_{a_{k-1}}^{b_{k-1}}C_{a_{k}}^{b_{k}} x^{b_{k}p^{k}+b_{k-1}p^{k-1}+…+b_{0}p^{0}}+其他项\\\\ (式①)&=C_{a_{0}}^{b_{0}}… C_{a_{k-1}}^{b_{k-1}}C_{a_{k}}^{b_{k}} x^{b}+其他项 \end{align} $$
∴有上式中等式左边$(1+x)^{a}$和右边累乘的$x^{b}$的系数分别为:
$$ C_{a}^{b} || C_{a_{0}}^{b_{0}}… C_{a_{k-1}}^{b_{k-1}}C_{a_{k}}^{b_{k}}(mod p) $$
结合①可知,
$$a_{0} = a \% p,b_{0} = b \% p \\\\ a_{1} = \frac{a}{p} \% p,b_{1} = \frac{b}{p} \% p\\\\ …\\\\ a_{k} = \frac{a}{p^{k}} \% p, b_{k} = \frac{b}{p^{k}} \% p $$
体现在代码中则只要$\frac{a}{p} 或者 \frac{b}{p} > p$就继续$lucas$递归下去直到$\frac{a}{p} 和 \frac{b}{p} < p$,递归的过程相当于自上向下将$C_{a_0}^{b_0}−>C_{a_k}^{b_k}$添加到乘式里，递归终点为$a_{k} < p {~}and{~}b_{k}<p$

#include<iostream>
#include<algorithm>

using namespace std;

typedef long long LL;

int qmi(int a,int k,int p)
{
    int res = 1;
    while(k)
    {
        if(k&1)res = (LL)res*a%p;
        a = (LL)a*a%p;
        k>>=1;
    }
    return res;
}

int C(int a,int b,int p)//自变量类型int
{
    if(b>a)return 0;//漏了边界条件
    int res = 1;
    // a!/(b!(a-b)!) = (a-b+1)*...*a / b! 分子有b项
    for(int i=1,j=a;i<=b;i++,j--)//i<=b而不是<
    {
        res = (LL)res*j%p;
        res = (LL)res*qmi(i,p-2,p)%p;
    }
    return res;
}
//对公式敲
int lucas(LL a,LL b,int p)
{
    if(a<p && b<p)return C(a,b,p);//lucas递归终点是C_{bk}^{ak}
    return (LL)C(a%p,b%p,p)*lucas(a/p,b/p,p)%p;//a%p后肯定是<p的,所以可以用C(),但a/p后不一定<p 所以用lucas继续递归
}

int main()
{
    int n;
    cin >> n;
    while(n--)
    {
        LL a,b;
        int p;
        cin >> a >> b >> p;
        cout << lucas(a,b,p) << endl;
    }
    return 0;
}