結果

問題 No.8030 ミラー・ラビン素数判定法のテスト
ユーザー wasd314
提出日時 2025-05-19 08:58:28
言語 C++23
(gcc 13.3.0 + boost 1.87.0)
結果
AC  
実行時間 56 ms / 9,973 ms
コード長 4,921 bytes
コンパイル時間 1,254 ms
コンパイル使用メモリ 102,784 KB
実行使用メモリ 6,272 KB
最終ジャッジ日時 2025-05-19 08:58:31
合計ジャッジ時間 2,538 ms
ジャッジサーバーID
（参考情報） 		judge2 / judge1
このコードへのチャレンジ
（要ログイン）
ファイルパターン 結果
other AC * 10
権限があれば一括ダウンロードができます

ソースコード

diff # 

#include <algorithm>
#include <bit>
#include <cassert>
#include <concepts>
#include <iostream>
#include <limits>
#include <type_traits>
#include <vector>
#include <cstdint>

namespace wasd314
{
    struct dynamic_modint1 {
        using U = uint64_t;
        using UU = __uint128_t;
        using I = std::make_signed_t<U>;
        using mint = dynamic_modint1;

        U rx;
        static inline U mod;
        static inline U R_1;
        static inline U R1;
        static inline U R2;
        static inline U N_;
        static constexpr int bits_u = std::numeric_limits<U>::digits;
        static constexpr int bits_half = bits_u / 2;

        static U get_N_()
        {
            U n_inv = mod;
            for (int bits = 2; bits < bits_u; bits <<= 1) {
                n_inv *= 2 - n_inv * mod;
            }
            return -n_inv;
        }
        static U get_R1() { return -mod % mod; }
        static U get_R2() { return -UU(mod) % mod; }
        static U get_R_1() { return (1 + UU(mod) * N_) >> bits_u; }

        static void set_mod(U new_mod)
        {
            assert(I(new_mod) > 0);
            assert(new_mod & 1);
            mod = new_mod;
            N_ = get_N_();
            R1 = get_R1();
            R2 = get_R2();
            R_1 = get_R_1();
        }

        static U safe_mod(I x)
        {
            x %= I(mod);
            if (x < 0) x += I(mod);
            return x;
        }
        static U reduce(const UU &x)
        {
            U y = (x + UU(U(x) * N_) * mod) >> bits_u;
            return y < mod ? y : y - mod;
        }

        dynamic_modint1(const I &x) : rx(reduce(UU(safe_mod(x)) * R2)) {}
        // for literal
        dynamic_modint1(I &&x) : rx(reduce(UU(safe_mod(x)) * R2)) {}

       private:
        dynamic_modint1(const U &x, auto) : rx(x) {}

       public:
        static mint raw(const U &x) { return mint(x, 0); }

        U val() const { return reduce(rx); }
        mint pow(uint64_t e) const
        {
            mint ans = raw(R1), b(*this);
            while (e) {
                if (e & 1) ans *= b;
                b *= b;
                e >>= 1;
            }
            return ans;
        }

        mint &operator+=(const mint &o)
        {
            rx += o.rx;
            if (rx >= mod) rx -= mod;
            return *this;
        }
        mint &operator-=(const mint &o)
        {
            if (__builtin_sub_overflow(rx, o.rx, &rx)) rx += mod;
            return *this;
        }
        mint &operator*=(const mint &o)
        {
            rx = reduce(UU(rx) * o.rx);
            return *this;
        }
        mint operator+() const { return mint(*this); }
        mint operator-() const { return mint::raw(0) -= *this; }
    };
    using mint = dynamic_modint1;
    bool operator==(const mint &x, const mint &y) { return x.rx == y.rx; }
    bool operator!=(const mint &x, const mint &y) { return !(x == y); }
    mint operator+(const mint &x, const mint &y) { return mint(x) += y; }
    mint operator-(const mint &x, const mint &y) { return mint(x) -= y; }
    mint operator*(const mint &x, const mint &y) { return mint(x) *= y; }

    using lint = long long;
    using u64 = uint64_t;
    bool is_prime(u64 n)
    {
        // using lint2 = __int128_t;
        if (n < 2) return false;
        for (u64 p : {2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37}) {
            if (n == p) return true;
            if (n % p == 0) return false;
        }
        if (n < 41 * 41) return true;
        mint::set_mod(n);

        const mint one = mint::raw(mint::R1), neg_one = -one;

        auto test_miller_rabin = [&](const std::vector<lint> &bases) {
            int e = std::countr_zero(n - 1);
            u64 o = n >> e;
            for (lint b : bases) {
                mint x = mint(b).pow(o);
                if (x == one) continue;
                for (int ei = 0; ei < e; ++ei) {
                    mint y = x * x;
                    if (y == one) {
                        if (x == neg_one) break;
                        return false;
                    }
                    x = y;
                    if (ei == e - 1) return false;
                }
            }
            return true;
        };
        if (n < 2047) return test_miller_rabin({2});
        if (n < 9080191) return test_miller_rabin({31, 73});
        if (n < 4759123141) return test_miller_rabin({2, 7, 61});
        if (n < 1122004669633) return test_miller_rabin({2, 13, 23, 1662803});
        if (n < 3770579582154547) return test_miller_rabin({2, 880937, 2570940, 610386380, 4130785767});
        return test_miller_rabin({2, 325, 9375, 28178, 450775, 9780504, 1795265022});
    }
}  // namespace wasd314


int main()
{
    using namespace std;
    using namespace wasd314;
    int q;
    cin >> q;
    for (int _ = 0; _ < q; ++_) {
        lint n;
        cin >> n;
        cout << n << ' ' << is_prime(n) << "\n";
    }
}
0