WA 30pts求调 有注释 码风优良
查看原帖
WA 30pts求调 有注释 码风优良
860817
Starstream楼主2023/11/16 21:57

Record

#include <iostream>
#define int long long

using namespace std;

const int N = 100010;

int n, m, p;
int primes[N], cnt;
int po[N]; // pow(p, po)

int exgcd(int a, int b, int &x, int &y) // exgcd: ax + by = gcd(a, b)
{
    if (!b)
    {
        x = 1, y = 0;
        return a;
    }
    int d = exgcd(b, a % b, y, x);
    y -= a / b * x;
    return d;
}

int inv(int t, int p) // inv of t in module p
{
    int x, y;
    exgcd(t, p, x, y);
    return (x % p + p) % p;
}

int qpow(int a, int b, int p) // quick_pow: a ^ b mod p
{
    int res = 1 % p;
    a %= p;
    while (b)
    {
        if (b & 1) res = res * a % p;
        b >>= 1;
        a = a * a % p;
    }
    return res;
}

int fact(int n, int p, int a) // n! mod pow(p, a)
{
    if (!n) return 1;
    int fr = 1, rm = 1; // circulate and remain
    for (int i = 1; i <= a; i ++ )
        if (i % p) fr = fr * i % a;
    fr = qpow(fr, n / a, a);
    for (int i = n / a * a; i <= n; i ++ )
        if (i % p) rm = rm * (i % a) % a;
    return fact(n / p, p, a) * fr % a * rm % a;
}

int get(int n, int p) // number of factor p in n!
{
    if (n < p) return 0;
    return get(n / p, p) + n / p;
}

int C(int n, int m, int p, int a) // C(n, m) mod pow(p, a)
{
    int fn = fact(n, p, a); // n!
    int ufm = inv(fact(m, p, a), a); // m!
    int ufnm = inv(fact(n - m, p, a), a); // (n - m)!
    int k = get(n, p) - get(m, p) - get(n - m, p); // C pow
    int res = qpow(p, k, a); // C
    return fn * ufm % a * ufnm % a * res % a;
}

int CRT(int p)
{
    int res = 0;
    for (int i = 1; i <= cnt; i ++ )
    {
        int mi = p / primes[i], ti = inv(mi, primes[i]);
        res = (res + po[i] * mi % p * ti % p) % p;
    }
    return res;
}

int exlucas(int n, int m, int p)
{
    int tmp = p; // divide p into primes
    for (int i = 2; i <= tmp / i; i ++ )
        if (tmp % i == 0)
        {
            int a = 1;
            while (tmp % i == 0) a ++ , tmp /= i;
            primes[ ++ cnt] = a;
            po[cnt] = C(n, m, i, a);
        }
    if (tmp > 1) primes[ ++ cnt] = tmp, po[cnt] = C(n, m, tmp, tmp);

    return CRT(p);
}

signed main()
{
    scanf("%lld%lld%lld", &n, &m, &p);
    printf("%lld\n", exlucas(n, m, p));
    return 0;
}
2023/11/16 21:57
加载中...