解题报告（二）多项式问题（多项式乘法及其各种运算）（ACM/ OI）超高质量题解

整理的算法模板合集： ACM模板

点我看算法全家桶系列！！！

实际上是一个全新的精炼模板整合计划

繁凡出品的全新系列：解题报告系列 —— 超高质量算法题单，配套我写的超高质量的题解和代码，题目难度不一定按照题号排序，我会在每道题后面加上题目难度指数（1∼51 \sim 51∼5），以模板题难度 111 为基准。

这样大家在学习算法的时候就可以执行这样的流程：

%
阅读【学习笔记】 / 【算法全家桶】学习算法 ⇒\Rightarrow⇒ 阅读相应算法的【解题报告】获得高质量题单 ⇒\Rightarrow⇒ 根据一句话题解的提示尝试自己解决问题 ⇒\Rightarrow⇒ 点开详细题解链接学习巩固（好耶）
%
要是26个英文字母用完了我就接上24个希腊字母，我就不信50道题不够我刷的hhh

%
解题报告系列合集：【解题报告系列】超高质量题单 + 题解（ICPC / CCPC / NOIP / NOI / CF / AT / NC / P / BZOJ）

本题单前置知识：【学习笔记】多项式全家桶（包含全套证明）

FFT备忘录：

两个多项式卷积之后的多项式的长度是 n+mn+mn+m 或者 2∗n2*n2∗n ！！！

最后除以 limitlimitlimit，以及 FFT 的四舍五入

生成函数相关的计数问题，尽量都开 long long！

A、【模板】多项式全家桶（1）
B、（P3338 [ZJOI2014]）力（多项式乘法FFT）（2）
C、（BZOJ 2194）快速傅立叶之二（FFT、卷积的概念、常用变换）
D、（BZOJ 3771）Triple（生成函数 + FFT）（3）
E、（BZOJ3513） [MUTC2013] idiots（FFT + 组合计数）
F、（CF438E） The Child and Binary Tree（生成函数、多项式求逆、开方）(3.5)
G、（BZOJ 3509 [CodeChef]）COUNTARI（FFT + 分块优化）

本专题题单主要收录多项式乘法相关的题目，以及多项式各种运算的应用如多项式除法，多项式快速幂，多项式对数指数等等。多项式除了会在一些公式题中直接只用以外，最常用的场景就是生成函数里，用于计算计数问题中得到的生成函数的卷积。根据推出的公式的不同会运用到不同的多项式运算。

这里简单复习一些重要的基础概念。

多项式乘法：

广义上的卷积为：

h(x)=∫−∞∞g(τ)⋅f(x−τ)dτ{h(x) = \int _{- \infty} ^{\infty}g(\tau) \cdot f(x - \tau)} \rm{d\tau} h(x)=∫−∞∞g(τ)⋅f(x−τ)dτ

我们这里仅讨论多项式域，对于两个多项式 / 数组 / 序列的卷积。

即给定两个多项式 f(x)f(x)f(x)，g(x)g(x)g(x) ：

f(x)=a0+a1x+⋯+anxnf(x)= a_0 + a_1x+ \dots +a_nx^nf(x)=a0+a1x+⋯+anxn

g(x)=b0+b1x+⋯+bmxmg(x)= b_0 + b_1x+ \dots +b_mx^mg(x)=b0+b1x+⋯+bmxm

计算多项式 Q(x)=f(x)⋅g(x):Q(x) = f(x)\ ·\ g(x):Q(x)=f(x) ⋅ g(x):

Q(x)=∑i=0n∑j=0maibjxi+j=c0+c1x+⋯+cn+mxn+m\boxed {Q(x) = \sum \limits_ {i = 0} ^ n \sum \limits_ {j = 0 } ^ m a_i b_j x ^ {i + j}} = c_0 + c_1 x + \dots + c_ {n + m} x ^ {n + m}Q(x)=i=0∑nj=0∑maibjxi+j=c0+c1x+⋯+cn+mxn+m

也可以写成这种形式：

C(x)=A(x)∗B(x)=∑k=02n−2(∑k=i+jaibj)xkC(x)=A(x)* B(x)=\sum_{k=0}^{2n-2}(\sum_{k=i+j}a_ib_j)x^k C(x)=A(x)∗B(x)=k=0∑2n−2(k=i+j∑aibj)xk

C(x)=A(x)∗B(x)=∑i=02×n(∑j=0min{i,n}ajbi−j)xj+i−j=∑i=02×n(∑j=0min{i,n}ajbi−j)xiC(x)=A(x)*B(x)=\sum_{i = 0}^{ 2\times n}(\sum_{j=0}^{min\{i,n\}}a_j b_{i-j})x^{j+i-j}=\sum_{i = 0}^{ 2\times n}(\sum_{j=0}^{min\{i,n\}}a_j b_{i-j})x^{i} C(x)=A(x)∗B(x)=i=0∑2×n(j=0∑min{i,n}ajbi−j)xj+i−j=i=0∑2×n(j=0∑min{i,n}ajbi−j)xi

也就是对于 C(x)C(x)C(x) 的第 iii 项的系数：

[i]C(x)=∑j=0iajbi−j[i]C(x)=\sum_{j=0}^{i}a_j b_{i-j} [i]C(x)=j=0∑iajbi−j

可以理解为只有小于 iii 的系数才有可能给 xix^ixi 做贡献，大的乘起来就大于 iii 了。

注意对于两个长度均为 nnn 的多项式 AAA 和 BBB ，卷积卷出来的多项式 CCC 的长度是 n+n=2×nn+n=2\times nn+n=2×n。切记，切记…

其中卷积的系数表示：

y3i=(∑j=02n−1ajxij)×(∑j=02n−1bjxij)=y1i×y2i{y_3}_i=(\sum_{j=0}^{2n-1}a_jx_i^j)\times(\sum_{j=0}^{2n-1}b_jx_i^j)={y_1}_i\times{y_2}_i y3i=(j=0∑2n−1ajxij)×(j=0∑2n−1bjxij)=y1i×y2i

A、【模板】多项式全家桶（1）

详细证明见：【学习笔记】多项式全家桶！！！（公式模板大全）

Code

//#pragma GCC optimize(2)
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;const int N = 3000007;
const int p = 998244353, gg = 3, ig = 332738118, img = 86583718;
const int mod = 998244353;template <typename T>void read(T &x)
{x = 0;register int f = 1;register char ch = getchar();while(ch < '0' || ch > '9') {if(ch == '-')f = -1;ch = getchar();}while(ch >= '0' && ch <= '9') {x = x * 10 + ch - '0';ch = getchar();}x *= f;
}int qpow(int a, int b)
{int res = 1;while(b) {if(b & 1) res = 1ll * res * a % mod;a = 1ll * a * a % mod;b >>= 1;}return res;
}namespace Poly
{#define mul(x, y) (1ll * x * y >= mod ? 1ll * x * y % mod : 1ll * x * y)#define minus(x, y) (1ll * x - y < 0 ? 1ll * x - y + mod : 1ll * x - y)#define plus(x, y) (1ll * x + y >= mod ? 1ll * x + y - mod : 1ll * x + y)#define ck(x) (x >= mod ? x - mod : x)//取模运算太慢了typedef vector<int> poly;const int G = 5;const int inv_G = qpow(G, mod - 2);int RR[N], deer[2][19][N], inv[N];void init(const int t) {//预处理出来NTT里需要的w和wn，砍掉了一个log的时间for(int p = 1; p <= t; ++ p) {int buf1 = qpow(G, (mod - 1) / (1 << p));int buf0 = qpow(inv_G, (mod - 1) / (1 << p));deer[0][p][0] = deer[1][p][0] = 1;for(int i = 1; i < (1 << p); ++ i) {deer[0][p][i] = 1ll * deer[0][p][i - 1] * buf0 % mod;//逆deer[1][p][i] = 1ll * deer[1][p][i - 1] * buf1 % mod;}}inv[1] = 1;for(int i = 2; i <= (1 << t); ++ i)inv[i] = 1ll * inv[mod % i] * (mod - mod / i) % mod;}int NTT_init(int n) {//快速数论变换预处理int limit = 1, L = 0;while(limit < n) limit <<= 1, L ++ ;for(int i = 0; i < limit; ++ i)RR[i] = (RR[i >> 1] >> 1) | ((i & 1) << (L - 1));return limit;}void NTT(poly &A, int type, int limit) {//快速数论变换A.resize(limit);for(int i = 0; i < limit; ++ i)if(i < RR[i])swap(A[i], A[RR[i]]);for(int mid = 2, j = 1; mid <= limit; mid <<= 1, ++ j) {int len = mid >> 1;for(int pos = 0; pos < limit; pos += mid) {int *wn = deer[type][j];for(int i = pos; i < pos + len; ++ i, ++ wn) {int tmp = 1ll * (*wn) * A[i + len] % mod;A[i + len] = ck(A[i] - tmp + mod);A[i] = ck(A[i] + tmp);}}}if(type == 0) {for(int i = 0; i < limit; ++ i)A[i] = 1ll * A[i] * inv[limit] % mod;}}poly poly_mul(poly A, poly B) {//多项式乘法int deg = A.size() + B.size() - 1;int limit = NTT_init(deg);poly C(limit);NTT(A, 1, limit);NTT(B, 1, limit);for(int i = 0; i < limit; ++ i)C[i] = 1ll * A[i] * B[i] % mod;NTT(C, 0, limit);C.resize(deg);return C;}poly poly_inv(poly &f, int deg) {//多项式求逆if(deg == 1)return poly(1, qpow(f[0], mod - 2));poly A(f.begin(), f.begin() + deg);poly B = poly_inv(f, (deg + 1) >> 1);int limit = NTT_init(deg << 1);NTT(A, 1, limit), NTT(B, 1, limit);for(int i = 0; i < limit; ++ i)A[i] = B[i] * (2 - 1ll * A[i] * B[i] % mod + mod) % mod;NTT(A, 0, limit);A.resize(deg);return A;}poly poly_dev(poly f) {//多项式求导int n = f.size();for(int i = 1; i < n; ++ i) f[i - 1] = 1ll * f[i] * i % mod;return f.resize(n - 1), f;//f[0] = 0，这里直接扔了,从1开始}poly poly_idev(poly f) {//多项式求积分int n = f.size();for(int i = n - 1; i ; -- i) f[i] = 1ll * f[i - 1] * inv[i] % mod;return f[0] = 0, f;}poly poly_ln(poly f, int deg) {//多项式求对数poly A = poly_idev(poly_mul(poly_dev(f), poly_inv(f, deg)));return A.resize(deg), A;}poly poly_exp(poly &f, int deg) {//多项式求指数if(deg == 1)return poly(1, 1);poly B = poly_exp(f, (deg + 1) >> 1);B.resize(deg);poly lnB = poly_ln(B, deg);for(int i = 0; i < deg; ++ i)lnB[i] = ck(f[i] - lnB[i] + mod);int limit = NTT_init(deg << 1);//n -> n^2NTT(B, 1, limit), NTT(lnB, 1, limit);for(int i = 0; i < limit; ++ i)B[i] = 1ll * B[i] * (1 + lnB[i]) % mod;NTT(B, 0, limit);B.resize(deg);return B;}poly poly_sqrt(poly &f, int deg) {//多项式开方if(deg == 1) return poly(1, 1);poly A(f.begin(), f.begin() + deg);poly B = poly_sqrt(f, (deg + 1) >> 1);poly IB = poly_inv(B, deg);int limit = NTT_init(deg << 1);NTT(A, 1, limit), NTT(IB, 1, limit);for(int i = 0; i < limit; ++ i)A[i] = 1ll * A[i] * IB[i] % mod;NTT(A, 0, limit);for(int i =0; i < deg; ++ i)A[i] = 1ll * (A[i] + B[i]) * inv[2] % mod;A.resize(deg);return A;}poly poly_pow(poly f, int k) {//多项式快速幂f = poly_ln(f, f.size());for(auto &x : f) x = 1ll * x * k % mod;return poly_exp(f, f.size());}poly poly_cos(poly f, int deg) {//多项式三角函数（cos）poly A(f.begin(), f.begin() + deg);poly B(deg), C(deg);for(int i = 0; i < deg; ++ i)A[i] = 1ll * A[i] * img % mod;B = poly_exp(A, deg);C = poly_inv(B, deg);int inv2 = qpow(2, mod - 2);for(int i = 0; i < deg; ++ i)A[i] = 1ll * (1ll * B[i] + C[i]) % mod * inv2 % mod;return A;}poly poly_sin(poly f, int deg) {//多项式三角函数（sin）poly A(f.begin(), f.begin() + deg);poly B(deg), C(deg);for(int i = 0; i < deg; ++ i)A[i] = 1ll * A[i] * img % mod;B = poly_exp(A, deg);C = poly_inv(B, deg);int inv2i = qpow(img << 1, mod - 2);for(int i = 0; i < deg; ++ i)A[i] = 1ll * (1ll * B[i] - C[i] + mod) % mod * inv2i % mod;return A;}poly poly_arcsin(poly f, int deg) {poly A(f.size()), B(f.size()), C(f.size());A = poly_dev(f);B = poly_mul(f, f);for(int i = 0; i < deg; ++ i)B[i] = minus(mod, B[i]);B[0] = plus(B[0], 1);C = poly_sqrt(B, deg);C = poly_inv(C, deg);C = poly_mul(A, C);C = poly_idev(C);return C;}poly poly_arctan(poly f, int deg) {poly A(f.size()), B(f.size()), C(f.size());A = poly_dev(f);B = poly_mul(f, f);B[0] = plus(B[0], 1);C = poly_inv(B, deg);C = poly_mul(A, C);C = poly_idev(C);return C;}
}using Poly::poly;
using Poly::poly_arcsin;
using Poly::poly_arctan;int n, m, x, k, type;
poly f, g;
char s[N];int main()
{Poly::init(18);//2^21 = 2,097,152,根据题目数据多项式项数的大小自由调整，注意大小需要跟deer数组同步(21+1=22)read(n), read(type);for(int i = 0; i < n; ++ i)read(x), f.push_back(x);if(type == 0) g = poly_arcsin(f, n);else g = poly_arctan(f, n);for(int i = 0; i < n; ++ i)printf("%d ", g[i]);return 0;
}

B、（P3338 [ZJOI2014]）力（多项式乘法FFT）（2）

Weblink

https://www.luogu.com.cn/problem/P3338

Problem

给出 nnn 个数 q1,q2,…qnq_1,q_2, \dots q_nq1,q2,…qn ，定义
Fj=∑i=1j−1qi×qj(i−j)2−∑i=j+1nqi×qj(i−j)2Ei=FiqiF_{j}=\sum_{i=1}^{j-1} \frac{q_{i} \times q_{j}}{(i-j)^{2}}-\sum_{i=j+1}^{n} \frac{q_{i} \times q_{j}}{(i-j)^{2}}\\E_{i}=\frac{F_{i}}{q_{i}} Fj=i=1∑j−1(i−j)2qi×qj−i=j+1∑n(i−j)2qi×qjEi=qiFi

对 1≤i≤n1 \leq i \leq n1≤i≤n，求 EiE_iEi的值。

1≤n≤105,0<qi<1091 \leq n \leq 10^{5}, 0<q_{i}<10^{9}1≤n≤105,0<qi<109

Hint

一句话题解： 化简得到 EEE 的式子以后，想办法设两个函数 f,gf,gf,g ，使得 f×gf\times gf×g 可以完全表示这个式子，这样我们，就可以用 FFT 解决了 ~

Solution

超高质量详细题解：点我

C、（BZOJ 2194）快速傅立叶之二（FFT、卷积的概念、常用变换）

Wbelink

https://darkbzoj.tk/problem/2194

Problem

请计算 C[k]=∑(a[i]∗b[i−k])\displaystyle C[k]=\sum(a[i]*b[i-k])C[k]=∑(a[i]∗b[i−k]) 其中 k<=i<nk < = i < nk<=i<n ，并且有 n<=105n < = 10 ^ 5n<=105。 a,ba,ba,b 中的元素均为小于等于 100100100 的非负整数。

n≤105n \leq 10^5n≤105

Solution

超高质量详细题解：点我

D、（BZOJ 3771）Triple（生成函数 + FFT）（3）

Link

https://darkbzoj.tk/problem/3771

Problem

我们讲一个悲伤的故事。
从前有一个贫穷的樵夫在河边砍柴。
这时候河里出现了一个水神，夺过了他的斧头，说：
“这把斧头，是不是你的？”
樵夫一看：“是啊是啊！”
水神把斧头扔在一边，又拿起一个东西问：
“这把斧头，是不是你的？”
樵夫看不清楚，但又怕真的是自己的斧头，只好又答：“是啊是啊！”
水神又把手上的东西扔在一边，拿起第三个东西问：
“这把斧头，是不是你的？”
樵夫还是看不清楚，但是他觉得再这样下去他就没法砍柴了。
于是他又一次答：“是啊是啊！真的是！”
水神看着他，哈哈大笑道：
“你看看你现在的样子，真是丑陋！”
之后就消失了。
樵夫觉得很坑爹，他今天不仅没有砍到柴，还丢了一把斧头给那个水神。
于是他准备回家换一把斧头。
回家之后他才发现真正坑爹的事情才刚开始。
水神拿着的的确是他的斧头。
但是不一定是他拿出去的那把，还有可能是水神不知道怎么偷偷从他家里拿走的。
换句话说，水神可能拿走了他的一把，两把或者三把斧头。
樵夫觉得今天真是倒霉透了，但不管怎么样日子还得过。
他想统计他的损失。
樵夫的每一把斧头都有一个价值，不同斧头的价值不同。总损失就是丢掉的斧头价值和。
他想对于每个可能的总损失，计算有几种可能的方案。
注意：如果水神拿走了两把斧头a和b，(a,b)和(b,a)视为一种方案。拿走三把斧头时，(a,b,c),(b,c,a),(c,a,b),(c,b,a),(b,a,c),(a,c,b)视为一种方案。

Input

第一行是整数N，表示有N把斧头。
接下来n行升序输入N个数字Ai，表示每把斧头的价值。

Output

若干行，按升序对于所有可能的总损失输出一行x y，x为损失值，y为方案数。

Solution

超高质量详细题解：点我

E、（BZOJ3513） [MUTC2013] idiots（FFT + 组合计数）

Weblink

https://darkbzoj.tk/problem/3513

Problem

给定 nnn 个长度分别为 aia_iai 的木棒，问随机选择 333 个木棒能够拼成三角形的概率。

input

第一行T( T≤100T\le 100T≤100 )，表示数据组数。
接下来若干行描述 T 组数据，每组数据第一行是 n ，接下来一行有 n 个数表示 aia_iai。
3≤N≤105,1≤ai≤1053≤N≤10^5,1≤a_i≤10^53≤N≤105,1≤ai≤105

output

T 行，每行一个整数，四舍五入保留7位小数。

Solution

超高质量详细题解：点我

F、（CF438E） The Child and Binary Tree（生成函数、多项式求逆、开方）(3.5)

Weblink

https://www.luogu.com.cn/problem/CF438E

Problem

首先我们发现模数为998244353998244353998244353，很自然的想到NTT！！！想到用多项式来求解本题。

具有同样效果的数字还有：1004535809、469762049…1004535809、469762049 \dots1004535809、469762049…

我们发现题目中所求的是所有的 sss 方案数，而1≤s≤m1\le s\le m1≤s≤m，我们发现这不就是多项式的xix^ixi嘛。

我们设fif_ifi 表示点权和为 iii 的符合条件的二叉树的个数，gig_igi 表示权值 iii 是否包含在 ccc 中

其中很明显 f0=1f_0=1f0=1

我们可以很轻松地根据题意得到这样一个式子：
fn=∑i=1ngi∑j=1n−ifjfn−i−j[n>0]f_n=\sum_{i=1}^ng_i\sum_{j=1}^{n-i}f_jf_{n-i-j}[n>0]fn=i=1∑ngij=1∑n−ifjfn−i−j[n>0]

其实就是一个递推式子，枚举每个权值是否存在（输入的时候存一下就行了）

权值和为 nnn 的方案数 === 权值和为 jjj 的方案数 +++ 权值和为 n−i−jn-i-jn−i−j 的方案数 +gi+\ g_i+ gi

我们只需要令 FFF 表示序列 fff 的生成函数，GGG 表示序列 ggg 的生成函数。那么可以得到：

F=G∗F2+1F=G*F^2+1F=G∗F2+1
由一元二次方程求根公式得方程的根：

F=1±1−4G2GF=\frac{1± \sqrt{1-4G}}{2G}F=2G1±1−4G

此时发现有两个解，分别讨论一下。

当取正号时lim⁡x→0为∞,舍去\text{当取正}\text{号时}\lim_{x→0}\text{为}\infty,\text{舍去}当取正号时limx→0为∞,舍去
当取负号时lim⁡x→0为1,满足\text{当取负}\text{号时}\lim_{x→0}\text{为1},\text{满足}当取负号时limx→0为1,满足

故F=1−1−4G2GF=\frac{1-\sqrt{1-4G}}{2G}F=2G1−1−4G

化简一下得到

F=21+1−4GF=\frac{2}{1+\sqrt{1-4G}}F=1+1−4G2

我们直接按照式子求一下即可，只需要用到一个多项式求逆和多项式开方。

//#pragma GCC optimize(2)
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;const int N = 1000007;
const int p = 998244353, gg = 3, ig = 332738118, img = 86583718;
const int mod = 998244353;template <typename T>void read(T &x)
{x = 0;register int f = 1;register char ch = getchar();while(ch < '0' || ch > '9') {if(ch == '-')f = -1;ch = getchar();}while(ch >= '0' && ch <= '9') {x = x * 10 + ch - '0';ch = getchar();}x *= f;
}int qpow(int a, int b)
{int res = 1;while(b) {if(b & 1) res = 1ll * res * a % mod;a = 1ll * a * a % mod;b >>= 1;}return res;
}namespace Poly
{#define mul(x, y) (1ll * x * y >= mod ? 1ll * x * y % mod : 1ll * x * y)#define minus(x, y) (1ll * x - y < 0 ? 1ll * x - y + mod : 1ll * x - y)#define plus(x, y) (1ll * x + y >= mod ? 1ll * x + y - mod : 1ll * x + y)#define ck(x) (x >= mod ? x - mod : x)//取模运算太慢了typedef vector<int> poly;const int G = 5;const int inv_G = qpow(G, mod - 2);int RR[N], deer[2][19][N], inv[N];void init(const int t) {//预处理出来NTT里需要的w和wn，砍掉了一个log的时间for(int p = 1; p <= t; ++ p) {int buf1 = qpow(G, (mod - 1) / (1 << p));int buf0 = qpow(inv_G, (mod - 1) / (1 << p));deer[0][p][0] = deer[1][p][0] = 1;for(int i = 1; i < (1 << p); ++ i) {deer[0][p][i] = 1ll * deer[0][p][i - 1] * buf0 % mod;//逆deer[1][p][i] = 1ll * deer[1][p][i - 1] * buf1 % mod;}}inv[1] = 1;for(int i = 2; i <= (1 << t); ++ i)inv[i] = 1ll * inv[mod % i] * (mod - mod / i) % mod;}int NTT_init(int n) {//快速数论变换预处理int limit = 1, L = 0;while(limit < n) limit <<= 1, L ++ ;for(int i = 0; i < limit; ++ i)RR[i] = (RR[i >> 1] >> 1) | ((i & 1) << (L - 1));return limit;}void NTT(poly &A, int type, int limit) {//快速数论变换A.resize(limit);for(int i = 0; i < limit; ++ i)if(i < RR[i])swap(A[i], A[RR[i]]);for(int mid = 2, j = 1; mid <= limit; mid <<= 1, ++ j) {int len = mid >> 1;for(int pos = 0; pos < limit; pos += mid) {int *wn = deer[type][j];for(int i = pos; i < pos + len; ++ i, ++ wn) {int tmp = 1ll * (*wn) * A[i + len] % mod;A[i + len] = ck(A[i] - tmp + mod);A[i] = ck(A[i] + tmp);}}}if(type == 0) {for(int i = 0; i < limit; ++ i)A[i] = 1ll * A[i] * inv[limit] % mod;}}poly poly_mul(poly A, poly B) {//多项式乘法int deg = A.size() + B.size() - 1;int limit = NTT_init(deg);poly C(limit);NTT(A, 1, limit);NTT(B, 1, limit);for(int i = 0; i < limit; ++ i)C[i] = 1ll * A[i] * B[i] % mod;NTT(C, 0, limit);C.resize(deg);return C;}poly poly_inv(poly &f, int deg) {//多项式求逆if(deg == 1)return poly(1, qpow(f[0], mod - 2));poly A(f.begin(), f.begin() + deg);poly B = poly_inv(f, (deg + 1) >> 1);int limit = NTT_init(deg << 1);NTT(A, 1, limit), NTT(B, 1, limit);for(int i = 0; i < limit; ++ i)A[i] = B[i] * (2 - 1ll * A[i] * B[i] % mod + mod) % mod;NTT(A, 0, limit);A.resize(deg);return A;}poly poly_sqrt(poly &f, int deg) {//多项式开方if(deg == 1) return poly(1, 1);poly A(f.begin(), f.begin() + deg);poly B = poly_sqrt(f, (deg + 1) >> 1);poly IB = poly_inv(B, deg);int limit = NTT_init(deg << 1);NTT(A, 1, limit), NTT(IB, 1, limit);for(int i = 0; i < limit; ++ i)A[i] = 1ll * A[i] * IB[i] % mod;NTT(A, 0, limit);for(int i =0; i < deg; ++ i)A[i] = 1ll * (A[i] + B[i]) * inv[2] % mod;A.resize(deg);return A;}
}using Poly::poly;
using Poly::poly_sqrt;
using Poly::poly_inv;int n, m, x, k, type;char s[N];
int cnt;int main()
{Poly::init(18);//2^21 = 2,097,152,根据题目数据多项式项数的大小自由调整，注意大小需要跟deer数组同步(21+1=22)read(n), read(m);poly f(N), g(N), tmp(N);for(int i = 0; i < n; ++ i)read(x), g[x] = 1;tmp[0] = 1;for(int i = 1; i <= m; ++ i) tmp[i] = mod - 4 * g[i] % mod;g = poly_sqrt(tmp, m + 1);++ g[0];f = poly_inv(g, m + 1);for(int i = 1; i <= m; ++ i) f[i] = f[i] * 2 % mod;for(int i = 1; i <= m; ++ i) printf("%d\n", f[i]);return 0;
}

G、（BZOJ 3509 [CodeChef]）COUNTARI（FFT + 分块优化）

Weblink

https://darkbzoj.tk/problem/3509

Problem

给定一个长度为 NNN 的数组 AAA，求有多少对 i,j,ki, j, ki,j,k（KaTeX parse error: Undefined control sequence: \lei at position 2: 1\̲l̲e̲i̲<j<k\leN）满足 A[k]−A[j]=A[j]−A[i]A[k]-A[j]=A[j]-A[i]A[k]−A[j]=A[j]−A[i]。

N≤105,A[i]≤30000N\le10^5,A[i]\le 30000N≤105,A[i]≤30000

Solution