#59729 (C++11) No.303 割れません

提出ソース
結果

問題	No.303 割れません
ユーザー	Min_25
提出日時	2015-11-19 22:42:15
言語	C++11 (gcc 11.4.0)
結果	AC
実行時間	400 ms / 10,000 ms
コード長	39,684 bytes
コンパイル時間	2,329 ms
コンパイル使用メモリ	103,744 KB
実行使用メモリ	7,640 KB
最終ジャッジ日時	2024-09-13 17:04:44
合計ジャッジ時間	6,760 ms
ジャッジサーバーID （参考情報）	judge4 / judge3
このコードへのチャレンジ
（要ログイン）
テストケース

テストケース表示
入力	結果	実行時間実行使用メモリ
testcase_00	AC	2 ms 5,248 KB
testcase_01	AC	2 ms 5,376 KB
testcase_02	AC	2 ms 5,376 KB
testcase_03	AC	141 ms 5,676 KB
testcase_04	AC	384 ms 7,588 KB
testcase_05	AC	382 ms 7,004 KB
testcase_06	AC	384 ms 7,516 KB
testcase_07	AC	383 ms 7,368 KB
testcase_08	AC	399 ms 7,640 KB
testcase_09	AC	395 ms 7,132 KB
testcase_10	AC	2 ms 5,376 KB
testcase_11	AC	384 ms 7,520 KB
testcase_12	AC	400 ms 7,484 KB
testcase_13	AC	217 ms 5,840 KB
testcase_14	AC	101 ms 5,376 KB
testcase_15	AC	64 ms 5,376 KB
testcase_16	AC	2 ms 5,376 KB
権限があれば一括ダウンロードができます
ソースコード

raw source code
#include <cstdio>
#include <iostream>

#include <algorithm>
#include <cstdint>
#include <cstdio>
#include <cmath>
#include <cassert>
#include <vector>
#include <string>

// M(n): O(n^1.59)

using uint32 = uint32_t;
using int64 = int64_t;
using uint64 = uint64_t;

#define HAVE_UINT128

#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 uint32 pop_count(uint64 n) {
  return __builtin_popcountll(n);
}

inline constexpr uint32 ctz(uint32 n) {
  return __builtin_ctz(n);
}

inline constexpr uint32 ctz(uint64 n) {
  return __builtin_ctzll(n);
}

inline constexpr uint32 clz(uint32 n) {
  return __builtin_clz(n);
}

inline constexpr uint32 clz(uint64 n) {
  return __builtin_clzll(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

}

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;
  }
  word_t magic_num;
  uint32 shamt;
  word_t mod;
};

inline word_t operator / (word_t n, FastDiv d) {
  return (dword_t(n) * d.magic_num) >> d.shamt;
}

inline word_t& operator /= (word_t& n, FastDiv d) {
  return n = n / d;
}

inline word_t operator % (word_t n, FastDiv d) {
  return n - n / d * d.mod;
}

inline word_t& operator %= (word_t& n, FastDiv d) {
  return n = n % d;
}

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_;
  }

  static inline 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_;
  }

  static inline word_t div(dword_t u, FastDiv21 fd) {
    word_t q, dummy;
    divmod(u, fd, q, dummy);
    return q;
  }

  static inline word_t mod(dword_t u, FastDiv21 fd) {
    word_t dummy, r;
    divmod(u, fd, dummy, r);
    return r;
  }

private:
  word_t reciprocal(word_t n) const {
    return bits::floor_div_small(~(dword_t(n) << (sizeof(word_t) * 8)), n);
  }
  word_t d_;
  uint32 shamt_;
  word_t v_;
};

static inline word_t operator / (dword_t n, FastDiv21 fd) {
  return FastDiv21::div(n, fd);
}

static inline dword_t& operator /= (dword_t& n, FastDiv21 fd) {
  return n = n / fd;
}

static inline word_t operator % (dword_t n, FastDiv21 fd) {
  return FastDiv21::mod(n, fd);
}

static inline dword_t& operator %= (dword_t& n, FastDiv21 fd) {
  return n = n % fd;
}

class BigInt {
public:
  struct BaseData {
    BaseData(uint32 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;
    uint32 e;
    bool lower;
private:
    BaseData() {}
  };

  using container_t = std::vector<word_t>;

  enum {
    W_BITS = sizeof(word_t) * 8,
#ifdef HAVE_UINT128
    MUL_TOOM22_THRESHOLD = 32,  // words
    SQ_TOOM22_THRESHOLD = 32,   // words
    STR_THRESHOLD = 1 << 13,    // bits
    DIVMOD_THRESHOLD = 1 << 13, // bits
#else
    MUL_TOOM22_THRESHOLD = 40,  // words
    SQ_TOOM22_THRESHOLD = 40,   // words
    STR_THRESHOLD = 1 << 13,    // bits
    DIVMOD_THRESHOLD = 1 << 13, // bits
#endif
  };

public:
  BigInt() : sign_(0) {
    n_ = container_t(1, 0);
  }
  BigInt(int n) : sign_(n < 0) {
    n_ = container_t(1, std::abs(n));
  }
  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, std::abs(n));
    (*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(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 uint32 size() const {
    return n_.size();
  }

  void push_back(word_t w) {
    n_.push_back(w);
  }

  void resize(uint32 size) {
    if (size == 0) {
      set_zero();
    } else {
      n_.resize(size, 0);
    }
  }

  inline word_t operator [] (uint32 idx) const {
    return n_[idx];
  }

  inline word_t& operator [] (uint32 idx) {
    return n_[idx];
  }

  inline word_t* data() const {
    return const_cast<word_t*>(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;
  }

  bool is_neg() const {
    return sign_ == 1;
  }

  bool is_zero() const {
    return size() == 1 && (*this)[0] == 0;
  }

  // unsigned
  bool is_one() const {
    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]) {
        resize(i + 1);
        return;
      }
    }
    resize(1);
    if (is_zero()) {
      clear_sign();
    }
    return;
  }

  uint32 bit_length() const {
    if (is_zero()) {
      return 0;
    }
    uint32 len = size();
    return W_BITS * len - bits::clz((*this)[len - 1]);
  }

  static uint32 bits_to_words(uint32 bits) {
    return (bits + W_BITS - 1) / W_BITS;
  }

  uint32 pop_count() const {
    uint32 ret = 0;
    for (uint32 i = 0; i < this->size(); ++i) {
      ret += bits::pop_count((*this)[i]);
    }
    return ret;
  }

  static BigInt mask(uint32 bit_length) {
    if (bit_length == 0) {
      return BigInt();
    }
    BigInt ret = BigInt();
    uint32 ret_size = 1 + (bit_length - 1) / W_BITS;
    ret.resize(ret_size);
    for (uint32 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 hex(bool lower=true) const {
    static const char* chars = lower 
      ? "0123456789abcdef" 
      : "0123456789ABCDEF";

    if (is_zero()) {
      return "0";
    }
    std::string ret;
    if (sign_) {
      ret += "-";
    }
    word_t h = (*this)[this->size() - 1];
    bool pad = false;
    for (uint32 j = 0; j < 2 * sizeof(word_t); ++j) {
      int d = (h >> (W_BITS - 4 * j - 4)) & 15;
      if (d == 0 && !pad) {
        continue;
      }
      pad = true;
      ret += chars[d];
    }
    for (int i = this->size() - 2; i >= 0; --i) {
      word_t w = (*this)[i];
      for (uint32 j = 0; j < 2 * sizeof(word_t); ++j) {
        int d = (w >> (W_BITS - 4 * j - 4)) & 15;
        ret += chars[d];
      }
    }
    return ret;
  }

  std::string str(uint32 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 (uint32 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 {
    return _add(*this, rhs, false);
  }

  /* slow */
  BigInt operator + (const int rhs) const {
    return _add(*this, BigInt(rhs), false);
  }

  BigInt& operator += (const BigInt& rhs) {
    return _add_assign(*this, rhs, false);
  }

  BigInt operator - (const BigInt& rhs) const {
    return _add(*this, rhs, true);
  }

  /* slow */
  BigInt operator - (const int rhs) const {
    return _add(*this, BigInt(rhs), true);
  }

  BigInt& operator -= (const BigInt& rhs) {
    return _add_assign(*this, rhs, true);
  }

  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();
    uint32 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 (rhs.size() > size()) {
      BigInt ret = BigInt(*this);
      _land(ret, rhs);
      return ret;
    } else {
      BigInt ret = BigInt(rhs);
      _land(ret, *this);
      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, rhs);
    return *this;
  }

  BigInt operator ^ (const BigInt& rhs) const {
    if (rhs.size() < size()) {
      BigInt ret = BigInt(*this);
      _lxor(ret, rhs);
      return ret;
    } else {
      BigInt ret = BigInt(rhs);
      _lxor(ret, *this);
      return ret;
    }
  }

  BigInt& operator ^= (const BigInt& rhs) {
    if (size() < rhs.size()) {
      resize(rhs.size());
    }
    _lxor(*this, rhs);
    return *this;
  }

  BigInt operator | (const BigInt& rhs) const {
    if (rhs.size() < size()) {
      BigInt ret = BigInt(*this);
      _lor(ret, rhs);
      return ret;
    } else {
      BigInt ret = BigInt(rhs);
      _lor(ret, *this);
      return ret;
    }
  }

  BigInt& operator |= (const BigInt& rhs) {
    if (size() < rhs.size()) {
      resize(rhs.size());
    }
    _lor(*this, rhs);
    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();
    uint32 bit_len = bit_length();
    ret.resize(bits_to_words(bit_len + shamt));
    _lshift(*this, this->size(), shamt, ret, ret.size());
    return ret;
  }

  BigInt operator << (uint32 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;
    }
    uint32 bit_len = bit_length();
    uint32 size = this->size();
    this->resize(bits_to_words(bit_len + shamt));
    _lshift(*this, size, shamt, *this, this->size());
    return *this;
  }

  BigInt& operator <<= (uint32 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();
    }
    uint32 bit_len = bit_length();
    if (bit_len <= uint32(shamt)) {
      return BigInt();
    }
    BigInt ret = BigInt();
    ret.resize(bits_to_words(bit_len - shamt));
    _rshift(*this, this->size(), shamt, ret, ret.size());
    return ret;
  }

  BigInt operator >> (uint32 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;
    }
    uint32 bit_len = bit_length();
    if (bit_len <= uint32(shamt)) {
      this->set_zero();
      return *this;
    }
    uint32 ret_size = bits_to_words(bit_len - shamt);
    _rshift(*this, this->size(), shamt, *this, ret_size);
    this->resize(ret_size);
    return *this;
  }

  BigInt& operator >>= (uint32 shamt) {
    return *this >>= int(shamt);
  }

  // subinteger ? strip ?
  BigInt subinteger(const uint32 bit_beg, uint32 bit_end) const {
    uint32 bit_len = bit_length();
    if (bit_end > bit_len) {
      bit_end = bit_len;
    }
    if (bit_beg >= bit_end) {
      return BigInt(0);
    }
    uint32 bit_size = bit_end - bit_beg;
    uint32 ret_size = bits_to_words(bit_size);
    uint32 in_size = bits_to_words(bit_end);
    BigInt ret = BigInt(ret_size, 0);
    _rshift(*this, in_size, bit_beg, ret, 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) << (W_BITS - 1 - bits::clz(e));
    while (mask) {
      if (mask & e) {
        ret *= *this;
      }
      mask >>= 1;
      if (!mask) {
        break;
      }
      ret *= ret;
    }
    return ret;
  }

  static BigInt fib(uint32 n) {
    BigInt a, b;
    _fib(n, a, b);
    return a;
  }

  static void fib(uint32 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:
  std::string _ustr(const BaseData& base, uint32 zpad=0) const {
    static const char* chars = base.lower
      ? "0123456789abcdefghijklmnopqrstuvwxyz"
      : "0123456789ABCDEFGHIJKLMNOPQRSTUVWXYZ";

    std::string ret;

    BigInt qs = BigInt(*this);

    word_t r = 0;
    while (qs >= base.base_l) {
      divmod_n1(qs, base.fd21, qs, r);
      for (uint32 i = 0; i < base.e; ++i) {
        word_t s = r / base.fd;
        ret += chars[r - s * base.fd.mod];
        r = s;
      }
    }
    r = qs[0];
    while (r > 0) {
      word_t s = r / base.fd;
      ret += chars[r - s * base.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, uint32 e, 
      const BaseData& base, const std::vector<BigInt>& large_bases, 
      std::string& ret, bool zpad=false) {

    if (e + 1 <= 11) {
      ret += n._ustr(base, zpad ? (1 << (e + 1)) : 0);
    } else {
      BigInt q, r;
      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 += "-";
    }

    uint32 bit_len = bit_length();
    if (bit_len <= STR_THRESHOLD) {
      ret += this->_ustr(base);
      return ret;
    }

    std::vector<BigInt> large_bases;
    BigInt base_pow = BigInt(base.base_s);

    uint32 e = 0;
    while (1) {
      large_bases.push_back(base_pow);
      if (base_pow.bit_length() * 2 - 1 > bit_len) {
        break;
      }
      e += 1;
      base_pow *= base_pow;
    }
    _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, uint32 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 BigInt& _add_assign(BigInt& a, const BigInt& b, bool is_sub) {
    if (a.is_neg() ^ b.is_neg() ^ is_sub) {
      if (_ucomp(a, b, true)) {
        // :(
        BigInt ret = BigInt(b);
        _usub(ret, a);
        if (is_sub) {
          ret.change_sign();
        }
        a = ret;
      } else {
        _usub(a, b);
      }
    } else {
      if (a.size() < b.size()) {
        a.resize(b.size());
      }
      _uadd(a, b);
    }
    return a;
  }

  static BigInt _add(const BigInt& a, const BigInt& b, bool is_sub) {
    if (a.is_neg() ^ b.is_neg() ^ is_sub) {
      if (_ucomp(a, b, true)) {
        BigInt ret = BigInt(b);
        _usub(ret, a);
        if (is_sub) {
          ret.change_sign();
        }
        return ret;
      } else {
        BigInt ret = BigInt(a);
        _usub(ret, b);
        return ret;
      }
    } else {
      if (a.size() < b.size()) {
        BigInt ret = BigInt(b);
        _uadd(ret, a);
        return ret;
      } else {
        BigInt ret = BigInt(a);
        _uadd(ret, b);
        return ret;
      }
    }
  }

  static void _uadd(BigInt& a, const BigInt& b) {
    // assert(a.size() >= b.size());
    word_t c = _uadd_core(a.data(), a.size(), b.data(), b.size(), a.data());
    if (c) {
      a.push_back(word_t(c));
    }
  }

  static word_t _uadd_core(word_t* a, uint32 a_size, const word_t* b, uint32 b_size, word_t* res) {
    // assert(a_size >= b_size);
    uint32 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;
    }
    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, uint32 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, uint32 a_size, const word_t* b, uint32 b_size, word_t* res) {
    uint32 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_core(word_t* a, uint32 a_size, const word_t* b, uint32 b_size) {
    // 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, a);
        return 1;
      }
    }
    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, a);
          return -1;
        } else {
          _usub_core(a, i + 1, b, i + 1, a);
          return 1;
        }
      }
      a[i] = 0;
    }
    return 1;
  }

  // mul

  static void _mul1(const BigInt& a, const word_t b, BigInt& res) {
    int sign = a.is_neg();
    _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_core(a.data(), a.size(), b, res.data());
    if (c) {
      res.push_back(c);
    }
  }

  static word_t _umul1_core(const word_t* a, uint32 a_size, const word_t b, word_t* res) {
    dword_t c = 0;
    for (uint32 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.is_neg() ^ b.is_neg();
      _umul(a, b, res);
      res.set_sign(sign);
    }
  }

  static void _umul(const BigInt& a, const BigInt& b, BigInt &res) {
    uint32 a_bit_len = a.bit_length();
    uint32 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;
    }

    uint32 a_size = a.size();
    uint32 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 {
      uint32 block_size = _umul_toom22_calc_block_size(a_size, b_size);
      uint32 nblock = (a_size + block_size - 1) / block_size;

      res.resize(block_size * (nblock + 1));

      word_t* work = new word_t[a.size() * 8 + 100]; // nonoptimal
      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 (uint32 i = 0; i < nblock - 1; ++i) {
          _umul_toom22(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_toom22(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_toom22(a.data(), c.data(), block_size, res.data(), work);
      }

      delete [] work;
    }
    res.normalize();
  }

  static uint32 _umul_toom22_calc_block_size(uint32 a_size, uint32 b_size) {
    double ratio = double(a_size) / b_size; // >= 1.0
    uint32 l = std::floor(ratio);
    uint32 h = std::ceil(ratio);

    uint32 block_size_a = (a_size + l - 1) / l;
    double cost1 = l * _umul_toom22_cost(block_size_a);
    double cost2 = h * _umul_toom22_cost(b_size);

    if (cost1 < cost2) {
      return block_size_a;
    } else {
      return b_size;
    }
  }

  static double _umul_toom22_cost(uint32 block_size) {
    return std::pow(double(block_size), 1.60);
  }

  static void _umul_toom22(
      const word_t* a, const word_t* b, uint32 size, word_t* res, word_t*& work) {

    if (size <= MUL_TOOM22_THRESHOLD) {
      return _umul_basecase(a, size, b, size, res);
    }

    const uint32 size_l = (size + 1) >> 1;
    const uint32 size_h = size - size_l;

    const uint32 work_size = size_l * 4 + 1;

    word_t* d1 = work;
    word_t* d2 = work + size_l;
    word_t* c2 = work + size_l * 2;
    work += work_size;

    _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

    std::copy(a, a + size_l, d1);
    int sign1 = _abs_sub_core(d1, size_l, a + size_l, size_h);
    std::copy(b, b + size_l, d2);
    int sign2 = _abs_sub_core(d2, size_l, b + size_l, size_h);

    _umul_toom22(d1, d2, size_l, c2, work);

    uint32 c2_size = size_l * 2;
    int sign;

    if (sign1 * sign2 > 0) {
      sign = _abs_sub_core(c2, c2_size, res, size_l * 2);
      if (sign < 0) {
        word_t c = _uadd_core(c2, c2_size, res + size_l * 2, size_h * 2, c2);
        if (c) {
          c2[c2_size++] = c;
        }
        sign = 1;
      } else {
        sign = -_abs_sub_core(c2, c2_size, res + size_l * 2, size_h * 2);
      }
    } else {
      word_t c = _uadd_core(c2, c2_size, res, size_l * 2, c2);
      if (c) {
        c2[c2_size++] = c;
      }
      c = _uadd_core(c2, c2_size, res + size_l * 2, size_h * 2, c2);
      if (c) {
        c2[c2_size++] = c;
      }
      sign = 1;
    }
    if (sign < 0) {
      _usub_core(res + size_l, size_l + size_h * 2, c2, c2_size, res + size_l);
    } else {
      _uadd_core(res + size_l, size_l + size_h * 2, c2, c2_size, res + size_l);
    }
    work -= work_size;
  }

  static void _umul_basecase(
      const word_t* a, uint32 a_size, const word_t* b, uint32 b_size, word_t* res) {

    std::fill(res, res + a_size + b_size, 0);
    for (uint32 i = 0; i < b_size; ++i) {
      if (b[i] == 0) {
        continue;
      }
      dword_t c = 0;
      word_t beta = b[i];
      for (uint32 j = 0; j < a_size; ++j) {
        dword_t s = dword_t(beta) * a[j] + res[i + j] + word_t(c);
        res[i + j] = s;
        c = s >> 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;
    }
    uint32 a_size = a.size();
    res.resize(a_size * 2);
    if (a_size <= SQ_TOOM22_THRESHOLD) {
      _square_basecase(a.data(), a_size, res.data());
    } else {
      word_t* work = new word_t[a_size * 6 + 100]; // nonoptimal
      _square_toom22(a.data(), a_size, res.data(), work);
      delete [] work;
    }
    res.normalize();
    return;
  }

  static void _square_toom22(
      const word_t* a, uint32 size, word_t* res, word_t*& work) {

    if (size <= SQ_TOOM22_THRESHOLD) {
      return _square_basecase(a, size, res);
    }

    const uint32 size_l = (size + 1) >> 1;
    const uint32 size_h = size - size_l;
    const uint32 work_size = size_l * 3 + 1;

    word_t* d = work;
    word_t* c2 = work + size_l;
    work += work_size;

    _square_toom22(a         , size_l, res             , work); // c0
    _square_toom22(a + size_l, size_h, res + size_l * 2, work); // c1

    std::copy(a, a + size_l, d);
    (void) _abs_sub_core(d, size_l, a + size_l, size_h);

    _square_toom22(d, size_l, c2, work);

    uint32 c2_size = size_l * 2;
    int sign = _abs_sub_core(c2, c2_size, res, size_l * 2);
    if (sign < 0) {
      word_t c = _uadd_core(c2, c2_size, res + size_l * 2, size_h * 2, c2);
      if (c) {
        c2[c2_size++] = c;
      }
      sign = 1;
    } else {
      sign = _abs_sub_core(c2, c2_size, res + size_l * 2, size_h * 2);
      sign = -sign;
    }
    if (sign < 0) {
      _usub_core(res + size_l, size_l + size_h * 2, c2, c2_size, res + size_l);
    } else {
      _uadd_core(res + size_l, size_l + size_h * 2, c2, c2_size, res + size_l);
    }
    work -= work_size;
  }

  static void _square_basecase(
      const word_t* a, uint32 a_size, word_t* res) {

    std::fill(res, res + 2 * a_size, word_t(0));
    for (uint32 i = 0; i < a_size; ++i) {
      word_t alpha = a[i];
      if (alpha == 0) {
        continue;
      }
      dword_t c = 0;
      for (uint32 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 (uint32 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 (a < b) {
      qs = BigInt();
      rs = BigInt(a);
      return;
    }
    if (b.is_zero()) {
      assert(0);
    }
    if (b.is_one()) {
      qs = BigInt(a);
      rs = BigInt();
      return;
    }
    BigInt denom;
    BigInt numer;

    uint32 b_bit_len = b.bit_length();
    uint32 shamt = (W_BITS - b_bit_len) & (W_BITS - 1);
    
    if (shamt > 0) {
      denom = b << shamt;
      numer = a << shamt;
    } else {
      denom = b;
      numer = a;
    }

    uint32 n_size = numer.size();
    uint32 d_size = denom.size();
    uint32 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];
    FastDiv21 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_core(denom.data(), denom.size(), q, tmp.data());
        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, uint32 n,
      BigInt& q, BigInt& r) {

    BigInt a2 = a21 >> n;
    if (a2 < b1) {
      _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, uint32 n, BigInt& q, BigInt& r) {
    if (a < b) {
      q = BigInt();
      r = BigInt(a);
      return;
    }
    if (n <= DIVMOD_THRESHOLD) {
      return _udivmod(a, b, q, r);
    }
    uint32 ofs = n & 1;
    uint32 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);
    BigInt q0;

    if (ofs) {
      b0 <<= 1;
      a0 <<= 1;
      b10 <<= 1;
    }

    _fast_div32(a32, a1, b10, b1, b0, nh, q, r);
    _fast_div32(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 (a < b) {
      q = BigInt();
      r = BigInt(a);
      return;
    }
    uint32 m = a.bit_length();
    uint32 n = b.bit_length();

    if (n <= DIVMOD_THRESHOLD) {
      return _udivmod(a, b, q, r);
    }

    uint32 q_bit_len = m - n + 1;
    uint32 nblock = 1 + m / n;

    q.resize(bits_to_words(q_bit_len));

    BigInt tmp_q;
    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, r);
      _lshifted_lor(q, tmp_q, i * n);
    }
    q.normalize();
  }

  // logical
  static void _land(BigInt& res, const BigInt& b) {
    // assert(res.size() <= b.size());
    for (uint32 i = 0; i < res.size(); ++i) {
      res[i] &= b[i];
    }
    res.normalize();
  }

  static void _lxor(BigInt& res, const BigInt& b) {
    // assert(res.size() >= b.size());
    for (uint32 i = 0; i < b.size(); ++i) {
      res[i] ^= b[i];
    }
    res.normalize();
  }

  static void _lor(BigInt& res, const BigInt& b) {
    // assert(res.size() >= b.size());
    _lor_core(res.data(), b.data(), b.size());
  }

  static void _lor_core(word_t* a, const word_t* b, uint32 b_size) {
    for (uint32 i = 0; i < b_size; ++i) {
      a[i] |= b[i];
    }
  }

  static void _lshifted_lor(BigInt& lhs, const BigInt& rhs, const uint32 shamt) {
    uint32 q = shamt / W_BITS;
    uint32 r = shamt % W_BITS;
    BigInt shifted_rhs = rhs << r;
    _lor_core(lhs.data() + q, shifted_rhs.data(), shifted_rhs.size());
  }

  static void _lshift(const BigInt& n, const uint32 size, 
      uint32 shamt, BigInt& res, const uint32 res_size) {

    uint32 q = shamt / W_BITS;
    uint32 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 (uint32 i = 0; i < q; ++i) {
      res[i] = 0;
    }
  }

  static void _rshift(const BigInt& n, const uint32 size, 
      uint32 shamt, BigInt& res, const uint32 res_size) {

    uint32 q = shamt / W_BITS;
    uint32 r = shamt % W_BITS;

    if (r) {
      word_t c = n[q] >> r;
      for (uint32 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 (uint32 i = 0; i < size - q; ++i) {
        res[i] = n[i + q];
      }
    }
  }

  // functions

  static void _fib(uint32 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;
    }

    uint32 n = a.bit_length();
    uint32 ofs = (n - 1) & 2;

    if (ofs > 0) {
      n += 2;
    }

    uint32 nq = (n + 1) >> 2;
    uint32 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;
}
yukicoder

結果

テストケース

ソースコード
