結果

問題 No.8030 ミラー・ラビン素数判定法のテスト
ユーザー kaoru murata
提出日時 2021-09-01 16:46:40
言語 C
(gcc 13.3.0)
結果
WA  
実行時間 -
コード長 2,013 bytes
コンパイル時間 1,554 ms
コンパイル使用メモリ 31,104 KB
実行使用メモリ 6,912 KB
最終ジャッジ日時 2024-11-27 21:58:44
合計ジャッジ時間 78,602 ms
ジャッジサーバーID
（参考情報） 		judge2 / judge3
このコードへのチャレンジ
（要ログイン）
ファイルパターン 結果
other WA * 3 TLE * 7
権限があれば一括ダウンロードができます

ソースコード

diff # 

#include <stdint.h>
#include <stdio.h>
#include <stdbool.h>

typedef int64_t i64;
typedef __int128_t i128;
typedef uint64_t u64;
typedef __uint128_t u128;

typedef u64 m64;

u64 mul_inv(u64 n) {
  u64 x = 1;
  for (int i = 6; i > 0; --i) {
    x = x * (2 - x * n);
  }
  return x;
}
m64 calc_mod(u64 m) {
  return m;
}
m64 calc_inv(u64 m) {
  return mul_inv(calc_mod(m));
}
m64 calc_r2(u64 m) {
  return -(u128)(calc_mod(m)) % calc_mod(m);
}
m64 MR(u128 x, m64 mod, m64 inv) {
  m64 y = (u64)(x >> 64) - (u64)(((u128)((u64)x * inv) * mod) >> 64);
  return (i64)y < 0 ? y + mod : y;
}
m64 toM64(u64 w, m64 mod, m64 inv, m64 r2) { return MR((u128)w * r2, mod, inv); }
u64 fromM64(m64 x, m64 mod, m64 inv) { return MR((u128)x, mod, inv); }
m64 addmod(m64 x, m64 y, m64 mod) {
  return (x + y >= mod) ? x + y - mod: x + y;
}
m64 submod(m64 x, m64 y, m64 mod) {
  return (i64)(x - y) < 0 ? x - y + mod : x - y;
}
m64 mulmod(m64 x, m64 y, m64 mod, m64 inv) {
  return MR((u128)x * y, mod, inv);
}
m64 powmod(m64 x, int i, m64 mod, m64 inv, m64 r2) {
  m64 ret = toM64(1, mod, inv, r2);
  for (; i; i >>= 1) {
    if (i & 1) ret = mulmod(ret, x, mod, inv);
    x = mulmod(x, x, mod, inv);
  }
  return ret;
}

bool is_prime(u64 N) {
  if (N < 2 || N % 6 % 4 != 1) return (N | 1) == 3;
  int s = __builtin_ctzll(N - 1);
  u64 d = (N - 1) >> s;
  u64 as[] = {2, 325, 9375, 28178, 450775, 9780504, 1795265022};
  m64 inv = calc_inv(N);
  m64 r2 = calc_r2(N);
  m64 e = toM64(1, N, inv, r2);
  m64 f = toM64(N - 1, N, inv, r2);
  for (int i = 0; i < 7; i++) {
    m64 a = toM64(as[i], N, inv, r2);
    m64 p = powmod(a, d, N, inv, r2);
    if (a == 0 || p == e || p == f) continue;
    for (int i = s - 1; p != f && i >= 0; --i) p = powmod(p, 2, N, inv, r2);
    if (p != f) return false;
  }
  return true;
}

int main() {
  int n; scanf("%d", &n);
  for (int i = 0; i < n; i++) {
    u64 x; scanf("%lu", &x);
    if (is_prime(x)) {
      printf("%d\n", 1);
    }
    else {
      printf("%d\n", 0);
    }
  }
  return 0;
}
0