#481438 (C++17) No.8030 ミラー・ラビン素数判定法のテスト

• 問題: No.8030 ミラー・ラビン素数判定法のテスト
• ユーザー: meowricator
• 提出日時: 2020-05-14 22:58:04
• 言語: C++17 (gcc 13.3.0 + boost 1.87.0)
• 結果: TLE
• 実行時間: -
• コード長: 21,959 bytes
• コンパイル時間: 3,519 ms
• コンパイル使用メモリ: 143,736 KB
• 最終ジャッジ日時: 2025-01-10 10:56:33
• ジャッジサーバーID （参考情報）: judge3 / judge2

ファイルパターン	結果
other	AC * 4 RE * 4 TLE * 2

ソースコード

// https://yukicoder.me/problems/no/3030

#include <iostream>
#include <string>
#include <cmath>
#include <type_traits>
#include <cassert>
#include <iomanip>
#include <vector>
#include <complex>
#include <random>

using i64 = int_fast64_t;
using u64 = uint_fast64_t;
using i32 = int_fast32_t;
using u32 = uint_fast32_t;

namespace mlib {

struct constexpr_string {
  constexpr constexpr_string(const char *str): str(str) {}

  const char *str;

  constexpr bool strcmp_begin(const constexpr_string st) const {
	return st.str[0] == '\x00' || (this->str[0] != '\x00' && this->str[0] == st.str[0] && constexpr_string(this->str + 1).strcmp_begin(st.str + 1));
  }
};

template <typename T>
constexpr bool is_type_in_root_namespace(const constexpr_string name) {
	return constexpr_string(__PRETTY_FUNCTION__ + (sizeof("constexpr bool mlib::is_type_in_root_namespace(mlib::constexpr_string) [with T = ") - 1)).strcmp_begin(name);
}

}

namespace mlib {

constexpr double pi = std::acos(-1);

using std::pow;
template<typename T, typename U, std::enable_if_t<is_type_in_root_namespace<T>("mlib"), std::nullptr_t> = nullptr>
T pow(T x, U n) {
	T res = 1;
	T a = x;
	while(n != 0) {
		if(n % 2 == 1) res *= a;
		a *= a;
		n /= 2;
	}
	return res;
}

template<typename T>
T mod_pow(const T &x, const T &n, const T &mod) {
	T res = 1;
	T a = x;
	T ll = n;
	while(ll != 0) {
		if(ll % 2 == 1) res = mod_mul(res, a, mod);
		a = mod_mul(a, a, mod);
		ll /= 2;
	}
	return res;
}

// @mod needs prime.
template <typename T>
T mod_inv(T x, T mod) {
	assert(mod > 0);
	if(mod == 1) return 0;
	return mod_pow(x, mod - 2, mod);
}

template <typename T>
T mod_add(T lhs, T rhs, T mod) {
	assert(mod > 0);
	lhs %= mod, rhs %= mod;
	lhs += rhs;
	if(lhs >= mod) lhs -= mod;
	return lhs;
}

template <typename T>
T mod_sub(T lhs, T rhs, T mod) {
	assert(mod > 0);
	lhs %= mod, rhs %= mod;
	if(lhs < rhs) lhs += mod;
	lhs -= rhs;
	return lhs;
}

template <typename T>
T mod_mul(T lhs, T rhs, T mod) {
	assert(mod > 0);
	return ((lhs % mod) * (rhs % mod)) % mod;
}

// @mod needs prime.
template <typename T>
T mod_div(T lhs, T rhs, T mod) {
	assert(mod > 0);
	assert(rhs > 0);
	lhs %= mod, rhs %= mod;
	return mod_mul(lhs, mod_inv(rhs, mod), mod);
}

}

namespace mlib {

class FFT {
public:
	static std::vector<std::complex<double>> fft(std::vector<std::complex<double>> &f);
	static std::vector<std::complex<double>> ifft(std::vector<std::complex<double>> &f);
};

}

using namespace mlib;

std::vector<std::complex<double>> FFT::fft(std::vector<std::complex<double>> &f) {
	u32 n = f.size();
	if(n <= 1) return f;
 
	std::vector<std::complex<double>> res = f;
	
	// Bit rev
	u32 s = 0;
	while((u32)(1 << s) < n) s++;
	for(u32 i = 0; i < n; i++) {
		u32 j = 0;
		for(u32 k = 0; k < s; k++) j |= (i >> k & 1) << (s - 1 - k);
		if(i < j) std::swap(res[i], res[j]);
	}
 
	double ang = -1.0 * pi / (2.0 * n);
	std::complex<double> zeta1, zeta3, x0, x1, x3;
	u32 m, mm, rev, k1, k2, k3;
 
	for(m = 4; m <= n; m <<= 1) {
		mm = m >> 2;
		for(u32 j = mm; j >= 1; j >>= 2) {
			for(u32 k = mm - j; k < mm - (j >> 1); k++) {
				k1 = k + mm;
				k2 = k1 + mm;
				k3 = k2 + mm;
				x1 = res[k] - res[k1];
				res[k] += res[k1];
				x3 = res[k3] - res[k2];
				res[k2] += res[k3];
				res[k1] = std::complex<double>(x1.real() - x3.imag(), x1.imag() + x3.real());	//res[k1] = x1 + x3 * ui;
				res[k3] = std::complex<double>(x1.real() + x3.imag(), x1.imag() - x3.real());	// res[k3] = x1 - x3 * ui;
			}
		}
		if(m == n) continue;
 
		rev = n >> 1;
		zeta1 = std::cos(ang * rev);
		for(u32 j = mm; j >= 1; j >>= 2) {
			for(u32 k = m + mm - j; k < m + mm - (j >> 1); k++) {
				k1 = k + mm;
				k2 = k1 + mm;
				k3 = k2 + mm;
				x1 = res[k] - res[k1];
				res[k] += res[k1];
				x3 = res[k3] - res[k2];
				res[k2] += res[k3];
				x0 = std::complex<double>(x1.real() - x3.imag(), x1.imag() + x3.real());	// x0 = x1 + x3 * ui; 
				x0 = std::complex<double>(x0.real() + x0.imag(), x0.imag() - x0.real());	// x0 *= 1.0 - ui;
				res[k1] = zeta1 * x0;
				x0 = std::complex<double>(x1.real() + x3.imag(), x1.imag() - x3.real());	// x0 = x1 - x3 * ui;
				x0 = std::complex<double>(-x0.real() + x0.imag(), -x0.imag() - x0.real());	// x0 *= -1.0 - ui;
				res[k3] = zeta1 * x0;
			}
		}
 
		for(u32 i = 2 * m; i < n; i += m) {
			for(u32 j = n >> 1; j > (rev ^= j); j >>= 1);
			zeta1 = std::polar(1.0, ang * rev);
			zeta3 = std::polar(1.0, 3.0 * ang * rev);
			for(u32 j = mm; j >= 1; j >>= 2) {
				for(u32 k = i + mm - j; k < i + mm - (j >> 1); k++) {
					k1 = k + mm;
					k2 = k1 + mm;
					k3 = k2 + mm;
					x1 = res[k] - res[k1];
					res[k] += res[k1];
					x3 = res[k3] - res[k2];
					res[k2] += res[k3];
					x0 = std::complex<double>(x1.real() - x3.imag(), x1.imag() + x3.real()); // x0 = x1 + x3 * ui;
					res[k1] = zeta1 * x0;
					x0 = std::complex<double>(x1.real() + x3.imag(), x1.imag() - x3.real()); // x0 = x1 - x3 * ui;
					res[k3] = zeta3 * x0;
				}
			}
		}
	}
 
	mm = n >> 1;
	for(u32 j = mm; j >= 1; j >>= 2) {
		for(u32 k = mm - j; k < mm - (j >> 1); k++) {
			k1 = mm + k;
			x0 = res[k] - res[k1];
			res[k] += res[k1];
			res[k1] = x0;
		}
	}
 
	return res;
}
 
std::vector<std::complex<double>> FFT::ifft(std::vector<std::complex<double>> &f) {
	u32 n = f.size();
	if(n <= 1) return f;
 
	std::vector<std::complex<double>> res = f;
	
	// Bit rev
	u32 s = 0;
	while((u32)(1 << s) < n) s++;
	for(u32 i = 0; i < n; i++) {
		u32 j = 0;
		for(u32 k = 0; k < s; k++) j |= (i >> k & 1) << (s - 1 - k);
		if(i < j) std::swap(res[i], res[j]);
	}
 
	double ang = 1.0 * pi / (2.0 * n);
	std::complex<double> zeta1, zeta3, x0, x1, x3;
	u32 m, mm, rev, k1, k2, k3;
 
	for(m = 4; m <= n; m <<= 1) {
		mm = m >> 2;
		for(u32 j = mm; j >= 1; j >>= 2) {
			for(u32 k = mm - j; k < mm - (j >> 1); k++) {
				k1 = k + mm;
				k2 = k1 + mm;
				k3 = k2 + mm;
				x1 = res[k] - res[k1];
				res[k] += res[k1];
				x3 = res[k3] - res[k2];
				res[k2] += res[k3];
				res[k1] = std::complex<double>(x1.real() + x3.imag(), x1.imag() - x3.real());	// res[k1] = x1 - x3 * ui;
				res[k3] = std::complex<double>(x1.real() - x3.imag(), x1.imag() + x3.real());	// res[k3] = x1 + x3 * ui;
			}
		}
		if(m == n) continue;
 
		rev = n >> 1;
		zeta1 = std::cos(ang * rev);
		for(u32 j = mm; j >= 1; j >>= 2) {
			for(u32 k = m + mm - j; k < m + mm - (j >> 1); k++) {
				k1 = k + mm;
				k2 = k1 + mm;
				k3 = k2 + mm;
				x1 = res[k] - res[k1];
				res[k] += res[k1];
				x3 = res[k3] - res[k2];
				res[k2] += res[k3];
				x0 = std::complex<double>(x1.real() + x3.imag(), x1.imag() - x3.real());	// x0 = x1 - x3 * ui;
				x0 = std::complex<double>(x0.real() - x0.imag(), x0.imag() + x0.real());	// x0 *= 1.0 + ui;
				res[k1] = zeta1 * x0;
				x0 = std::complex<double>(x1.real() - x3.imag(), x1.imag() + x3.real());	// x0 = x1 + x3 * ui;
				x0 = std::complex<double>(-x0.real() - x0.imag(), x0.real() + x0.imag());	// x0 *= -1.0 + ui;
				res[k3] = zeta1 * x0;
			}
		}
 
		for(u32 i = 2 * m; i < n; i += m) {
			for(u32 j = n >> 1; j > (rev ^= j); j >>= 1);
			zeta1 = std::polar(1.0, ang * rev);
			zeta3 = std::polar(1.0, 3.0 * ang * rev);
			for(u32 j = mm; j >= 1; j >>= 2) {
				for(u32 k = i + mm - j; k < i + mm - (j >> 1); k++) {
					k1 = k + mm;
					k2 = k1 + mm;
					k3 = k2 + mm;
					x1 = res[k] - res[k1];
					res[k] += res[k1];
					x3 = res[k3] - res[k2];
					res[k2] += res[k3];
					x0 = std::complex<double>(x1.real() + x3.imag(), x1.imag() - x3.real());	// x0 = x1 - x3 * ui;
					res[k1] = zeta1 * x0;
					x0 = std::complex<double>(x1.real() - x3.imag(), x1.imag() + x3.real());	// x0 = x1 + x3 * ui;
					res[k3] = zeta3 * x0;
				}
			}
		}
	}
 
	mm = n >> 1;
	for(u32 j = mm; j >= 1; j >>= 2) {
		for(u32 k = mm - j; k < mm - (j >> 1); k++) {
			k1 = k + mm;
			x0 = res[k] - res[k1];
			res[k] += res[k1];
			res[k1] = x0;
		}
	}
 
	return res;
}

namespace mlib {

// number of base digits
constexpr u64 MAX_DIGITNUM = 5;
// radix base
constexpr u64 RADIX_BASE = 10;
// radix
constexpr u64 RADIX = std::pow(RADIX_BASE, MAX_DIGITNUM);

// sign
using SIGN = bool;
constexpr bool PLUS = true;
constexpr bool MINUS = false;

using DIGITS = std::vector<u64>;

class Int {
private:
	// 12345678 => digits = {78, 56, 34, 12}
	DIGITS digits;
	SIGN sign;

	void assignment(const Int &x) noexcept;
	void to_zero() noexcept;
	void to_radix_power(u64 n) noexcept;
	void push_up(u64 num) noexcept;	// link as upper digit

	// unsigned operations	
	// absolute peer arithmetics
	// The argument does not have to be non-negative.
	void add(const Int &x) noexcept;	
	SIGN sub(const Int &x) noexcept;	// automatically subtracted from the larger one
	void mul(const Int &x) noexcept;	
	void div_short(const Int &x, bool rem = false) noexcept;
	void div(const Int &x, bool rem = false) noexcept;
	bool is_equal(const Int &x) const noexcept;	
	bool less_than(const Int &x, bool equal = false) const noexcept;	
	bool greater_than(const Int &x, bool equal = false) const noexcept;

	// absolute
	Int abs() const noexcept;	

	// signed operations (arithmetic)
	void a_add(const Int &x) noexcept;	
	void a_sub(const Int &x) noexcept; 
	void a_mul(const Int &x) noexcept;	
	void a_div(const Int &x) noexcept;
	void a_mod(const Int &x) noexcept;
	bool a_is_equal(const Int &x) const noexcept;	
	bool a_less_than(const Int &x, bool equal = false) const noexcept;	
	bool a_greater_than(const Int &x, bool equal = false) const noexcept;

public:
	// default constructor
	Int();

	// integer constructor
	template<typename T, std::enable_if_t<std::is_integral<T>::value, std::nullptr_t> = nullptr>
	Int(T n): Int(std::to_string(n)) {
	}

	Int(const DIGITS &digits);
	
	// std::string constructor
	Int(const std::string &x);

	size_t size() const noexcept;
	u64 at(i64 n) const noexcept;
	u64 most_sig() const noexcept;

	void resize(size_t size) noexcept;
	void normalize() noexcept;	// remove meaningless 0 from the top digit
	void align() noexcept;		// adjusted so that each digit does not exceed RADIX

	/* operators */
	// io
	friend std::ostream &operator <<(std::ostream &os, Int &x) noexcept {
		if(x.sign == MINUS && x != 0) os << "-";
		os << x.digits.at(x.size() - 1);
		for(size_t i = x.size() - 1; i --> 0;) os << std::setw(MAX_DIGITNUM) << std::setfill('0') << x.digits.at(i);
		return os;
	}

	friend std::ostream &operator <<(std::ostream &os, Int &&x) noexcept {
		if(x.sign == MINUS && x != 0) os << "-";
		os << x.digits.at(x.size() - 1);
		for(size_t i = x.size() - 1; i --> 0;) os << std::setw(MAX_DIGITNUM) << std::setfill('0') << x.digits.at(i);
		return os;
	}

	friend std::istream &operator >>(std::istream &is, Int &x) noexcept {
		std::string s;
		is >> s;
		x = Int(s);
		return is;
	}

	// arithmetic
	Int &operator +=(const Int &x) noexcept {
		a_add(x);
		return *this;
	}

	Int &operator -=(const Int &x) noexcept {
		a_sub(x);
		return *this;
	}

	Int &operator *=(const Int &x) noexcept {
		a_mul(x);
		return *this;
	}

	Int &operator /=(const Int &x) noexcept {
		a_div(x);
		return *this;
	}

	Int &operator %=(const Int &x) noexcept {
		a_mod(x);
		return *this;
	}

	friend Int operator +(Int lhs, const Int &rhs) noexcept {
		lhs += rhs;
		return lhs;
	}

	friend Int operator -(Int lhs, const Int &rhs) noexcept {
		lhs -= rhs;
		return lhs;
	}

	friend Int operator *(Int lhs, const Int &rhs) noexcept {
		lhs *= rhs;
		return lhs;
	}

	friend Int operator /(Int lhs, const Int &rhs) noexcept {
		lhs /= rhs;
		return lhs;
	}

	friend Int operator %(Int lhs, const Int &rhs) noexcept {
		lhs %= rhs;
		return lhs;
	}

	// unary
	Int &operator ++() noexcept {
		operator +=(1);
		return *this;
	}

	Int operator ++(int) noexcept {
		Int x = *this;
		operator +=(1);
		return x; 
	}

	Int &operator --() noexcept {
		operator -=(1);
		return *this;
	}

	Int operator --(int) noexcept {
		Int x = *this;
		operator -=(1);
		return x;
	}

	Int operator +() const noexcept {
		return *this;
	}

	Int operator -() const noexcept {
		Int x = *this;
		x.sign = !x.sign;
		return x;
	}

	// comparison
	friend bool operator ==(const Int &lhs, const Int &rhs) noexcept {
		return lhs.a_is_equal(rhs);
	}

	friend bool operator !=(const Int &lhs, const Int &rhs) noexcept {
		return !lhs.a_is_equal(rhs);
	}

	friend bool operator <(const Int &lhs, const Int &rhs) noexcept {
		return lhs.a_less_than(rhs);
	}

	friend bool operator <=(const Int &lhs, const Int &rhs) noexcept {
		return lhs.a_less_than(rhs, true);
	}

	friend bool operator >(const Int &lhs, const Int &rhs) noexcept {
		return lhs.a_greater_than(rhs);
	}

	friend bool operator >=(const Int &lhs, const Int &rhs) noexcept {
		return lhs.a_greater_than(rhs, true);
	}
};

inline Int Int::abs() const noexcept {
	if(this->sign == PLUS) return *this;

	Int res = *this;
	res.sign = PLUS;
	return res;
}

inline void Int::resize(size_t size) noexcept {
	this->digits.resize(size, 0);
}

inline void Int::normalize() noexcept {
	while(size() > 0) {
		if(this->digits.back() == 0) this->digits.pop_back();
		else break;
	}

	if(size() == 0) to_zero();
}

inline u64 Int::at(i64 n) const noexcept {
	return this->digits.at(n);
}

inline void Int::assignment(const Int &x) noexcept {
	this->digits = x.digits;
	this->sign = x.sign;
}

inline void Int::to_zero() noexcept {
	std::vector<u64>().swap(this->digits);
	push_up(0);
}

inline void Int::to_radix_power(u64 n) noexcept {
	if(n == 0) {
		assignment(1);
		return;
	}

	to_zero();
	this->digits.reserve(n);
	for(u64 i = 0; i < n -1; i++) push_up(0);
	push_up(1);
}

inline void Int::push_up(u64 num) noexcept {
	this->digits.emplace_back(num);
}

inline size_t Int::size() const noexcept {
	return this->digits.size();
}

inline u64 Int::most_sig() const noexcept {
	return this->digits.back();
}

inline void Int::align() noexcept {
	for(size_t i = 0; i < this->digits.size() - 1; i++) {
		if(this->digits.at(i) >= RADIX) {
			this->digits.at(i + 1) += this->digits.at(i) / RADIX;
			this->digits.at(i) %= RADIX;
		}
	}

	if(this->digits.back() >= RADIX) {
		this->push_up(this->digits.back() / RADIX);
		this->digits.at(this->digits.size() - 2) %= RADIX;
	}
}

}

using namespace mlib;

Int::Int(): Int("0") {
}

Int::Int(const DIGITS &digits) {
	this->sign = PLUS;
	this->digits = digits;
}

Int::Int(const std::string &x) {
	std::string num = x;

	this->sign = PLUS;
	if(num.size() != 0) {
		if(num.front() == '+') num.erase(num.begin());
		else if(num.front() == '-') {
			num.erase(num.begin());
			this->sign = MINUS;
		}
	}

	if(num.size() == 0) num = "0";

	while(num.size() > MAX_DIGITNUM) {
		push_up(std::stoull(num.substr(num.size() - MAX_DIGITNUM, MAX_DIGITNUM)));
		for(size_t i = 0; i < MAX_DIGITNUM; i++) num.pop_back();
	}
	push_up(std::stoull(num));
}

void Int::add(const Int &x) noexcept {
	if(x == 0) return;

	if(size() < x.size()) resize(x.size());

	u64 buf, carry = 0;
	u64 rd;
	for(size_t i = 0; i < size(); i++) {
		rd = (i < x.size())? x.digits.at(i) : 0;
		buf = this->digits.at(i) + rd + carry;
		carry = buf / RADIX;
		this->digits.at(i) = buf % RADIX;
	}

	if(carry != 0) push_up(carry);
	normalize();
}

SIGN Int::sub(const Int &x) noexcept {
	if(is_equal(x)) {
		to_zero();
		return PLUS;
	}

	u32 borrow = 0;
	u64 rd;

	// *this - x
	if(greater_than(x)) {
		for(size_t i = 0; i < size(); i++) {
			rd = (i < x.size())? x.digits.at(i) : 0;
			if(this->digits.at(i) >= rd + borrow) {
				this->digits.at(i) = this->digits.at(i) - borrow - rd;
				borrow = 0;
			} else {
				this->digits.at(i) = this->digits.at(i) + RADIX - borrow - rd;
				borrow = 1;
			}
		}
		normalize();
		return PLUS;
	}

	// x - *this
	if(size() < x.size()) resize(x.size());
	for(size_t i = 0; i < x.size(); i ++) {
		rd = (i < size())? this->digits.at(i) : 0;
		if(x.digits.at(i) >= rd + borrow) {
			this->digits.at(i) = x.digits.at(i) - borrow - rd;
			borrow = 0;
		} else {
			this->digits.at(i) = x.digits.at(i) + RADIX - borrow - rd;
			borrow = 1;
		}
	}
	normalize();
	return MINUS;
}

void Int::mul(const Int &x) noexcept {
	if(x == 0) {
		to_zero();
		return;
	}

	// adjust the number of digits to 2^n
	size_t n = 1;
	while (n < size() + x.size()) n <<= 1;

	std::vector<std::complex<double>> a(n);
	for(size_t i = 0; i < n; i++) {
		a.at(i).real((i >= size()) ? 0 : this->digits.at(i));
		a.at(i).imag((i >= x.size()) ? 0 : x.digits.at(i));
	}

	a = FFT::fft(a);
	std::vector<std::complex<double>> res(n);
	for(size_t i = 0; i < n; i++) {
		res.at(i) = (a.at(i) + std::conj(a.at((n - i) % n))) * 0.5 * (a.at(i) - std::conj(a.at((n - i) % n))) * std::complex<double>(0, -0.5);
	}
		
	res = FFT::ifft(res);
	resize(n);
	for(size_t i = 0; i < n; i++) this->digits.at(i) = std::round(res.at(i).real() / static_cast<double>(n));

	normalize();
	align();
}

// short divisor
// x is absolute value
void Int::div_short(const Int &x, bool rem) noexcept {
	if(x < RADIX) {
		u64 dividend = this->digits.back();
		for(u64 i = this->digits.size(); i --> 0;) {
			this->digits.at(i) = dividend / x.at(0);
			if(i == 0) {
				if(rem) assignment(dividend % x.at(0));
				else normalize();
				return;
			}
			dividend = (dividend % x.at(0)) * RADIX + this->digits.at(i - 1);
		}
	}
	
	u64 radix_per_10 = RADIX / 10;
	Int dividend(this->digits.back() / radix_per_10);

	u64 qq, xx;
	Int rr;
	for(u64 i = this->digits.size(); i --> 0;) {
		xx = 0;
		for(u64 j = MAX_DIGITNUM; j --> 0;) {
			qq = 0;
			rr = dividend;
			while(rr >= x) {
				rr -= x;
				qq++;
			}
			xx += qq * std::pow(RADIX_BASE, j);
			
			if(i == 0 && j == 0) {
				this->digits.at(i) = xx;
				if(!rem) {
					this->digits.at(i) = xx;
					normalize();
				} else {
					assignment(rr);
				}
				return;
			}
			
			if(qq != 0) dividend = rr;
			dividend *= RADIX_BASE;
			if(j == 0) dividend += at(i - 1) / radix_per_10;
			else dividend += at(i) / static_cast<u64>(std::pow(RADIX_BASE, j - 1)) % 10;
		}
		this->digits.at(i) = xx;
	}
}


// Knuth D
void Int::div(const Int &x, bool rem) noexcept {
	assert(x != 0);
	if(less_than(x)) {
		if(!rem) to_zero();
		return;
	}

	if(is_equal(x)) {
		if(!rem) assignment(1);
		else assignment(0);
		return;
	}

	i64 s, t;
	s = size() - 1, t = x.size() - 1;
	if(t <= 1) {
		div_short(x.abs(), rem);
		return;
	}

	u64 k = RADIX / (x.digits.back() + 1);
	Int a = abs() * k;
	Int b = x.abs() * k;

	DIGITS q(a.size(), 0);

	i64 u;
	u64 qp, rp;
	u64 fi = 1;
	while(true) {
		s = a.size() - 1 + fi;
		a.push_up(0);
		u = (a.at(s) < b.at(t))? s - t - 1 : s - t;
		qp = (a.at(s) * RADIX + a.at(s - 1)) / b.at(t);
		rp = (a.at(s) * RADIX + a.at(s - 1)) % b.at(t);
		while(rp < RADIX && qp * b.at(t - 1) > RADIX * rp + a.at(s - 2)) {
			qp--;
			rp += b.at(t);
		}

		q.at(u) = qp;
		Int xx;
		xx.to_radix_power(u);
		xx *= b;
		a -= qp * xx;
		if(a < 0) {
			--q.at(u);
			a += xx;
		}

		fi = 0;
		if(u == 0) {
			if(!rem) {
				this->digits = q;
				normalize();
			} else {
				assignment(a / k);
			}
			return;
		}
	}
}

bool Int::is_equal(const Int &x) const noexcept {
	if(size() != x.size()) return false;

	for(size_t i = 0; i < size(); i++) {
		if(this->digits.at(i) != x.digits.at(i)) return false;
	}

	return true;
}

bool Int::less_than(const Int &x, bool equal) const noexcept {
	if(size() < x.size()) return true;
	if(size() > x.size()) return false;

	for(size_t i = size(); i --> 0;) {
		if(this->digits.at(i) < x.digits.at(i)) return true;
		if(this->digits.at(i) > x.digits.at(i)) return false;
	}

	// equal
	return equal;
}

bool Int::greater_than(const Int &x, bool equal) const noexcept {
	return !less_than(x, !equal);
}

void Int::a_add(const Int &x) noexcept {
	if(this->sign == x.sign) {
		add(x);
	} else {
		if(sub(x) == MINUS) this->sign = x.sign;
	}
}

void Int::a_sub(const Int &x) noexcept {
	a_add(-x);
}

void Int::a_mul(const Int &x) noexcept {
	mul(x);
	this->sign = !(this->sign ^ x.sign);
}

void Int::a_div(const Int &x) noexcept {
	div(x);
	this->sign = !(this->sign ^ x.sign);
}

void Int::a_mod(const Int &x) noexcept {
	div(x, true);
}

bool Int::a_is_equal(const Int &x) const noexcept {
	if(is_equal(0) && x.is_equal(0)) return true;
	if(this->sign != x.sign) return false;
	return is_equal(x);
}

bool Int::a_less_than(const Int &x, bool equal) const noexcept {
	// +, - zero
	if(is_equal(0)) {
		if(x.sign == PLUS) return less_than(x, equal);
		else return greater_than(x, equal);
	}
	
	// different sign
	if(!x.is_equal(0)) {
		if(this->sign == MINUS && x.sign == PLUS) return true;
		if(this->sign == PLUS && x.sign == MINUS) return false;
	}

	// same sign
	if(this->sign == PLUS) return less_than(x, equal);
	return greater_than(x, equal);
}

bool Int::a_greater_than(const Int &x, bool equal) const noexcept {
	return !a_less_than(x, !equal);
}

namespace mlib {
namespace random {

// return random Int value in [0, x)
Int randint(const Int &x);

}
}

mlib::Int mlib::random::randint(const mlib::Int &x) {
	assert(x.size() != 0);

	std::random_device seed;
	std::mt19937 mt(seed());
	std::uniform_int_distribution<> rand(0, mlib::RADIX - 1);
	mlib::DIGITS dig(x.size());
	
	for(u64 i = 0; i < dig.size() - 1; i++) dig[i] = rand(mt);	

	std::uniform_int_distribution<> rand2(0, x.most_sig() - 1);
	dig[dig.size() - 1] = rand2(mt);

	mlib::Int ret(dig);
	ret.normalize();
	return ret;
}

namespace mlib {

bool is_prime(const Int &n);

}

// Miller-Rabin algorithm
bool mlib::is_prime(const mlib::Int &n) {
	if(n == 2) return true;
	if(n == 1) return false;
	if(n % 2 == 0) return false;

	mlib::Int d = (n - 1) / 2;
	while(d % 2 == 0) d /= 2;

	for(u64 i = 0; i < 100; i++) {
		auto a = mlib::random::randint(n - 2) + 1;
		auto x = mlib::mod_pow(a, d, n);
		auto t = d;

		while(t != n - 1 && x != 1 && x != n - 1) {
			x = mlib::mod_pow(x, mlib::Int(2), n);
			t *= 2;
		}

		if(x != n - 1 && x % 2 == 0) return false;
	}

	return true;
}

int main() {
	u64 n;
	std::cin >> n;
	for(u64 i = 0; i < n; i++) {
		mlib::Int x;
		std::cin >> x;
		std::cout << x << " ";
		if(mlib::is_prime(x)) std::cout << 1 << std::endl;
		else std::cout << 0 << std::endl;
	}
	return 0;
}