exLucas

/*
    exLucas : 求 C(n,m)%P (P不保证是质数&&P<=1e6) 的问题

    算法
    设 P= p1^a1 * p2^a2 * ... pr^ar <-唯一分解定理
    求出
    {
        C(n,m) % p1^a1 
        C(n,m) % p2^a2
        ...
        C(n,m) % pr^ar
    }
    最后用中国剩余定理合并

    分解成要求 C(n,m) % p^k(p为质数)
                        n! 
             so =   ———————————— % p^k
                      m!(n-m)!
    


    由于式子是在模 p^k 意义下,所以分母要算乘法逆元.

    我们无法求得m!,(n−m)!的逆元, why? m!可能含有 p 因子.

    同余方程 a*x ≡ 1(%P)(即乘法逆元)有解 的充要条件为 <-裴蜀定理：gcd⁡(a,b)=1

    so我们换一个式子：

         n!
        —————
         p^x
    —————————————— * p ^ (x-y-z) % p^k
     m!    (n-m)!
    ——— * ————————
    p^y      p^z
    x 表示 n! 中包含多少个p因子,y,z同理.

    现在问题转化为问题等价于求
      n!
    ————— % p^k
     p^x

    我们对n!进行变形:

    n!=1*2*3*...*n
    =(P*2P*3P...)(1*2*...)
    左边是1∼n中所有 <=n的 P 的倍数,右边反之.

    因为1∼n中有 ⌊Pn​⌋ 个 P 的倍数,

    so n! = p^(⌊n/P​⌋) * (1*2*3*...) * (1*2*3*...)
                                  (     n        )
          = p^(⌊n/P​⌋) * (⌊n/P​⌋)! *(     ∏ i      )
                                  ( i=1,(i%p!=0) )
    显然后面这个 % P 是有循环节的.
                                    (    p^k       )^⌊n/(p^k)​⌋    (           n              )
          = p^(⌊n/P​⌋) * (⌊n/P​⌋)! *  (     ∏ i      )              (           ∏ i            )
                                    ( i=1,(i%p!=0) )              ( i=p^k*⌊n/(p^k)​⌋,(i%p!=0) )
    其中

        (    p^k       )^⌊n/(p^k)​⌋
        (     ∏ i      )
        ( i=1,(i%p!=0) )

    是循环节1∼P中所有无 P 因子数的乘积.

        (           n              )
        (           ∏ i            )
        ( i=p^k*⌊n/(p^k)​⌋,(i%p!=0) )

    是余项的乘积.


    发现前面这个 P^⌊n/P⌋ 是要除掉的.

    然而(⌊n/P⌋)!里可能还包含 P .

    所以,我们定义函数：
             n!
    f(n) = ———————
             P^x

  so 
  
                                     (    p^k       )^⌊n/(p^k)​⌋    (           n              )
        f(n) = f(⌊n/P​⌋) * (⌊n/P​⌋) *  (     ∏ i      )              (           ∏ i            )
                                     ( i=1,(i%p!=0) )              ( i=p^k*⌊n/(p^k)​⌋,(i%p!=0) )
    这样就好了.时间复杂度是O(log⁡Pn).
    边界f(0)=1
    看看原式


         n!
        —————
         p^x
    —————————————— * p ^ (x-y-z) % p^k
     m!    (n-m)!
    ——— * ————————
    p^y      p^z

           f(n)
======= ——————————— * p ^ (x-y-z) % p^k
        f(m)*f(n-m)
    
    由于f(m)一定与P^k互质,所以我们可以直接用exgcd求逆元.

    How to solve P ^ (x-y-z)?
    
    easy!

    For example : we want to get f(n) = n! / (p^x) 中的x
    设g(n)=x(above)
    再看看阶乘的式子
                                    (     n        )
        n!  = p^(⌊n/P​⌋) * (⌊n/P​⌋)! *(     ∏ i      )
                                    ( i=1,(i%p!=0) )
    显然 we want p^(⌊n/P​⌋)

    而且(⌊nP⌋)!可能还有P因子。

    所以我们可以得到以下递推式:
        g(n)=⌊n​/P⌋+g(⌊n​/P⌋)
    
    时间复杂度是O(log⁡(n/p))

    边界g(n)=0(n<P)

    所以答案就是
                   n!/(p^x)
    =    ———————————————————————————— * P^{g(n)-g(m)-g(n-m)} % (p^k) 
            m!/(p^y) * (n-m)!/(p^z)

    最后用中国剩余定理合并答案即可得到C(n,m)

*/
#include<bits/stdc++.h>
#define Bessie moo~
#define int long long
using namespace std;
int read(){
    int a=0,fl=1;
    char c;
    c=getchar();
    while(c<'0'||c>'9')fl= c=='-' ? -1 : 1, c=getchar();
    while('0'<=c&&c<='9')a=(a<<3)+(a<<1)+(c^48),c=getchar();
    return a*fl;
}
void exgcd(int a,int b,int &x,int &y){
    if(!b){
        x=1,y=0;
        return ;
    }
    exgcd(b,a%b,x,y);
    int t=x;
    x=y,y=t-a/b*y;
}
int INV(int a,int MOD){
    int x,y;
    exgcd(a,MOD,x,y);
    return (x+MOD)%MOD;
}
int qpow(int a,int b,int MOD){
    int base=1;
    a%=MOD;
    while(b){
        if(b&1)base=(base*a)%MOD;
        b>>=1;
        a=(a*a)%MOD;
    }
    return base;
}
int F(int n,int MOD,int PK){
    if(n==0)return 1;
    int rou=1,//循环节
        rem=1;//余项
    for(int i=1;i<=PK;i++){
        if(i%MOD)rou=rou*i%PK;
    }
    rou=qpow(rou,n/PK,PK);

    for(int i=PK*(n/PK);i<=n;i++){
        if(i%MOD)rem=rem*(i%PK)%PK;
    }
    return F(n/MOD,MOD,PK)*rou%PK*rem%PK;
}
int G(int n,int MOD){
    if(n<MOD)return 0;
    return G(n/MOD,MOD)+(n/MOD);
}
int C_PK(int n,int m,int MOD,int PK){
    int fz=F(n,MOD,PK),fm1=INV(F(m,MOD,PK),PK),fm2=INV(F(n-m,MOD,PK),PK);
    int mi=qpow(MOD,G(n,MOD)-G(m,MOD)-G(n-m,MOD),PK);
    return fz*fm1%PK*fm2%PK*mi%PK;
}
int A[1005],B[1005];
//x≡B(MOD A)
int exLucas(int n,int m,int MOD){
    int ljc=MOD,tot=0;
    for(int tmp=2;tmp*tmp<=MOD;tmp++){
        if(!(ljc%tmp)){
            int PK=1;//p^k
            while(!(ljc%tmp)){
                PK*=tmp;
                ljc/=tmp;
            }
            A[++tot]=PK;
            B[tot]=C_PK(n,m,tmp,PK);
        }
    }
    if(ljc!=1){
        A[++tot]=ljc;
        B[tot]=C_PK(n,m,ljc,ljc);
    }
    int ans=0;
    for(int i=1;i<=tot;i++){//CRT
        int M=MOD/A[i],T=INV(M,A[i]);
        ans=(ans+B[i]*M%MOD*T%MOD)%MOD;
    }
    return ans;
}

signed main(){
    int n=read(),m=read(),MOD=read();
    printf("%lld",exLucas(n,m,MOD));
    return 0;
}