CP-Library

This documentation is automatically generated by online-judge-tools/verification-helper

View the Project on GitHub joji0/CP-Library

:heavy_check_mark: test/1_library_checker/math/convolution/Convolution107.test.cpp

Depends on

Code

// 我思う、故に我在り
#define PROBLEM "https://judge.yosupo.jp/problem/convolution_mod_1000000007"
#include <bits/stdc++.h>
#include "math/convolution/Convolution.hpp"

void solve() {
        int N, M; std::cin >> N >> M;
        std::vector<modint107> A(N), B(M);
        for (auto &x : A) std::cin >> x;
        for (auto &x : B) std::cin >> x;
        std::vector<modint107> C = Convolution::convolution(A, B);
        for (auto &x : C) std::cout << x << " ";
        std::cout << '\n';
}

int main() {
        std::ios::sync_with_stdio(false); std::cin.tie(nullptr);
        int tc = 1;
        // std::cin >> tc;
        while (tc--) solve();
        return 0;
}
#line 1 "test/1_library_checker/math/convolution/Convolution107.test.cpp"
// 我思う、故に我在り
#define PROBLEM "https://judge.yosupo.jp/problem/convolution_mod_1000000007"
#include <bits/stdc++.h>
#line 4 "math/convolution/Convolution.hpp"
#include <type_traits>
#line 4 "math/Barrett.hpp"

// reference: https://maspypy.github.io/library/mod/barrett.hpp
struct Barrett {
        uint32_t mod;
        uint64_t inv_mod;
        explicit Barrett(uint32_t m = 1) : mod(m), inv_mod((uint64_t)(-1) / m + 1) {}
        uint32_t umod() const { return mod; }
        uint32_t modulo(uint64_t z) const {
                if (mod == 1) return 0;
                uint64_t x = (uint64_t)(((__uint128_t)(z)*inv_mod) >> 64);
                uint64_t y = x * mod;
                return (uint32_t)(z - y + (z < y ? mod : 0));
        }
        uint64_t floor(uint64_t z) const {
                if (mod == 1) return z;
                uint64_t x = (uint64_t)(((__uint128_t)(z)*inv_mod) >> 64);
                uint64_t y = x * mod;
                return (z < y ? x - 1 : x);
        }
        std::pair<uint64_t, uint32_t> divmod(uint64_t z) const {
                if (mod == 1) return {z, 0};
                uint64_t x = (uint64_t)(((__uint128_t)(z)*inv_mod) >> 64);
                uint64_t y = x * mod;
                if (z < y) return {x - 1, z - y + mod};
                return {x, z - y};
        }
        uint32_t mul(uint32_t a, uint32_t b) const { return modulo((uint64_t)(a)*b); }
};
#line 8 "math/ModInt.hpp"

template <uint32_t Mod>
class StaticModInt {
        static_assert(Mod > 0, "Modulus must be positive");
        static_assert(Mod <= std::numeric_limits<uint32_t>::max() / 2,
                        "Modulus * 2 must be <= uint32_t::max()");
private:
        uint32_t val;
        static constexpr uint32_t inv1000000007[] = {
                0, 1, 500000004, 333333336, 250000002, 400000003,
                166666668, 142857144, 125000001, 111111112, 700000005
        };
        static constexpr uint32_t inv998244353[] = {
                0, 1, 499122177, 332748118, 748683265, 598946612,
                166374059, 855638017, 873463809, 443664157, 299473306
        };
        static constexpr uint32_t mod_inv(int64_t a) {
                int64_t b = Mod, x = 1, u = 0, t = 0, tmp = 0;
                while (b) {
                        t = a / b;
                        tmp = a - t * b;
                        a = b;
                        b = tmp;
                        tmp = x - t * u;
                        x = u;
                        u = tmp;
                }
                if (x < 0) x += Mod;
                return (uint32_t)x;
        }
        static constexpr uint32_t mod_reduce(int64_t v) {
                int64_t x = v % (int64_t)Mod;
                if (x < 0) x += Mod;
                return (uint32_t)x;
        }
public:
        constexpr StaticModInt() : val(0) {}
        template <class U, std::enable_if_t<std::is_integral_v<U> && std::is_signed_v<U>, int> = 0>
        constexpr StaticModInt(U v) : val(mod_reduce((int64_t)v)) {}
        template <class U, std::enable_if_t<std::is_integral_v<U> && std::is_unsigned_v<U>, int> = 0>
        constexpr StaticModInt(U v) : val((uint32_t)(v % Mod)) {}
        constexpr uint32_t get() const { return val; }
        static constexpr uint32_t get_mod() { return Mod; }
        static constexpr StaticModInt raw(uint32_t v) {
                StaticModInt res;
                res.val = v;
                return res;
        }
        constexpr StaticModInt inv() const {
                if constexpr (Mod == 1000000007) {
                        if (val <= 10) return raw(inv1000000007[val]);
                } else if constexpr (Mod == 998244353) {
                        if (val <= 10) return raw(inv998244353[val]);
                }
                return raw(mod_inv(val));
        }
        constexpr StaticModInt& operator++() {
                if (++val == Mod) val = 0;
                return *this;
        }
        constexpr StaticModInt operator++(int) {
                StaticModInt res = *this;
                ++*this;
                return res;
        }
        constexpr StaticModInt& operator--() {
                if (val == 0) val = Mod;
                --val;
                return *this;
        }
        constexpr StaticModInt operator--(int) {
                StaticModInt res = *this;
                --*this;
                return res;
        }
        constexpr StaticModInt& operator+=(const StaticModInt& other) {
                val += other.val;
                if (val >= Mod) val -= Mod;
                return *this;
        }
        constexpr StaticModInt& operator-=(const StaticModInt& other) {
                if (val < other.val) val += Mod;
                val -= other.val;
                return *this;
        }
        constexpr StaticModInt& operator*=(const StaticModInt& other) {
                val = (uint32_t)((uint64_t)val * other.val % Mod);
                return *this;
        }
        constexpr StaticModInt& operator/=(const StaticModInt& other) {
                *this *= other.inv();
                return *this;
        }
        friend constexpr StaticModInt operator+(const StaticModInt& lhs, const StaticModInt& rhs) { return StaticModInt(lhs) += rhs; }
        friend constexpr StaticModInt operator-(const StaticModInt& lhs, const StaticModInt& rhs) { return StaticModInt(lhs) -= rhs; }
        friend constexpr StaticModInt operator*(const StaticModInt& lhs, const StaticModInt& rhs) { return StaticModInt(lhs) *= rhs; }
        friend constexpr StaticModInt operator/(const StaticModInt& lhs, const StaticModInt& rhs) { return StaticModInt(lhs) /= rhs; }
        constexpr StaticModInt operator+() const { return *this; }
        constexpr StaticModInt operator-() const { return StaticModInt() - *this; }
        friend constexpr bool operator==(const StaticModInt& lhs, const StaticModInt& rhs) { return lhs.val == rhs.val; }
        friend constexpr bool operator!=(const StaticModInt& lhs, const StaticModInt& rhs) { return lhs.val != rhs.val; }
        constexpr StaticModInt pow(int64_t a) const {
                StaticModInt v = *this, res = 1;
                while (a) {
                        if (a & 1) res *= v;
                        a >>= 1;
                        v *= v;
                }
                return res;
        }
        friend std::ostream& operator<<(std::ostream& os, const StaticModInt& m) { return os << m.val; }
        friend std::istream& operator>>(std::istream& is, StaticModInt& m) {
                int64_t v; is >> v;
                m = StaticModInt(v);
                return is;
        }
};

using modint107 = StaticModInt<1000000007>;
using modint998 = StaticModInt<998244353>;

template <int ID>
class DynamicModInt {
private:
        uint32_t val;
        static Barrett bt;
        static uint32_t mod_inv(int64_t a) {
                int64_t b = bt.mod, x = 1, u = 0, t = 0, tmp = 0;
                while (b) {
                        t = a / b;
                        tmp = a - t * b;
                        a = b;
                        b = tmp;
                        tmp = x - t * u;
                        x = u;
                        u = tmp;
                }
                if (x < 0) x += bt.mod;
                return (uint32_t)x;
        }
public:
        DynamicModInt() : val(0) {}
        template <class U, std::enable_if_t<std::is_integral_v<U> && std::is_signed_v<U>, int> = 0>
        DynamicModInt(U v) {
                int64_t x = (int64_t)(v % (int64_t)bt.mod);
                if (x < 0) x += bt.mod;
                val = (uint32_t)x;
        }
        template <class U, std::enable_if_t<std::is_integral_v<U> && std::is_unsigned_v<U>, int> = 0>
        DynamicModInt(U v) : val((uint32_t)(v % bt.mod)) {}
        uint32_t get() const { return val; }
        static uint32_t get_mod() { return bt.mod; }
        static void set_mod(uint32_t m) {
                assert(m > 0 && m <= std::numeric_limits<uint32_t>::max() / 2);
                bt = Barrett(m);
        }
        static DynamicModInt raw(uint32_t v) {
                DynamicModInt res;
                res.val = v;
                return res;
        }
        DynamicModInt inv() const { return raw(mod_inv(val)); }
        DynamicModInt& operator++() { if (++val == bt.mod) val = 0; return *this; }
        DynamicModInt operator++(int) { DynamicModInt res = *this; ++*this; return res; }
        DynamicModInt& operator--() { if (val == 0) val = bt.mod; --val; return *this; }
        DynamicModInt operator--(int) { DynamicModInt res = *this; --*this; return res; }
        DynamicModInt& operator+=(const DynamicModInt& other) { val += other.val; if (val >= bt.mod) val -= bt.mod; return *this; }
        DynamicModInt& operator-=(const DynamicModInt& other) { if (val < other.val) val += bt.mod; val -= other.val; return *this; }
        DynamicModInt& operator*=(const DynamicModInt& other) { val = bt.mul(val, other.val); return *this; }
        DynamicModInt& operator/=(const DynamicModInt& other) { *this *= other.inv(); return *this; }
        friend DynamicModInt operator+(const DynamicModInt& lhs, const DynamicModInt& rhs) { return DynamicModInt(lhs) += rhs; }
        friend DynamicModInt operator-(const DynamicModInt& lhs, const DynamicModInt& rhs) { return DynamicModInt(lhs) -= rhs; }
        friend DynamicModInt operator*(const DynamicModInt& lhs, const DynamicModInt& rhs) { return DynamicModInt(lhs) *= rhs; }
        friend DynamicModInt operator/(const DynamicModInt& lhs, const DynamicModInt& rhs) { return DynamicModInt(lhs) /= rhs; }
        DynamicModInt operator+() const { return *this; }
        DynamicModInt operator-() const { return DynamicModInt() - *this; }
        friend bool operator==(const DynamicModInt& lhs, const DynamicModInt& rhs) { return lhs.val == rhs.val; }
        friend bool operator!=(const DynamicModInt& lhs, const DynamicModInt& rhs) { return lhs.val != rhs.val; }
        DynamicModInt pow(int64_t a) const {
                DynamicModInt v = *this, res = 1;
                while (a) {
                        if (a & 1) res *= v;
                        a >>= 1;
                        v *= v;
                }
                return res;
        }
        friend std::ostream& operator<<(std::ostream& os, const DynamicModInt& m) { return os << m.val; }
        friend std::istream& operator>>(std::istream& is, DynamicModInt& m) {
                int64_t v; is >> v;
                m = DynamicModInt(v);
                return is;
        }
};

template <int ID>
Barrett DynamicModInt<ID>::bt = Barrett(998244353);

using dmint = DynamicModInt<0>;
#line 8 "math/convolution/Convolution.hpp"

namespace Convolution {

constexpr uint64_t primitive_root_for_convolution(uint64_t p) {
        if (p == 2) return 1;
        if (p == 998244353) return 3;
        if (p == 469762049) return 3;
        if (p == 1811939329) return 11;
        if (p == 2013265921) return 11;
        for (uint64_t g = 2; g < p; ++g) {
                uint64_t base = g, exp = (p - 1) >> 1, res = 1;
                while (exp > 0) {
                        if (exp & 1) res = (res * base) % p;
                        base = (base * base) % p;
                        exp >>= 1;
                }
                if (res != 1) return g;
        }
        return -1;
}

namespace internal {

constexpr int msb(uint32_t x) { return x == 0 ? -1 : 31 - __builtin_clz(x); }
constexpr int ceil_log2(uint32_t x) { return x <= 1 ? 0 : 32 - __builtin_clz(x - 1); }

template <class T>
class NthRoot {
private:
        static constexpr unsigned int lg = msb((T::get_mod() - 1) & (1 - T::get_mod()));
        T root[lg + 1];
        T inv_root[lg + 1];
        T rate[lg + 1];
        T inv_rate[lg + 1];

public:
        constexpr NthRoot() : root{}, inv_root{}, rate{}, inv_rate{} {
                root[lg] = T(primitive_root_for_convolution(T::get_mod())).pow((T::get_mod() - 1) >> lg);
                inv_root[lg] = root[lg].inv();
                for (int i = lg - 1; i >= 0; i--) {
                        root[i] = root[i + 1] * root[i + 1];
                        inv_root[i] = inv_root[i + 1] * inv_root[i + 1];
                }
                T r = 1;
                for (int i = 2; i <= (int)lg; i++) {
                        rate[i - 2] = r * root[i];
                        r = r * inv_root[i];
                }
                r = 1;
                for (int i = 2; i <= (int)lg; i++) {
                        inv_rate[i - 2] = r * inv_root[i];
                        r = r * root[i];
                }
        }
        static constexpr unsigned int get_lg() { return lg; }
        constexpr T get(int n) const { return root[n]; }
        constexpr T inv(int n) const { return inv_root[n]; }
        constexpr T get_rate(int n) const { return rate[n]; }
        constexpr T get_inv_rate(int n) const { return inv_rate[n]; }
};

template <class T>
void number_theoretic_transform(std::vector<T>& a) {
        static constexpr NthRoot<T> nth_root;
        int n = (int)a.size();
        for (int i = n >> 1; i > 0; i >>= 1) {
                T z = T::raw(1);
                for (int j = 0; j < n; j += (i << 1)) {
                        for (int k = 0; k < i; k++) {
                                const T x = a[j + k];
                                const T y = a[j + i + k] * z;
                                a[j + k] = x + y;
                                a[j + i + k] = x - y;
                        }
                        z *= nth_root.get_rate(__builtin_popcount(j & ~(j + (i << 1))));
                }
        }
}

template <class T>
void inverse_number_theoretic_transform(std::vector<T>& a) {
        static constexpr NthRoot<T> nth_root;
        int n = (int)a.size();
        for (int i = 1; i < n; i <<= 1) {
                T z = T::raw(1);
                for (int j = 0; j < n; j += (i << 1)) {
                        for (int k = 0; k < i; k++) {
                                const T x = a[j + k];
                                const T y = a[j + i + k];
                                a[j + k] = x + y;
                                a[j + i + k] = (x - y) * z;
                        }
                        z *= nth_root.get_inv_rate(__builtin_popcount(j & ~(j + (i << 1))));
                }
        }
        T inv_n = T(1) / T(n);
        for (auto& x : a) x *= inv_n;
}

template <class T>
std::vector<T> convolution_naive(const std::vector<T>& a, const std::vector<T>& b) {
        int n = (int)a.size(), m = (int)b.size();
        std::vector<T> c(n + m - 1);
        for (int i = 0; i < n; i++) {
                for (int j = 0; j < m; j++) c[i + j] += a[i] * b[j];
        }
        return c;
}

template <class T>
std::vector<T> convolution_pow2(std::vector<T> a) {
        int n = (int)a.size() * 2 - 1;
        if (n <= 1) return convolution_naive(a, a);
        int lg = msb(n - 1) + 1;
        if (n - (1 << (lg - 1)) <= 5) {
                --lg;
                int m = (int)a.size() - (1 << (lg - 1));
                std::vector<T> a1(a.begin(), a.begin() + m);
                std::vector<T> a2(a.begin() + m, a.end());
                std::vector<T> c(n);
                std::vector<T> c1 = convolution_naive(a1, a1);
                std::vector<T> c2 = convolution_naive(a1, a2);
                std::vector<T> c3 = convolution_pow2(a2);
                for (int i = 0; i < (int)c1.size(); i++) c[i] += c1[i];
                for (int i = 0; i < (int)c2.size(); i++) c[i + m] += c2[i] * 2;
                for (int i = 0; i < (int)c3.size(); i++) c[i + m * 2] += c3[i];
                return c;
        }
        int m = 1 << lg;
        a.resize(m);
        number_theoretic_transform(a);
        for (int i = 0; i < m; i++) a[i] *= a[i];
        inverse_number_theoretic_transform(a);
        a.resize(n);
        return a;
}

template <class T>
std::vector<T> convolution(std::vector<T> a, std::vector<T> b) {
        int n = (int)a.size() + (int)b.size() - 1;
        if (n <= 1) return convolution_naive(a, b);
        int lg = ceil_log2(n);
        int m = 1 << lg;
        if (n - (1 << (lg - 1)) <= 5) {
                --lg;
                if (a.size() < b.size()) std::swap(a, b);
                int half_m = n - (1 << lg);
                std::vector<T> a1(a.begin(), a.begin() + half_m);
                std::vector<T> a2(a.begin() + half_m, a.end());
                std::vector<T> c(n);
                std::vector<T> c1 = convolution_naive(a1, b);
                std::vector<T> c2 = convolution(a2, b);
                for (int i = 0; i < (int)c1.size(); i++) c[i] += c1[i];
                for (int i = 0; i < (int)c2.size(); i++) c[i + half_m] += c2[i];
                return c;
        }
        a.resize(m);
        b.resize(m);
        number_theoretic_transform(a);
        number_theoretic_transform(b);
        for (int i = 0; i < m; i++) a[i] *= b[i];
        inverse_number_theoretic_transform(a);
        a.resize(n);
        return a;
}

} // namespace internal

template <uint32_t p>
struct is_ntt_friendly : std::false_type {};
template <> struct is_ntt_friendly<998244353> : std::true_type {};
template <> struct is_ntt_friendly<469762049> : std::true_type {};
template <> struct is_ntt_friendly<1811939329> : std::true_type {};
template <> struct is_ntt_friendly<2013265921> : std::true_type {};

template <uint32_t p>
std::vector<StaticModInt<p>> convolution_for_any_mod(const std::vector<StaticModInt<p>>& a, const std::vector<StaticModInt<p>>& b);

template <uint32_t p>
std::vector<StaticModInt<p>> convolution_large(const std::vector<StaticModInt<p>>& a, const std::vector<StaticModInt<p>>& b);

template <uint32_t p>
std::vector<StaticModInt<p>> convolution(const std::vector<StaticModInt<p>>& a, const std::vector<StaticModInt<p>>& b) {
        int n = (int)a.size(), m = (int)b.size();
        if (n == 0 || m == 0) return {};
        if (n <= 60 || m <= 60) return internal::convolution_naive(a, b);
        if (n + m - 1 <= (int)((1 - p) & (p - 1))) {
                if (n == m && a == b) return internal::convolution_pow2(a);
                return internal::convolution(a, b);
        }
        if constexpr (is_ntt_friendly<p>::value) {
                return convolution_large<p>(a, b);
        }
        return convolution_for_any_mod(a, b);
}

template <uint32_t p>
std::vector<StaticModInt<p>> convolution_for_any_mod(const std::vector<StaticModInt<p>>& a, const std::vector<StaticModInt<p>>& b) {
        int n = (int)a.size(), m = (int)b.size();
        assert(n + m - 1 <= (1 << 26));
        static constexpr uint32_t MOD1 = 469762049;
        static constexpr uint32_t MOD2 = 1811939329;
        static constexpr uint32_t MOD3 = 2013265921;
        using mint1 = StaticModInt<MOD1>;
        using mint2 = StaticModInt<MOD2>;
        using mint3 = StaticModInt<MOD3>;
        std::vector<mint1> a1(n), b1(m);
        std::vector<mint2> a2(n), b2(m);
        std::vector<mint3> a3(n), b3(m);
        for (int i = 0; i < n; i++) {
                uint32_t v = a[i].get();
                a1[i] = v; a2[i] = v; a3[i] = v;
        }
        for (int i = 0; i < m; i++) {
                uint32_t v = b[i].get();
                b1[i] = v; b2[i] = v; b3[i] = v;
        }
        auto c1 = internal::convolution(a1, b1);
        auto c2 = internal::convolution(a2, b2);
        auto c3 = internal::convolution(a3, b3);
        static constexpr uint32_t INV1_2 = mint2(MOD1).inv().get();
        static constexpr uint32_t INV1_3 = mint3(MOD1).inv().get();
        static constexpr uint32_t INV2_3 = mint3(MOD2).inv().get();
        std::vector<StaticModInt<p>> res(n + m - 1);
        for (int i = 0; i < n + m - 1; i++) {
                int64_t t1 = c1[i].get();
                int64_t t2 = (c2[i].get() - t1 + MOD2) * INV1_2 % MOD2;
                if (t2 < 0) t2 += MOD2;
                int64_t t3 = ((c3[i].get() - t1 + MOD3) * INV1_3 % MOD3 - t2 + MOD3) * INV2_3 % MOD3;
                if (t3 < 0) t3 += MOD3;
                res[i] = StaticModInt<p>(t1 + (t2 + t3 * MOD2) % p * (MOD1 % p));
        }
        return res;
}

template <uint32_t p>
std::vector<StaticModInt<p>> convolution_large(const std::vector<StaticModInt<p>>& a, const std::vector<StaticModInt<p>>& b) {
        int n = (int)a.size(), m = (int)b.size();
        int limit = (int)((1 - p) & (p - 1));
        int L = limit >> 1;
        int k = (n + L - 1) / L;
        int l = (m + L - 1) / L;
        std::vector<std::vector<StaticModInt<p>>> A_hat(k, std::vector<StaticModInt<p>>(limit));
        std::vector<std::vector<StaticModInt<p>>> B_hat(l, std::vector<StaticModInt<p>>(limit));
        for (int i = 0; i < k; ++i) {
                int start = i * L;
                int end = std::min(n, start + L);
                for (int j = start; j < end; ++j) {
                        A_hat[i][j - start] = a[j];
                }
                internal::number_theoretic_transform(A_hat[i]);
        }
        for (int i = 0; i < l; ++i) {
                int start = i * L;
                int end = std::min(m, start + L);
                for (int j = start; j < end; ++j) {
                        B_hat[i][j - start] = b[j];
                }
                internal::number_theoretic_transform(B_hat[i]);
        }
        std::vector<StaticModInt<p>> res(n + m - 1);
        for (int c = 0; c < k + l - 1; ++c) {
                std::vector<StaticModInt<p>> C_chunk(limit);
                for (int i = 0; i < k; ++i) {
                        int j = c - i;
                        if (0 <= j && j < l) {
                                for (int x = 0; x < limit; ++x) {
                                        C_chunk[x] += A_hat[i][x] * B_hat[j][x];
                                }
                        }
                }
                internal::inverse_number_theoretic_transform(C_chunk);
                int start = c * L;
                int end = std::min(n + m - 1, start + limit);
                for (int j = start; j < end; ++j) {
                        res[j] += C_chunk[j - start];
                }
        }
        return res;
}

template <uint32_t p>
std::vector<int64_t> convolution(const std::vector<int64_t>& a, const std::vector<int64_t>& b) {
        int n = (int)a.size(), m = (int)b.size();
        std::vector<StaticModInt<p>> a2(n), b2(m);
        for (int i = 0; i < n; i++) a2[i] = a[i];
        for (int i = 0; i < m; i++) b2[i] = b[i];
        auto c2 = convolution<p>(a2, b2);
        std::vector<int64_t> c(c2.size());
        for (int i = 0; i < (int)c2.size(); i++) c[i] = c2[i].get();
        return c;
}

inline std::vector<int64_t> convolution_ll(const std::vector<int64_t>& a, const std::vector<int64_t>& b) {
        int n = (int)a.size(), m = (int)b.size();
        if (n == 0 || m == 0) return {};
        static constexpr uint32_t MOD1 = 469762049;
        static constexpr uint32_t MOD2 = 1811939329;
        static constexpr uint32_t MOD3 = 2013265921;
        using mint1 = StaticModInt<MOD1>;
        using mint2 = StaticModInt<MOD2>;
        using mint3 = StaticModInt<MOD3>;
        std::vector<mint1> a1(n), b1(m);
        std::vector<mint2> a2(n), b2(m);
        std::vector<mint3> a3(n), b3(m);
        for (int i = 0; i < n; i++) {
                int64_t v = a[i] % MOD1; if (v < 0) v += MOD1; a1[i] = v;
                v = a[i] % MOD2; if (v < 0) v += MOD2; a2[i] = v;
                v = a[i] % MOD3; if (v < 0) v += MOD3; a3[i] = v;
        }
        for (int i = 0; i < m; i++) {
                int64_t v = b[i] % MOD1; if (v < 0) v += MOD1; b1[i] = v;
                v = b[i] % MOD2; if (v < 0) v += MOD2; b2[i] = v;
                v = b[i] % MOD3; if (v < 0) v += MOD3; b3[i] = v;
        }
        auto c1 = internal::convolution(a1, b1);
        auto c2 = internal::convolution(a2, b2);
        auto c3 = internal::convolution(a3, b3);
        static constexpr uint32_t INV1_2 = mint2(MOD1).inv().get();
        static constexpr uint32_t INV1_3 = mint3(MOD1).inv().get();
        static constexpr uint32_t INV2_3 = mint3(MOD2).inv().get();
        std::vector<int64_t> res(n + m - 1);
        for (int i = 0; i < n + m - 1; i++) {
                int64_t t1 = c1[i].get();
                int64_t t2 = c2[i].get() - t1;
                t2 %= MOD2;
                if (t2 < 0) t2 += MOD2;
                t2 = (t2 * INV1_2) % MOD2;
                int64_t t3 = c3[i].get() - t1;
                t3 %= MOD3;
                if (t3 < 0) t3 += MOD3;
                t3 = (t3 * INV1_3) % MOD3;
                t3 = (t3 - t2) % MOD3;
                if (t3 < 0) t3 += MOD3;
                t3 = (t3 * INV2_3) % MOD3;
                __int128_t result = t1;
                result += (__int128_t)t2 * MOD1;
                result += (__int128_t)t3 * MOD1 * MOD2;
                res[i] = (int64_t)result;
        }
        return res;
}

// https://judge.yosupo.jp/problem/convolution_mod_2_64
inline std::vector<uint64_t> convolution_64(const std::vector<uint64_t>& a, const std::vector<uint64_t>& b) {
        int n = (int)a.size(), m = (int)b.size();
        if (n == 0 || m == 0) return {};
        std::vector<int64_t> a_lo(n), a_hi(n), b_lo(m), b_hi(m);
        for (int i = 0; i < n; i++) {
                a_lo[i] = a[i] & 0xFFFFFFFF;
                a_hi[i] = a[i] >> 32;
        }
        for (int i = 0; i < m; i++) {
                b_lo[i] = b[i] & 0xFFFFFFFF;
                b_hi[i] = b[i] >> 32;
        }
        auto c_lolo = convolution_ll(a_lo, b_lo);
        auto c_hilo = convolution_ll(a_hi, b_lo);
        auto c_lohi = convolution_ll(a_lo, b_hi);
        std::vector<uint64_t> res(n + m - 1);
        for (int i = 0; i < n + m - 1; i++) {
                uint64_t lolo = (uint64_t)c_lolo[i];
                uint64_t hilo = (uint64_t)c_hilo[i];
                uint64_t lohi = (uint64_t)c_lohi[i];
                res[i] = lolo + (hilo << 32) + (lohi << 32);
        }
        return res;
}

template <class T>
void ntt_doubling_(std::vector<T>& a, std::vector<T> b) {
        static constexpr internal::NthRoot<T> nth_root;
        int n = (int)a.size();
        const T z = nth_root.get(internal::msb(n) + 1);
        T r = 1;
        for (int i = 0; i < n; i++) {
                b[i] *= r;
                r *= z;
        }
        internal::number_theoretic_transform(b);
        a.reserve(2 * n);
        a.insert(a.end(), b.begin(), b.end());
}

template <class T>
void ntt_doubling_(std::vector<T>& a) {
        auto b = a;
        internal::inverse_number_theoretic_transform(b);
        ntt_doubling_(a, std::move(b));
}

template <int ID>
std::vector<DynamicModInt<ID>> convolution(const std::vector<DynamicModInt<ID>>& a, const std::vector<DynamicModInt<ID>>& b) {
        int n = (int)a.size(), m = (int)b.size();
        if (n == 0 || m == 0) return {};
        if (n <= 60 || m <= 60) return internal::convolution_naive(a, b);
        uint32_t mod = DynamicModInt<ID>::get_mod();
        if (n + m - 1 <= (int)((1 - mod) & (mod - 1))) {
                if (n == m && a == b) return internal::convolution_pow2(a);
                return internal::convolution(a, b);
        }
        static constexpr uint32_t MOD1 = 469762049;
        static constexpr uint32_t MOD2 = 1811939329;
        static constexpr uint32_t MOD3 = 2013265921;
        using mint1 = StaticModInt<MOD1>;
        using mint2 = StaticModInt<MOD2>;
        using mint3 = StaticModInt<MOD3>;
        std::vector<mint1> a1(n), b1(m);
        std::vector<mint2> a2(n), b2(m);
        std::vector<mint3> a3(n), b3(m);
        for (int i = 0; i < n; i++) {
                uint32_t v = a[i].get();
                a1[i] = v; a2[i] = v; a3[i] = v;
        }
        for (int i = 0; i < m; i++) {
                uint32_t v = b[i].get();
                b1[i] = v; b2[i] = v; b3[i] = v;
        }
        auto c1 = internal::convolution(a1, b1);
        auto c2 = internal::convolution(a2, b2);
        auto c3 = internal::convolution(a3, b3);
        static constexpr uint32_t INV1_2 = mint2(MOD1).inv().get();
        static constexpr uint32_t INV1_3 = mint3(MOD1).inv().get();
        static constexpr uint32_t INV2_3 = mint3(MOD2).inv().get();
        std::vector<DynamicModInt<ID>> res(n + m - 1);
        for (int i = 0; i < n + m - 1; i++) {
                int64_t t1 = c1[i].get();
                int64_t t2 = (c2[i].get() - t1 + MOD2) * INV1_2 % MOD2;
                if (t2 < 0) t2 += MOD2;
                int64_t t3 = ((c3[i].get() - t1 + MOD3) * INV1_3 % MOD3 - t2 + MOD3) * INV2_3 % MOD3;
                if (t3 < 0) t3 += MOD3;
                res[i] = DynamicModInt<ID>(t1 + (t2 + t3 * MOD2) % mod * MOD1);
        }
        return res;
}

} // namespace Convolution
#line 5 "test/1_library_checker/math/convolution/Convolution107.test.cpp"

void solve() {
        int N, M; std::cin >> N >> M;
        std::vector<modint107> A(N), B(M);
        for (auto &x : A) std::cin >> x;
        for (auto &x : B) std::cin >> x;
        std::vector<modint107> C = Convolution::convolution(A, B);
        for (auto &x : C) std::cout << x << " ";
        std::cout << '\n';
}

int main() {
        std::ios::sync_with_stdio(false); std::cin.tie(nullptr);
        int tc = 1;
        // std::cin >> tc;
        while (tc--) solve();
        return 0;
}
Back to top page