結果

問題 No.3030 ミラー・ラビン素数判定法のテスト
ユーザー nonamaenonamae
提出日時 2022-07-21 01:06:54
言語 C++17
(gcc 12.3.0 + boost 1.83.0)
結果
WA  
実行時間 -
コード長 6,062 bytes
コンパイル時間 2,035 ms
コンパイル使用メモリ 200,360 KB
実行使用メモリ 4,380 KB
最終ジャッジ日時 2023-09-15 10:24:44
合計ジャッジ時間 2,706 ms
ジャッジサーバーID
（参考情報） 		judge14 / judge15
このコードへのチャレンジ
（要ログイン）

テストケース

テストケース表示
入力 結果 実行時間
実行使用メモリ
testcase_00 AC 2 ms
4,376 KB
testcase_01 AC 1 ms
4,380 KB
testcase_02 AC 1 ms
4,380 KB
testcase_03 AC 1 ms
4,376 KB
testcase_04 WA -
testcase_05 WA -
testcase_06 WA -
testcase_07 WA -
testcase_08 WA -
testcase_09 WA -
権限があれば一括ダウンロードができます

ソースコード

diff # 

#include <bits/stdc++.h>

using i32 = std::int32_t;
using i64 = std::int64_t;
using i128 = __int128_t;
using u32 = std::uint32_t;
using u64 = std::uint64_t;
using u128 = __uint128_t;

int jacobi_symbol(i64 a, u64 n) {
    u64 t;
    int j = 1;
    while (a) {
        if (a < 0) {
            a = -a;
            if ((n & 3) == 3) j = -j;
        }
        int s = __builtin_ctzll(a);
        a >>= s;
        if (((n & 7) == 3 || (n & 7) == 5) && (s & 1)) j = -j;
        if ((a & n & 3) == 3) j = -j;
        t = a, a = n, n = t;
        a %= n;
        if (u64(a) > n / 2) a -= n;
    }
    return n == 1 ? j : 0;
}

struct Runtime_modint64 {
    using m64 = u64;

    inline static m64 one, r2, n, mod;

    m64 x;

    static void set_mod(u64 m) {
        mod = m;
        one = (u64)-1ull % m + 1;
        r2 = (u128)(i128)-1 % m + 1;
        u64 nn = m;
        for (int _ = 0; _ < 5; ++_) nn *= 2 - nn * m;
        n = nn;
    }

    static m64 reduce(u128 a) {
        u64 y = (u64(a >> 64)) - (u64((u128(u64(a) * n) * mod) >> 64));
        return i64(y) < 0 ? y + mod : y;
    }

    Runtime_modint64() : x(0) { }
    Runtime_modint64(u64 x) : x(reduce(u128(x) * r2)) { }
    Runtime_modint64(u64 x, bool is_montgomery) : x(is_montgomery ? x : reduce(u128(x) * r2)) { }

    u64 get_val() const { return reduce(u128(x)); }
    u64 get_raw() const { return x; }

    Runtime_modint64 &operator+=(Runtime_modint64 y) {
        x += y.x - mod;
        if (i64(x) < 0) x += mod;
        return *this;
    }
    Runtime_modint64 &operator-=(Runtime_modint64 y) {
        if (i64(x -= y.x) < 0) x += 2 * mod;
        return *this;
    }
    Runtime_modint64 &operator*=(Runtime_modint64 y) {
        x = reduce(u128(x) * y.x);
        return *this;
    }
    Runtime_modint64 &operator/=(Runtime_modint64 y) {
        return *this *= y.inv();
    }
    Runtime_modint64 &operator++() {
        if (x + one >= mod) x = x + one - mod;
        else x += one;
        return *this;
    }
    Runtime_modint64 &operator--() {
        if (x < one) x = mod + x - one;
        else x -= one;
        return *this;
    }
    Runtime_modint64 &operator>>=(int k) {
        x >>= k;
        return *this;
    }
    Runtime_modint64 &operator<<=(int k) {
        x <<= k;
        return *this;
    }
    Runtime_modint64 operator+(Runtime_modint64 y) const { return Runtime_modint64(*this) += y; }
    Runtime_modint64 operator-(Runtime_modint64 y) const { return Runtime_modint64(*this) -= y; }
    Runtime_modint64 operator*(Runtime_modint64 y) const { return Runtime_modint64(*this) *= y; }
    Runtime_modint64 operator/(Runtime_modint64 y) const { return Runtime_modint64(*this) /= y; }
    Runtime_modint64 operator-() const { return Runtime_modint64() - Runtime_modint64(*this); }
    Runtime_modint64 operator++(int) {
        Runtime_modint64 temp = *this;
        ++*this;
        return temp;
    }
    Runtime_modint64 operator--(int) {
        Runtime_modint64 temp = *this;
        --*this;
        return temp;
    }
    Runtime_modint64 operator>>(int k) const { return Runtime_modint64(*this) >>= k; }
    Runtime_modint64 operator<<(int k) const { return Runtime_modint64(*this) <<= k; }
    Runtime_modint64 pow(u64 k) {
        Runtime_modint64 y = 1, z = *this;
        for ( ; k; k >>= 1, z *= z) if (k & 1) y *= z;
        return y;
    }
    Runtime_modint64 inv() { return (*this).pow(mod - 2); }

    bool operator==(Runtime_modint64 y) const { return (x >= mod ? x - mod : x) == (y.x >= mod ? y.x - mod : y.x); }
    bool operator!=(Runtime_modint64 y) const { return not operator==(y); }
    bool operator<(const Runtime_modint64& other) { return (*this).get_val() < other.get_val(); }
    bool operator<=(const Runtime_modint64& other) { return (*this).get_val() <= other.get_val(); }
    bool operator>(const Runtime_modint64& other) { return (*this).get_val() > other.get_val(); }
    bool operator>=(const Runtime_modint64& other) { return (*this).get_val() >= other.get_val(); }

};

int is_prime(u64 n) {
    
    {
        if (n <= 1) return 0;
        if (n <= 3) return 1;
        if (!(n & 1)) return 0;
    }
    
    Runtime_modint64::set_mod(n);

    Runtime_modint64 o(1);
    {
        u64 d = (n - 1) << __builtin_clzll(n - 1);
        Runtime_modint64 t = o << 1;
        for (d <<= 1; d; d <<= 1) {
            t *= t;
            if (d >> 63) t <<= 1;
        }
        if (t != o) {
            u64 x = (n - 1) & -(n - 1);
            Runtime_modint64 mone = -o;
            for (x >>= 1; t != mone; x >>= 1) {
                if (x == 0) return 0;
                t *= t;
            }
        }
    }

    {
        i64 D_ = 5;
        for (int i = 0; jacobi_symbol(D_, n) != -1 && i < 64; i++) {
            if (i == 32) {
                u32 k = round(sqrtl(n));
                if (k * k == n) return 0;
            }
            if (i & 1) D_ -= 2;
            else D_ += 2;
            D_ = -D_;
        }

        Runtime_modint64 Q(D_ < 0 ? (1 - D_) / 4 % n : n - (D_ - 1) / 4 % n);

        Runtime_modint64 u, v, Qn;
        u64 k = (n + 1) << __builtin_clzll(n + 1);
        u = o;
        v = o;
        Qn = Q;
        D_ %= (i64)n;
        Runtime_modint64 D(D_ < 0 ? n + D_ : D_);

        for (k <<= 1; k; k <<= 1) {
            u *= v;
            v = (v * v) - (Qn + Qn);
            Qn *= Qn;
            if (k >> 63) {
                u64 uu = u.get_raw() + v.get_raw();
                if (uu & 1) uu += n;
                uu >>= 1;
                v = D * u + v;
                v >>= 1;
                u = uu;
                Qn *= Q;
            }
        }

        if (u.get_val() == 0 || v.get_val() == 0) return 1;
        u64 x = (n + 1) & ~n;
        for (x >>= 1; x; x >>= 1) {
            u *= v;
            v = (v * v) - (Qn + Qn);
            if (v.get_val() == 0) return 1;
            Qn *= Qn;
        }
    }

    return 0;
}


int main() {
    i32 T; scanf("%d", &T);
    while (T--) {
        u64 x;
        scanf("%llu", &x);
        printf("%llu %d\n", x, is_prime(x));
    }
}
0