提出 #64335602


ソースコード 拡げる

#include <math.h>
#include <stdint.h>
#include <algorithm>
#include <array>
#include <bitset>
#include <cassert>
#include <cmath>
#include <cstring>
#include <deque>
#include <fstream>
#include <functional>
#include <iomanip>
#include <iostream>
#include <list>
#include <map>
#include <numeric>
#include <optional>
#include <queue>
#include <random>
#include <set>
#include <stack>
#include <string>
#include <unordered_map>
#include <unordered_set>
#include <vector>

#include <algorithm>
#include <array>
#include <bitset>
#include <cassert>
#include <cstring>
#include <deque>
#include <fstream>
#include <functional>
#include <iomanip>
#include <iostream>
#include <list>
#include <map>
#include <random>
#include <set>
#include <stack>
#include <string>
#include <unordered_map>
#include <unordered_set>
#include <vector>

struct xor_space {
    std::vector<uint64_t> basis;
    bool                  in(uint64_t x) {
        for (uint64_t b : basis)
            if ((b ^ x) < x)
                x = b ^ x;
        return x == 0;
    }
    void add(uint64_t x) {
        for (uint64_t b : basis)
            if ((b ^ x) < x)
                x = b ^ x;
        if (x != 0)
            basis.push_back(x);
    }
};

template<typename T>
struct matrix {
    std::vector<std::vector<T>> dat;
    int                         r;
    int                         c;

    matrix() {
    }
    matrix(int r, int c) : dat(r, std::vector<T>(c)), r(r), c(c) {
    }
    matrix(int r, int c, T val) : dat(r, std::vector<T>(c, val)), r(r), c(c) {
    }
    matrix(const matrix<T>& src) : dat(src.dat), r(src.r), c(src.c) {
    }

    const matrix<T>& operator=(const matrix<T>& src) {
        dat = src.dat;
        r = src.r;
        c = src.c;
        return *this;
    }

    std::vector<T>& operator[](int ind) {
        return dat[ind];
    }
    const std::vector<T>& operator[](int ind) const {
        return dat[ind];
    }

    matrix<T> operator*(const matrix<T>& m) {
        matrix<T> ans{ r, m.c };

        int K = std::max<int>(c, m.r);

        for (int i = 0; i < r; ++i) {
            for (int j = 0; j < m.c; ++j) {
                ans[i][j] = dat[i][0] * m[0][j];
                for (int k = 1; k < K; ++k)
                    ans[i][j] += dat[i][k] * m[k][j];
            }
        }
        return ans;
    }

    void init(int r, int c, T v) {
        dat.resize(r, std::vector<T>(c, v));
        this->r = r;
        this->c = c;
    }
};

template<int col>
struct bitmatrix {
    const int                     c = col;
    std::vector<std::bitset<col>> dat;
    int                           r;

    bitmatrix() : r(0) {
    }
    bitmatrix(int r) : dat(r), r(r) {
    }
    bitmatrix(const bitmatrix<col>& src) : dat(src.dat), r(src.r) {
    }
    const bitmatrix<col>& operator=(const bitmatrix<col>& src) {
        dat = src.dat;
        r = src.r;
        return *this;
    }

    std::bitset<col>& operator[](int ind) {
        return dat[ind];
    }
    const std::bitset<col>& operator[](int ind) const {
        return dat[ind];
    }

    void init(int r) {
        dat.resize(r);
        this->r = r;
    }
};

template<typename T>
int gauss_jordan(const matrix<T>& src, matrix<T>& ans, bool is_extended = false) {
    ans = src;

    int nxtr = 0;
    for (int c = 0; c < ans.c; ++c) {
        if (is_extended && ans.c - 1 == c)
            break;

        int r = nxtr;
        for (; r < ans.r; ++r)
            if (ans[r][c] != 0)
                break;

        if (r == ans.r)
            continue;

        std::swap(ans[r], ans[nxtr]);
        r = nxtr++;

        for (int i = r + 1; i < ans.r; ++i) {
            if (ans[i][c] == 0)
                continue;
            for (int j = c; j < ans.c; ++j)
                ans[i][j] = ans[i][j] / ans[i][c] - ans[r][j];
        }
    }
    return nxtr;
}

int gauss_jordan(const matrix<double>& src, matrix<double>& ans, bool is_extended = false, double eps = -1e9) {
    ans = src;

    int nxtr = 0;
    for (int c = 0; c < ans.c; ++c) {
        if (is_extended && ans.c - 1 == c)
            break;

        int r = nxtr;
        for (; r < ans.r; ++r)
            if (abs(ans[r][c]) < eps)
                break;

        if (r == ans.r)
            continue;

        std::swap(ans[r], ans[nxtr]);
        r = nxtr++;

        for (int i = r + 1; i < ans.r; ++i) {
            for (int j = c; j < ans.c; ++j) {
                if (abs(ans[i][c]) < eps)
                    continue;
                ans[i][j] = ans[i][j] / ans[i][c] - ans[r][j];
            }
        }
    }
    return nxtr;
}

template<unsigned int col>
int gauss_jordan(const bitmatrix<col>& src, bitmatrix<col>& ans, bool is_extended = false) {
    ans = src;

    int nxtr = 0;
    for (int c = 0; c < ans.c; ++c) {
        if (is_extended && ans.c - 1 == c)
            break;

        int r = nxtr;
        for (; r < ans.r; ++r)
            if (ans[r][c])
                break;

        if (r == ans.r)
            continue;

        std::swap(ans[r], ans[nxtr]);
        r = nxtr++;

        for (int i = r + 1; i < ans.r; ++i)
            if (ans[i][c])
                ans[i] ^= ans[r];
    }
    return nxtr;
}

template<typename T>
class BIT {
public:
    long long      N;
    std::vector<T> seg;

    BIT() {
    }
    BIT(long long N) : N(N) {
        seg.resize(N + 1);
        for (int i = 0; i < N + 1; ++i)
            seg[i] = 0;
    }
    void operator=(const BIT<T>& source) {
        N = source.N;
        seg = source.seg;
    }

    //i番目にxを足す
    void add(int i, T x) {
        ++i;
        while (i <= N) {
            seg[i] += x;
            i += i & (-i);
        }
    }

    //[0,i)の累積
    T cumulative(int i) {
        T s = 0;
        while (i > 0) {
            s += seg[i];
            i -= i & (-i);
        }
        return s;
    }

    //[a,b)の累積
    T sum(int a, int b) {
        return cumulative(b) - cumulative(a);
    }

    //i番目の要素
    T elem(int i) {
        return sum(i, i + 1);
    }

    //i番目にxを代入
    void assign(int i, T x) {
        add(i, -elem(i));
        add(i, x);
    }
};

template<typename T>
class BIT2D {
public:
    long long                   H, W;
    std::vector<std::vector<T>> v;
    BIT2D() {
    }
    BIT2D(long long H_, long long W_) {
        H = H_ + 1;
        W = W_ + 1;
        v.assign(H, std::vector<T>(W, 0));
    }
    void add(int h, int w, const T& x) {
        ++h;
        ++w;
        for (int i = h; i < H; i += (i & -i)) {
            for (int j = w; j < W; j += (j & -j)) {
                v[i][j] += x;
            }
        }
    }
    //0<=y<=h, 0<=x<=w???
    T sum(int h, int w) {
        ++h;
        ++w;
        T ret{ 0 };
        for (int i = h; i > 0; i -= (i & -i)) {
            for (int j = w; j > 0; j -= (j & -j)) {
                ret += v[i][j];
            }
        }
        return ret;
    }
    //h1<=y<=h2,w1<=x<=w2???
    T sum(int h1, int w1, int h2, int w2) {
        return sum(h2, w2) - sum(h2, w1 - 1) - sum(h1 - 1, w2) + sum(h1 - 1, w1 - 1);
    }

    T elem(int h, int w) {
        return sum(h, w, h, w);
    }

    void assign(int h, int w, const T& x) {
        add(h, w, -elem(h, w));
        add(h, w, x);
    }
};

int next_subset_descending(int sub, int super) {
    return (sub - 1) & super;
}

int next_subset_ascending(int sub, int super) {
    return (sub - super) & super;
}

int next_subset_including_x(int x, int sub) {
    return (sub + 1) | x;
}

int next_combination(int sub) {
    int x = sub & -sub, y = sub + x;
    return (((sub & ~y) / x) >> 1) | y;
}

long long fill_lower(long long b) {
    b |= b >> 1;
    b |= b >> 2;
    b |= b >> 4;
    b |= b >> 8;
    b |= b >> 16;
    b |= b >> 32;
    return b;
}

long long fill_upper(long long b) {
    b |= b << 1;
    b |= b << 2;
    b |= b << 4;
    b |= b << 8;
    b |= b << 16;
    b |= b << 32;
    return b;
}

long long msb(long long b) {
    long long tmp = fill_lower(b);
    return b & ~(tmp >> 1);
}

long long lsb(long long b) {
    return b & (-b);
}

template<typename T>
class com_table {
public:
    std::vector<T> frac;
    std::vector<T> finv;

    void init(long long size) {
        frac.resize(size + 1);
        finv.resize(size + 1);
        frac[0] = 1;
        finv[0] = 1;
        for (int i = 1; i <= size; ++i)
            frac[i] = frac[i - 1] * T(i);
        finv[size] = T(1) / frac[size];
        for (int i = size; i >= 2; --i)
            finv[i - 1] = finv[i] * T(i);
    }
    T operator()(int a, int b) {
        if (a < b)
            return 0;
        return frac[a] * finv[b] * finv[a - b];
    }

    com_table() {
    }
    com_table(long long size) {
        init(size);
    }
};

template<typename T>
class frac_table {
public:
    std::vector<T> frac;

    void init(long long size) {
        frac.resize(size + 1);
        frac[0] = 1;
        for (int i = 1; i <= size; ++i)
            frac[i] = frac[i - 1] * T(i);
    }
    T operator()(int a) {
        return frac[a];
    }

    frac_table() {
    }
    frac_table(long long size) {
        init(size);
    }
};

/*
segment tree
LAZY         = true -> 遅延セグ木.  注意 : applyは左にかける
LAZY = BEATS = true -> SGT Beats
BEATS = trueのとき必要なこと
(1) T::fail : bool がある
(2) mappingが失敗したときのみfail = true, 区間幅1の要素へのmappingは失敗しない
*/
template<typename T, bool LAZY = false, typename M = int, bool BEATS = false>
struct segtree {
private:
    std::function<M(M, M)> composition;
    std::function<T(M, T)> mapping;
    std::function<T(T, T)> op;
    std::function<M()>     id;
    std::function<T()>     e;
    std::vector<M>         lazy;
    std::vector<T>         dat;
    long long              size;
    long long              log;

    void update(int ind) {
        dat[ind] = op(dat[ind << 1 | 0], dat[ind << 1 | 1]);
    }
    void node_apply(int ind, const M& m) {
        dat[ind] = mapping(m, dat[ind]);
        if (ind < size) {
            lazy[ind] = composition(m, lazy[ind]);
            if constexpr (BEATS) {
                if (dat[ind].fail)
                    prop(ind), update(ind);
            }
        }
    }
    void prop(int ind) {
        node_apply(ind << 1 | 0, lazy[ind]);
        node_apply(ind << 1 | 1, lazy[ind]);
        lazy[ind] = id();
    }
    void prop_topdown(int ind) {
        for (int i = log; i > 0; --i)
            prop(ind >> i);
    }

public:
    segtree() : size{ 0 }, log{ 0 } {
    }
    segtree(long long _n, std::function<T(T, T)> _op, std::function<T()> _e) : op{ _op }, e{ _e } {
        static_assert(!LAZY);
        size = 1, log = 0;
        while (size < _n)
            size <<= 1, ++log;
        dat.resize(2 * size, e());
    }
    segtree(const std::vector<T>& _src, std::function<T(T, T)> _op, std::function<T()> _e) : op{ _op }, e{ _e } {
        static_assert(!LAZY);
        size = 1, log = 0;
        while (size < _src.size())
            size <<= 1, ++log;
        dat.resize(2 * size);
        for (int i = 0; i < _src.size(); ++i)
            dat[i + size] = _src[i];
        for (int i = _src.size(); i < size; ++i)
            dat[i + size] = e();
        for (int i = size - 1; i > 0; --i)
            update(i);
    }
    segtree(long long              _n,            //
        std::function<T(T, T)> _op,           //
        std::function<T()>     _e,            //
        std::function<T(M, T)> _mapping,      //
        std::function<M(M, M)> _composition,  //
        std::function<M()>     _id)               //
        : op{ _op }
        ,  //
        e{ _e }
        ,  //
        mapping{ _mapping }
        ,  //
        composition{ _composition }
        ,           //
        id{ _id }  //
    {
        static_assert(LAZY);
        size = 1, log = 0;
        while (size < _n)
            size <<= 1, ++log;
        dat.resize(2 * size, e());
        lazy.resize(2 * size, id());
    }
    segtree(const std::vector<T>& _src,          //
        std::function<T(T, T)> _op,           //
        std::function<T()>     _e,            //
        std::function<T(M, T)> _mapping,      //
        std::function<M(M, M)> _composition,  //
        std::function<M()>     _id)               //
        : op{ _op }
        ,  //
        e{ _e }
        ,  //
        mapping{ _mapping }
        ,  //
        composition{ _composition }
        ,           //
        id{ _id }  //
    {
        static_assert(LAZY);
        size = 1, log = 0;
        while (size < _src.size())
            size <<= 1, ++log;
        dat.resize(2 * size);
        lazy.resize(2 * size, id());
        for (int i = 0; i < _src.size(); ++i)
            dat[i + size] = _src[i];
        for (int i = _src.size(); i < size; ++i)
            dat[i + size] = e();
        for (int i = size - 1; i > 0; --i)
            update(i);
    }
    //O(1)
    T all_prod() {
        return dat[1];
    }

    //[l,r) mは左にかける
    void apply(long long l, long long r, const M& m) {
        static_assert(LAZY);
        assert(0 <= l && l <= r && r <= size);
        if (l == r)
            return;
        l += size;
        r += size;
        for (int i = log; i >= 1; i--) {
            if (((l >> i) << i) != l)
                prop(l >> i);
            if (((r >> i) << i) != r)
                prop((r - 1) >> i);
        }
        long long a = l, b = r;
        while (a < b) {
            if (a & 1)
                node_apply(a++, m);
            if (b & 1)
                node_apply(--b, m);
            a >>= 1, b >>= 1;
        }
        for (int i = 1; i <= log; ++i) {
            if ((l >> i) << i != l)
                update(l >> i);
            if ((r >> i) << i != r)
                update((r - 1) >> i);
        }
    }
    //[l,r)   l = r -> e
    T query(long long l, long long r) {
        static_assert(!BEATS || LAZY);
        assert(0 <= l && l <= r && r <= size);
        l += size, r += size;
        if constexpr (LAZY) {
            for (int i = log; i >= 1; i--) {
                if (((l >> i) << i) != l)
                    prop(l >> i);
                if (((r >> i) << i) != r)
                    prop((r - 1) >> i);
            }
        }
        T tmp_l = e();
        T tmp_r = e();
        while (l < r) {
            if (l & 1)
                tmp_l = op(tmp_l, dat[l++]);
            if (r & 1)
                tmp_r = op(dat[--r], tmp_r);
            l >>= 1, r >>= 1;
        }
        return op(tmp_l, tmp_r);
    }

    T elem(long long i) {
        assert(0 <= i && i < size);
        i += size;
        if constexpr (LAZY)
            prop_topdown(i);
        return dat[i];
    }

    void set(long long i, const T& value) {
        assert(0 <= i && i < size);
        i += size;
        if constexpr (LAZY)
            prop_topdown(i);
        dat[i] = value;
        for (int j = 1; j <= log; ++j)
            update(i >> j);
    }
    //cond(op[x,r))は左に単調であるとする
    //cond(e) = trueとする
    //x <= rで、cond[x-1, r) = falseとなる最大のx   全区間trueなら0
    //secondはop[x,r)
    //ただしr = 0 -> (0, e)
    std::pair<int, T> search_left(long long r, std::function<bool(T)> cond) {
        assert(0 <= r && r <= size);
        assert(cond(e()));
        if (r == 0)
            return { 0, e() };
        r += size;
        T tmp = e();
        if constexpr (LAZY)
            prop_topdown(r - 1);
        do {
            --r;
            while (r % 2 && r > 1)
                r >>= 1;
            T next = op(dat[r], tmp);
            if (cond(next)) {
                tmp = next;
                continue;
            }
            while (r < size) {
                if constexpr (LAZY)
                    prop(r);
                r <<= 1, ++r;
                T next = op(dat[r], tmp);
                if (cond(next)) {
                    --r;
                    tmp = next;
                }
            }
            return { r - size + 1, tmp };
        } while ((r & -r) != r);
        return { 0, tmp };
    }
    //cond(op[l, x))は右に単調であるとする
    //cond(e) = trueとする
    //l <= x で、cond(op[l,x+1))=falseとなる最小のx  全区間trueならsize
    //secondはop[l,x)
    //ただしl = size -> (size, e)
    std::pair<int, T> search_right(long long l, std::function<bool(T)> cond) {
        assert(0 <= l && l <= size);
        assert(cond(e()));
        if (l == size)
            return { size, e() };
        l += size;
        T tmp = e();
        if constexpr (LAZY)
            prop_topdown(l);
        do {
            while (l % 2 == 0)
                l >>= 1;
            T next = op(tmp, dat[l]);
            if (cond(next)) {
                ++l;
                tmp = next;
                continue;
            }
            while (l < size) {
                if constexpr (LAZY)
                    prop(l);
                l <<= 1;
                next = op(tmp, dat[l]);
                if (cond(next)) {
                    ++l;
                    tmp = next;
                }
            }
            return { l - size, tmp };
        } while ((l & -l) != l);
        return { size, tmp };
    }

    auto get_size() const {
        return size;
    }
};
template<typename T, typename M>
using lazy_segtree = segtree<T, true, M>;
template<typename T, typename M>
using segtree_beats = segtree<T, true, M, true>;



/*
Lucasの定理
pを素数として
mCn == m1Cn1 * m2Cn2 * ... (mod p)
ただしmi,ni は m, nのp進数表記のi桁目
Wilsonの定理
p > 2について
(p-1)! == p-1 (mod p)   <=>   pは素数
Primeswingみる
a == b → fabs(a - b) < EPS
a < b → a + EPS < b
a <= b → a < b + EPS
sum k = 0 to n   pCk * qC(n-k) = (p+q)Cn    (a < b -> aCb=0)
sum k = 0 to a   kCb = (a+1)C(b+1)
*/

//エラトステネスの篩。Mまで
std::vector<int> get_prime(int M) {
    std::vector<int> ans;
    std::vector<int> F(M + 1, 0);
    for (int i = 2; i <= M; ++i) {
        if (F[i])
            continue;
        ans.push_back(i);
        for (int j = i + i; j <= M; j += i)
            F[j] = 1;
    }
    return ans;
}
//素因数分解。√n以下の素数をPで渡す
std::vector<std::pair<long long, int>> factorize(long long n, const std::vector<int>& P) {
    std::vector<std::pair<long long, int>> ans;
    for (int p : P) {
        if (n < p)
            break;
        if (n % p)
            continue;
        ans.push_back({ p, 0 });
        do {
            ++ans.back().second;
            n /= p;
        } while (n % p == 0);
    }
    if (n != 1)
        ans.push_back({ n, 1 });
    return ans;
}
//mまでの各数の最小の素因数を格納O(m)
std::vector<int> min_factor(int m) {
    std::vector<int> ans(m + 1);
    std::iota(ans.begin(), ans.end(), 0);
    for (int i = 2; i * i <= m; ++i) {
        if (ans[i] < i)
            continue;
        for (int j = i * i; j <= m; j += i)
            if (ans[j] == j)
                ans[j] = i;
    }
    return ans;
}
//素因数分解。nまでの最小の素因数をmin_factorで渡す。O(logn)
std::vector<std::pair<long long, int>> factorize2(long long n, const std::vector<int>& min_fact) {
    std::vector<std::pair<long long, int>> ans;
    while (n > 1) {
        if (ans.empty() || min_fact[n] != ans.back().first)
            ans.push_back({ min_fact[n], 0 });
        ++ans.back().second;
        n /= min_fact[n];
    }
    return ans;
}
//素因数から約数列挙。
std::vector<long long> divisor(const std::vector<std::pair<long long, int>>& fact) {
    std::vector<long long> ans(1, 1);
    for (std::pair<long long, int> p : fact) {
        int last = ans.size();
        for (int i = 0; i < last; ++i) {
            long long tmp = p.first;
            for (int j = 0; j < p.second; ++j) {
                ans.push_back(tmp * ans[i]);
                tmp *= p.first;
            }
        }
    }
    return ans;
}

template<typename T>
T power(T r, long long n) {
    T ans = 1;
    T tmp = r;
    while (n > 0) {
        if ((n & 1) > 0) {
            ans = ans * tmp;
        }
        tmp = tmp * tmp;
        n >>= 1;
    }
    return ans;
}

template<typename T>
T power(T r, long long n, T e) {
    T ans = e;
    T tmp = r;
    while (n > 0) {
        if ((n & 1) > 0) {
            ans = ans * tmp;
        }
        tmp = tmp * tmp;
        n >>= 1;
    }
    return ans;
}

//切り捨て除算
template<typename T>
T floor(T a, T b) {
    if ((a >= 0 && b > 0) || (a <= 0 && b < 0))
        return abs(a) / abs(b);
    else {
        return -((abs(a) + abs(b) - 1) / abs(b));
    }
}

//切り上げ除算
template<typename T>
T ceil(T a, T b) {
    if ((a >= 0 && b > 0) || (a <= 0 && b < 0))
        return (abs(a) + abs(b) - 1) / abs(b);
    else {
        return (-abs(a)) / abs(b);
    }
}

//最大公約数
template<typename T>
T gcd(T a, T b) {
    if (b == 0)
        return a;
    return gcd(b, a % b);
}
//最小公倍数
template<typename T>
T lcm(T a, T b) {
    return a * b / gcd(a, b);
}

//拡張ユークリッドの互除法。戻り値は最大公約数、x,yはax+by=gcd(a,b)を満たす組の一つ
template<typename T>
T gcdext(T a, T b, T& x, T& y) {
    if (b == 0) {
        x = 1;
        y = 0;
        return a;
    }
    T g = gcdext(b, a % b, y, x);
    y -= a / b * x;
    return g;
}

template<typename T>
T modinv(T v, T mod) {
    T x, y;
    gcdext(v, mod, x, y);
    return (x % mod + mod) % mod;
}

template<typename T>
T modpow(T r, T n, T mod) {
    T ans = 1;
    T tmp = r % mod;
    while (n > 0) {
        if ((n & 1) > 0) {
            ans = ans * tmp % mod;
        }
        tmp = tmp * tmp % mod;
        n >>= 1;
    }
    return ans;
}

//素数判定O(log n)
bool is_prime(unsigned int n) {
    if (n == 1 || n % 2 == 0)
        return n == 2;
    unsigned int d = n - 1, s = 0;
    while (d % 2 == 0) {
        d /= 2;
        ++s;
    }
    constexpr std::array<unsigned int, 7> test = { 2, 325, 9375, 28178, 450775, 9780504, 1795265022ULL };
    for (auto a : test) {
        if (a % n == 0)
            return true;
        unsigned long long tmp = modpow<unsigned long long>(a, d, n);
        if (tmp == 1)
            continue;
        bool no_prime = true;
        for (int r = 0; r < (int)s; ++r) {
            if (tmp == n - 1) {
                no_prime = false;
                break;
            }
            tmp = tmp * tmp % n;
        }
        if (no_prime)
            return false;
    }
    return true;
}
/*
bool is_prime_ull(unsigned long long n) {
    if (n == 1 || n % 2 == 0) return n == 2;
    unsigned long long d = n - 1, s = 0;
    while (d % 2 == 0) {
        d /= 2;
        ++s;
    }
    constexpr std::array<unsigned long long, 7> test = {2, 325, 9375, 28178, 450775, 9780504, 1795265022ULL};
    for (auto a : test) {
        if (a % n == 0) return true;
        __uint128_t tmp = modpow<__uint128_t>(a, d, n);
        if (tmp == 1) continue;
        bool no_prime = true;
        for (int r = 0; r < s; ++r) {
            if (tmp == n - 1) {
                no_prime = false;
                break;
            }
            tmp = tmp * tmp % n;
        }
        if (no_prime) return false;
    }
    return true;
}*/

//x == v[i].first (mod v[i].second) を満たす最小のxを求める。mは互いに素
template<typename T>
T garner(const std::vector<std::pair<T, T>>& v) {
    int n = v.size();
    T   m_prod = 1;
    T   x = v[0].first % v[0].second;
    for (int i = 1; i < n; ++i) {
        m_prod *= v[i - 1].second;
        T t = ((v[i].first - x) * modinv<T>(m_prod, v[i].second)) % v[i].second;
        if (t < 0)
            t += v[i].second;
        x += t * m_prod;
    }
    return x;
}

//x == v[i].first (mod v[i].second) を、互に素にする。そもそも解がなければfalse
template<typename T>
bool garner_normalization(std::vector<std::pair<T, T>>& v) {
    using ll = long long;
    for (int i = 0; i < v.size(); ++i) {
        for (int j = i + 1; j < v.size(); ++j) {
            ll g = gcd(v[i].second, v[j].second);
            if ((v[i].first - v[j].first) % g != 0) {
                return false;
            }
            ll afilter = gcd(g, v[i].second / g);
            ll bfilter = gcd(g, v[j].second / g);
            ll gcur = g;
            ll a = 1, b = 1;
            while (true) {
                ll tmp = gcd(afilter, gcur);
                if (tmp == 1)
                    break;
                a *= tmp;
                gcur /= tmp;
            }
            while (true) {
                ll tmp = gcd(bfilter, gcur);
                if (tmp == 1)
                    break;
                b *= tmp;
                gcur /= tmp;
            }
            ll na = v[i].second / lcm(a, g) * a * gcur;
            ll nb = v[j].second / lcm(b, g) * b;
            v[i].second = na;
            v[j].second = nb;
            v[i].first %= na;
            v[j].first %= nb;
        }
    }
    return true;
}

// 等比数列 r^i (i = 0, 1,.., n)の和を計算, O(log(n))、公式が使えないときに
template<typename T>
T geometric_seq_sum(T r, long long n) {
    using mat22 = std::array<std::array<T, 2>, 2>;
    auto mul = [](mat22 a, mat22 b) {
        mat22 ans;
        for (int i = 0; i < 2; ++i)
            for (int j = 0; j < 2; ++j) {
                ans[i][j] = T(0);
                for (int k = 0; k < 2; ++k)
                    ans[i][j] += a[i][k] * b[k][j];
            }
        return ans;
        };

    mat22 mat, tmp;
    mat[0] = { T(1), T(0) };
    mat[1] = { T(0), T(1) };
    tmp[0] = { T(1), T(1) };
    tmp[1] = { T(0), r };

    while (n > 0) {
        if ((n & 1) > 0)
            mat = mul(tmp, mat);
        tmp = mul(tmp, tmp);
        n >>= 1;
    }
    return (mat[0][0] + mat[0][1] * r);
}

template<long long mod>
class modint {
private:
    long long gcdext(long long a, long long b, long long& x, long long& y) const {
        if (b == 0) {
            x = 1;
            y = 0;
            return a;
        }
        long long g = gcdext(b, a % b, y, x);
        y -= a / b * x;
        return g;
    }

public:
    long long v;

    modint() : v(0) {
    }
    modint(long long val) : v((val% mod + mod) % mod) {
    }

    operator long long() const {
        return v;
    }
    modint<mod> operator-() {
        return modint<mod>(v - mod);
    }

    void operator=(long long val) {
        v = (val % mod + mod) % mod;
    }
    void operator=(modint<mod> src) {
        v = src.v;
    }
    bool operator==(long long a) const {
        return v == a;
    }

    modint<mod> operator+(modint<mod> a) const {
        return modint<mod>(v + a.v);
    }
    modint<mod> operator-(modint<mod> a) const {
        return modint<mod>(v - a.v);
    }
    modint<mod> operator*(modint<mod> a) const {
        return modint<mod>(v * a.v);
    }
    modint<mod> operator+(long long a) const {
        return modint<mod>(v + a % mod);
    }
    modint<mod> operator-(long long a) const {
        return modint<mod>(v - a % mod);
    }
    modint<mod> operator*(long long a) const {
        return modint<mod>(v * (a % mod));
    }

    void operator+=(modint<mod> a) {
        v = (v + a.v) % mod;
    }
    void operator-=(modint<mod> a) {
        v = ((v - a.v) % mod + mod) % mod;
    }
    void operator*=(modint<mod> a) {
        v = v * a.v % mod;
    }
    void operator+=(long long a) {
        v = (v + modint<mod>(a).v) % mod;
    }
    void operator-=(long long a) {
        v = ((v - modint<mod>(a).v) % mod + mod) % mod;
    }
    void operator*=(long long a) {
        v = v * modint<mod>(a).v % mod;
    }

    // モジュラ逆数。vとmodが互いに素でないとき未定義
    modint<mod> inv() const {
        long long x, y;
        gcdext(v, mod, x, y);
        return modint<mod>(x);
    }

    //aとmodが互いに素でないとき未定義
    modint<mod> operator/(modint<mod> a) const {
        return v * a.inv();
    }
    modint<mod> operator/(long long a) const {
        return v * modint<mod>(a).inv();
    }

    //aとmodが互いに素でないとき未定義
    void operator/=(modint<mod> a) {
        v = v * a.inv() % mod;
    }
    void operator/=(long long a) {
        v = v * modint<mod>(a).inv() % mod;
    }

    modint<mod>& operator++() {
        v = (v + 1) % mod;
        return *this;
    }
    modint<mod>& operator--() {
        v = ((v - 1) % mod + mod) % mod;
        return *this;
    }
    modint<mod> operator++(int) {
        modint<mod> ret = *this;
        v = (v + 1) % mod;
        return ret;
    }
    modint<mod> operator--(int) {
        modint<mod> ret = *this;
        v = ((v - 1) % mod + mod) % mod;
        return ret;
    }
};
using mint7 = modint<1000000007>;
using mint998 = modint<998244353>;

struct unionfind {
    int              N;
    int              set_num;
    std::vector<int> par;
    std::vector<int> size;
    std::vector<int> buf;
    unionfind() : N{ 0 }, set_num{ 0 } {
    }
    unionfind(long long N_) : N((int)N_), set_num((int)N_), par(N), size(N), buf(N) {
        for (int i = 0; i < N; ++i) {
            par[i] = i;
            size[i] = 1;
        }
    }
    int get(int ind) {
        int cur = ind;
        buf.clear();
        while (cur != par[cur]) {
            buf.push_back(cur);
            cur = par[cur];
        }
        for (int i : buf)
            par[i] = cur;
        return cur;
    }
    int getsize(int ind) {
        int p = get(ind);
        return size[p];
    }
    void merge(int i, int j) {
        int pi = get(i);
        int pj = get(j);
        if (pi == pj)
            return;
        --set_num;
        if (getsize(pi) < getsize(pj)) {
            par[pi] = pj;
            size[pj] += size[pi];
        }
        else {
            par[pj] = pi;
            size[pi] += size[pj];
        }
    }
    bool same(int i, int j) {
        return get(i) == get(j);
    }
};

class Rational {

    long long _bunsi;
    long long _bunbo;

    void normalize() {
        if (_bunbo == 0) {
            std::cerr << "Rational: bunbo is 0" << std::endl;
            exit(1);
        }
        if (_bunbo < 0) {
            _bunsi *= -1;
            _bunbo *= -1;
        }
        long long g = gcd(abs(_bunsi), _bunbo);
        _bunsi /= g;
        _bunbo /= g;
    }

public:
    Rational() : _bunsi(0), _bunbo(1) {};
    Rational(long long a) : _bunsi(a), _bunbo(1) {};
    Rational(long long a, long long b) : _bunsi(a), _bunbo(b) {
        normalize();
    };
    Rational(const Rational& r) : _bunsi(r._bunsi), _bunbo(r._bunbo) {};
    Rational(Rational&& r) : _bunsi(r._bunsi), _bunbo(r._bunbo) {};

    long long bunsi() const {
        return _bunsi;
    }
    long long bunbo() const {
        return _bunbo;
    }

    Rational& operator=(const Rational& r) = default;
    Rational& operator=(Rational&& r) = default;

    Rational operator-() const {
        return Rational(-_bunsi, _bunbo);
    }

    Rational operator+(const Rational& r) const {
        return Rational(_bunsi * r._bunbo + r._bunsi * _bunbo, _bunbo * r._bunbo);
    }
    Rational operator-(const Rational& r) const {
        return Rational(_bunsi * r._bunbo - r._bunsi * _bunbo, _bunbo * r._bunbo);
    }
    Rational operator*(const Rational& r) const {
        return Rational(_bunsi * r._bunsi, _bunbo * r._bunbo);
    }
    Rational operator/(const Rational& r) const {
        return Rational(_bunsi * r._bunbo, _bunbo * r._bunsi);
    }
    template<typename I>
    Rational operator+(I a) const {
        return Rational(_bunsi + a * _bunbo, _bunbo);
    }
    template<typename I>
    Rational operator-(I a) const {
        return Rational(_bunsi - a * _bunbo, _bunbo);
    }
    template<typename I>
    Rational operator*(I a) const {
        return Rational(_bunsi * a, _bunbo);
    }
    template<typename I>
    Rational operator/(I a) const {
        return Rational(_bunsi, _bunbo * a);
    }

    Rational& operator+=(const Rational& r) {
        *this = *this + r;
        return *this;
    }
    Rational& operator-=(const Rational& r) {
        *this = *this - r;
        return *this;
    }
    Rational& operator*=(const Rational& r) {
        *this = *this * r;
        return *this;
    }
    Rational& operator/=(const Rational& r) {
        *this = *this / r;
        return *this;
    }
    template<typename I>
    Rational& operator+=(I a) {
        *this = *this + a;
        return *this;
    }
    template<typename I>
    Rational& operator-=(I a) {
        *this = *this - a;
        return *this;
    }
    template<typename I>
    Rational& operator*=(I a) {
        *this = *this * a;
        return *this;
    }
    template<typename I>
    Rational& operator/=(I a) {
        *this = *this / a;
        return *this;
    }

    bool operator==(const Rational& r) const {
        return _bunsi == r._bunsi && _bunbo == r._bunbo;
    }
    bool operator!=(const Rational& r) const {
        return _bunsi != r._bunsi || _bunbo != r._bunbo;
    }
    bool operator<(const Rational& r) const {
        return _bunsi * r._bunbo < r._bunsi * _bunbo;
    }
    bool operator>(const Rational& r) const {
        return _bunsi * r._bunbo > r._bunsi * _bunbo;
    }
    bool operator<=(const Rational& r) const {
        return _bunsi * r._bunbo <= r._bunsi * _bunbo;
    }
    bool operator>=(const Rational& r) const {
        return _bunsi * r._bunbo >= r._bunsi * _bunbo;
    }
};

namespace Geometry {

    template<typename T>
    struct NumberTrait {
        static bool isZero(const T& a) {
            return a == 0;
        }
    };

    template<>
    struct NumberTrait<Rational> {
        static bool isZero(const Rational& a) {
            return a == Rational(0);
        }
    };

    template<typename T>
    struct Point {
        T x;
        T y;
        Point() : x(0), y(0) {}
        Point(T x, T y) : x(x), y(y) {}
    };

    //ax + by + c = 0
    template<typename T>
    using Line = std::tuple<T, T, T>;

    template<typename T>
    std::optional<Point<T>> crossPoint(const Line<T>& l1, const Line<T>& l2) {
        auto [a1, b1, c1] = l1;
        auto [a2, b2, c2] = l2;

        T d = a1 * b2 - a2 * b1;
        if (NumberTrait<T>::isZero(d)) return std::nullopt;

        T x = (b1 * c2 - b2 * c1) / d;
        T y = (a2 * c1 - a1 * c2) / d;
        return Point<T>(x, y);
    }

    template<typename T>
    Line<T> lineThrough(const Point<T>& p1, const Point<T>& p2) {
        T a = p2.y - p1.y;
        T b = p1.x - p2.x;
        T c = -a * p1.x - b * p1.y;
        return Line<T>(a, b, c);
    }

    template<typename T>
    bool onLine(const Point<T>& p, const Line<T>& l) {
        auto [a, b, c] = l;
        return NumberTrait<T>::isZero(a * p.x + b * p.y + c);
    }
}

using namespace std;
using uint = unsigned int;
using ll = long long;
using ull = unsigned long long;
using lpair = std::pair<ll, ll>;
#define REP(i, a, b)    for(ll i = a; i < b; ++i)
#define REPREV(i, a, b) for(ll i = a; i > b; --i)

const int _ = []() {
    std::cout << std::fixed << std::setprecision(10);
    return 0;
    }();

template<typename A, typename B>
void chmin(A& a, const B& b) {
    a = min(a, static_cast<A>(b));
};
template<typename A, typename B>
void chmax(A& a, const B& b) {
    a = max(a, static_cast<A>(b));
};


int main() {

    vector<mint998> A[11];

    ll N, K; cin >> N >> K;
    A[0] = vector<mint998>(N);
    REP(i, 0, N) {
        A[0][i] = 1LL;
    }
    A[1] = vector<mint998>(N);
    REP(i, 0, N) {
        ll x; cin >> x;
        A[1][i] = x;
    }
    REP(k, 2, 11) {
        A[k] = vector<mint998>(N);
        REP(i, 0, N) {
            A[k][i] = A[k - 1][i] * A[1][i];
        }
    }

    com_table<mint998> com(100);

    vector<mint998> B[11];
    REP(k, 0, 11) {
        B[k] = vector<mint998>(N + 1, mint998(0));
        REP(i, 0, N) {
            B[k][i + 1] = A[k][i];
            REP(m, 0, k + 1) {
                B[k][i + 1] += com(k, m) * A[m][i] * B[k - m][i];
            }
        }
    }

    mint998 ans = 0LL;
    REP(i, 0, N) {
        ans += B[K][i + 1];
    }
    cout << ans << endl;

    return 0;
}

提出情報

提出日時
問題 F - Range Power Sum
ユーザ tarakojo1019
言語 C++ 20 (gcc 12.2)
得点 550
コード長 37605 Byte
結果 AC
実行時間 247 ms
メモリ 37540 KiB

コンパイルエラー

Main.cpp: In function ‘int main()’:
Main.cpp:1431:43: warning: implicitly-declared ‘constexpr modint<998244353>::modint(const modint<998244353>&)’ is deprecated [-Wdeprecated-copy]
 1431 |             A[k][i] = A[k - 1][i] * A[1][i];
      |                                           ^
Main.cpp:1068:10: note: because ‘modint<998244353>’ has user-provided ‘void modint<mod>::operator=(modint<mod>) [with long long int mod = 998244353]’
 1068 |     void operator=(modint<mod> src) {
      |          ^~~~~~~~
Main.cpp:1081:39: note:   initializing argument 1 of ‘modint<mod> modint<mod>::operator*(modint<mod>) const [with long long int mod = 998244353]’
 1081 |     modint<mod> operator*(modint<mod> a) const {
      |                           ~~~~~~~~~~~~^
Main.cpp:1441:33: warning: implicitly-declared ‘constexpr modint<998244353>::modint(const modint<998244353>&)’ is deprecated [-Wdeprecated-copy]
 1441 |             B[k][i + 1] = A[k][i];
      |                                 ^
Main.cpp:1068:10: note: because ‘modint<998244353>’ has user-provided ‘void modint<mod>::operator=(modint<mod>) [with long long int mod = 998244353]’
 1068 |     void operator=(modint<mod> src) {
      |          ^~~~~~~~
Main.cpp:1068:32: note:   initializing argument 1 of ‘void modint<mod>::operator=(modint<mod>) [with long long int mod = 998244353]’
 1068 |     void operator=(modint<mod> src) {
      |                    ~~~~~~~~~~~~^~~
Main.cpp:1443:50: warning: implicitly-declared ‘constexpr modint<998244353>::modint(const modint<998244353>&)’ is deprecated [-Wdeprecated-copy]
 1443 |                 B[k][i + 1] += com(k, m) * A[m][i] * B[k - m][i];
      |                                                  ^
Main.cpp:1068:10: note: because ‘modint<998244353>’ has user-provided ‘void modint<mod>::operator=(modint<mod>) [with long long int mod = 998244353]’
 1068 |     void operator=(modint<mod> src) {
      |          ^~~~~~~~
Main.cpp:1081:39: note:   initializing argument 1 of ‘modint<mod> modint<mod>::operator*(modint<mod>) const [with long long int mod = 998244353]’
 1081 |     modint<mod> operator*(modint<mod> a) const {
      |                           ~~~~~~~~~~~~^
Main.cpp:1443:64: warning: implicitly-declared ‘constexpr modint<998244353>::modint(const modint<998244353>&)’ is deprecated [-Wdeprecated-copy]
 1443 |                 B[k][i + 1] += com(k, m) * A[m][i] * B[k - m][i];
      |                                                                ^
Main.cpp:1068:10: note: because ‘modint<998244353>’ has user-provided ‘void modint<mod>::operator=(modint<mod>) [with long long int mod = 998244353]’
 1068 |     void operator=(modint<mod> src) {
      |          ^~~~~~~~
Main.cpp:1081:39: note:   initializing argument 1 of ‘modint<mod> modint<mod>::operator*(modint<mod>) const [with long long int mod = 998244353]’
 1081 |     modint<mod> operator*(modint<mod> a) const {
      |                           ~~~~~~~~~~~~^
Main.cpp:1450:26: warning: implicitly-declared ‘constexpr modint<998244353>::modint(const modint<998244353>&)’ is deprecated [-Wdeprecated-copy]
 1450 |         ans += B[K][i + 1];
      |                          ^
Main.cpp:1068:10: note: because ‘modint<998244353>’ has user-provided ‘void modint<mod>::operator=(modint<mod>) [with long long int mod = 998244353]’
 1068 |     void operator=(modint<mod> src) {
      |          ^~~~~~~~
Main.cpp:1094:33: note:   initializing argument 1 of ‘void modint<mod>::operator+=(modint<mod>) [with long long int mod = 998244353]’
 1094 |     void operator+=(modint<mod> a) {
      |                     ~~~~~~~~~~~~^
Main.cpp: In instantiation of ‘T com_table<T>::operator()(int, int) [with T = modint<998244353>]’:
Main.cpp:1443:35:   required from here
Main.cpp:401:24: warning: implicitly-declared ‘constexpr modint<998244353>::modint(const modint<998244353>&)’ is deprecated [-Wdeprecated-copy]
  401 |         return frac[a] * finv[b] * finv[a - b];
      |                ~~~~~~~~^~~~~~~
Main.cpp:1068:10: note: because ‘modint<998244353>’ has user-provided ‘void modint<mod>::operator=(modint<mod>) [with long long int mod = 998244353]’
 1068 |     void operator=(modint<mod> src) {
      |          ^~~~~~~~
Main.cpp:1081:39: note:   initializing argument 1 of ‘modint<mod> modint<mod>::operator*(modint<mod>) const [with long long int mod = 998244353]’
 1081 |     modint<mod> operator*(modint<mod> a) const {
      |                           ~~~~~~~~~~~~^
Main.cpp:401:34: warning: implicitly-declared ‘constexpr modint<998244353>::modint(const modint<998244353>&)’ is deprecated [-Wdeprecated-copy]
  401 |         return frac[a] * finv[b] * finv[a - b];
      |                ~~~~~~~~~~~~~~~~~~^~~~~~~
Main.cpp:1068:10: note: because ‘modint<998244353>’ has user-provided ‘void modint<mod>::operator=(modint<mod>) [with long long int mod = 998244353]’
 1068 |     void operator=(modint<mod> src) {
      |          ^~~~~~~~
Main.cpp:1081:39: note:   initializing argument 1 of ‘modint<mod> modint<mod>::operator*(modint<mod>) const [with long long int mod = 998244353]’
 1081 |     modint<mod> operator*(modint<mod> a) const {
      |                           ~~~~~~~~~~~~^
Main.cpp: In instantiation of ‘void com_table<T>::init(long long int) [with T = modint<998244353>]’:
Main.cpp:407:9:   required from ‘com_table<T>::com_table(long long int) [with T = modint<998244353>]’
Main.cpp:1435:31:   required from here
Main.cpp:394:27: warning: implicitly-declared ‘constexpr modint<998244353>::modint(const modint<998244353>&)’ is deprecated [-Wdeprecated-copy]
  394 |         finv[size] = T(1) / frac[size];
      |                      ~~~~~^~~~~~~
Main.cpp:1068:10: note: because ‘modint<998244353>’ has user-provided ‘void modint<mod>::operator=(modint<mod>) [with long long int mod = 998244353]’
 1068 |     void operator=(modint<mod> src) {
      |          ^~~~~~~~
Main.cpp:1121:39: note:   initializing argument 1 of ‘modint<mod> modint<mod>::operator/(modint<mod>) const [with long long int mod = 998244353]’
 1121 |     modint<mod> operator/(modint<mod> a) const {
      |                           ~~~~~~~~~~~~^

ジャッジ結果

セット名 Sample All
得点 / 配点 0 / 0 550 / 550
結果
AC × 3
AC × 44
セット名 テストケース
Sample 00_sample_00.txt, 00_sample_01.txt, 00_sample_02.txt
All 00_sample_00.txt, 00_sample_01.txt, 00_sample_02.txt, 01_random_00.txt, 01_random_01.txt, 01_random_02.txt, 01_random_03.txt, 01_random_04.txt, 01_random_05.txt, 01_random_06.txt, 01_random_07.txt, 01_random_08.txt, 01_random_09.txt, 01_random_10.txt, 01_random_11.txt, 01_random_12.txt, 01_random_13.txt, 01_random_14.txt, 01_random_15.txt, 01_random_16.txt, 01_random_17.txt, 01_random_18.txt, 01_random_19.txt, 01_random_20.txt, 01_random_21.txt, 01_random_22.txt, 01_random_23.txt, 01_random_24.txt, 01_random_25.txt, 01_random_26.txt, 01_random_27.txt, 01_random_28.txt, 01_random_29.txt, 02_random2_00.txt, 02_random2_01.txt, 02_random2_02.txt, 03_minmax_00.txt, 03_minmax_01.txt, 03_minmax_02.txt, 03_minmax_03.txt, 03_minmax_04.txt, 03_minmax_05.txt, 03_minmax_06.txt, 03_minmax_07.txt
ケース名 結果 実行時間 メモリ
00_sample_00.txt AC 1 ms 3524 KiB
00_sample_01.txt AC 1 ms 3464 KiB
00_sample_02.txt AC 1 ms 3388 KiB
01_random_00.txt AC 10 ms 4800 KiB
01_random_01.txt AC 247 ms 37396 KiB
01_random_02.txt AC 245 ms 37464 KiB
01_random_03.txt AC 181 ms 28212 KiB
01_random_04.txt AC 244 ms 37448 KiB
01_random_05.txt AC 245 ms 37396 KiB
01_random_06.txt AC 169 ms 26648 KiB
01_random_07.txt AC 244 ms 37304 KiB
01_random_08.txt AC 244 ms 37312 KiB
01_random_09.txt AC 165 ms 26116 KiB
01_random_10.txt AC 244 ms 37312 KiB
01_random_11.txt AC 245 ms 37384 KiB
01_random_12.txt AC 31 ms 7260 KiB
01_random_13.txt AC 246 ms 37444 KiB
01_random_14.txt AC 246 ms 37472 KiB
01_random_15.txt AC 115 ms 19176 KiB
01_random_16.txt AC 245 ms 37416 KiB
01_random_17.txt AC 246 ms 37380 KiB
01_random_18.txt AC 137 ms 22088 KiB
01_random_19.txt AC 246 ms 37540 KiB
01_random_20.txt AC 245 ms 37308 KiB
01_random_21.txt AC 121 ms 19956 KiB
01_random_22.txt AC 245 ms 37532 KiB
01_random_23.txt AC 246 ms 37532 KiB
01_random_24.txt AC 22 ms 6160 KiB
01_random_25.txt AC 245 ms 37448 KiB
01_random_26.txt AC 245 ms 37376 KiB
01_random_27.txt AC 170 ms 26864 KiB
01_random_28.txt AC 245 ms 37308 KiB
01_random_29.txt AC 245 ms 37404 KiB
02_random2_00.txt AC 246 ms 37392 KiB
02_random2_01.txt AC 245 ms 37308 KiB
02_random2_02.txt AC 245 ms 37336 KiB
03_minmax_00.txt AC 1 ms 3464 KiB
03_minmax_01.txt AC 1 ms 3456 KiB
03_minmax_02.txt AC 1 ms 3460 KiB
03_minmax_03.txt AC 1 ms 3524 KiB
03_minmax_04.txt AC 215 ms 37460 KiB
03_minmax_05.txt AC 245 ms 37400 KiB
03_minmax_06.txt AC 213 ms 37396 KiB
03_minmax_07.txt AC 246 ms 37380 KiB