Poly

poly（多项式）

一般的， poly 问题就是将题目转化为离散卷积的形式然后用 NTT 或者 FFT 求出卷积的系数

离散卷积：给定两个多项式

\[A(x)=\sum_{i=0}^{n}a_ix^i \]

\[B(x)=\sum_{i=0}^{m}b_ix^i \]

设

\[C(x)=A(x)B(x)=\sum_{k=0}^{n+m}c_kx^k \]

则有

\[c_k=\sum_{i+j=k}a_ib_j \]

称序列 \(c_k\) 为序列 \(a_i,b_i\) 的离散卷积。

复数：我们定义 \(i^2=-1\) 表示虚数单位，则有复数 \(z=a+bi\)，\(a\) 为实部，\(b\) 为虚部。

复平面：在一个二维平面中，横轴表示实部，纵轴表示虚部，即 \((a,b)\)。

欧拉公式：

\[e^{i\theta}=\cos\theta+i\sin\theta \]

证明：

首先可以发现证明这个东西就等价于证明 \(\frac{\cos\theta+i\sin\theta}{e^{i\theta}}=1\) 恒成立，显然有一个 \(1\) 即要证明 \(f(\theta)=\frac{\cos\theta+i\sin\theta}{e^{i\theta}}\) 的 \(f'(\theta)=0\) 恒成立。

\[\begin{aligned}(e^{i\theta})'&=ie^{i\theta}\\f'(\theta)&=\frac{-\sin\theta e^{i\theta}+i\sin\theta\cdot ie^{i\theta}}{e^{2i\theta}}\\&=0\end{aligned} \]

证毕。

单位根：

若复数 \(\omega_n=1\) 则称 \(\omega\) 是 \(n\) 次单位根。

这个 \(\omega_n\) 可以形成一个类似于环的东西即 \(\omega_n^0=\omega_n^n\) （有没有点像欧拉定理 \(\varphi\) 那个）

一般的在 FFT 中 \(\omega_n=e^{\frac{2\pi i}{n}}\) 而在 NTT 中一般为 \(3\)（模数为 \(998244353\)）

点值表示法：对于一个多项式存系数 \((a_0,a_1,...,a_n)\) 叫做系数表示法，若我们对于多项式 \(F(x)\) 对于每个不同的 \(x\) 都直接求出他的 \(F(x)\) 的真实值，那么这个多项式可以表示为 \((x_0,F(x_0)),(x_1,F(x_1)),...,(x_n,F(x_n))\) 那么这个就叫做点值表示法。

FFT

模板

套用点值表示法，如果把两个多项式 \(F(x),G(x)\) 都换成点值表示法，那么 \(H(x)=F(x)G(x)\) 的点值表示法就是他们两个点值表示法乘起来，然后再转化为点值表示法，这个就是 FFT 和 NTT 的基本思想。

首先观察多项式 \(F(x)=a_0+a_1x^1+...+a_{n-1}x^{n-1}\)

如果做可爱的因式分解多的话第一眼就能想到想去给他奇偶划分，然后再套用单位根去思考，具体的，我们划分为 \((a_0+a_2x^2+...)+x(a_1+a_3x^2+...)\)

即 \(F(x)=F_0(x^2)+xF_1(x^2)\)

那么我们发现原问题分裂成两个原问题了，所以我们要这么做，但是原问题得计算量还是没有缩小。

将 \(\omega_n^k\) 套入得

\[\begin{aligned}F(\omega_n^k)&=F_0(\omega_n^{2k})+\omega_n^kF_1(\omega_n^{2k})\\&=F_0(\omega_{\frac{n}{2}}^{k})+\omega_n^kF_1(\omega_{\frac{n}{2}}^{k})\end{aligned} \]

将 \(\omega_n^{k+\frac{n}{2}}\) 套入得

\[\begin{aligned}F(\omega_n^{k+\frac{n}{2}})&= F_0(\omega_n^{2k+n})+\omega_n^{k+\frac{n}{2}}F_1(\omega_n^{2k+n})\\ &=F_0(\omega_{\frac{n}{2}}^{k})-\omega_n^kF_1(\omega_{\frac{n}{2}}^{k})\end{aligned}\]

可以发现系数不同，所以我们的问题就缩小了一半了（对应：\(k\rightarrow k+\frac{n}{2}\)）

到这里我们对这个由系数表示法转化为点值表示法的问题已经能 \(O(n\log n)\) 做了，但是对于点值表示法转化为系数表示法还是不知道。

我们构造一个方法重塑，设 \(F(x)\) 的点值表示为 \((x_0,x_1,x_2,...,x_{n-1})\)（\(x\) 分别等于 \((\omega_n^0,\omega_n^1,...\omega_n^{n-1})\) 一一对应）

定义 \(c_k=\sum\limits_{i=0}^{n-1}x_i(\omega_n^{-k})^i\)

有 \(c_k=\sum\limits_{i=0}^{n-1}\sum\limits_{j=0}^{n-1}a_j(\omega_n^{j-k})^i=\sum\limits_{j=0}^{n-1}a_j(\sum\limits_{i=0}^{n-1}(\omega_n^{j-k})^i)\)

对于 \(\sum(\omega_n^k)^i\) 且 \(k\not=0\) 显然有，这个东西等于 \(0\) （等比数列求和和性质可得）。

对于 \(k=0\) 这个式子等于 \(n\)

所以只有 \(k=j\) 的东西被保留下来了。

所以我们只需要给他做一个负的单位根然后再除以 \(n\) 就可以了。

我们可以递归写，但是递归常数太大，至于递归代码：略。

这里将给出非递归写法：

\(\mathscr{Code}\)

const int N=1<<21;
const double pi=acos(-1); // 并不懂（没学过几何），知道表示 pi 就行，似乎是 cos 的反向吧
struct cd{ // 手写复数可以减少很多常数似乎，为什么这么写自己推一推吧
    double x,y;
    cd(double _x=0,double _y=0){x=_x,y=_y;}
    cd operator-(const cd&b)const{return {x-b.x,y-b.y};}
    cd operator+(const cd&b)const{return {x+b.x,y+b.y};}
    cd operator*(const cd&b)const{return {x*b.x-y*b.y,x*b.y+y*b.x};}
}f[N],g[N];
int n,m,rev[N],len;
void fft(cd *x,int op){
    for(int i=0;i<len;++i)if(i<rev[i])swap(x[i],x[rev[i]]); // 交换操作
    for(int h=2;h<=len;h<<=1){ // 表示长度即 n 
        cd wmi(cos(op*2*pi/h),sin(op*2*pi/h)); // cos 的 op 可以省略
        for(int i=0;i<len;i+=h){
            cd w(1,0);
            for(int j=0;j<h/2;++j){
                cd u=x[i+j],t=x[i+j+h/2]*w;
                x[i+j]=u+t,x[i+j+h/2]=u-t;
                w=w*wmi;
            }
        }
    }
    if(op==-1){
        for(int i=0;i<len;++i)
            x[i].x/=len;
    }
}
void Solve(){
    rd(n,m);
    for(int i=0;i<=n;++i)scanf("%lf",&f[i].x),f[i].y=0;
    for(int i=0;i<=m;++i)scanf("%lf",&g[i].x),g[i].y=0;
    len=1<<int(__lg(n+m))+1;
    for(int i=0;i<len;++i)rev[i]=(rev[i>>1]>>1)|((i&1)*len>>1); // 这里是对其操作，可以看成二进制，好理解的，不过这个码风不优良后面会给出好的
    fft(f,1),fft(g,1);
    for(int i=0;i<len;++i)f[i]=f[i]*g[i];
    fft(f,-1);
    for(int i=0;i<=n+m;++i)printf("%d ",(int)round(f[i].x)); // 一定要四舍五入，原因忘了
}

至此 FFT 终于讲完了 /ll

NTT

和 FFT 差不多思路一毛一样就是换了一个单位根，这俩最大区别就是一个是不取模，一个是取模的了。

模板

首先我们可以用原根去代替单位根，原根是啥？先鸽着。

先说阶： \(a^k\equiv1\pmod m\)

那么满足这个条件的最小正整数 \(k\)，叫做 \(a\) 模 \(m\) 的阶。

记作 \(k=\operatorname{ord}_m(a)\)

原根：首先我们知道欧拉定理 \(a^{\varphi(m)}\equiv 1\pmod m\) （这玩意群论可以成环，是不是很像单位根了）。

如果一个数 \(g\) 的阶等于 \(\varphi(m)\) 那么 \(g\) 是模 \(m\) 的原根。

在这里 \(\omega_n^k=g^{(p-1)\cdot \frac{k}{n}}\bmod p\)

那么我们考虑一下他能替换单位根的必要条件是啥（也就是 FFT 中使用到了单位根的啥）

这 \(n\) 个数互不相等
\(\omega_{2n}^{2k}=\omega_n^k\)
\(\omega_n^{k+\frac{n}{2}}=-\omega_n^k\)

\[\begin{aligned}\omega_n^{k+\frac{n}{2}}=-\omega_n^k\Rightarrow \omega_n^{\frac{n}{2}}=-1\Rightarrow g^{\frac{1}{2}(p-1)}=-1\end{aligned} \]

那么对于 \(g^\frac{p-1}{2}\) 可能有两种解，但是发现如果等于 \(1\) 是谬论，所以得证。

\(\sum\limits_{k\not=0}(\omega_n^k)^i=0\)

对于 NTT 如果模数是 \(998244353=7\times 17\times 2^{23}+1\) 的话才能用，因为 \(2\) 真的很多，这个东西的原根是 \(3\)。

至此终于讲完 NTT 了，下面给出代码？

其实这里有两种写法，就是对于点指表示转化为系数表示的时候，这里讲一下常数小的做法，其实就是小技巧吧，给出代码。

\(\mathscr{Code:}\)

#define ll long long
const int mod=998244353;
const int N=1<<20;
int n,m,len=1,rev[N];
ll f[N],g[N];
ll qpow(ll a,int b){
    ll res=1;
    while(b){
        if(b&1)res=res*a%mod;
        a=a*a%mod;
        b>>=1;
    }
    return res;
}
void ntt(ll *x,int on){
    for(int i=0;i<len;++i)if(i<rev[i])swap(x[i],x[rev[i]]);
    for(int h=2;h<=len;h<<=1){
        ll gn=qpow(3,(mod-1)/h); // 这里先补药疑惑，看后面
        for(int i=0;i<len;i+=h){
            ll g=1;
            for(int j=0;j<h/2;++j){
                int u=x[i+j];
                int t=x[i+j+h/2]*g%mod;
                x[i+j]=(u+t)%mod;
                x[i+j+h/2]=(u-t)%mod;
                g=g*gn%mod;
            }
        }
    }
    if(on==-1){
        reverse(x+1,x+len);
        // 这么倒过来就是因为对称性质，注意都是整体式子的点值表示，成立的原因就是因为有 g^{-ik}=g^{(p-i)k} 注意是交换哪里啊
        ll inv=qpow(len,mod-2);
        for(int i=0;i<len;++i)
            x[i]=x[i]*inv%mod;
    }
}
inline void Solve(){
    n=rd(),m=rd();
    for(int i=0;i<n;++i)f[i]=rd();
    for(int i=0;i<m;++i)g[i]=rd();
    while(len<=n+m-2)len<<=1;
    for(int i=0;i<len;++i)rev[i]=rev[i>>1]>>1|(i&1)*(len>>1);
    ntt(f,1),ntt(g,1);
    for(int i=0;i<len;++i)f[i]=f[i]*g[i]%mod;
    ntt(f,-1);
    for(int i=0;i<=n+m-2;++i)printf("%lld ",(f[i]+mod)%mod);
}

至此，终于把这两坨答辩写完了啊 /ll ~~虽然我的代码更答辩就是了~~

多项式乘法逆

泰勒展开（证明）：

泰勒展开的核心思想是：用函数在某一点的各阶导数信息，构造一个多项式，使它和原函数在该点附近完全一致。

首先泰勒展开：\(f(x)=\sum\limits_{k=0}^\infty\frac{f^{(k)}(a)}{k!}(x-a)^k\)

设展开式：\(f(x)=\sum\limits_{k=0}^\infty c_k(x-a)^k\) （注意 \(c\) 和 \(a\) 都是常量）

求 \(c_0\) 显然为 \(c_0=f(a)\)
求 \(c_1\) 导数一次然后再带入得 \(c_1=f'(a)\)
求 \(c_2\) 导数两次然后再带入得 \(c_1=\frac{f''(a)}{2!}\)
求 \(c_k\) 导数 \(k\) 次然后再带入得 \(c_k=\frac{f^{(k)}(a)}{k!}\)

即得到原式：\(f(x)=\sum\limits_{k=0}^\infty c_k(x-a)^k\)

严格数学证明还需要证明余项趋于 \(0\)，怎么证明？不会

复合多项式：\((f\circ g)(x)=f_0+\sum\limits_{k=1}^\infty f_kg^k(x)\)（这里 \(f\) 和 \(g\) 都是幂级数，幂级数就是一种特殊的无限多项式）

复合多项式有 \(n\log^2n\) 的做法，自己看 wiki 吧

乘法逆元：

例题

首先我们有

\[\frac{1}{1-x}=1+x+x^2+... \]

证明显然。

所以我们可以发现多项式的倒数是是能展开成无穷个的，存在倒数当且仅当常数项不为 \(0\)，倒数也满足常数项不为 \(0\)。

我们设多项式 \(F\) 的倒数为 \(G=F^{-1}\) 即 \(FG=1\)

首先我们对于这个问题他的常数项就是 \(F_0G_0=1\)

对于第 \(x^n\) 项的系数我们必然是 \(\sum\limits_{k=0}^nF_kG_{n-k}=0\) 即 \(\sum\limits_{k=}^nF_kG_{n-k}+F_0G_n=0\)

那么我们便有递推式：\(F^{-1}_0=\frac{1}{F_0},F_n^{-1}=\frac{-1}{F_0}\sum\limits_{k=0}^{n-1}F_k^{-1}F_{n-k}\)

有两种推导方式，先说牛顿迭代：

首先我们考虑在 \(x_n\) 处的泰勒展开有

\[f(x)=f(x_n)+f'(x_n)(x-x_n)+O((x-x_n)^2) \]

先忽略 \(O((x-x_n)^2)\) 带入 \(f(x)=0\) 得

\[f(x_n)+f'(x_n)(x-x_n)=0 \]

即：

\[x=x_n-\frac{f(x)}{f'(x_n)} \]

那么显然得我们可以由 \(x_n\) 推到 \(x_{n+1}\) 然后当 \(n\rightarrow \infty\) 的时候这个东西就是对的，即 \(\lim\limits_{n\rightarrow\infty}f(x_n)\rightarrow0\)。

注意到我们每个 \(x_n\) 基本上都是估计值啊。

这个东西这里写的并不严谨与详细，可以自己搜搜。

其实牛顿迭代每一次都是增加正确位数。

我们要求 \(F(x)\equiv\frac{1}{G(x)}\) 这个东西

设 \(H(t)\equiv F(x)-\frac{1}{t}\)

那么有 \(H(G(x))\equiv0\)

套入公式得

\[G'\equiv G-\frac{H(G)}{H'(G)} \]

有

\[G'\equiv G-\frac{F-\frac{1}{G}}{\frac{1}{G^2}} \]

算一下得

\[G'\equiv G(2-GF) \]

到这里我们已经会了枚举 \(m\) 次的做法。

但是我们根据同余的性质可以发现，设模 \(x^n\) 的逆元为 \(B\)，模 \(x^{2n}\) 的逆元为 \(B'\)

则有 \(AB'\equiv1\pmod{x^n}\) 和 \(AB\equiv1\pmod {x^{n}}\)

那么必然有 \(B\equiv B'\pmod{x^n}\)

所以当我们进行一次迭代的时候都可以看作乘了 \(2\) 那么我们就只进行 \(\log n\) 次操作就可以实现。

倍增法

由 \(B'-B\equiv0\pmod{x^n}\)

推出 \((B'-B)^2\equiv 0\pmod{x^{2n}}\)

推一下式子得：

\[\begin{aligned}B'^2-2B'B+B^2&\equiv0\pmod{x^{2n}}\\AB'^2-2AB'B+AB^2&\equiv0\pmod{x^{2n}}\\ B'-2B+AB^2&\equiv0\pmod{x^{2n}}\\ B'&\equiv B(2-AB)\pmod{x^{2n}}\end{aligned}\]

那么剩下的和牛顿迭代一个样子了。

\(\mathscr{Code:}\)

static inline Poly inv(const Poly& a) {
	int n = a.p.size();
	if (n == 0) return Poly();
	if (n == 1) return Poly(vi{qpow(a.p[0], mod - 2)});
	int mid = (n + 1) >> 1;
	Poly b = inv(Poly(vi(a.p.begin(), a.p.begin() + mid)));
	Poly two(vi{2}), tmp = a * b; tmp.p.resize(n, 0);
	Poly diff = two - tmp; diff.p.resize(n, 0);
	Poly res = b * diff; res.p.resize(n, 0);
	return res;
}

poly 优美模板

~~上面的没看懂？没关系把模板抄下来吧~~

\(\mathscr{Code:}\)

const int N = 1 << 20!!!!!!!!;
int rev[N];
struct Poly {
	static constexpr int mod = 998244353;
	ll& operator[](int idx) { return p[idx]; }
	const ll& operator[](int idx) const { return p[idx]; }
	using ll = long long;
	using vi = vector<ll>;
	vi p;
	Poly(int n) {p.resize(n); } // 注意这里 {n} 会进来，注意啊！！！！！！！！！！ 
	Poly(vi _p = {}) {p = _p; }
	static inline ll qpow(ll a, int b) {
		ll res = 1;
		while (b) {
			if (b & 1) res = res * a % mod;
			a = a * a % mod;
			b >>= 1;
		}
		return res;
	}
	static inline void NTT(Poly& x, int len, int on) {
		int L = 0; while ((1 << L) < len) L++;
		for(int i = 0; i < len; ++i) {
			rev[i] = (rev[i >> 1] >> 1) | ((i & 1) << (L - 1));
			if(i < rev[i]) swap(x.p[i], x.p[rev[i]]);
		}
		for(int h = 2; h <= len; h <<= 1){
			int gn = qpow(3, (mod - 1) / h);
			for(int i = 0; i < len; i += h) {
				int g = 1;
				for(int j = 0; j < h / 2; ++j){
					int u = x.p[i + j], t = 1LL * x.p[i + j + h / 2] * g % mod;
					(x.p[i + j] = u + t) >= mod ? x.p[i + j] -= mod : 0;
					(x.p[i + j + h / 2] = u - t) < 0 ? x.p[i + j + h / 2] += mod : 0;
					g = 1LL * g * gn % mod;
				}
			}
		}
		if(on == -1) {
			reverse(x.p.begin() + 1, x.p.begin() + len);
			int inv = qpow(len, mod - 2);
			for(int i = 0; i < len; ++i) x.p[i] = 1LL * x.p[i] * inv % mod;
		}
	}
	Poly operator+ (const Poly& b) const {
		Poly c(p); c.p.resize(max(p.size(), b.p.size()));
		for (int i = 0; i < b.p.size(); ++i) (c.p[i] += b.p[i]) >= mod ? c.p[i] -= mod : 0;
		return c;
	}
	Poly operator- (const Poly& b) const {
	    Poly c(p); c.p.resize(max(p.size(), b.p.size()));
	    for (int i = 0; i < b.p.size(); ++i) (c.p[i] -= b.p[i]) < 0 ? c.p[i] += mod : 0;
	    return c;
	}
	Poly operator* (const Poly& b) const {
		if (p.empty() || b.p.empty()) return Poly({0});
		int n = p.size(), m = b.p.size();
		int len = 1;
		while (len < n + m - 1) len <<= 1;
		Poly A(p), B(b.p);
		A.p.resize(len, 0); NTT(A, len, 1);
		B.p.resize(len, 0); NTT(B, len, 1);
		Poly C; C.p.resize(len, 0);
		for (int i = 0; i < len; ++i) C[i] = 1LL * A[i] * B[i] % mod;
		NTT(C, len, -1);
		C.p.resize(n + m - 1);
		return C;
	}
	int size() {return p.size();}
	static inline Poly inv(const Poly& a) {
		int n = a.p.size();
		if (n == 0) return Poly();
		if (n == 1) return Poly(vi{qpow(a.p[0], mod - 2)});
		int mid = (n + 1) >> 1;
		Poly b = inv(Poly(vi(a.p.begin(), a.p.begin() + mid)));
		Poly two(vi{2}), tmp = a * b; tmp.p.resize(n, 0);
		Poly diff = two - tmp; diff.p.resize(n, 0);
		Poly res = b * diff; res.p.resize(n, 0);
		return res;
	}
};

~~暂时扔个这个半成品，还有一堆没写~~

生成函数

又称母函数，是一种形式幂级数，其每一项的系数可以提供关于这个序列的信息。

生成函数一般表示成这个样子：

\[F(x)=\sum a_nk_n(x) \]

普通生成函数：\(k_n(x)=x^n\)
指数生成函数：\(k_n(x)=\frac{x^n}{n!}\)
狄利克雷生成函数：\(k_n(x)=\frac{1}{n^x}\)

posted @ 2026-08-02 11:49 static_inline 阅读(6) 评论(0) 收藏举报

刷新页面返回顶部

state

Poly

poly（多项式）

FFT

NTT

多项式乘法逆

poly 优美模板

生成函数

公告