mint.cpp

#include <iostream>
#include <cstdint>
#include <chrono>
#include <vector>
const long long MOD = (long long)1e9+7;
using namespace std;

// 剰余を考慮した整数のクラス
template <uint_fast64_t Mod> class Mint {
    using u64 = uint64_t;
    public:
    u64 a;
    constexpr Mint(const u64 x = 0) noexcept : a(x % Mod){}
    constexpr void set(const u64 x = 0) noexcept { a = x % Mod; }
    constexpr u64 &value() noexcept { return a; }
    constexpr const u64 &value() const noexcept { return a; }
    constexpr Mint operator+(const Mint rhs) const noexcept {
        return Mint(*this) += rhs;
    }
    constexpr Mint operator-(const Mint rhs) const noexcept {
        return Mint(*this) -= rhs;
    }
    constexpr Mint operator*(const Mint rhs) const noexcept {
        return Mint(*this) *= rhs;
    }
    constexpr Mint operator/(const Mint rhs) const noexcept {
        return Mint(*this) /= rhs;
    }
    constexpr Mint &operator+=(const Mint rhs) noexcept {
        a += rhs.a;
        if (a >= Mod) a -= Mod;
        return *this;
    }
    constexpr Mint &operator++() noexcept {
        a += 1;
        if (a >= Mod) a -= Mod;
        return *this;
    }
    constexpr Mint operator++(int) noexcept {
        Mint m = *this;
        a += 1;
        if (a >= Mod) a -= Mod;
        return m;
    }
    constexpr Mint &operator-=(const Mint rhs) noexcept {
        if (a < rhs.a) a += Mod;
        a -= rhs.a;
        return *this;
    }
    constexpr Mint &operator--() noexcept {
        if (a < 0) a += Mod;
        a -= 1;
        return *this;
    }
    constexpr Mint operator--(int) noexcept {
        Mint m = *this;
        if (a < 0) a += Mod;
        a -= 1;
        return m;
    }
    constexpr Mint &operator*=(const Mint rhs) noexcept {
        a = a * rhs.a % Mod;
        return *this;
    }
    constexpr Mint &operator/=(const Mint rhs) noexcept {
        u64 exp = Mod - 2;
        while(exp) {
            if(exp % 2) {
                *this *= rhs;
            }
            rhs *= rhs;
            exp /= 2;
        }
        return *this;
    }
    constexpr operator long long() const {return a;}
    friend ostream& operator<< (ostream& os, const Mint<Mod> &m) {
        os << m.a;
        return os;
    };
    friend istream& operator>> (istream& is, Mint<Mod> &m) {
        long long l;
        is >> l;
        m.set(l);
        return is;
    };
};
typedef Mint<MOD> mint;
typedef vector<Mint<MOD>> vmint;

// 階乗の剰余
vmint factional(const int m) {
    vmint ret;
    mint n = 1;
    ret.push_back(n);
    for(mint i = 1; i.a <= m; i++)
    {
        n *= i;
        ret.push_back(n);
    }
    return ret;
}

// 繰り返し二乗和
mint RepSqMod(long long N, long long P) {
    if(P==0) return mint(1);
    if(P % 2 == 0) {
        mint t = RepSqMod(N, P / 2);
        return t * t;
    }
    return mint(N * RepSqMod(N, P - 1)).value();
}

// 階乗の剰余の逆元
vmint invfactional(const int m) {
    vmint ret;
    mint n = 1, tmp;
    ret.push_back(n);
    for(mint i = 1; i.a <= m; i++)
    {
        n *= i;
        tmp = RepSqMod(n, MOD - 2);
        ret.push_back(tmp);
    }
    return ret;
}

// 二項係数の剰余
mint ModComb(long long n, long long r, vmint fac, vmint inv) {
    return fac[n] * inv[r] * inv[n - r];
}

int main() {
    vmint fac, inv;
    int N;
    long long n, r;
    cin >> N >> n >> r;
    auto start = chrono::system_clock::now();
    fac = factional(N);
    auto end = chrono::system_clock::now();
    cout << "階乗の剰余の計算　　　 :" << (double)chrono::duration_cast<chrono::milliseconds>(end - start).count() << "ms" << endl;
    start = chrono::system_clock::now();
    inv = invfactional(N);
    end = chrono::system_clock::now();
    cout << "階乗の剰余の逆元の計算 :" << (double)chrono::duration_cast<chrono::milliseconds>(end - start).count() << "ms" << endl;

    
    start = chrono::system_clock::now();
    mint m = ModComb(n, r, fac, inv);
    end = chrono::system_clock::now();
    cout << "組み合わせの剰余計算　 :" << (double)chrono::duration_cast<chrono::milliseconds>(end - start).count() << "ms" << endl;
    cout << n << "C" << r << " mod " << MOD << " = " << m << endl;
}