P12416 的多项式推法
不知名用户
·
·
题解
题意
已知非负整数数列 b_1,b_2,\cdots,b_n 的值，另有非负整数数列 a 满足 a_1+a_2+\cdots+a_n = m，a_1\le a_2\le\cdots\le a_n。请求出对于所有满足条件的数列 a，a_1b_1+a_2b_2+\cdots+a_nb_n 的和。由于答案可能很大，你只需要输出答案对 10^9+7 取模的结果。
### 做法
提供一个无脑多项式推法，不过最后可以避开多项式操作。
感谢 @[Konata28](https://www.luogu.com.cn/user/177029) 提出本做法，@[Milmon](https://www.luogu.com.cn/user/234641) 编写了该做法题解。本人将该题解做了很多补充，写成了下面的题解。
首先做转化，问题变为把 $m$ 拆分为若干个 $\leq n$ 的正整数的和，正整数 $i$ 的权值为 $s_i=b_{n-i+1}+b_{n-i+2}+\cdots+b_n$，需要求出所有拆分的权值之和。意思是把 $a$ 视作一些后缀叠起来的。
考虑写出其生成函数，但是发现权值不太好记录到系数里面，因为系数是做乘法而权值是做加法，而在值域做乘法下定义域做加法的函数是指数函数！所以再引入变量 $y$，用 $y$ 的幂表示权值。记二元生成函数 $F(x,y)=\prod\limits_{i=1}^n \dfrac{1}{1-x^iy^{s_i}}$，则只需要求所有 $x^m$ 项的 $y$ 指数与系数的积总和，等价于对 $y$ 求导之后带入 $y=1$ 的 $x^n$ 项系数，这个可以写出 $x,y$ 都是非负指数的形式后看出。
形如 $\dfrac{1}{1-x}$ 的式子的乘积可以转化成先求 $\ln$ 做加法再 $\exp$ 回去。显然不能一个一个做 $\ln\dfrac{1}{1-x}$，但是可以先对其进行一些推导加速运算。设 $f(x)=\ln\dfrac{1}{1-x}$，对 $f(x)$ 求导再积分：
$$
\begin{aligned}
\ln\dfrac{1}{1-x}&=f(x)\\
(-(-\frac{1}{(1-x)^2}))\cdot\frac{1}{\frac{1}{1-x}}&=f'(x)\\
\frac{1}{1-x}&=f'(x)\\
\sum_{i=0}^{\infty}x^i&=f'(x)\\
\sum_{i=1}^{\infty}\frac{1}{i}x^i&=f(x)
\end{aligned}
$$
（为了方便大家理解，展开了复合函数求导法则）
将上式的 $x$ 改成 $x^iy^{s_i}$ 代入：
$$
\exp\sum_{i=1}^n\ln\dfrac{1}{1-x^iy^{s_i}}
=\exp\sum_{i=1}^n\sum_{j=1}^{\infty}\dfrac{x^{ij}y^{s_ij}}{j}
$$
根据复合函数求导法则和 $(\exp x)'=\exp x$，$(\exp f(x))'=f'(x)\exp f(x)$ 以及函数加法导数加法原则，对 $y$ 求导，再代入 $y=1$ 的结果为：
$$
\left(\sum_{i=1}^n\sum_{j=1}^{\infty}s_ix^{ij}\right)\exp\left(\sum_{i=1}^n\sum_{j=1}^{\infty}\frac{1}{j}x^{ij}\right)
$$
只需求出 $x^m$ 的系数即可。因为只关心 $x^m$ 的系数，所以不需要额外的多项式乘法，直接枚举前面的 $s_ix^{ij}$ 然后和后半部分 $x^{m-ij}$ 的系数相乘即可。
后面的式子处理出来是调和级数复杂度 $\Theta(m\log n)$ 的，$\exp$ 复杂度是 $\Theta(m\log m)$。总时间复杂度 $\Theta(m(\log m+\log n))$。
需自备模数为 $10^9+7$ 的多项式 exp 板子。虽然复杂度很低，但 MTT exp 常数很大，可能需要比较优的板子。不过出这题本意就是不想让多项式过的。省去多项式板子的代码：
```cpp
#include<bits/stdc++.h>
using namespace std;
const int N = 1.1e5 + 10, mod = 1e9 + 7;
int a[N], s[N], inv[N];
int main()
{
int n, m, i, j, ans = 0;
scanf("%d%d", &n, &m);
for(i=1;i<=n;i++) scanf("%d", &a[n-i+1]);
for(i=1;i<=n;i++) s[i] = (s[i-1] + a[i]) % mod;
poly p;
p.resize(m+1);
inv[1] = 1;
for(i=2;i<=m;i++) inv[i] = 1ll * inv[mod%i] * (mod - mod / i) % mod;
for(i=1;i<=n;i++)
for(j=i;j<=m;j+=i)
p[j] = (p[j] + inv[j/i]) % mod;
exp(p);
for(i=1;i<=n;i++)
for(j=i;j<=m;j+=i)
ans = (1ll * s[i] * p[m-j] + ans) % mod;
printf("%d", ans);
return 0;
}
```
### 如何避开大常数 MTT exp
观察 $\displaystyle\exp\left(\sum_{i=1}^n\sum_{j=1}^{\infty}\frac{1}{j}x^{ij}\right)$ 括号里面的式子，这就是 $\displaystyle\sum\limits_{i=1}^{\infty}\ln\frac{1}{1-x^i}$！这个式子就是商品体积为 $1,2,\cdots,n$，凑得一个价格的方案数的生成函数（[P4389](https://www.luogu.com.cn/problem/P4389)，没理解可以去翻那题的题解）。在 $m=n$ 的时候这就是[拆分数](https://www.luogu.com.cn/problem/P6189)，$m>n$ 的时候可以先跑 $m$ 的拆分数再把 $n+1,n+2,\cdots,m$ 从计数背包撤销，这样就用 $\Theta(m\sqrt m+(m-n)^2)$ 的复杂度代替了多项式 $\exp$ 操作。总复杂度 $\Theta(m\sqrt m+(m-n)^2+m\log n)$。
```cpp
#include<bits/stdc++.h>
using namespace std;
const int N = 1.1e5 + 10, sqN = 350, mod = 1e9 + 7;
int f[sqN+5][N], p[N], s[N];
signed main()
{
int n, m, i, j, ans = 0;
scanf("%d%d", &n, &m);
for(i=1;i<=n;i++) scanf("%d", &s[n-i+1]);
for(i=1;i<=n;i++) s[i] = (s[i] + s[i-1]) % mod;
int sqm = (int)sqrt(m);
f[0][0] = p[0] = 1;
for(i=1;i<=m/sqm;i++) for(j=0;j<=m;j++)
{
if(j>=i) f[i][j] = f[i][j-i];
if(j>=sqm) f[i][j] = (f[i][j] + f[i-1][j-sqm]) % mod;
p[j] = (p[j] + f[i][j]) % mod;
}
for(i=1;i<sqm;i++) for(j=i;j<=m;j++) p[j] = (p[j] + p[j-i]) % mod;
for(i=m;i>n;i--) for(j=m;j>=i;j--) p[j] = (p[j] - p[j-i] + mod) % mod;
for(i=1;i<=n;i++)
for(j=i;j<=m;j+=i)
ans = (ans + 1ll * s[i] * p[m-j] % mod) % mod;
printf("%d", ans);
return 0;
}
```
### 花絮
对每个 $s_i$ 算贡献不使用多项式也可以做，代码和上面是一样的。但是如果要用 MTT exp 的话对多项式板子要求非常严格，也许没经过任何卡常的板子跑不过去。
其实最开始出出来还是暴力 $\Theta(nm)$ 背包做的，然后我发现物品是独立的，可以算贡献，但是要先求出来 $1,2,\cdots,n$ 去凑出 $1,2,\cdots,m$ 的方案数，具体参见[这篇题解](https://www.luogu.com.cn/article/2eztwi30)，和他做法一样。这可以拆分数+撤销也可以直接套用 P4389 的 poly 做法，推出来的式子跟上面的式子是一样的。验题人提出了本篇题解的做法。