結果
| 問題 |
No.303 割れません
|
| ユーザー |
Min_25
|
| 提出日時 | 2015-11-22 18:44:18 |
| 言語 | C++11(廃止可能性あり) (gcc 13.3.0) |
| 結果 |
AC
|
| 実行時間 | 323 ms / 10,000 ms |
| コード長 | 49,902 bytes |
| コンパイル時間 | 2,057 ms |
| コンパイル使用メモリ | 107,016 KB |
| 実行使用メモリ | 8,156 KB |
| 最終ジャッジ日時 | 2024-09-13 17:11:31 |
| 合計ジャッジ時間 | 5,823 ms |
|
ジャッジサーバーID (参考情報) |
judge5 / judge2 |
(要ログイン)
| ファイルパターン | 結果 |
|---|---|
| sample | AC * 3 |
| other | AC * 14 |
ソースコード
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <cmath>
#include <cstdint>
#include <cassert>
#include <iostream>
#include <vector>
#include <algorithm>
#include <string>
// M(n): O(n^1.47)
#define HAVE_UINT128
using uint32 = uint32_t;
using int64 = int64_t;
using uint64 = uint64_t;
using uint32 = uint32_t;
using int64 = int64_t;
using uint64 = uint64_t;
#ifdef HAVE_UINT128
using uint128 = __uint128_t;
using word_t = uint64;
using dword_t = uint128;
#else
using word_t = uint32;
using dword_t = uint64;
#endif
namespace bits {
inline constexpr uint32 pop_count(uint32 n) {
return __builtin_popcount(n);
}
inline constexpr uint64 pop_count(uint64 n) {
return __builtin_popcountll(n);
}
inline constexpr uint32 ctz(uint32 n) {
return __builtin_ctz(n);
}
inline constexpr uint64 ctz(uint64 n) {
return __builtin_ctzll(n);
}
inline constexpr uint32 clz(uint32 n) {
return __builtin_clz(n);
}
inline constexpr uint64 clz(uint64 n) {
return __builtin_clzll(n);
}
template <typename T>
inline constexpr T ilog2(T n) {
return n == 0 ? 0 : (sizeof(n) * 8 - 1) - clz(n);
}
#ifdef HAVE_UINT128
inline uint64 floor_div_small(uint128 a, uint64 b) {
uint64 q, r;
__asm__ (
"divq\t%4"
: "=a"(q), "=d"(r)
: "0"(uint64(a)), "1"(uint64(a >> 64)), "rm"(b)
);
return q;
}
#else
inline uint32 floor_div_small(uint64 a, uint32 b) {
uint32 q, r;
__asm__ (
"divl\t%4"
: "=a"(q), "=d"(r)
: "0"(uint32(a)), "1"(uint32(a >> 32)), "rm"(b)
);
return q;
}
#endif
inline word_t mul_inv(word_t n) {
word_t x = n;
while (1) {
word_t t = n * x;
if (t == 1) {
break;
}
x *= 2 - t;
}
return x;
}
}
struct FastDiv {
FastDiv() {}
FastDiv(word_t n) {
if (n == 1) {
shamt = 0;
magic_num = 1;
} else {
shamt = (8 * sizeof(dword_t)) - 1 - bits::clz(n - 1);
magic_num = bits::floor_div_small((dword_t(1) << shamt) + n - 1, n);
}
mod = n;
}
friend word_t operator / (word_t n, FastDiv d) {
return (dword_t(n) * d.magic_num) >> d.shamt;
}
friend word_t& operator /= (word_t& n, FastDiv d) {
return n = n / d;
}
friend word_t operator % (word_t n, FastDiv d) {
return n - n / d * d.mod;
}
friend word_t& operator %= (word_t& n, FastDiv d) {
return n = n % d;
}
word_t magic_num;
word_t shamt;
word_t mod;
};
class FastDiv21 {
public:
FastDiv21(word_t n) {
shamt_ = bits::clz(n);
d_ = n << shamt_;
v_ = reciprocal(d_);
};
FastDiv21() {}
word_t divisor() const {
return d_ >> shamt_;
}
// divl in x86.
static void divmod(
dword_t u, FastDiv21 fd, word_t& ret_q, word_t& ret_r) {
u <<= fd.shamt_;
word_t u_hi = u >> (8 * sizeof(word_t));
word_t u_lo = word_t(u);
dword_t q = dword_t(u_hi) * fd.v_ + u;
word_t q_hi = (q >> (8 * sizeof(word_t))) + 1;
word_t r = u_lo - q_hi * fd.d_;
if (r > word_t(q)) {
q_hi -= 1;
r += fd.d_;
}
if (r >= fd.d_) {
q_hi += 1;
r -= fd.d_;
}
ret_q = q_hi;
ret_r = r >> fd.shamt_;
}
friend word_t operator / (dword_t n, FastDiv21 fd) {
return FastDiv21::div(n, fd);
}
friend dword_t& operator /= (dword_t& n, FastDiv21 fd) {
return n = n / fd;
}
friend word_t operator % (dword_t n, FastDiv21 fd) {
return FastDiv21::mod(n, fd);
}
friend dword_t& operator %= (dword_t& n, FastDiv21 fd) {
return n = n % fd;
}
private:
static word_t div(dword_t u, FastDiv21 fd) {
word_t q, dummy;
divmod(u, fd, q, dummy);
return q;
}
static word_t mod(dword_t u, FastDiv21 fd) {
word_t dummy, r;
divmod(u, fd, dummy, r);
return r;
}
word_t reciprocal(word_t n) const {
return bits::floor_div_small(~(dword_t(n) << (sizeof(word_t) * 8)), n);
}
word_t d_;
word_t shamt_;
word_t v_;
};
class BigInt {
public:
enum {
W_BITS = sizeof(word_t) * 8,
#ifdef HAVE_UINT128
MUL_TOOM22_THRESHOLD = 24, // words
MUL_TOOM33_THRESHOLD = 120, // >= 4 words
SQ_TOOM22_THRESHOLD = 32, // words
SQ_TOOM33_THRESHOLD = 160, // >= 4 words
INT_THRESHOLD = 64, // words
STR_THRESHOLD = 1 << 12, // bits
DIVMOD_THRESHOLD = 1 << 13, // bits
#else
MUL_TOOM22_THRESHOLD = 16, // words
MUL_TOOM33_THRESHOLD = 80, // >= 4 words
SQ_TOOM22_THRESHOLD = 28, // words
SQ_TOOM33_THRESHOLD = 140, // >= 4 words
INT_THRESHOLD = 64, // words
STR_THRESHOLD = 1 << 9, // bits
DIVMOD_THRESHOLD = 1 << 12, // bits
#endif
};
static_assert(MUL_TOOM33_THRESHOLD >= 4, "mul_toom33");
static_assert(SQ_TOOM33_THRESHOLD >= 4, "sq_toom33");
public:
struct BarrettReduction {
~BarrettReduction() {
delete divisor;
delete reciprocal;
}
BarrettReduction() : divisor(nullptr), reciprocal(nullptr), shamt(0) {}
BarrettReduction(BigInt* d, BigInt* r, word_t s) :
divisor(d), reciprocal(r), shamt(s) {}
static void divmod(const BigInt& n, const BarrettReduction& d, BigInt& q, BigInt& r) {
if (d.reciprocal == nullptr) {
return BigInt::divmod(n, *d.divisor, q, r);
}
q = ((n >> (d.shamt - 1)) * *d.reciprocal) >> (d.shamt + 1);
r = n - q * *d.divisor;
while (r >= *d.divisor) {
r -= *d.divisor;
q += 1;
}
}
BigInt* divisor;
BigInt* reciprocal;
word_t shamt;
};
struct BaseData {
BaseData(word_t base, bool lower) : base_s(base), lower(lower) {
fd = FastDiv(base);
word_t lim = (word_t(1) << (W_BITS - 1)) - 1;
base_l = 1;
e = 0;
while (lim >= base_s) {
base_l *= base_s;
lim /= fd;
e += 1;
}
fd21 = FastDiv21(base_l);
}
FastDiv fd;
FastDiv21 fd21;
word_t base_s;
word_t base_l;
word_t e;
bool lower;
private:
BaseData() {}
};
using container_t = std::vector<word_t>;
private:
struct _bigint_t {
_bigint_t(int s, const word_t* d, word_t sz) : sign(s), data(const_cast<word_t*>(d)), size(sz) {}
_bigint_t(word_t* d) : sign(0), data(d), size(0) {}
_bigint_t() {}
int sign;
word_t* data;
word_t size;
};
public:
BigInt() : sign_(0) {
n_ = container_t(1, 0);
}
BigInt(int n) : sign_(n < 0) {
n_ = container_t(1, uint32(std::abs(n))); // fixed
}
BigInt(uint32 n) : sign_(0) {
n_ = container_t(1, n);
}
#ifdef HAVE_UINT128
BigInt(int64 n) : sign_(n < 0) {
n_ = container_t(1, std::abs(n));
}
BigInt(uint64 n) : sign_(0) {
n_ = container_t(1, n);
}
inline operator uint64() const {
uint64 ret = (*this)[0];
if (is_neg()) {
ret = -ret;
}
return ret;
}
#else
BigInt(int64 n) : sign_(n < 0) {
n_ = container_t(2);
(*this)[0] = n;
(*this)[1] = std::abs(int(n >> W_BITS));
}
BigInt(uint64 n) : sign_(0) {
n_ = container_t(2);
(*this)[0] = n;
(*this)[1] = n >> W_BITS;
}
inline operator uint64() const {
uint64 ret = (size() >= 2 ? ((uint64((*this)[1]) << 32) | (*this)[0]) : (*this)[0]);
if (is_neg()) {
ret = -ret;
}
return ret;
}
#endif
BigInt(const char* cstr, word_t base=10) {
_int(cstr, std::strlen(cstr), base, *this);
}
BigInt(const std::string& str, word_t base=10) {
_int(str.c_str(), str.size(), base, *this);
}
BigInt(word_t size, int n) : sign_(n < 0) {
assert(size > 0);
n_ = container_t(size, 0);
n_[0] = n;
}
BigInt(const BigInt& n) {
this->sign_ = n.sign_;
this->n_ = n.n_;
}
// container
inline word_t size() const {
return n_.size();
}
void push_back(word_t w) {
n_.push_back(w);
}
void resize(word_t size) {
if (size == 0) {
set_zero();
} else {
n_.resize(size, 0);
}
}
inline word_t operator [] (word_t idx) const {
return n_[idx];
}
inline word_t& operator [] (word_t idx) {
return n_[idx];
}
inline word_t* data() {
return n_.data();
}
inline const word_t* data() const {
return n_.data();
}
// basic
void change_sign() {
if (!is_zero()) {
sign_ ^= 1;
}
}
void set_sign(int s) {
if (!is_zero()) {
sign_ = s;
}
}
void clear_sign() {
sign_ = 0;
}
int sign() const {
return sign_;
}
bool is_neg() const {
return sign_ == 1;
}
bool is_zero() const {
return size() == 1 && (*this)[0] == 0;
}
bool is_one() const {
// unsigned
return size() == 1 && (*this)[0] == 1;
}
BigInt abs() const {
BigInt ret = BigInt(*this);
ret.clear_sign();
return ret;
}
void set_zero() {
resize(1);
clear_sign();
(*this)[0] = 0;
}
void normalize() {
for (int i = size() - 1; i >= 0; --i) {
if ((*this)[i]) {
return resize(i + 1);
}
}
return set_zero();
}
word_t bit_length() const {
if (is_zero()) {
return 0;
}
word_t len = size();
return W_BITS * len - bits::clz((*this)[len - 1]);
}
static word_t bits_to_words(word_t bits) {
return (bits + W_BITS - 1) / W_BITS;
}
word_t pop_count() const {
word_t ret = 0;
for (word_t i = 0; i < this->size(); ++i) {
ret += bits::pop_count((*this)[i]);
}
return ret;
}
static BigInt mask(word_t bit_length) {
if (bit_length == 0) {
return BigInt();
}
BigInt ret = BigInt();
word_t ret_size = 1 + (bit_length - 1) / W_BITS;
ret.resize(ret_size);
for (word_t i = 0; i < ret_size; ++i) {
ret[i] = ~word_t(0);
}
ret[ret_size - 1] &= ~word_t(0) >> ((W_BITS - bit_length) & (W_BITS - 1));
return ret;
}
// string
std::string str(word_t base=10, bool lower=true) const {
if (is_zero()) {
return "0";
}
BaseData bdata = BaseData(base, lower);
return _fast_str(bdata);
}
// relational
bool operator == (const BigInt& rhs) const {
if (this->is_neg() ^ rhs.is_neg()) {
return false;
}
if (size() != rhs.size()) {
return false;
}
for (word_t i = 0; i < rhs.size(); ++i) {
if ((*this)[i] != rhs[i]) {
return false;
}
}
return true;
}
bool operator != (const BigInt& rhs) const {
return !(*this == rhs);
}
bool operator < (const BigInt& rhs) const {
if (this->is_neg()) {
if (!rhs.is_neg()) {
return true;
} else {
return _ucomp(rhs, *this, false);
}
} else {
if (rhs.is_neg()) {
return false;
} else {
return _ucomp(*this, rhs, false);
}
}
}
/* slow */
bool operator < (const int rhs) const {
return *this < BigInt(rhs);
}
#ifndef HAVE_UINT128
/* tmp */
bool operator < (const uint64 rhs) const {
return *this < BigInt(rhs);
}
#endif
bool operator < (const word_t rhs) const {
if (is_neg()) {
return true;
}
return !(size() >= 2) && (*this)[0] < rhs;
}
bool operator <= (const BigInt& rhs) const {
if (this->is_neg()) {
if (!rhs.is_neg()) {
return true;
} else {
return _ucomp(rhs, *this, true);
}
} else {
if (rhs.is_neg()) {
return false;
} else {
return _ucomp(*this, rhs, true);
}
}
}
bool operator <= (const word_t rhs) const {
if (is_neg()) {
return true;
}
return !(size() >= 2) && (*this)[0] <= rhs;
}
bool operator >= (const BigInt& rhs) const {
return rhs <= *this;
}
bool operator >= (const word_t rhs) const {
return !(*this < rhs);
}
bool operator > (const BigInt& rhs) const {
return rhs < *this;
}
bool operator > (const word_t rhs) const {
return !(*this <= rhs);
}
// unary operator
bool operator ! () const {
return is_zero();
}
BigInt operator + () const {
return BigInt(*this);
}
BigInt operator - () const {
BigInt ret = BigInt(*this);
ret.change_sign();
return ret;
}
// arith
BigInt operator + (const BigInt& rhs) const {
BigInt res;
_add(*this, rhs, false, res);
return res;
}
/* slow */
BigInt operator + (const int rhs) const {
BigInt res;
_add(*this, BigInt(rhs), false, res);
return res;
}
BigInt& operator += (const BigInt& rhs) {
_add(*this, rhs, false, *this);
return *this;
}
BigInt operator - (const BigInt& rhs) const {
BigInt res;
_add(*this, rhs, true, res);
return res;
}
/* slow */
BigInt operator - (const int rhs) const {
BigInt res;
_add(*this, BigInt(rhs), true, res);
return res;
}
BigInt& operator -= (const BigInt& rhs) {
_add(*this, rhs, true, *this);
return *this;
}
BigInt operator * (const BigInt& rhs) const {
BigInt ret;
_mul(*this, rhs, ret);
return ret;
}
BigInt operator * (const word_t rhs) const {
BigInt ret;
_mul1(*this, rhs, ret);
return ret;
}
BigInt& operator *= (const BigInt& rhs) {
return *this = *this * rhs;
}
BigInt& operator *= (const word_t rhs) {
_mul1(*this, rhs, *this);
return *this;
}
BigInt operator / (const BigInt& rhs) const {
BigInt q, r;
divmod(*this, rhs, q, r);
return q;
}
BigInt operator / (const word_t rhs) const {
BigInt ret;
word_t dummy;
divmod_n1(*this, rhs, ret, dummy);
return ret;
}
BigInt operator % (const BigInt& rhs) const {
BigInt q, r;
divmod(*this, rhs, q, r);
return r;
}
word_t operator % (const word_t rhs) const {
BigInt dummy;
word_t ret;
divmod_n1(*this, rhs, dummy, ret);
return ret;
}
static void divmod(const BigInt& a, const BigInt& b, BigInt& q, BigInt& r) {
/* sign */
_fast_udivmod(a, b, q, r);
}
static void divmod_n1(const BigInt& n, const word_t d, BigInt& qs, word_t& r) {
/* sign */
auto fd = FastDiv21(d);
return divmod_n1(n, fd, qs, r);
}
static void divmod_n1(const BigInt& n, const FastDiv21 fd, BigInt& qs, word_t& r) {
const word_t d = fd.divisor();
word_t size = n.size();
r = n[size - 1];
word_t q = 0;
if (r >= d) {
FastDiv21::divmod(r, fd, q, r);
qs.resize(size);
qs[size - 1] = q;
} else {
qs.resize(size - 1);
}
for (int i = size - 2; i >= 0; --i) {
FastDiv21::divmod((dword_t(r) << W_BITS) | n[i], fd, q, r);
qs[i] = q;
}
}
// logical [unsigned]
BigInt operator & (const BigInt& rhs) const {
if (size() > rhs.size()) {
return rhs & *this;
}
BigInt ret = BigInt(*this);
_land(ret.data(), rhs.data(), ret.size());
ret.normalize();
return ret;
}
BigInt operator & (const int rhs) const {
return BigInt((*this)[0] & rhs);
}
BigInt& operator &= (const BigInt& rhs) {
if (size() > rhs.size()) {
resize(rhs.size());
}
_land(this->data(), rhs.data(), this->size());
this->normalize();
return *this;
}
BigInt operator ^ (const BigInt& rhs) const {
if (size() < rhs.size()) {
return rhs ^ *this;
}
BigInt ret = BigInt(*this);
_lxor(ret.data(), rhs.data(), rhs.size());
ret.normalize();
return ret;
}
BigInt& operator ^= (const BigInt& rhs) {
if (size() < rhs.size()) {
resize(rhs.size());
}
_lxor(this->data(), rhs.data(), this->size());
this->normalize();
return *this;
}
BigInt operator | (const BigInt& rhs) const {
if (size() < rhs.size()) {
return rhs | *this;
}
BigInt ret = BigInt(*this);
_lor(ret.data(), rhs.data(), rhs.size());
return ret;
}
BigInt& operator |= (const BigInt& rhs) {
if (size() < rhs.size()) {
resize(rhs.size());
}
_lor(this->data(), rhs.data(), rhs.size());
return *this;
}
// shift [unsigned]
BigInt operator << (int shamt) const {
if (shamt < 0) {
return *this >> -shamt;
}
if (shamt == 0) {
return BigInt(*this);
}
if (is_zero()) {
return BigInt();
}
BigInt ret = BigInt();
word_t bit_len = bit_length();
ret.resize(bits_to_words(bit_len + shamt));
_lshift(this->data(), this->size(), shamt, ret.data(), ret.size());
return ret;
}
BigInt operator << (word_t shamt) const {
return *this << int(shamt);
}
BigInt& operator <<= (int shamt) {
if (shamt < 0) {
return *this >>= -shamt;
}
if (shamt == 0) {
return *this;
}
if (is_zero()) {
return *this;
}
word_t bit_len = bit_length();
word_t size = this->size();
this->resize(bits_to_words(bit_len + shamt));
_lshift(this->data(), size, shamt, this->data(), this->size());
return *this;
}
BigInt& operator <<= (word_t shamt) {
return *this <<= int(shamt);
}
BigInt operator >> (int shamt) const {
if (shamt < 0) {
return *this << -shamt;
}
if (shamt == 0) {
return BigInt(*this);
}
if (is_zero()) {
return BigInt();
}
word_t bit_len = bit_length();
if (bit_len <= word_t(shamt)) {
return BigInt();
}
BigInt ret = BigInt();
ret.resize(bits_to_words(bit_len - shamt));
_rshift(this->data(), this->size(), shamt, ret.data(), ret.size());
return ret;
}
BigInt operator >> (word_t shamt) const {
return *this >> int(shamt);
}
BigInt& operator >>= (int shamt) {
if (shamt < 0) {
return *this <<= -shamt;
}
if (shamt == 0) {
return *this;
}
if (is_zero()) {
return *this;
}
word_t bit_len = bit_length();
if (bit_len <= word_t(shamt)) {
this->set_zero();
return *this;
}
word_t ret_size = bits_to_words(bit_len - shamt);
_rshift(this->data(), this->size(), shamt, this->data(), ret_size);
this->resize(ret_size);
return *this;
}
BigInt& operator >>= (word_t shamt) {
return *this >>= int(shamt);
}
// subinteger ? strip ?
BigInt subinteger(const word_t bit_beg, word_t bit_end) const {
word_t bit_len = bit_length();
if (bit_end > bit_len) {
bit_end = bit_len;
}
if (bit_beg >= bit_end) {
return BigInt(0);
}
word_t bit_size = bit_end - bit_beg;
word_t ret_size = bits_to_words(bit_size);
word_t in_size = bits_to_words(bit_end);
BigInt ret = BigInt(ret_size, 0);
_rshift(this->data(), in_size, bit_beg, ret.data(), ret_size);
word_t mask = ~word_t(0) >> ((W_BITS - bit_size) & (W_BITS - 1));
ret[ret_size - 1] &= mask;
ret.normalize();
return ret;
}
// functions
BigInt isqrt() const {
if (is_neg()) {
assert(0);
}
BigInt s, r;
_isqrt(*this, s, r);
return s;
}
BigInt sqrt() const {
return isqrt();
}
BigInt pow(word_t e) const {
BigInt ret = BigInt(1);
if (e == 0) {
return ret;
}
word_t mask = word_t(1) << bits::ilog2(e);
while (mask) {
if (mask & e) {
ret *= *this;
}
mask >>= 1;
if (!mask) {
break;
}
ret *= ret;
}
return ret;
}
static BigInt fib(word_t n) {
BigInt a, b;
_fib(n, a, b);
return a;
}
static void fib(word_t n, BigInt& a, BigInt& b) {
_fib(n, a, b);
}
// output
friend std::ostream & operator << (std::ostream& os, const BigInt& n) {
return os << n.str();
}
private:
static BigInt _int_small(word_t size, const word_t* smalls, const word_t base) {
BigInt ret;
for (int i = size - 1; i >= 0; --i) {
word_t c = _umul1_add(ret.data(), ret.size(), base, ret.data(), smalls[i]);
if (c) {
ret.push_back(c);
}
}
return ret;
}
static BigInt _int_rec(word_t size, word_t e, word_t base_l,
const word_t* smalls, const std::vector<BigInt>& large_bases) {
if (size <= INT_THRESHOLD) {
return _int_small(size, smalls, base_l);
} else {
word_t d = word_t(1) << e;
if (size > d) {
return _int_rec(d, e - 1, base_l, smalls , large_bases) +
_int_rec(size - d, e - 1, base_l, smalls + d, large_bases) *
large_bases[e];
} else {
return _int_rec(size, e - 1, base_l, smalls , large_bases);
}
}
}
static void _int(const char* cstr, word_t size, word_t base, BigInt& res) {
static const word_t char_to_word[128] = {
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 0, 0, 0, 0, 0, 0,
0, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24,
25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 0, 0, 0, 0, 0,
0, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24,
25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 0, 0, 0, 0, 0
};
res.set_zero();
if (size == 0) {
return;
}
bool neg = false;
word_t pos = 0;
if (cstr[pos] == '-') {
neg = true;
++pos;
}
word_t lim = ~word_t(0);
word_t base_l = 1;
auto fd = FastDiv(base);
word_t e = 0;
while (lim >= base) {
lim /= fd;
e += 1;
base_l *= base;
}
// cstr => word_size
word_t size_s = (size - pos + e - 1) / e;
auto smalls = std::vector<word_t>(size_s, 0);
for (word_t i = 0; i < size_s; ++i) {
int end = size - i * e;
int beg = std::max(int(pos), end - int(e));
word_t t = 0;
for (int j = beg; j < end; ++j) {
t = t * base + char_to_word[cstr[j] & 0x7f];
}
smalls[i] = t;
}
if (size_s == 1) {
res = BigInt(smalls[0]);
return;
}
word_t large_e = bits::ilog2(size_s - 1);
auto larges = std::vector<BigInt>(large_e + 1);
BigInt large_base = BigInt(base_l);
for (word_t i = 0; i < large_e; ++i) {
larges[i] = large_base;
large_base *= large_base;
}
larges[large_e] = large_base;
res = _int_rec(size_s, large_e, base_l, smalls.data(), larges);
if (neg) {
res.change_sign();
}
}
std::string _ustr(const BaseData& base, word_t zpad=0) const {
static const char* chars = base.lower
? "0123456789abcdefghijklmnopqrstuvwxyz"
: "0123456789ABCDEFGHIJKLMNOPQRSTUVWXYZ";
std::string ret;
BigInt qs = BigInt(*this); qs.clear_sign();
auto fd = base.fd;
word_t r = 0;
while (qs >= base.base_l) {
divmod_n1(qs, base.fd21, qs, r);
for (word_t i = 0; i < base.e; ++i) {
word_t s = r / fd;
ret += chars[r - s * fd.mod];
r = s;
}
}
r = qs[0];
while (r > 0) {
word_t s = r / fd;
ret += chars[r - s * fd.mod];
r = s;
}
if (zpad && ret.size() < zpad) {
ret.append(zpad - ret.size(), '0');
}
std::reverse(ret.begin(), ret.end());
return ret;
}
static void _fast_str_rec(
const BigInt& n, word_t e,
const BaseData& base, const std::vector<BarrettReduction>& large_bases,
std::string& ret, bool zpad=false) {
if (e + 1 <= 9) {
ret += n._ustr(base, zpad ? (1 << (e + 1)) : 0);
} else {
BigInt q, r;
BarrettReduction::divmod(n, large_bases[e], q, r);
if (zpad || !q.is_zero()) {
_fast_str_rec(q, e - 1, base, large_bases, ret, zpad);
_fast_str_rec(r, e - 1, base, large_bases, ret, true);
} else {
_fast_str_rec(r, e - 1, base, large_bases, ret, zpad);
}
}
}
std::string _fast_str(const BaseData& base) const {
std::string ret;
if (is_neg()) {
ret += "-";
}
word_t bit_len = bit_length();
if (bit_len <= STR_THRESHOLD) {
ret += this->_ustr(base);
return ret;
}
auto large_bases = std::vector<BarrettReduction>(50);
BigInt base_pow = BigInt(base.base_s);
word_t e = 0;
while (1) {
large_bases[e].divisor = new BigInt(base_pow);
if (base_pow.bit_length() * 2 - 1 > bit_len) {
break;
}
e += 1;
base_pow *= base_pow;
}
BigInt two = BigInt(1);
word_t prev_k = 0;
for (int i = 0; i < int(e) - 2; ++i) {
word_t k = large_bases[i].divisor->bit_length();
two <<= 2 * (k - prev_k);
large_bases[i].reciprocal = new BigInt(two / *large_bases[i].divisor);
large_bases[i].shamt = k;
prev_k = k;
}
_fast_str_rec(*this, e, base, large_bases, ret);
return ret;
}
// relational
static bool _ucomp(const BigInt& a, const BigInt& b, bool case_eq) {
if (a.size() < b.size()) {
return true;
} else if (a.size() > b.size()) {
return false;
}
for (int i = a.size() - 1; i >= 0; --i) {
if (a[i] != b[i]) {
return a[i] < b[i];
}
}
return case_eq;
}
static bool _shifted_ucomp(const BigInt& a, word_t ofs, const BigInt& b, bool case_eq) {
// assert(a >= b);
for (int i = b.size() - 1; i >= 0; --i) {
if (a[i + ofs] != b[i]) {
return a[i + ofs] < b[i];
}
}
return case_eq;
}
// add sub
static void _add(const BigInt& a, const BigInt& b, bool is_sub, BigInt& res) {
word_t a_size = a.size();
word_t b_size = b.size();
res.resize(std::max(a_size, b_size) + 1); // a.data() (and b.data()) might be changed.
_bigint_t a_ = _bigint_t(a.sign(), a.data(), a_size);
_bigint_t b_ = _bigint_t(b.sign(), b.data(), b_size);
_bigint_t r_ = _bigint_t(res.data());
_add_ar(a_, b_, is_sub, r_);
res.set_sign(r_.sign);
res.normalize();
}
static void _add_ar(const _bigint_t& a, const _bigint_t& b, bool is_sub, _bigint_t& res) {
// when a_size < b_size and is_sub == true, b must be normalized.
word_t b_size = b.size;
if (is_sub) {
while (b_size > a.size && b.data[b_size - 1] == 0) {
b_size -= 1;
}
}
if (a.size < b_size) {
_add_ar(b, a, is_sub, res);
if (is_sub) {
res.sign ^= 1;
}
} else {
int res_sign = a.sign;
res.size = a.size;
if (a.sign ^ b.sign ^ is_sub) {
res_sign ^= _abs_sub(a.data, a.size, b.data, b_size, res.data);
} else {
word_t c = _uadd_core(a.data, a.size, b.data, b_size, res.data);
if (c) {
res.data[res.size++] = c;
}
}
res.sign = res_sign;
}
}
static word_t _uadd_core(const word_t* a, word_t a_size, const word_t* b, word_t b_size, word_t* res) {
// assert(a_size >= b_size);
word_t i = 0;
dword_t c = 0;
for (; i < b_size; ++i) {
c = dword_t(a[i]) + b[i] + word_t(c);
res[i] = word_t(c);
c >>= W_BITS;
}
while (c && i < a_size) {
c += a[i];
res[i++] = word_t(c);
c >>= W_BITS;
}
if (res != a) {
for ( ; i < a_size; ++i) {
res[i] = a[i];
}
}
return c;
}
static void _usub(BigInt& a, const BigInt& b) {
// assert(a >= b);
_usub_core(a.data(), a.size(), b.data(), b.size(), a.data());
a.normalize();
}
static void _usub_unnorm(BigInt& a, const BigInt& b) {
_usub_core(a.data(), a.size(), b.data(), b.size(), a.data());
}
static void _shifted_usub_unnorm(BigInt& a, word_t ofs, const BigInt& b) {
_usub_core(a.data() + ofs, a.size() - ofs, b.data(), b.size(), a.data() + ofs);
}
static void _usub_core(const word_t* a, word_t a_size, const word_t* b, word_t b_size, word_t* res) {
word_t i = 0;
dword_t c = 0;
for (; i < b_size; ++i) {
c = dword_t(a[i]) - b[i] - word_t(c);
res[i] = word_t(c);
c = ((c >> W_BITS) >> (W_BITS - 1));
}
while (c && i < a_size) {
c = dword_t(a[i]) - word_t(c);
res[i++] = word_t(c);
c = ((c >> W_BITS) >> (W_BITS - 1));
}
if (res != a) {
for (; i < a_size; ++i) {
res[i] = a[i];
}
}
}
static int _abs_sub(const word_t* a, word_t a_size, const word_t* b, word_t b_size, word_t* res) {
// assert(a_size >= b_size);
for (int i = a_size - 1; i >= int(b_size); --i) {
if (a[i]) {
_usub_core(a, i + 1, b, b_size, res);
return 0;
}
res[i] = 0;
}
for (int i = b_size - 1; i >= 0; --i) {
if (a[i] != b[i]) {
if (a[i] < b[i]) {
_usub_core(b, i + 1, a, i + 1, res);
return 1;
} else {
_usub_core(a, i + 1, b, i + 1, res);
return 0;
}
}
res[i] = 0;
}
return 0;
}
// mul
static void _mul1(const BigInt& a, const word_t b, BigInt& res) {
int sign = a.sign();
_umul1(a, b, res);
res.set_sign(sign);
}
static void _umul1(const BigInt& a, const word_t b, BigInt& res) {
if (a.is_zero() || b == 0) {
res = BigInt();
return;
}
if (b == 1) {
res = BigInt(a);
return;
}
res.resize(a.size());
word_t c = _umul1_add(a.data(), a.size(), b, res.data(), 0);
if (c) {
res.push_back(c);
}
}
static word_t _umul1_add(const word_t* a, word_t a_size, const word_t b,
word_t* res, const word_t gamma) {
dword_t c = gamma;
for (word_t i = 0; i < a_size; ++i) {
c += dword_t(a[i]) * b;
res[i] = c;
c >>= W_BITS;
}
return c;
}
static void _mul(const BigInt& a, const BigInt& b, BigInt &res) {
if (&a == &b) {
_square(a, res);
} else {
int sign = a.sign() ^ b.sign();
_umul(a, b, res);
res.set_sign(sign);
}
}
static void _umul(const BigInt& a, const BigInt& b, BigInt &res) {
word_t a_bit_len = a.bit_length();
word_t b_bit_len = b.bit_length();
if (a_bit_len < b_bit_len) {
return _umul(b, a, res);
}
if (b.size() == 1) {
return _umul1(a, b[0], res);
}
if (a.is_zero() || b.is_zero()) {
res = BigInt();
return;
}
if (a.is_one()) {
res = BigInt(b);
return;
}
if (b.is_one()) {
res = BigInt(a);
return;
}
word_t a_size = a.size();
word_t b_size = b.size();
if (b.size() <= MUL_TOOM22_THRESHOLD) {
res.resize(a_size + b_size);
_umul_basecase(a.data(), a_size, b.data(), b_size, res.data());
} else {
word_t block_size = _umul_toom33_calc_block_size(a_size, b_size);
word_t nblock = (a_size + block_size - 1) / block_size;
res.resize(block_size * (nblock + 1));
word_t* work = new word_t[(block_size * 16 / 3) + 3 * 128];
BigInt c = BigInt(b);
if (block_size > b_size) {
c.resize(block_size);
}
if (nblock > 1) {
BigInt res_s = BigInt();
res_s.resize(2 * block_size);
for (word_t i = 0; i < nblock - 1; ++i) {
_umul_toom33(a.data() + i * block_size, c.data(), block_size,
res_s.data(), work);
(void) _uadd_core(
res.data() + i * block_size, 2 * block_size,
res_s.data(), 2 * block_size, res.data() + i * block_size);
}
// ... :(
BigInt last_a = a.subinteger((nblock - 1) * block_size * W_BITS, a.bit_length());
last_a.resize(block_size);
_umul_toom33(last_a.data(), c.data(), block_size, res_s.data(), work);
(void) _uadd_core(
res.data() + (nblock - 1) * block_size, res.size() - (nblock - 1) * block_size,
res_s.data(), 2 * block_size, res.data() + (nblock - 1) * block_size);
} else {
_umul_toom33(a.data(), c.data(), block_size, res.data(), work);
}
delete [] work;
}
res.normalize();
}
static word_t _umul_toom33_calc_block_size(word_t a_size, word_t b_size) {
double ratio = double(a_size) / b_size; // >= 1.0
word_t l = std::floor(ratio);
word_t h = std::ceil(ratio);
word_t block_size_a = (a_size + l - 1) / l;
double cost1 = l * _umul_toom33_cost(block_size_a);
double cost2 = h * _umul_toom33_cost(b_size);
if (cost1 < cost2) {
return block_size_a;
} else {
return b_size;
}
}
static double _umul_toom33_cost(word_t block_size) {
return std::pow(double(block_size), 1.47);
}
static int _umul_toom33_d1_d3(
const word_t* a0, const word_t* a1, const word_t* a2,
word_t size_l, word_t size_h, word_t* d1, word_t* d3) {
// d1 = a2 + a1 + a1
// d3 = a2 - a1 + a0 (sign)
word_t d1_size = size_l;
word_t c = _uadd_core(a0, size_l, a2, size_h, d1);
if (c) {
d1[d1_size++] = c;
}
word_t d3_size = d1_size;
int sign = _abs_sub(d1, d1_size, a1, size_l, d3);
c = _uadd_core(d1, d1_size, a1, size_l, d1);
if (c) {
d1[d1_size++] = c;
}
if (d3_size == size_l) {
if (d1_size == size_l) {
d1[d1_size++] = 0;
}
d3[d3_size++] = 0;
}
return sign;
}
static void _umul_toom33_d2(
const word_t* a0, const word_t* a1, const word_t* a2,
word_t size_l, word_t size_h, word_t* d2) {
// d2 = 4 * a2 + 2 * a1 + a0
word_t i = 0;
dword_t c = 0;
for ( ; i < size_h; ++i) {
c = (((dword_t(a2[i]) << 1) + a1[i]) << 1) + a0[i] + word_t(c);
d2[i] = c;
c >>= W_BITS;
}
for ( ; i < size_l; ++i) {
c = (dword_t(a1[i]) << 1) + a0[i] + word_t(c);
d2[i] = c;
c >>= W_BITS;
}
d2[i++] = c;
}
static void _umul_toom33(
const word_t* a, const word_t* b, word_t size, word_t* res, word_t*& work) {
while (size > 1 && (a[size - 1] | b[size - 1]) == 0) {
size -= 1;
res[size * 2] = res[size * 2 + 1] = 0;
}
if (size < MUL_TOOM33_THRESHOLD) {
return _umul_toom22(a, b, size, res, work);
}
const word_t size_l = (size + 2) / 3;
const word_t size_h = size - size_l * 2;
const word_t work_size = 8 * (size_l + 1);
auto* a0 = a;
auto* a1 = a0 + size_l;
auto* a2 = a1 + size_l;
auto* b0 = b;
auto* b1 = b0 + size_l;
auto* b2 = b1 + size_l;
auto* ad1 = work;
auto* bd1 = ad1 + (size_l + 1);
auto* ad3 = bd1 + (size_l + 1);
auto* bd3 = ad3 + (size_l + 1);
auto* ad2 = ad1;
auto* bd2 = bd1;
auto* s0 = res;
auto* s1 = work + 6 * (size_l + 1);
auto* s3 = work + 4 * (size_l + 1);
auto* s2 = work + 2 * (size_l + 1);
auto* s4 = res + 4 * size_l;
work += work_size;
int sign_ad3 = _umul_toom33_d1_d3(a0, a1, a2, size_l, size_h, ad1, ad3);
int sign_bd3 = _umul_toom33_d1_d3(b0, b1, b2, size_l, size_h, bd1, bd3);
_umul_toom33(ad1, bd1, size_l + 1, s1, work);
_umul_toom33(ad3, bd3, size_l + 1, s3, work);
_umul_toom33_d2(a0, a1, a2, size_l, size_h, ad2);
_umul_toom33_d2(b0, b1, b2, size_l, size_h, bd2);
_umul_toom33(ad2, bd2, size_l + 1, s2, work);
_umul_toom33(a0, b0, size_l, s0, work);
_umul_toom33(a2, b2, size_h, s4, work);
auto s0_ = _bigint_t(0, s0, 2 * size_l);
auto s1_ = _bigint_t(0, s1, 2 * (size_l + 1));
auto s2_ = _bigint_t(0, s2, 2 * (size_l + 1));
auto s3_ = _bigint_t(sign_ad3 ^ sign_bd3, s3, 2 * (size_l + 1));
auto s4_ = _bigint_t(0, s4, 2 * size_h);
auto r_ = _bigint_t(0, res + size_l, size * 2 - size_l);
_add_ar(s2_, s3_, 1, s2_);
_exact_div(s2, s2_.size, 3);
_add_ar(s1_, s3_, 1, s3_);
_rshift(s3, s3_.size, 1, s3, s3_.size);
_add_ar(s1_, s0_, 1, s1_);
_add_ar(s2_, s1_, 1, s2_);
_rshift(s2, s2_.size, 1, s2, s2_.size);
_add_ar(s1_, s3_, 1, s1_);
_add_ar(s1_, s4_, 1, s1_);
_add_ar(s2_, s4_, 1, s2_);
_add_ar(s2_, s4_, 1, s2_);
_add_ar(s3_, s2_, 1, s3_);
std::fill(res + size_l * 2, res + size_l * 4, 0);
s2_.size = std::min(s2_.size, r_.size - size_l * 2);
_add_ar(r_, s3_, 0, r_); r_.data += size_l; r_.size -= size_l;
_add_ar(r_, s1_, 0, r_); r_.data += size_l; r_.size -= size_l;
_add_ar(r_, s2_, 0, r_);
work -= work_size;
}
static void _umul_toom22(
const word_t* a, const word_t* b, word_t size, word_t* res, word_t*& work) {
if (size <= MUL_TOOM22_THRESHOLD) {
return _umul_basecase(a, size, b, size, res);
}
const word_t size_l = (size + 1) >> 1;
const word_t size_h = size - size_l;
const word_t work_size = size_l * 2 + 1;
word_t* d1 = res;
word_t* d2 = res + size_l;
word_t* c2 = work;
work += work_size;
int sign1 = _abs_sub(a, size_l, a + size_l, size_h, d1);
int sign2 = _abs_sub(b, size_l, b + size_l, size_h, d2);
_umul_toom22(d1, d2, size_l, c2, work);
_umul_toom22(a, b, size_l, res, work); // c0
_umul_toom22(a + size_l, b + size_l, size_h, res + size_l * 2, work); // c1
word_t c2_size = size_l * 2;
_bigint_t c2_ = _bigint_t(sign1 ^ sign2 ^ 1, c2, c2_size);
_bigint_t c1_ = _bigint_t(0, res + size_l * 2, size_h * 2);
_bigint_t c0_ = _bigint_t(0, res, size_l * 2);
_bigint_t r_ = _bigint_t(0, res + size_l , size_l + size_h * 2);
_add_ar(c2_, c0_, 0, c2_);
_add_ar(c2_, c1_, 0, c2_);
_add_ar(r_, c2_, 0, r_);
work -= work_size;
}
static void _umul_basecase(
const word_t* a, word_t a_size, const word_t* b, word_t b_size, word_t* res) {
std::fill(res, res + a_size + b_size, 0);
for (word_t i = 0; i < b_size; ++i) {
if (b[i] == 0) {
continue;
}
dword_t c = 0;
word_t beta = b[i];
for (word_t j = 0; j < a_size; ++j) {
c = dword_t(beta) * a[j] + res[i + j] + word_t(c);
res[i + j] = c;
c >>= W_BITS;
}
if (c) {
res[i + a_size] = c;
}
}
}
// square
static void _square(const BigInt& a, BigInt& res) {
if (a.is_zero()) {
res = BigInt();
return;
}
if (a.is_one()) {
res = BigInt(1);
return;
}
word_t a_size = a.size();
res.resize(a_size * 2);
if (a_size <= SQ_TOOM22_THRESHOLD) {
_square_basecase(a.data(), a_size, res.data());
} else if (a_size <= SQ_TOOM33_THRESHOLD) {
word_t* work = new word_t[a_size * 2 + 3 * 32];
_square_toom22(a.data(), a_size, res.data(), work);
delete [] work;
} else {
word_t* work = new word_t[a_size * 4 + 6 * 96];
_square_toom33(a.data(), a_size, res.data(), work);
delete [] work;
}
res.normalize();
return;
}
static void _square_toom33(
const word_t* a, word_t size, word_t* res, word_t*& work) {
while (size > 1 && a[size - 1] == 0) {
size -= 1;
res[size * 2] = res[size * 2 + 1] = 0;
}
if (size <= SQ_TOOM33_THRESHOLD) {
return _square_toom22(a, size, res, work);
}
const word_t size_l = (size + 2) / 3;
const word_t size_h = size - size_l * 2;
const word_t work_size = 6 * (size_l + 1);
auto* a0 = a;
auto* a1 = a0 + size_l;
auto* a2 = a1 + size_l;
auto* d1 = work;
auto* d2 = d1 + size_l + 1;
auto* s0 = res;
auto* s1 = d1 + (size_l + 1) * 2;
auto* s2 = s1 + (size_l + 1) * 2;
auto* s3 = d1;
auto* s4 = res + size_l * 4;
work += work_size;
// d1 = a2 + a1 + a0
// d2 = a2 - a1 + a0
(void) _umul_toom33_d1_d3(a0, a1, a2, size_l, size_h, d1, d2);
std::copy(a2, a2 + size_h, res);
std::fill(res + size_h, res + size_l, 0);
auto* a2_ = res;
_square_toom33(d2, size_l + 1, s2, work); // s2 = d2 ** 2 [w4, w6)
_square_toom33(d1, size_l + 1, s1, work); // s1 = d1 ** 2 [w2, w4)
_umul_toom33(a1, a2_, size_l, s3, work); // s3 = 2 * a1 * a2 [w0, w2)
_lshift(s3, size_l * 2, 1, s3, size_l * 2 + 1);
_square_toom33(a0, size_l, s0, work); // s0 = a0 ** 2 [r0, r2)
_square_toom33(a2, size_h, s4, work); // s4 = a2 ** 2 [r4, r6)
auto s0_ = _bigint_t(0, s0, 2 * size_l);
auto s1_ = _bigint_t(0, s1, 2 * (size_l + 1));
auto s2_ = _bigint_t(0, s2, 2 * (size_l + 1));
auto s3_ = _bigint_t(0, s3, 2 * size_l + 1);
auto s4_ = _bigint_t(0, s4, 2 * size_h);
auto r_ = _bigint_t(0, res + size_l, size * 2 - size_l);
_add_ar(s2_, s1_, 0, s2_);
_rshift(s2, s2_.size, 1, s2, s2_.size);
_add_ar(s1_, s2_, 1, s1_);
_add_ar(s1_, s3_, 1, s1_);
_add_ar(s2_, s4_, 1, s2_);
_add_ar(s2_, s0_, 1, s2_);
std::fill(res + size_l * 2, res + size_l * 4, 0);
s3_.size = std::min(s3_.size, r_.size - size_l * 2);
_add_ar(r_, s1_, 0, r_); r_.data += size_l; r_.size -= size_l;
_add_ar(r_, s2_, 0, r_); r_.data += size_l; r_.size -= size_l;
_add_ar(r_, s3_, 0, r_);
work -= work_size;
}
static void _square_toom22(
const word_t* a, word_t size, word_t* res, word_t*& work) {
if (size <= SQ_TOOM22_THRESHOLD) {
return _square_basecase(a, size, res);
}
const word_t size_l = (size + 1) >> 1;
const word_t size_h = size - size_l;
const word_t work_size = size_l * 2 + 1;
word_t* d = res;
word_t* c2 = work;
work += work_size;
(void) _abs_sub(a, size_l, a + size_l, size_h, d);
_square_toom22(d, size_l, c2, work);
_square_toom22(a, size_l, res, work); // c0
_square_toom22(a + size_l, size_h, res + size_l * 2, work); // c1
word_t c2_size = size_l * 2;
_bigint_t c2_ = _bigint_t(1, c2, c2_size);
_bigint_t c1_ = _bigint_t(0, res + size_l * 2, size_h * 2);
_bigint_t c0_ = _bigint_t(0, res, size_l * 2);
_bigint_t r_ = _bigint_t(0, res + size_l, size_l + size_h * 2);
_add_ar(c2_, c0_, 0, c2_);
_add_ar(c2_, c1_, 0, c2_);
_add_ar(r_, c2_, 0, r_);
work -= work_size;
}
static void _square_basecase(
const word_t* a, word_t a_size, word_t* res) {
std::fill(res, res + 2 * a_size, word_t(0));
for (word_t i = 0; i < a_size; ++i) {
word_t alpha = a[i];
if (alpha == 0) {
continue;
}
dword_t c = 0;
for (word_t j = i + 1; j < a_size; ++j) {
dword_t s = dword_t(alpha) * a[j] + res[i + j] + word_t(c);
res[i + j] = s;
c = s >> W_BITS;
}
if (c) {
res[i + a_size] = c;
}
}
word_t c = 0;
for (word_t i = 0; i < a_size; ++i) {
dword_t a2 = dword_t(a[i]) * a[i] + c;
dword_t t = (dword_t(res[2 * i]) << 1) + word_t(a2);
res[2 * i] = t;
t = (dword_t(res[2 * i + 1]) << 1) + word_t(t >> W_BITS) + word_t(a2 >> W_BITS);
res[2 * i + 1] = t;
c = t >> W_BITS;
}
}
// div-mod
static void _udivmod(const BigInt& a, const BigInt& b, BigInt& qs, BigInt& rs) {
// Modern computer arithmetic
if (_ucomp(a, b, false)) { // fixed
qs = BigInt();
rs = BigInt(a);
return;
}
if (b.is_zero()) {
assert(0);
}
if (b.is_one()) {
qs = BigInt(a);
rs = BigInt();
return;
}
word_t b_bit_len = b.bit_length();
word_t shamt = (W_BITS - b_bit_len) & (W_BITS - 1);
BigInt numer = a << shamt;
BigInt denom = b << shamt;
word_t n_size = numer.size();
word_t d_size = denom.size();
word_t ofs = n_size - d_size;
if (!_shifted_ucomp(numer, ofs, denom, false)) {
qs.resize(ofs + 1);
qs[ofs] = 1;
_shifted_usub_unnorm(numer, ofs, denom);
} else {
qs.resize(ofs);
}
word_t d = denom[d_size - 1];
auto fd = FastDiv21(d);
BigInt tmp = BigInt(d_size + 1, 0);
for (int i = ofs - 1; i >= 0; --i) {
word_t q;
if (numer[i + d_size] >= d) {
q = ~word_t(0);
} else {
q = ((dword_t(numer[i + d_size]) << W_BITS) | numer[i + d_size - 1]) / fd;
}
if (q > 0) {
tmp[d_size] = _umul1_add(denom.data(), denom.size(), q, tmp.data(), 0);
while (_shifted_ucomp(numer, i, tmp, false)) {
_usub_unnorm(tmp, denom);
q -= 1;
}
_shifted_usub_unnorm(numer, i, tmp);
}
qs[i] = q;
numer.resize(numer.size() - 1);
}
numer.normalize();
rs = BigInt(numer);
if (shamt) {
rs >>= shamt;
}
}
static void _fast_div32(
const BigInt& a21, const BigInt& a0,
const BigInt& b10, const BigInt& b1, const BigInt& b0, word_t n,
BigInt& q, BigInt& r) {
BigInt a2 = a21 >> n;
if (_ucomp(a2, b1, false)) { // fixed
_fast_div21(a21, b1, n, q, r);
} else {
q = mask(n);
r = a21; r += b1; r -= b1 << n;
}
r <<= n; r |= a0; r -= b0 * q; // bottleneck
while (r < 0) {
q -= 1;
r += b10;
}
}
static void _fast_div21(const BigInt& a, const BigInt& b, word_t n, BigInt& q, BigInt& r) {
if (_ucomp(a, b, false)) {
q = BigInt();
r = BigInt(a);
return;
}
if (n <= DIVMOD_THRESHOLD) {
return _udivmod(a, b, q, r);
}
word_t ofs = n & 1;
word_t nh = (n + ofs) >> 1;
BigInt b10 = BigInt(b);
BigInt b1 = b >> (nh - ofs);
BigInt b0 = b.subinteger(0, nh - ofs);
BigInt a32 = a >> n;
BigInt a1 = a.subinteger(nh - ofs, n);
BigInt a0 = a.subinteger(0, nh - ofs);
if (ofs) {
b0 <<= 1;
a0 <<= 1;
b10 <<= 1;
}
BigInt tmp_r;
BigInt q0;
_fast_div32(a32, a1, b10, b1, b0, nh, q, tmp_r); // fixed
_fast_div32(tmp_r, a0, b10, b1, b0, nh, q0, r);
q <<= nh;
q |= q0;
if (ofs) {
r >>= 1;
}
}
static void _fast_udivmod(const BigInt& a, const BigInt& b, BigInt& q, BigInt& r) {
if (_ucomp(a, b, false)) {
q = BigInt();
r = BigInt(a);
return;
}
word_t m = a.bit_length();
word_t n = b.bit_length();
if (n <= DIVMOD_THRESHOLD) {
return _udivmod(a, b, q, r);
}
word_t q_bit_len = m - n + 1;
word_t nblock = 1 + m / n;
q.set_zero(); // fixed
q.resize(bits_to_words(q_bit_len));
BigInt tmp_q, tmp_r;
r = a.subinteger((nblock - 1) * n, nblock * n);
for (int i = nblock - 2; i >= 0; --i) {
r <<= n;
r |= a.subinteger(i * n, (i + 1) * n);
_fast_div21(r, b, n, tmp_q, tmp_r); // fixed
_lshifted_lor(q, tmp_q, i * n);
r = tmp_r;
}
q.normalize();
}
static void _exact_div(word_t* a, word_t size, word_t d) {
assert(d & 1);
word_t inv = bits::mul_inv(d);
word_t c = 0;
for (word_t i = 0; i < size; ++i) {
a[i] = inv * (a[i] - c);
c = (dword_t(a[i]) * d) >> W_BITS;
}
}
// logical
static void _land(word_t* a, const word_t* b, word_t b_size) {
// assert(a_size == b_size);
for (word_t i = 0; i < b_size; ++i) {
a[i] &= b[i];
}
}
static void _lxor(word_t* a, const word_t* b, word_t b_size) {
// assert(a_size >= b.size);
for (word_t i = 0; i < b_size; ++i) {
a[i] ^= b[i];
}
}
static void _lor(word_t* a, const word_t* b, word_t b_size) {
// assert(a_size >= b_size);
for (word_t i = 0; i < b_size; ++i) {
a[i] |= b[i];
}
}
static void _lshifted_lor(BigInt& lhs, const BigInt& rhs, const word_t shamt) {
word_t q = shamt / W_BITS;
word_t r = shamt % W_BITS;
BigInt shifted_rhs = rhs << r;
_lor(lhs.data() + q, shifted_rhs.data(), shifted_rhs.size());
}
static void _lshift(const word_t* n, const word_t size,
word_t shamt, word_t* res, const word_t res_size) {
word_t q = shamt / W_BITS;
word_t r = shamt % W_BITS;
if (r) {
word_t c = n[size - 1];
if (res_size > size + q) {
res[res_size - 1] = c >> (W_BITS - r);
}
c <<= r;
for (int i = size - 1; i >= 1; --i) {
word_t t = n[i - 1];
res[i + q] = c | (t >> (W_BITS - r));
c = t << r;
}
res[q] = c;
} else {
for (int i = size - 1; i >= 0; --i) {
res[i + q] = n[i];
}
}
for (word_t i = 0; i < q; ++i) {
res[i] = 0;
}
}
static void _rshift(const word_t* n, const word_t size,
word_t shamt, word_t* res, const word_t res_size) {
word_t q = shamt / W_BITS;
word_t r = shamt % W_BITS;
if (r) {
word_t c = n[q] >> r;
for (word_t i = 0; i < res_size - 1; ++i) {
word_t t = n[i + q + 1];
res[i] = c | (t << (W_BITS - r));
c = t >> r;
}
if (res_size + q < size) {
c |= n[size - 1] << (W_BITS - r);
}
res[res_size - 1] = c;
} else {
for (word_t i = 0; i < size - q; ++i) {
res[i] = n[i + q];
}
}
}
// functions
static void _fib(word_t n, BigInt& a, BigInt& b) {
if (n <= 2) {
a = BigInt((n + 1) >> 1);
b = BigInt((n + 2) >> 1);
return;
}
_fib((n >> 1) - 1, a, b);
a *= a;
b *= b;
BigInt c = b + a;
BigInt d = (b << 2) - a + ((n >> 1) & 1 ? -2 : 2);
if (n & 1) {
a = d;
b = (d << 1) - c;
} else {
a = d - c;
b = d;
}
}
static void _isqrt(const BigInt& a, BigInt& s, BigInt& r) {
if (a < (uint64(1) << 52)) {
s = BigInt(word_t(std::sqrt(uint64(a))));
r = a - s * s;
return;
}
word_t n = a.bit_length();
word_t ofs = (n - 1) & 2;
if (ofs > 0) {
n += 2;
}
word_t nq = (n + 1) >> 2;
word_t nh = nq << 1;
BigInt ah = a >> (nh - ofs);
BigInt al = a.subinteger(0, nh - ofs);
BigInt a1 = al >> (nq - ofs);
BigInt a0 = al.subinteger(0, nq - ofs);
if (ofs > 0) {
a0 <<= ofs;
}
BigInt s1, r1;
_isqrt(ah, s1, r1);
BigInt q, u;
divmod((r1 << nq) + a1, s1 << 1, q, u);
s = (s1 << nq) + q;
r = (u << nq) + a0 - q * q;
if (r < 0) {
r += (s << 1);
r -= 1;
s -= 1;
}
if (ofs > 0) {
r >>= 2;
if (s & 1) {
s >>= 1;
r += s;
r += 1;
} else {
s >>= 1;
}
}
}
int sign_;
container_t n_;
};
int main() {
using namespace std;
uint32 n;
while (~scanf("%u", &n)) {
if (n == 2) {
cout << "3\nINF" << endl;
} else {
BigInt ans;
if (n & 1) {
ans = BigInt::fib(n);
} else {
BigInt a, b;
BigInt::fib(n / 2 - 1, a, b);
ans = a * b;
ans <<= 1;
}
cout << n << endl;
cout << ans << endl;
}
}
return 0;
}
Min_25