#include #include using namespace std; using i64 = long long; // {{{ mod_mul, extgcd, mod_inv, mod_pow inline i64 mod_mul(const i64& x, const i64& y, const i64& m) { // yとmの最大ビット長はともに48ビット。 // http://codeforces.com/blog/entry/15884 // 12 < floor(63 - log2(MAX_VALUE)) = 15 // 0xfff = (1<<12) - 1 i64 b, res = y * (x & 0xfff); res += (b=(y<<12) % m) * ((x>>12) & 0xfff); res += (b=(b<<12) % m) * ((x>>24) & 0xfff); res += (b<<12) % m * ((x>>36) & 0xfff); res %= m; return res; } inline tuple extgcd(const i64& x, const i64& y) { if(y == 0) { return make_tuple(1LL, 0LL, x); } i64 a, b, d; tie(a, b, d) = extgcd(y, x%y); return make_tuple(b, a-x/y*b, d); } inline i64 mod_inv(const i64 a, const i64 m) { i64 x, _, d; tie(x, _, d) = extgcd(a, m); return (x % m + m) % m; } inline i64 mod_pow(const i64& a, const i64& r, const i64& n, const i64& m) { if(n == 0) { return r; } if(n & 1) { return mod_pow(a, mod_mul(a, r, m), n-1, m); } return mod_pow(mod_mul(a, a, m), r, n>>1, m); } inline i64 mod_pow(const i64 a, const i64 n, const i64 m) { return mod_pow(a, 1, n, m); } // }}} struct FMT { static void fmt(const int, vector&); static void ifmt(const int, vector&); static void convol(const int, vector&, vector&); private: static constexpr i64 N = 5LL << 37; static constexpr i64 O = 2003LL; static constexpr i64 M = N * 211 + 1; static void myfmt(const int, vector&, const bool); }; // {{{ FMT 実装 void FMT::myfmt(const int logn, vector &a, const bool inv) { const int n = 1 << logn; for(int H=1, W=n>>1; H>=1) { vector y = a; i64 r = mod_pow(O, N/(H*2), M); if(inv) { r = mod_inv(r, M); } i64 w = 1; for(int k=0; k &a) { myfmt(logn, a, false); } void FMT::ifmt(const int logn, vector &a) { myfmt(logn, a, true); int n = 1 << logn; i64 inv = mod_inv(n, M); for(int i=0; i &a, vector &b) { fmt(logn, a); fmt(logn, b); for(int i=0; i<1< v; static bigint init(const int size); inline bool absGreater(const bigint&, const bigint&) const; inline bigint absAdd(const bigint&, const bigint&) const; inline bigint sub(const bigint&, const bigint&) const; inline bigint fftmul(const bigint&, const bigint&) const; inline vector karatsuba(const vector&, const vector&) const; inline bigint karatsuba(const bigint&, const bigint&) const; static bigint pow(const bigint&, const bigint&, const int); }; // {{{ bigint 実装 bigint::bigint(void) {} bigint bigint::init(const int size) { bigint res; res.n = size; res.v.assign(size, 0LL); return res; } bigint::bigint(const string s) { minus = s[0] == '-'; int last = int(s.size()); int m = last - minus; // 桁数 n = (m + BASE_DIGITS - 1) / BASE_DIGITS; v.assign(n, 0LL); for(int i=last-1, k=0; i>=minus; i-=BASE_DIGITS, ++k) { for(int j=max(minus, i-BASE_DIGITS+1); j<=i; ++j) { v[k] = v[k] * 10 + s[j] - '0'; } } if(v[n-1] == 0) { minus = false; } // -0にはしない } bigint::bigint(const i64 x) { constexpr int N = 5; *this = bigint::init(N); minus = x < 0; i64 y = x; for(int i=0; i 1) { --n; } v.resize(n); if(v[n-1] == 0) { minus = false; } // -0にはしない } bigint& bigint::operator = (const bigint &other) { minus = other.minus; n = other.n; v = other.v; return *this; } bool bigint::absGreater(const bigint& a, const bigint& b) const { // abs(a) < abs(b)かどうか。minusフラグを見ずにチェックすればいい。 if(a.n != b.n) { return a.n < b.n; } for(int i=a.n-1; i>=0; --i) { if(a.v[i] != b.v[i]) { return a.v[i] < b.v[i]; } } return false; // abs(a) == abs(b)のときはfalseを返す } bool bigint::operator < (bigint& other) const { if(minus != other.minus) { return minus; } return absGreater(*this, other) ^ minus; } bigint bigint::absAdd(const bigint& a, const bigint& b) const { // abs(a) + abs(b) を計算する int size = max(a.n, b.n) + 1; bigint res = bigint::init(size); for(int i=0; i= BASE) { // 繰り上がり ++res.v[i+1]; res.v[i] -= BASE; } } while(res.v[res.n-1] == 0 && res.n > 1) { --res.n; } res.v.resize(res.n); return res; } bigint bigint::abs() { bigint res = *this; res.minus = false; return res; } bigint bigint::sub(const bigint& a, const bigint& b) const { // 必ず abs(a) >= abs(b) になるように引数は来る // abs(a) - abs(b) を計算する int size = max(a.n, b.n) + 1; bigint res = bigint::init(size); for(int i=0; i 1) { --res.n; } res.v.resize(res.n); return res; } bigint bigint::operator + (const bigint& other) const { bigint res; if(minus == other.minus) { // +, + or -, - res = absAdd(*this, other); res.minus = minus; } else if(absGreater(*this, other)) { // +, - or -, + res = sub(other, *this); // sub()を簡単に実装するために、こことelseの部分の二つはまとめない方がいいと思う res.minus = other.minus; } else { res = sub(*this, other); res.minus = minus; } if(res.n == 1 && res.v[0] == 0) { res.minus = false; } return res; } bigint operator + (const i64& other, const bigint& x) { return x + other; } bigint bigint::operator - () const { bigint obj = *this; obj.minus = !minus; return obj; } bigint bigint::operator - (const bigint& other) const { return *this + (-other); } bigint bigint::fftmul(const bigint& a, const bigint& b) const { int size = max(a.n, b.n) * 2; int m; for(m=0; 1< c(1< 1) { --res.n; } res.v.resize(res.n); res.minus = a.minus ^ b.minus; return res; } vector bigint::karatsuba(const vector& a, const vector& b) const { int n = int(a.size()); vector res(n + n); if(n <= 32) { for(int i=0; i> 1; vector a1(begin(a), begin(a) + k), a2(begin(a) + k, end(a)), b1(begin(b), begin(b) + k), b2(begin(b) + k, end(b)), a1b1 = karatsuba(a1, b1), a2b2 = karatsuba(a2, b2); for(int i=0; i r = karatsuba(a2, b2); for(int i=0, m=int(a1b1.size()); i a(begin(ba.v), end(ba.v)), b(begin(bb.v), end(bb.v)); int size = max(ba.n, bb.n); int m; for(m=0 ; 1< c = karatsuba(a, b); int cn = int(c.size()); bigint res = bigint::init(cn); for(int i=0; i 1) { --res.n; } res.v.resize(res.n); res.minus = ba.minus ^ bb.minus; return res; } bigint bigint::operator * (const bigint& other) const { return max(n, other.n) < 1000 ? karatsuba(*this, other) : fftmul(*this, other); } bigint bigint::operator * (const i64& other) const { return karatsuba(*this, bigint(other)); } bigint& bigint::operator += (const bigint& other) { *this = *this + other; return *this; } bigint& bigint::operator -= (const bigint& other) { *this = *this - other; return *this; } bigint& bigint::operator *= (const bigint& other) { *this = *this * other; return *this; } bigint& bigint::operator *= (const i64& other) { *this = *this * bigint(other); return *this; } bigint bigint::pow(const bigint &x, const bigint &r, const int n) { if(n == 0) { return r; } if(n & 1) { return pow(x, r*x, n-1); } return pow(x*x, r, n>>1); } bigint bigint::pow(const int n) const { return bigint::pow(*this, 1, n); } void bigint::print() const { if(minus) { putchar('-'); } printf("%lld", v[n-1]); for(int i=n-2; i>=0; --i) { printf("%04lld", v[i]); } } void bigint::println() const { print(); putchar('\n'); } inline bigint recfact(const int start, const int n) { if(n <= 2) { bigint r = bigint(start); for(int i=start+1; i struct BIT { int N; vector data; BIT(int n); T sum(int k); void add(int k, T x); }; // {{{ BIT 実装 template BIT::BIT(int n) { N = n + 1; data.assign(N, 0); } template T BIT::sum(int k) { T res = 0; for(++k; k>0; k-=k&-k) { res += data[k]; } return res; } template void BIT::add(int k, T x) { for(++k; k v; pair f(int lo, int hi) { if(hi - lo == 1) { return {v[lo], lo+1}; } int md = (lo + hi) / 2; bigint x0, x1, y0, y1; tie(x0, x1) = f(lo, md); tie(y0, y1) = f(md, hi); return {y0*x1+x0, x1*y1}; } int main(void) { int n; scanf("%d", &n); vector a(n); v.resize(n); BIT bit(n+1); for(int i=0; i