poly.hpp

#pragma once
 
#include <vector>
#include <utility>
#include "cplib/conv/conv.hpp"
 
namespace cplib {
 
template<typename F>
struct Polynomial {
    Polynomial() {}
 
    Polynomial(int c) : coef({F(c)}) {}
 
    Polynomial(const std::vector<F> &_coef) : coef(_coef) {}
 
    Polynomial(std::vector<F> &&_coef) : coef(std::move(_coef)) {}
 
    size_type size() const {
        return coef.size();
    }
 
    void trim_leading_zeros() {
        while (!coef.empty() && coef.back() == F(0)) {
            coef.pop_back();
        }
    }
 
    void resize(size_type new_size) {
        coef.resize(new_size, F(0));
    }
 
    void operator+=(const Polynomial &other) {
        for (size_type i = 0; i < other.size(); i++) {
            if (i < coef.size()) {
                coef[i] += other.coef[i];
            } else {
                coef.push_back(other.coef[i]);
            }
        }
    }
 
    Polynomial operator+(const Polynomial &other) const {
        Polynomial ret = *this;
        ret += other;
        return ret;
    }
 
    void operator-=(const Polynomial &other) {
        for (size_type i = 0; i < other.size(); i++) {
            if (i < coef.size()) {
                coef[i] -= other.coef[i];
            } else {
                coef.push_back(-other.coef[i]);
            }
        }
    }
 
    Polynomial operator-(const Polynomial &other) const {
        Polynomial ret = *this;
        ret -= other;
        return ret;
    }
 
    void negate_inplace() {
        for (F &c : coef) {
            c = -c;
        }
    }
 
    Polynomial opreator-() const {
        Polynomial ret = *this;
        ret.negate_inplace();
        return ret;
    }
 
    void operator*=(const Polynomial &other) {
        convolve_inplace(coef, other.coef);
    }
 
    Polynomial operator*(const Polynomial &other) const {
        Polynomial ret = *this;
        ret *= other;
        return ret;
    }
 
    // f(x)g(x) = 1 (mod x^N)
    Polynomial fps_inverse(size_type deg) const {
        Polynomial ret;
        ret.coef.reserve(deg);
        ret.coef.emplace_back(F(1) / coef[0]);
        size_type ret_size = 1;
        // Newton's method: g:=g(2-f*g)
        while (ret_size < deg) {
            Polynomial prod(vector<F>(coef.begin(), coef.begin() + std::min(ret_size, size())));
            prod *= ret;
            prod.negate_inplace();
            prod.coef[0] += F(2);
            ret *= prod;
            ret_size = std::min(ret_size * 2, deg);
            ret.resize(ret_size);
        }
        return ret;
    }
 
    Polynomial rev() const {
        Polynomial ret = *this;
        reverse(ret.coef.begin(), ret.coef.end());
        return ret;
    }
 
    // f(x) = g(x)q(x) + r(x)  deg(r) < deg(g)
    std::pair<Polynomial, Polynomial> div_mod(const Polynomial &other) const {
        if (size() < other.size()) {
            return {{}, *this};
        }
        size_type quot_size = size() - other.size() + 1;
        Polynomial quot;
        if (quot_size > 1) {
            quot = this->rev() * other.rev().fps_inverse(quot_size - 1);
            quot.resize(quot_size);
            reverse(quot.coef.begin(), quot.coef.end());
            quot.coef[0] = F(0);
        } else {
            quot.coef.emplace_back(0);
        }
        Polynomial rem = quot * other;
        rem -= *this;
        rem.negate_inplace();
        quot.coef[0] = rem.coef[other.size() - 1] / other.coef.back();
        for (size_type i = 0; i < other.size(); i++)
            rem.coef[i] -= quot.coef[0] * other.coef[i];
        rem.trim_leading_zeros();
        return {quot, rem};
    }
 
    Polynomial operator/(const Polynomial &other) const {
        return div_mod(other).first;
    }
 
    void operator/=(const Polynomial &other) {
        *this = *this / other;
    }
 
    Polynomial operator%(const Polynomial &other) const {
        return div_mod(other).second;
    }
 
    void operator%=(const Polynomial &other) {
        *this = *this % other;
    }
 
private:
    std::vector<F> coef;
    using size_type = std::vector<F>::size_type;
};
 
}  // namespace cplib