zoukankan      html  css  js  c++  java
  • 拓展卢卡斯定理模板完整注释

    #define _CRT_SECURE_NO_WARNINGS
    #include<bits/stdc++.h>
    using namespace std;
    //#define ll long long
    typedef long long ll;
    //拓展gcd,求解ax+by=gcd(a,b)
    void ex_gcd(ll a, ll b, ll & x, ll &y) {
        if (!b) {
            x = 1; y = 0;
        }
        else {
            ex_gcd(b, a%b, y, x);
            y -= x * (a / b);
        }
    }
    //求解a关于p的逆元,其中gcd(a,p)==1
    ll inverse(ll a, ll p) {
        if (!a)return 0ll;
        ll x, y;
        ex_gcd(a, p, x, y);
        x = (x%p + p) % p;
        return x;
    }
    //快速幂a^b
    ll qpow(ll a, ll b, ll p) {
        ll res = 1;
        while (b) {
            if (b & 1)res = res * a%p;
            a = a * a%p;
            b >>= 1;
        }
        return res;
    }
    //中国剩余定理,求解x=a[i]%m[i]方程组
    ll CRT(int n, ll *a, ll *m) {
        ll M = 1, res = 0;
        for (int i = 0; i < n; i++)M = M * m[i];
        for (int i = 0; i < n; i++) {
            ll x, y;
            ll tm = M / m[i];
            ex_gcd(tm, m[i], x, y);
            res = (res + tm * x*a[i] % M) % M;
        }
        return (res + M) % M;
    }
    //x表示质数,P表示分解的x^t
    //calc计算n!%P,P不一定质数,复杂度P*log(x)(n)
    ll calc(ll n, ll x, ll P) {
        if (!n)return 1;
        ll s = 1;
        //线性时间内计算1~P中不是x的倍数的乘积
        //由于P不是质数,不能用威尔逊定理
        //如果多组,这里可以提前打表优化,不过好像也不太方便,P总数太多
        for (ll i = 2; i <= P; i++) {
            if (i%x)s = s * i%P;
        }
        s = qpow(s, n / P, P);//前面一部分的数关于P下面的循环
        //这里是计算多余的,也就是n%P剩下的
        for (ll i =2; i <=n%P; i++)
            if (i%x)s = s*i%P;
        return s * calc(n / x, x, P) % P;
    }
    //x表示质数,P表示分解的x^t
    //计算C(n,m)%P,P为一个质数的幂
    ll multilucas(ll n, ll m, ll x, ll P) {
        ll cnt = 0;
        //下面分别用log的复杂度三个阶乘中含x的指数数目
        //这一块不懂得话可以了解n!%p的相关性质
        for (register ll i = n; i>0; i /= x)
            cnt += i / x;
        for ( ll i = m; i>0; i /= x)cnt -= i / x;
        for (ll i = n - m; i>0; i /= x)cnt -= i / x;
        return qpow(x, cnt, P) % P*calc(n, x, P) % P*inverse(calc(m, x, P), P) % P*inverse(calc(n - m, x, P), P) % P;
    }
    //求解C(n,m)%P的最终答案
    ll ex_lucas(ll n, ll m, ll P) {
        ll cnt = 0;
        //a数组保存C(n,m)对每一种分解之后质因数求余的结果
        //由于每一种质因数的t次方是互质的,所以用CRT即可求解
        //p数组保存每一类质因数分解之后的答案
        ll p[20], a[20];//最多质因数分解20个
        for (ll i = 2; i*i <= P; i++) {
            //提前打素数表可以优化一点
            if (P%i == 0) {
                p[cnt] = 1;
                while (P%i == 0) { p[cnt] = p[cnt] * i; P /= i; }
                a[cnt] = multilucas(n, m, i, p[cnt]);
                cnt++;
            }
        }
        //最后是一个质数
        if (P > 1)p[cnt] = P, a[cnt++] = multilucas(n, m, P, P);
        return CRT(cnt, a, p);
    }
    int main()
    {
        ll n, m, p;
        scanf("%lld%lld%lld", &n, &m, &p);
        printf("%lld
    ", ex_lucas(n, m, p));
        cin >> n;
        return 0;
    }
    不疯魔不成活
  • 相关阅读:
    Netty应用
    原生JDK网络编程- NIO之Reactor模式
    Kafka入门教程
    Java CAS
    Java读写锁
    遍历map的四种方法
    java selector
    Range Sum Query
    Increasing Triplet Subsequence
    Count Numbers with Unique Digits
  • 原文地址:https://www.cnblogs.com/gzr2018/p/11607182.html
Copyright © 2011-2022 走看看