polynomial/test/pow_of_formal_power_series_acl.test.cpp

• View this file on GitHub
• Last update: 2022-08-25 20:03:29+09:00
• Problem: https://judge.yosupo.jp/problem/pow_of_formal_power_series

Depends on

• polynomial/fps_exp.hpp
• polynomial/fps_inv.hpp
• polynomial/fps_log.hpp
• polynomial/fps_pow.hpp
• polynomial/polynomial.hpp
• polynomial/with_acl.hpp
• template/template.hpp

Code

#define PROBLEM "https://judge.yosupo.jp/problem/pow_of_formal_power_series"

#define FAST_IO

#include "../../template/template.hpp"
#include "../../polynomial/with_acl.hpp"
#include "../../polynomial/fps_pow.hpp"

using Mint = atcoder::modint998244353;

istream &operator>>(istream &is, Mint &x) {
    i32 val;
    is >> val;
    x = Mint(val);
    return is;
}
ostream &operator<<(ostream &os, Mint x) {
    os << x.val();
    return os;
}

int main() {
    using FPS = Polynomial<Mint, ACLNTT<998244353>>;
    
    i32 n;
    i64 m;
    cin >> n >> m;
    FPS f(n);
    REP(i, n) {
        cin >> f[i];
    }
    FPS g = fps_pow(f, m);
    REP(i, n) {
        cout << g[i] << " \n"[i + 1 == n];
    }
}

#line 1 "polynomial/test/pow_of_formal_power_series_acl.test.cpp"
#define PROBLEM "https://judge.yosupo.jp/problem/pow_of_formal_power_series"

#define FAST_IO

#line 1 "template/template.hpp"
#include <algorithm>
#include <array>
#include <bitset>
#include <cassert>
#include <cmath>
#include <iomanip>
#include <iostream>
#include <list>
#include <map>
#include <numeric>
#include <queue>
#include <random>
#include <set>
#include <stack>
#include <string>
#include <tuple>
#include <unordered_map>
#include <unordered_set>
#include <utility>
#include <vector>

#define OVERRIDE(a, b, c, d, ...) d
#define REP2(i, n) for (i32 i = 0; i < (i32) (n); ++i)
#define REP3(i, m, n) for (i32 i = (i32) (m); i < (i32) (n); ++i)
#define REP(...) OVERRIDE(__VA_ARGS__, REP3, REP2)(__VA_ARGS__)
#define PER(i, n) for (i32 i = (i32) (n) - 1; i >= 0; --i)
#define ALL(x) begin(x), end(x)

using namespace std;

using u32 = unsigned int;
using u64 = unsigned long long;
using u128 = __uint128_t;
using i32 = signed int;
using i64 = signed long long;
using i128 = __int128_t;
using f64 = double;
using f80 = long double;

template <typename T>
using Vec = vector<T>;

template <typename T>
bool chmin(T &x, const T &y) {
    if (x > y) {
        x = y;
        return true;
    }
    return false;
}
template <typename T>
bool chmax(T &x, const T &y) {
    if (x < y) {
        x = y;
        return true;
    }
    return false;
}

istream &operator>>(istream &is, i128 &x) {
    i64 v;
    is >> v;
    x = v;
    return is;
}
ostream &operator<<(ostream &os, i128 x) {
    os << (i64) x;
    return os;
}
istream &operator>>(istream &is, u128 &x) {
    u64 v;
    is >> v;
    x = v;
    return is;
}
ostream &operator<<(ostream &os, u128 x) {
    os << (u64) x;
    return os;
}

[[maybe_unused]] constexpr i32 INF = 1000000100;
[[maybe_unused]] constexpr i64 INF64 = 3000000000000000100;
struct SetUpIO {
    SetUpIO() {
#ifdef FAST_IO
        ios::sync_with_stdio(false);
        cin.tie(nullptr);
#endif
        cout << fixed << setprecision(15);
    }
} set_up_io;
#line 2 "polynomial/with_acl.hpp"

#include <atcoder/convolution>

template <int mod>
class ACLNTT {
public:
    using Mint = atcoder::static_modint<mod>;
    
    static void dft(std::vector<Mint> &a) {
        atcoder::internal::butterfly(a);
    }
    static void idft(std::vector<Mint> &a) {
        atcoder::internal::butterfly_inv(a);
        Mint inv = Mint::raw(a.size()).inv();
        for (Mint &ele : a) {
            ele *= inv;
        }
    }
    static std::vector<Mint> mul(std::vector<Mint> a, std::vector<Mint> b) {
        return atcoder::convolution(a, b);
    }
};
#line 2 "polynomial/fps_pow.hpp"

#line 2 "polynomial/fps_log.hpp"

#line 2 "polynomial/fps_inv.hpp"

#line 2 "polynomial/polynomial.hpp"

#line 7 "polynomial/polynomial.hpp"

template <typename T, typename Mul>
class Polynomial {
    std::vector<T> coeff;
    
public:
    using This = Polynomial<T, Mul>;
    
    Polynomial() : coeff() {}
    Polynomial(int n) : coeff(n, T(0)) {}
    Polynomial(std::vector<T> c) : coeff(std::move(c)) {}
    
    const std::vector<T> &vec() const {
        return coeff;
    }
    
    int size() const {
        return (int) coeff.size();
    }
    
    const T &operator[](int idx) const {
        return coeff[idx];
    }
    T &operator[](int idx) {
        return coeff[idx];
    }
    
    T at(int idx) const {
        if (idx < size()) {
            return coeff[idx];
        } else {
            return T(0);
        }
    }
    
    void pre_(int n) {
        assert(n >= 0);
        coeff.resize(n, T(0));
    }
    This pre(int n) const {
        This tmp(*this);
        tmp.pre_(n);
        return tmp;
    }
    
    T operator()(const T &x) const {
        T p(1), sum(0);
        for (const T &ele : coeff) {
            sum += p * ele;
            p *= x;
        }
        return sum;
    }
    
    This &operator+=(const This &rhs) {
        if (coeff.size() < rhs.coeff.size()) {
            coeff.resize(rhs.coeff.size(), T(0));
        }
        for (int i = 0; i < (int) rhs.coeff.size(); ++i) {
            coeff[i] += rhs.coeff[i];
        }
        return *this;
    }
    friend This operator+(This lhs, const This &rhs) {
        lhs += rhs;
        return lhs;
    }
    This &operator-=(const This &rhs) {
        if (coeff.size() < rhs.coeff.size()) {
            coeff.resize(rhs.coeff.size(), T(0));
        }
        for (int i = 0; i < (int) rhs.coeff.size(); ++i) {
            coeff[i] -= rhs.coeff[i];
        }
        return *this;
    }
    friend This operator-(This lhs, const This &rhs) {
        lhs -= rhs;
        return lhs;
    }
    
    This &operator*=(This rhs) {
        coeff = Mul::mul(std::move(coeff), std::move(rhs.coeff));
        return *this;
    }
    friend This operator*(This lhs, This rhs) {
        return This(Mul::mul(std::move(lhs.coeff), std::move(rhs.coeff)));
    }
    
    This diff() const {
        if (coeff.empty()) {
            return This();
        }
        std::vector<T> c(coeff.size() - 1);
        for (int i = 0; i < (int) c.size(); ++i) {
            c[i] = T(i + 1) * coeff[i + 1];
        }
        return This(c);
    }
    This integ() const {
        std::vector<T> c(coeff.size() + 1, T(0));
        for (int i = 0; i < (int) coeff.size(); ++i) {
            c[i + 1] = coeff[i] / T(i + 1);
        }
        return This(c);
    }
};
#line 4 "polynomial/fps_inv.hpp"

template <typename T, typename Mul>
Polynomial<T, Mul> fps_inv(const Polynomial<T, Mul> &f, int sz = -1) {
    const std::vector<T> &coeff = f.vec();
    assert(!coeff.empty() && coeff[0] != T(0));
    if (sz == -1) {
        sz = (int) coeff.size();
    }
    assert(sz >= 0);
    std::vector<T> g({T(1) / coeff[0]});
    while ((int) g.size() < sz) {
        std::vector<T> fg;
        if (2 * g.size() <= coeff.size()) {
            fg = std::vector<T>(coeff.begin(), coeff.begin() + 2 * g.size());
        } else {
            fg = coeff;
            fg.resize(2 * g.size());
        }
        Mul::dft(fg);
        std::vector<T> dft_g = g;
        dft_g.resize(2 * g.size(), T(0));
        Mul::dft(dft_g);
        for (int i = 0; i < 2 * (int) g.size(); ++i) {
            fg[i] *= dft_g[i];
        }
        Mul::idft(fg);
        std::fill(fg.begin(), fg.begin() + g.size(), T(0));
        Mul::dft(fg);
        for (int i = 0; i < 2 * (int) g.size(); ++i) {
            fg[i] *= dft_g[i];
        }
        Mul::idft(fg);
        g.resize(2 * g.size());
        for (int i = (int) g.size() / 2; i < (int) g.size(); ++i) {
            g[i] = -fg[i];
        }
    }
    g.resize(sz);
    return Polynomial<T, Mul>(g);
}
#line 4 "polynomial/fps_log.hpp"

template <typename T, typename Mul>
Polynomial<T, Mul> fps_log(const Polynomial<T, Mul> &f, int sz = -1) {
    const std::vector<T> &coeff = f.vec();
    assert(!coeff.empty() && coeff[0] == T(1));
    if (sz == -1) {
        sz = (int) coeff.size();
    }
    assert(sz >= 0);
    if (sz == 0) {
        return Polynomial<T, Mul>();
    }
    return (f.diff().pre(sz - 1) * fps_inv(f, sz - 1)).pre(sz - 1).integ();
}
#line 2 "polynomial/fps_exp.hpp"

#line 4 "polynomial/fps_exp.hpp"

template <typename T, typename Mul>
Polynomial<T, Mul> fps_exp(const Polynomial<T, Mul> &h, int sz = -1) {
    const std::vector<T> &coeff = h.vec();
    assert(!coeff.empty() && coeff[0] == T(0));
    if (sz == -1) {
        sz = (int) coeff.size();
    }
    assert(sz >= 0);
    std::vector<T> f({T(1)});
    std::vector<T> g({T(1)});
    std::vector<T> dft_f_({T(1), T(1)});
    
    while ((int) f.size() < sz) {
        int n = (int) f.size();
        
        // F_{2n}(g_0)
        std::vector<T> dft_g_2 = g;
        dft_g_2.resize(2 * n, T(0));
        Mul::dft(dft_g_2);
        
        // \delta
        std::vector<T> delta(n, T(0));
        for (int i = 0; i < n; ++i) {
            delta[i] = dft_f_[i] * dft_g_2[i];
        }
        Mul::idft(delta);
        delta.resize(2 * n);
        for (int i = 0; i < n; ++i) {
            std::swap(delta[i], delta[n + i]);
        }
        delta[n] -= T(1);
        
        // F_n(D(f_0))
        std::vector<T> dft_d_f(n, T(0));
        for (int i = 0; i < n - 1; ++i) {
            dft_d_f[i] = T(i + 1) * f[i + 1];
        }
        Mul::dft(dft_d_f);
        
        // D(f_0) g_0
        std::vector<T> d_f_g(n, T(0));
        for (int i = 0; i < n; ++i) {
            d_f_g[i] = dft_d_f[i] * dft_g_2[i];
        }
        Mul::idft(d_f_g);
        d_f_g.resize(2 * n, T(0));
        for (int i = 0; i < n - 1; ++i) {
            T tmp = T(i + 1) * h.at(i + 1);
            d_f_g[n + i] = d_f_g[i] - tmp;
            d_f_g[i] = tmp;
        }
        
        // \delta D(h_0)
        std::vector<T> dft_delta = delta;
        Mul::dft(dft_delta);
        std::vector<T> delta_d_h(2 * n);
        for (int i = 0; i < n - 1; ++i) {
            delta_d_h[i] = T(i + 1) * h.at(i + 1);
        }
        Mul::dft(delta_d_h);
        for (int i = 0; i < 2 * n; ++i) {
            delta_d_h[i] *= dft_delta[i];
        }
        Mul::idft(delta_d_h);
        std::fill(delta_d_h.begin(), delta_d_h.begin() + n, T(0));
        
        // \epsilon
        std::vector<T> eps = std::move(d_f_g);
        for (int i = 0; i < 2 * n; ++i) {
            eps[i] -= T(i + 1) * h.at(i + 1) + delta_d_h[i];
        }
        for (int i = 2 * n - 1; i > 0; --i) {
            eps[i] = eps[i - 1] / T(i);
        }
        eps[0] = T(0);
        
        // \epsilon f_0
        std::vector<T> dft_eps = eps;
        Mul::dft(dft_eps);
        std::vector<T> eps_f(2 * n);
        for (int i = 0; i < 2 * n; ++i) {
            eps_f[i] = dft_eps[i] * dft_f_[i];
        }
        Mul::idft(eps_f);
        std::fill(eps_f.begin(), eps_f.begin() + n - 1, T(0));
        
        // update f
        f.resize(2 * n, T(0));
        for (int i = 0; i < 2 * n; ++i) {
            f[i] -= eps_f[i];
        }
        
        if ((int) f.size() >= sz) {
            break;
        }
        
        // update F_{2n}(f)
        dft_f_ = f;
        dft_f_.resize(4 * n);
        Mul::dft(dft_f_);
        
        // update g
        std::vector<T> fg(dft_f_.begin(), dft_f_.begin() + 2 * n);
        for (int i = 0; i < 2 * n; ++i) {
            fg[i] *= dft_g_2[i];
        }
        Mul::idft(fg);
        std::fill(fg.begin(), fg.begin() + n, T(0));
        Mul::dft(fg);
        for (int i = 0; i < 2 * n; ++i) {
            fg[i] *= dft_g_2[i];
        }
        Mul::idft(fg);
        g.resize(2 * n);
        for (int i = n; i < 2 * n; ++i) {
            g[i] = -fg[i];
        }
    }
    
    f.resize(sz);
    return Polynomial<T, Mul>(f);
}
#line 5 "polynomial/fps_pow.hpp"

template <typename T, typename Mul>
Polynomial<T, Mul> fps_pow(const Polynomial<T, Mul> &h, unsigned long long m, int sz = -1) {
    const std::vector<T> &coeff = h.vec();
    if (sz == -1) {
        sz = (int) coeff.size();
    }
    assert(sz >= 0);
    
    if (m == 0) {
        std::vector<T> a(sz, T(0));
        if (sz >= 1) {
            a[0] = T(1);
        }
        return Polynomial<T, Mul>(a);
    }
    
    int ord = (int) coeff.size();
    for (int i = 0; i < (int) coeff.size(); ++i) {
        if (coeff[i] != T(0)) {
            ord = i;
            break;
        }
    }
    if (ord == (int) coeff.size() || (unsigned long long) ord >= ((unsigned long long) sz + m - 1) / m) {
        return Polynomial<T, Mul>(sz);
    }
    int zero = (int) (ord * m);
    int nonzero = sz - zero;
    
    Polynomial<T, Mul> f(std::vector<T>(coeff.begin() + ord, coeff.end()));
    T cf = f[0];
    T cf_inv = cf.inv();
    for (int i = 0; i < f.size(); ++i) {
        f[i] *= cf_inv;
    }
    f = fps_log(f, nonzero); //f.log(nonzero);
    T tm = T(m);
    for (int i = 0; i < f.size(); ++i) {
        f[i] *= tm;
    }
    f = fps_exp(f);
    T cf_m = cf.pow(m);
    for (int i = 0; i < f.size(); ++i) {
        f[i] *= cf_m;
    }
    std::vector<T> ans = f.vec();
    ans.insert(ans.begin(), zero, T(0));
    return Polynomial<T, Mul>(ans);
}
#line 8 "polynomial/test/pow_of_formal_power_series_acl.test.cpp"

using Mint = atcoder::modint998244353;

istream &operator>>(istream &is, Mint &x) {
    i32 val;
    is >> val;
    x = Mint(val);
    return is;
}
ostream &operator<<(ostream &os, Mint x) {
    os << x.val();
    return os;
}

int main() {
    using FPS = Polynomial<Mint, ACLNTT<998244353>>;
    
    i32 n;
    i64 m;
    cin >> n >> m;
    FPS f(n);
    REP(i, n) {
        cin >> f[i];
    }
    FPS g = fps_pow(f, m);
    REP(i, n) {
        cout << g[i] << " \n"[i + 1 == n];
    }
}

Back to top page