QOJ.ac

QOJ

IDProblemSubmitterResultTimeMemoryLanguageFile sizeSubmit timeJudge time
#371479#618. 多项式乘法IsrothyCompile Error//C++2320.0kb2024-03-30 13:10:282024-03-30 13:10:29

Judging History

你现在查看的是最新测评结果

  • [2024-03-30 13:10:29]
  • 评测
  • [2024-03-30 13:10:28]
  • 提交

answer

#include <bit>
#include <cassert>
#include <cmath>
#include <cstdio>
#include <functional>
#include <numeric>
#include <optional>
#include <queue>
#include <random>
#include <span>
#include <unordered_map>
#include <vector>


constexpr int mod = 998244353;
constexpr int g = 3;

int32_t constexpr mul_mod(int32_t a, int32_t b, int32_t mod) {
    return static_cast<int>(static_cast<int64_t>(a) * b % mod);
}
int64_t constexpr mul_mod(int64_t a, int64_t b, int64_t mod) {
#ifdef __SIZEOF_INT128__
    return static_cast<int64_t>(static_cast<__int128>(a) * b % mod);
#else
    int64_t res = 0;
    for (b = (b % mod + mod) % mod; b; b >>= 1) {
        if (b & 1) {
            res = (res + a) % mod;
        }
        a = (a + a) % mod;
    }
    return res;
#endif
}
template<typename T>
constexpr std::tuple<T, T, T> ex_gcd(const T &a, const T &b) {
    if (b == 0) {
        return {a, T(1), T(0)};
    }
    auto [d, x, y] = ex_gcd(b, a % b);
    return {d, y, x - a / b * y};
}
template<typename T, typename U>
constexpr T power(T x, U k, const std::function<T(T, T)> &multiply) {
    T res{1};
    for (; k; k >>= 1) {
        if (k & 1) {
            res = multiply(res, x);
        }
        x = multiply(x, x);
    }
    return res;
}
template<typename T, typename U>
constexpr T power(T x, U k, const T &mod) {
    T res{1};
    for (; k; k >>= 1) {
        if (k & 1) {
            res = mul_mod(res, x, mod);
        }
        x = mul_mod(x, x, mod);
    }
    return res;
}
template<typename T>
constexpr T inverse(const T &a, const T &mod) {
    auto [d, x, y] = ex_gcd(a, mod);
    return d * x;
}
template<typename T>
std::optional<T> bsgs(const T &a, T b, const T &mod) {
    if (mod == 1) {
        return 0;
    }
    T w{1}, x{1}, s{static_cast<T>(std::sqrt(mod)) + 1};
    std::unordered_map<T, T> map;
    map.reserve(s);
    for (T k = 1; k <= s; ++k) {
        b = mul_mod(b, a, mod);
        w = mul_mod(w, a, mod);
        map[b] = k;
    }
    for (T k = 1; k <= s; ++k) {
        x = mul_mod(x, w, mod);
        if (map.contains(x)) {
            return (k * s - map[x]) % (mod - 1);
        }
    }
    return std::nullopt;
}
template<typename T>
std::optional<T> ex_bsgs(T a, T b, const T &mod) {
    a = (a % mod + mod) % mod;
    b = (b % mod + mod) % mod;
    if (b == 1 || mod == 1) {
        return 0;
    }
    auto d = gcd(a, mod);
    if (b % d) {
        return std::nullopt;
    }
    if (d == 1) {
        return bsgs(a, b, mod);
    }
    auto g = inverse(a / d, mod / d);
    auto x = ex_bsgs(a, b / d * g, mod / d);
    if (!x.has_value()) {
        return std::nullopt;
    }
    return x.value() + 1;
}
template<typename T>
struct Crt {
    std::vector<T> mt;
    T m{};
    Crt() = default;
    explicit Crt(std::span<T> a) : mt(a.size()) {
        m = std::accumulate(a.begin(), a.end(), T{1}, std::multiplies<>());
        for (int i = 0; i < a.size(); ++i) {
            auto mi = m / a[i];
            mt[i] = mi * inverse(mi, a[i]) % m;
        }
    }
    T query(std::span<T> b) {
        assert(b.size() == mt.size());
        T res = 0;
        for (int i = 0; i < mt.size(); ++i) {
            res = (res + b[i] * mt[i]) % m;
        }
        return res;
    }
};
template<typename T>
auto ex_crt(T a1, T m1, T a2, T m2) -> std::optional<std::pair<T, T>> {
    auto [d, x, y] = ex_gcd(m1, m2);
    if ((a2 - a1) % d) {
        return std::nullopt;
    }
    auto m = m1 / d * m2;
    auto t = ((a2 - a1) / d * x % (m2 / d)) * m1 % m;
    auto a = (a1 + t) % m;
    if (a < 0) {
        a += m;
    }
    return std::pair{a, m};
}
auto sieve_of_euler(std::span<bool> is_composite) {
    auto n = is_composite.size();
    std::vector<int> primes;
    primes.reserve(static_cast<int>(static_cast<double>(n) / std::log(n)));
    primes.push_back(0);
    for (int i = 2; i < n; ++i) {
        if (!is_composite[i]) {
            primes.push_back(i);
        }
        for (int j = 1; j < primes.size() && i * primes[j] < n; ++j) {
            is_composite[i * primes[j]] = true;
            if (i % primes[j] == 0) {
                break;
            }
        }
    }
    return primes;
}
template<typename T>
std::optional<T> primitive_root(T n, std::span<int> primes) {
    if (n == 2 || n == 4) {
        return n - 1;
    }
    if (n == 1 || (n & 3) == 0) {
        return std::nullopt;
    }
    auto a = prime_factors(n, primes);
    if (2 < a.size() || (a.size() == 2 && a[0] != 2)) {
        return std::nullopt;
    }
    T m = a.size() == 2 ? n / 2 / a[1] * (a[1] - 1) : n / a[0] * (a[0] - 1);
    auto b = prime_factors(m, primes);
    for (T g{2}; g < n; ++g) {
        if (power(g, m, n) == 1
            && std::all_of(b.begin(), b.end(), [&](auto p) { return power(g, m / p, n) != 1; })) {
            return g;
        }
    }
    return std::nullopt;
}
template<typename T>
bool is_prime(const T &n) {
    if (n < 2) {
        return false;
    }
    if (~n & 1) {
        return n == 2;
    }
    auto d = n - 1, s = 0;
    for (; ~d & 1; d >>= 1) {
        ++s;
    }
    static constexpr auto p = {2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37};
    return std::none_of(p.begin(), p.end(), [=](auto a) {
        if (a == n) {
            return false;
        }
        T x = power<T, T>(a, d, n);
        if (x == 1 || x == n - 1) {
            return false;
        }
        for (int i = 1; i < s; ++i) {
            x = mul_mod(x, x, n);
            if (x == n - 1) {
                return false;
            }
        }
        return true;
    });
}
template<typename T>
T pollard_rho(const T &n) {
    if (is_prime(n)) {
        return n;
    }
    for (auto p: {2, 3, 5, 7, 11, 13, 17, 19, 23, 29}) {
        if (n % p == 0) {
            return p;
        }
    }
    std::uniform_int_distribution<T> dist(1, n - 1);
    while (true) {
        static std::mt19937 mt_rand(std::random_device{}());
        auto c = dist(mt_rand);
        auto f = [&](const T &x) {
            return (mul_mod(x, x, n) + c) % n;
        };
        auto t = f(0), r = f(t);
        int steps = 1;
        while (t != r) {
            T prod{1};
            for (int i = 0; i < steps; ++i) {
                if (auto tmp = mul_mod(prod, std::abs(t - r), n)) {
                    prod = tmp;
                    t = f(t);
                    r = f(f(r));
                } else {
                    break;
                }
            }
            if (auto d = std::gcd(prod, n); d != 1) {
                return d;
            }
            steps = std::min(128, steps << 1);
        }
    }
}
template<typename T>
auto prime_factors(const T &n) {
    std::queue<T> q;
    std::vector<T> res;
    for (q.push(n); !q.empty(); q.pop()) {
        if (auto x = q.front(); is_prime(x)) {
            res.push_back(x);
        } else {
            auto d = pollard_rho(x);
            q.push(d);
            q.push(x / d);
        }
    }
    std::sort(res.begin(), res.end());
    res.erase(std::unique(res.begin(), res.end()), res.end());
    return res;
}
namespace polynomial {
    auto congruence_equation(int64_t a, int64_t b, int64_t p) {
        b = (b % p + p) % p;
        auto d = std::gcd(a, p);
        std::tie(a, b, p) = std::make_tuple(a / d, b / d, p / d);
        return (b * inverse(a, p) % p + p) % p;
    }
    auto quadratic_residue(int64_t x, int64_t mod, int64_t g) {
        if (!x) {
            return 0ll;
        }
        x = bsgs(g, x, mod).value();
        x = power(g, congruence_equation(2, x, mod - 1), mod);
        return std::min(x, mod - x);
    }
    template<int Mod, int G>
    struct polynomial : private std::vector<int> {
        using std::vector<int>::vector;
        using std::vector<int>::operator[];
        using std::vector<int>::begin;
        using std::vector<int>::end;
        using std::vector<int>::size;
        using std::vector<int>::resize;
        auto &operator*=(int64_t k) {
            for (auto &x: *this) {
                x = x * k % Mod;
            }
            return *this;
        }
        auto &operator+=(const polynomial<Mod, G> &rhs) {
            if (rhs.size() > size()) {
                resize(rhs.size());
            }
            for (int i = 0; i < rhs.size(); ++i) {
                (*this)[i] = ((*this)[i] + rhs[i]) % Mod;
            }
            return *this;
        }
        auto &operator-=(const polynomial<Mod, G> &rhs) {
            if (rhs.size() > size()) {
                resize(rhs.size());
            }
            for (int i = 0; i < rhs.size(); ++i) {
                (*this)[i] = ((*this)[i] - rhs[i]) % Mod;
            }
            return *this;
        }
        auto &operator+=(int64_t x) {
            (*this)[0] = ((*this)[0] + x) % Mod;
            return *this;
        }
        auto &operator-=(int64_t x) {
            (*this)[0] = ((*this)[0] - x) % Mod;
            return *this;
        }
        auto &operator*=(polynomial<Mod, G> rhs) {
            auto m = size() + rhs.size() - 1;
            auto n = std::bit_ceil(m);
            *this = dft(modXN(std::move(*this), n), 1);
            rhs = dft(modXN(std::move(rhs), n), 1);
            for (int i = 0; i < n; ++i) {
                (*this)[i] = (int64_t) (*this)[i] * rhs[i] % Mod;
            }
            *this = modXN(dft(std::move(*this), -1), m);
            return *this;
        }
    };
    template<int Mod, int G>
    auto dft(polynomial<Mod, G> a, int f) {
        static constexpr auto wn{[]() constexpr {
            constexpr auto len = std::countr_zero(static_cast<uint64_t>(Mod) - 1);
            std::array<std::array<int, len>, 2> wn{};
            for (int i = 0; i < len; ++i) {
                wn[0][i] = power(G, (Mod - 1) >> (i + 1), Mod);
                wn[1][i] = inverse(wn[0][i], Mod);
            }
            return wn;
        }()};
        int n = a.size();
        std::vector<int> w(n);
        for (int i = 0, j = 0; i < n; ++i) {
            if (i < j) {
                std::swap(a[i], a[j]);
            }
            for (int l = n >> 1; (j ^= l) < l; l >>= 1)
                ;
        }
        w[0] = 1;
        for (int i = 0; 1 << i < n; ++i) {
            for (int j = (1 << (i + 1)) - 1; j; --j) {
                w[j] = j & 1 ? (int64_t) w[j >> 1] * wn[(1 - f) / 2][i] % Mod : w[j >> 1];
            }
            for (int j = 0; j < n; j += 1 << (i + 1)) {
                auto *p = &a[j], *q = &a[j | 1 << i], *r = &w[0];
                for (int k = 0; k < 1 << i; ++k) {
                    auto t = (int64_t) q[k] * r[k];
                    q[k] = (p[k] - t) % Mod;
                    p[k] = (p[k] + t) % Mod;
                }
            }
        }
        if (f == -1) {
            int64_t in = ::inverse(n, Mod);
            for (auto &x: a) {
                x = x * in % Mod;
            }
        }
        return a;
    }
    template<int Mod, int G>
    auto modXN(polynomial<Mod, G> &&p, int n) {
        p.resize(n);
        return p;
    }
    template<int Mod, int G>
    auto modXN(const polynomial<Mod, G> &p, int n) {
        polynomial<Mod, G> res(n);
        std::copy(p.begin(), p.begin() + std::min(n, p.size()), res.begin());
        return res;
    }
    template<int Mod, int G>
    auto divXN(polynomial<Mod, G> &&p, int n) {
        std::copy(p.begin() + n, p.end(), p.begin());
        p.resize(p.size() - n);
        return p;
    }
    template<int Mod, int G>
    auto divXN(const polynomial<Mod, G> &p, int n) {
        polynomial res(p.size() - n);
        std::copy(p.begin() + n, p.end(), res.begin());
        return res;
    }
    template<int Mod, int G>
    auto reverse(polynomial<Mod, G> p) {
        std::reverse(p.begin(), p.end());
        return p;
    }
    template<int Mod, int G>
    auto operator+(polynomial<Mod, G> lhs, const polynomial<Mod, G> &rhs) {
        return lhs += rhs;
    }
    template<int Mod, int G>
    auto operator-(polynomial<Mod, G> lhs, const polynomial<Mod, G> &rhs) {
        return lhs -= rhs;
    }
    template<int Mod, int G>
    auto operator+(int64_t x, polynomial<Mod, G> p) {
        return p += x;
    }
    template<int Mod, int G>
    auto operator+(polynomial<Mod, G> p, int64_t x) {
        return p += x;
    }
    template<int Mod, int G>
    auto operator-(polynomial<Mod, G> p, int64_t x) {
        return p -= x;
    }
    template<int Mod, int G>
    auto operator-(int64_t x, polynomial<Mod, G> p) {
        return p -= x;
    }
    template<int Mod, int G>
    auto operator*(int64_t x, polynomial<Mod, G> p) {
        return p *= x;
    }
    template<int Mod, int G>
    auto operator*(polynomial<Mod, G> p, int64_t x) {
        return p *= x;
    }
    template<int Mod, int G>
    auto operator*(polynomial<Mod, G> lhs, const polynomial<Mod, G> &rhs) {
        return lhs *= rhs;
    }
    template<int Mod, int G>
    auto inverse(const polynomial<Mod, G> &p) {
        polynomial<Mod, G> res = {static_cast<int>(inv(p[0], Mod))};
        auto n = std::bit_ceil(p.size());
        for (int i = 2; i <= n; i <<= 1) {
            auto a = dft(modXN(modXN(p, i), i << 1), 1);
            auto b = dft(modXN(std::move(res), i << 1), 1);
            for (int j = 0; j < i << 1; ++j) {
                b[j] = b[j] * (2 - (int64_t) a[j] * b[j] % Mod) % Mod;
            }
            res = modXN(dft(std::move(b), -1), i);
        }
        return modXN(std::move(res), p.size());
    }
    template<int Mod, int G>
    auto derivative(polynomial<Mod, G> p) {
        for (int i = 1; i < p.size(); ++i) {
            p[i - 1] = (int64_t) i * p[i] % Mod;
        }
        p.resize(p.size() - 1);
        return p;
    }
    template<int Mod, int G>
    auto integral(polynomial<Mod, G> p) {
        p.resize(p.size() + 1);
        for (int i = (int) p.size() - 1; i >= 0; --i) {
            p[i] = ::inverse(i, Mod) * p[i - 1] % Mod;
        }
        p[0] = 0;
        return p;
    }
    template<int Mod, int G>
    auto log(const polynomial<Mod, G> &p) {
        return modXN(integral(derivative(p) * inverse(p)), p.size());
    }
    template<int Mod, int G>
    auto exp(const polynomial<Mod, G> &p) {
        polynomial<Mod, G> res = {1};
        auto n = std::bit_ceil(p.size());
        for (int i = 2; i <= n; i <<= 1) {
            auto a = dft(modXN(modXN(p, i), i << 1), 1);
            auto b = dft(modXN(res, i << 1), 1);
            auto c = dft(modXN(log(modXN(std::move(res), i)), i << 1), 1);
            for (int j = 0; j < i << 1; ++j) {
                b[j] = (int64_t) b[j] * (1 + a[j] - c[j]) % Mod;
            }
            res = modXN(dft(std::move(b), -1), i);
        }
        return modXN(std::move(res), p.size());
    }
    template<int Mod, int G>
    auto pow(const polynomial<Mod, G> &p, int64_t k) {
        return exp(log(p) * k);
    }
    template<int Mod, int G>
    auto sqrt(const polynomial<Mod, G> &p) {
        polynomial<Mod, G> res = {static_cast<int>(quadratic_residue(p[0], Mod, G))};
        constexpr auto inv2 = ::inverse(2, Mod);
        auto n = std::bit_ceil(p.size());
        for (int i = 2; i <= n; i <<= 1) {
            auto a = dft(modXN(modXN(p, i), i << 1), 1);
            auto b = dft(modXN(res, i << 1), 1);
            auto c = dft(modXN(inverse(modXN(std::move(res), i)), i << 1), 1);
            for (int j = 0; j < i << 1; ++j) {
                b[j] = (b[j] + (int64_t) a[j] * c[j]) % Mod * inv2 % Mod;
            }
            res = modXN(dft(std::move(b), -1), i);
        }
        return modXN(std::move(res), p.size());
    }
    template<int Mod, int G>
    auto operator/(const polynomial<Mod, G> &lhs, const polynomial<Mod, G> &rhs) {
        auto n = lhs.size();
        auto m = rhs.size();
        if (n < m) {
            return polynomial<Mod, G>{0};
        }
        auto a = modXN(reverse(lhs), n - m + 1);
        auto b = modXN(reverse(rhs), n - m + 1);
        return reverse(modXN(a * inverse(b), n - m + 1));
    }
    template<int Mod, int G>
    auto operator%(const polynomial<Mod, G> &lhs, const polynomial<Mod, G> &rhs) {
        return modXN(lhs - lhs / rhs * rhs, rhs.size() - 1);
    }
    template<int Mod, int G>
    auto operator/=(polynomial<Mod, G> &lhs, const polynomial<Mod, G> &rhs) {
        return lhs = lhs / rhs;
    }
    template<int Mod, int G>
    auto operator%=(polynomial<Mod, G> &lhs, const polynomial<Mod, G> &rhs) {
        return lhs = lhs % rhs;
    }
    template<int Mod, int G>
    auto
    eva_build(int p, int l, int r, const std::vector<int> &x, std::vector<polynomial<Mod, G>> &a) {
        if (l == r) {
            a[p] = {1, l < x.size() ? -x[l] : 0};
            return;
        }
        auto mid = (l + r) >> 1;
        eva_build(p << 1, l, mid, x, a);
        eva_build(p << 1 | 1, mid + 1, r, x, a);
        a[p] = a[p << 1] * a[p << 1 | 1];
    }
    template<int Mod, int G>
    auto eva_work(
        int p,
        int l,
        int r,
        const polynomial<Mod, G> &f,
        std::vector<polynomial<Mod, G>> &a,
        std::vector<int> &res
    ) {
        if (l == r) {
            if (l < res.size()) {
                res[l] = f[0];
            }
            return;
        }
        int mid = (l + r) >> 1;
        auto fsize = f.size();
        auto n = std::bit_ceil(fsize);
        auto x = dft(modXN(f, n), 1);
        auto helper = [n, fsize](polynomial<Mod, G> x, const polynomial<Mod, G> &g) {
            auto b = dft(modXN(g, n), 1);
            for (int i = 0; i < n; ++i) {
                x[i] = (int64_t) x[i] * b[i] % Mod;
            }
            return divXN(modXN(dft(std::move(x), -1), fsize), g.size() - 1);
        };
        auto lf = helper(x, a[p << 1 | 1]);
        auto rf = helper(x, a[p << 1]);
        eva_work(p << 1, l, mid, lf, a, res);
        eva_work(p << 1 | 1, mid + 1, r, rf, a, res);
    }
    template<int Mod, int G>
    auto evaluation(const polynomial<Mod, G> &p, const std::vector<int> &x) {
        int m = std::max(x.size(), p.size() - 1);
        std::vector<polynomial<Mod, G>> a(m << 2);
        std::vector<int> res(x.size());
        eva_build(1, 0, m - 1, x, a);
        auto f = modXN(reverse(modXN(p, m + 1)) * inverse(a[1]), m + 1);
        eva_work(1, 0, m - 1, f, a, res);
        for (int i = 0; i < x.size(); ++i) {
            res[i] = (p[0] + (int64_t) res[i] * x[i]) % Mod;
        }
        return res;
    }
    template<int Mod, int G>
    polynomial<Mod, G> interpolation_work(
        int p,
        int l,
        int r,
        const std::vector<int> &y,
        std::vector<polynomial<Mod, G>> &a,
        const std::vector<int> &b
    ) {
        if (l == r) {
            return {(int) (y[l] * ::inverse<int64_t>(b[l], Mod) % Mod)};
        }
        auto mid = (l + r) >> 1;
        auto lf = interpolation_work(p << 1, l, mid, y, a, b);
        auto rf = interpolation_work(p << 1 | 1, mid + 1, r, y, a, b);
        return lf * reverse(a[p << 1 | 1]) + rf * reverse(a[p << 1]);
    }
    template<int Mod, int G>
    auto interpolation(const std::vector<int> &x, const std::vector<int> &y) {
        auto n = x.size();
        std::vector<polynomial<Mod, G>> a(n << 2);
        std::vector<int> b(n);
        eva_build(1, 0, n - 1, x, a);
        auto f = derivative(reverse(a[1]));
        auto g = modXN(reverse(modXN(f, n + 1)) * inverse(a[1]), n + 1);
        eva_work(1, 0, n - 1, g, a, b);
        for (int i = 0; i < n; ++i) {
            b[i] = (f[0] + (int64_t) b[i] * x[i]) % Mod;
        }
        return interpolation_work(1, 0, n - 1, y, a, b);
    }
}// namespace polynomial

int main() {
    int n, m;
    scanf("%d%d", &n, &m);
    polynomial::polynomial<mod, g> a(n + 1), b(m + 1);
    for (int i = 0; i <= n; ++i) {
        scanf("%d", &a[i]);
    }
    for (int i = 0; i <= m; ++i) {
        scanf("%d", &b[i]);
    }
    auto c = a * b;
    for (auto x: c) {
        printf("%d ", (x + mod) % mod);
    }
    puts("");
    return 0;
}

Details

answer.code: In function ‘auto polynomial::quadratic_residue(int64_t, int64_t, int64_t)’:
answer.code:284:24: error: inconsistent deduction for auto return type: ‘long long int’ and then ‘long int’
  284 |         return std::min(x, mod - x);
      |                ~~~~~~~~^~~~~~~~~~~~
answer.code: In function ‘int main()’:
answer.code:628:10: warning: ignoring return value of ‘int scanf(const char*, ...)’ declared with attribute ‘warn_unused_result’ [-Wunused-result]
  628 |     scanf("%d%d", &n, &m);
      |     ~~~~~^~~~~~~~~~~~~~~~
answer.code:631:14: warning: ignoring return value of ‘int scanf(const char*, ...)’ declared with attribute ‘warn_unused_result’ [-Wunused-result]
  631 |         scanf("%d", &a[i]);
      |         ~~~~~^~~~~~~~~~~~~
answer.code:634:14: warning: ignoring return value of ‘int scanf(const char*, ...)’ declared with attribute ‘warn_unused_result’ [-Wunused-result]
  634 |         scanf("%d", &b[i]);
      |         ~~~~~^~~~~~~~~~~~~