形式的冪級数 LIB::FormalPowerSeries

LIB::FormalPowerSeries

modの畳み込みでn次式fとm次式gの乗除算を高速に計算する
https://zenn.dev/koboshi/books/600753d917cdbc/viewer/3a148e
https://web.archive.org/web/20220813112322/https://opt-cp.com/fps-implementation/
https://judge.yosupo.jp/submission/42413
https://web.archive.org/web/20220813113030/https://opt-cp.com/fps-fast-algorithms/
https://maspypy.com/category/%e5%bd%a2%e5%bc%8f%e7%9a%84%e3%81%b9%e3%81%8d%e7%b4%9a%e6%95%b0%e8%a7%a3%e8%aa%ac

使い方

と に対して
積 を計算する

#include <bits/stdc++.h>
#include <atcoder/all>
using namespace std;
using namespace atcoder;
const long long mint1=1000000007;
const long long mint9=998244353; 

using mint=static_modint<mint9>;
using fps=FormalPowerSeries<mint>;
using sfps=vector<pair<int,mint>>;

int main(void){
	fps f={9,8,7},g1={9,0,8,7};
	f.resize(6);
	fps ans1=f*g1; // 畳み込みを使った演算
	cout<<ans1[5].val()<<endl;
	sfps g2={{0,9},{1,0},{2,8},{3,7}};
	fps ans2=f*g2; // 疎な多項式の演算
	cout<<ans2[5].val()<<endl;
	return 0;
}

49
49

3行の using を書いてからメイン処理を書く。fは3次式で積が5次式なので、fを5次にリサイズする必要がある（リサイズしなかった場合、3次以上の情報は捨てられる）。畳み込みを使う場合はfps型、疎な（項の少ない）多項式にはsfps型を使うと高速に計算できる。なお atcoder ライブラリを利用するので<atcoder/all>が必要。

逆元を使った高速化テク

より
また （ただし は定数項を持たない多項式）

とすると
となり が答え
https://atcoder.jp/contests/tdpc/submissions/45972918

メンバ関数

• F &operator*=(const F &g): 畳み込みで高速化したn次式とm次式の積
• F &operator*=(vector<pair<int,T>> g): 疎な多項式の積 （n,mは項数）
• F &operator/=(const F &g): 畳み込みで高速化したn次式とm次式の商
• F &operator/=(vector<pair<int,T>> g): 疎な多項式の商 （n,mは項数）

#include <bits/stdc++.h>
#include <atcoder/all>
using namespace std;
using namespace atcoder;

namespace LIB{
#define rep2(i,m,n) for(int i=(m);i<(n);++i)
#define rep1(i,n) rep2(i,0,n)
#define drep2(i,m,n) for(int i=(m)-1;i>=(n);--i)
#define drep1(i,n) drep2(i,n,0)

	template<class T>struct FormalPowerSeries:vector<T>{
		using vector<T>::vector;
		using vector<T>::operator=;
		using F=FormalPowerSeries;

		F operator-() const {
			F res(*this);
			for (auto &e : res) e = -e;
			return res;
		}
		F &operator*=(const T &g) {
			for (auto &e : *this) e *= g;
			return *this;
		}
		F &operator/=(const T &g) {
			assert(g != T(0));
			*this *= g.inv();
			return *this;
		}
		F &operator+=(const F &g) {
			int n = (*this).size(), m = g.size();
			rep1(i, min(n, m)) (*this)[i] += g[i];
			return *this;
		}
		F &operator-=(const F &g) {
			int n = (*this).size(), m = g.size();
			rep1(i, min(n, m)) (*this)[i] -= g[i];
			return *this;
		}
		F &operator<<=(const int d) {
			int n = (*this).size();
			(*this).insert((*this).begin(), d, 0);
			(*this).resize(n);
			return *this;
		}
		F &operator>>=(const int d) {
			int n = (*this).size();
			(*this).erase((*this).begin(), (*this).begin() + min(n, d));
			(*this).resize(n);
			return *this;
		}
		F inv(int d = -1) const {
			int n = (*this).size();
			assert(n != 0 && (*this)[0] != 0);
			if (d == -1) d = n;
			assert(d > 0);
			F res{(*this)[0].inv()};
			while (res.size() < d) {
				int m = size(res);
				F f(begin(*this), begin(*this) + min(n, 2*m));
				F r(res);
				f.resize(2*m), internal::butterfly(f);
				r.resize(2*m), internal::butterfly(r);
				rep1(i, 2*m) f[i] *= r[i];
				internal::butterfly_inv(f);
				f.erase(f.begin(), f.begin() + m);
				f.resize(2*m), internal::butterfly(f);
				rep1(i, 2*m) f[i] *= r[i];
				internal::butterfly_inv(f);
				T iz = T(2*m).inv(); iz *= -iz;
				rep1(i, m) f[i] *= iz;
				res.insert(res.end(), f.begin(), f.begin() + m);
			}
			return {res.begin(), res.begin() + d};
		}
		F &operator*=(const F &g) {
			int n = (int)(*this).size();
			*this = convolution(*this, g);
			(*this).resize(n);
			return *this;
		}
		F &operator/=(const F &g) {
		   int n = (*this).size();
		   *this = convolution(*this, g.inv(n));
		   (*this).resize(n);
		   return *this;
		}
		F &operator*=(vector<pair<int, T>> g) {
			int n = (int)(*this).size();
			auto [d, c] = g.front();
			if (d == 0) g.erase(g.begin());
			else c = 0;
			drep1(i, n) {
				(*this)[i] *= c;
				for (auto &[j, b] : g) {
					if (j > i) break;
					(*this)[i] += (*this)[i-j] * b;
				}
			}
			return *this;
		}
		F &operator/=(vector<pair<int, T>> g) {
			int n = (*this).size();
			auto [d, c] = g.front();
			assert(d == 0 && c != T(0));
			T ic = c.inv();
			g.erase(g.begin());
			rep1(i, n) {
				for (auto &[j, b] : g) {
					if (j > i) break;
					(*this)[i] -= (*this)[i-j] * b;
				}
				(*this)[i] *= ic;
			}
			return *this;
		}
		void multiply(const int d, const T c) { 
			int n = (*this).size();
			if (c == T(1)) drep1(i, n-d) (*this)[i+d] += (*this)[i];
			else if (c == T(-1)) drep1(i, n-d) (*this)[i+d] -= (*this)[i];
			else drep1(i, n-d) (*this)[i+d] += (*this)[i] * c;
		}
		void divide(const int d, const T c) {
			int n = (*this).size();
			if (c == T(1)) rep1(i, n-d) (*this)[i+d] -= (*this)[i];
			else if (c == T(-1)) rep1(i, n-d) (*this)[i+d] += (*this)[i];
			else rep1(i, n-d) (*this)[i+d] -= (*this)[i] * c;
		}
		T eval(const T &a) const {
			T x(1), res(0);
			for (auto e : *this) res += e * x, x *= a;
			return res;
		}
		F operator*(const T &g) const { return F(*this) *= g; }
		F operator/(const T &g) const { return F(*this) /= g; }
		F operator+(const F &g) const { return F(*this) += g; }
		F operator-(const F &g) const { return F(*this) -= g; }
		F operator<<(const int d) const { return F(*this) <<= d; }
		F operator>>(const int d) const { return F(*this) >>= d; }
		F operator*(const F &g) const { return F(*this) *= g; }
		F operator/(const F &g) const { return F(*this) /= g; }
		F operator*(vector<pair<int, T>> g) const { return F(*this) *= g; }
		F operator/(vector<pair<int, T>> g) const { return F(*this) /= g; }
	};

#undef rep2
#undef rep1
#undef drep2
#undef drep1
}

実装例

https://atcoder.jp/contests/tdpc/submissions/45969929