結果

問題 No.8030 ミラー・ラビン素数判定法のテスト
ユーザー nonamae
提出日時 2022-07-12 16:42:48
言語 C++17
(gcc 13.3.0 + boost 1.87.0)
結果
RE  
実行時間 -
コード長 2,689 bytes
コンパイル時間 1,994 ms
コンパイル使用メモリ 194,744 KB
最終ジャッジ日時 2025-01-30 06:47:18
ジャッジサーバーID
（参考情報） 		judge5 / judge1
このコードへのチャレンジ
（要ログイン）
ファイルパターン 結果
other AC * 4 RE * 6
権限があれば一括ダウンロードができます
コンパイルメッセージ
main.cpp: In function ‘bool isPrime(uint64_t)’:
main.cpp:85:1: warning: control reaches end of non-void function [-Wreturn-type]
   85 | }
      | ^

ソースコード

diff # 

#include <bits/stdc++.h>

struct Runtime_Montgomery_Int64 {
    
    using m64 = Runtime_Montgomery_Int64;

    using i64 = std::int64_t;
    using u64 = std::uint64_t;
    using u128 = __uint128_t;

    inline static u64 m, r, n2;

    static void set_mod(u64 m) {
        assert(m < (1ull << 62));
        assert((m & 1) == 1);
        m64::m = m;
        n2 = -u128(m) % m;
        r = m;
        for (int _ = 0; _ < 5; _++) r *= 2 - m * r;
        r = -r;
        assert(r * m == -1ull);
    }

    static u64 reduce(u128 b) { return (b + u128(u64(b) * r) * m) >> 64; }

    u64 x;

    Runtime_Montgomery_Int64() : x(0) { }
    Runtime_Montgomery_Int64(u64 x) : x(reduce(u128(x) * n2)) { }
    u64 val() const {
        u64 y = reduce(x);
        return y >= m ? y - m : y;
    }
    m64 &operator+=(m64 y) {
        x += y.x - (m << 1);
        x = (i64(x) < 0 ? x + (m << 1) : x);
        return *this;
    }
    m64 &operator-=(m64 y) {
        x -= y.x;
        x = (i64(x) < 0 ? x + (m << 1) : x);
        return *this;
    }
    m64 &operator*=(m64 y) {
        x = reduce(u128(x) * y.x);
        return *this;
    }
    m64 operator+(m64 y) const { return m64(*this) += y; }
    m64 operator-(m64 y) const { return m64(*this) -= y; }
    m64 operator*(m64 y) const { return m64(*this) *= y; }
    bool operator==(m64 y) const { return (x >= m ? x - m : x) == (y.x >= m ? y.x - m : y.x); }
    bool operator!=(m64 y) const { return not operator==(y); }
    m64 pow(u64 n) const {
        m64 y = 1, z = *this;
        for ( ; n; n >>= 1, z *= z) if (n & 1) y *= z;
        return y;
    }
};

bool isPrime(const std::uint64_t x) {
    using u64 = std::uint64_t;
    using m64 = Runtime_Montgomery_Int64;
    if (x == 2 or x == 3 or x == 5 or x == 7) return true;
    if (x % 2 == 0 or x % 3 == 0 or x % 5 == 0 or x % 7 == 0) return false;
    if (x < 121) return x > 1;
    const u64 d = (x - 1) >> __builtin_ctzll(x - 1);
    m64::set_mod(x);
    const m64 one{1}, rev{x - 1};
    auto ok = [&](u64 a) {
        auto y = m64(a).pow(d);
        u64 t = d;
        while (y != one and y != rev and t != x - 1) y *= y, t <<= 1;
        if (y != rev and t % 2 == 0) return false;
        return true;
    };
    if (x < (1ull << 32)) {
        for (u64 a : { 2, 7, 61 }) if (not ok(a)) return false;
    }
    else {
        for (u64 a : { 2, 325, 9375, 28178, 450775, 9780504, 1795265022 }) {
            if (x <= a) return true;
            if (not ok(a)) return false;
        }
    }
}

int main() {
    std::uint64_t a, b;
    std::cin >> a;
    while (a--) {
        std::cin >> b;
        std::cout << b;
        std::cout << " " << (isPrime(b) ? "1" : "0") << "\n";
    }
}
0