[ACNOI2022]Patience

题目

题目背景
哪怕旁边的人说着：“我们是同伴，要将真正的同伴留在心里！” 哪怕旁边的人说着：“尽管这力量很微小，但如果它能左右世界，那我们更应该谨小慎微。” 然而到最后，我也只是作了炮灰。

$\mathsf{O}\mathsf{U}\mathsf{Y}\mathsf{E}$ 挺身而出，挑起大梁；$\mathsf{D}\mathsf{D}(\mathsf{X}\mathsf{Y}\mathsf{X})$ 一往无前，渴望孤独。那么，我呢？失去了 $\mathsf{O}\mathsf{U}\mathsf{Y}\mathsf{E}$ 的「卷坷垃」保护的我，发挥过哪怕一点点的作用吗？

题目描述
面前或许有很多数字，但最终有用的寥寥，就像 $\mathsf{D}\mathsf{D}(\mathsf{X}\mathsf{Y}\mathsf{X})$ 和 $\mathsf{O}\mathsf{U}\mathsf{Y}\mathsf{E}$ 那样。

给出自然数的集合 $A={\bigcup }_{i=1}^n[l_i,r_i]$，求有多少种方式将 $s$ 表示为 $A$ 中 $4$ 个不同数字的和。表示方式不同，当且仅当选出的 $4$ 个数字的集合不同。

数据范围与提示
$s⩽8\times 10^8$ 但 $n⩽800$ 。

思路

说白了是个背包问题，也就是个生成函数问题。写出 $F(x)=\sum \frac{x^{l_i}-x^{r_i+1}}{1-x}$ 是容易的，但是 $F(x{)}^4$ 并非答案，因为其可能选择了相同的元素。

有一个略去思考过程的结论，见此文最后一节。若不知道，那我们可以从头开始，求出每种个数分配对应的生成函数。例如，$G_{1,2}(x)$ 表示某个数选择 $1$ 次、另一个数选择 $2$ 次，对应的生成函数。当然，选择次数相同的数字是本质相同的（没有顺序关系）。那么我们可以写出
$$\begin{array}{rl}{G}_{1,1}(x)& =2^{-1}[F(x{)}^2-F(x^2)]\\ G_{1,1,1}(x)& =3^{-1}[G_{1,1}(x)F(x)-G_{1,2}(x)]\\ G_{1,1,1,1}(x)& =4^{-1}[G_{1,1,1}(x)F(x)-G_{1,1,2}(x)]\end{array}$$

形式是相仿的，应该容易理解，$G_{\cdots ,2}(x)$ 就是出现重复的情况。除以 $2/3/4$ 是取消选择次数相同的数字之间的顺序关系。

类似地，我们可以把需要用到的值都补齐
$$\begin{array}{rl}{G}_{1,2}(x)& =F(x^2)F(x)-F(x^3)\\ G_{1,1,2}(x)& =G_{1,1}(x)F(x^2)-G_{1,3}(x)\\ G_{1,3}(x)& =F(x)F(x^3)-F(x^4)\end{array}$$

然后按照定义，层层展开。过程较长且无技术含量，这里略去。
$$24G_{1,1,1,1}(x)=F(x{)}^4-6F(x^2)F(x{)}^2+3F(x^2{)}^2+8F(x^3)F(x)-6F(x^4)$$

答案就是 $[x^s]G_{1,1,1,1}(x)$ 。对于 $F(x^2{)}^2$ 和 $F(x^3)F(x)$ 和 $F(x^4)$，它们的分子都不超过 $\mathcal{O}(n^2)$ 项，可以暴力卷积。事实上，直接按照定义去计算，复杂度就能做到 $\mathcal{O}(n^2)$ 了；只是没必要。把分子乘出来，枚举每一项，就知道它除以分母后对 $x^s$ 的贡献。

可惜的是，指数过于离散，需要最后再排序。所以逃不过 $\mathcal{O}(n^2\log n)$ 了。

而 $F(x{)}^4$ 呢？可以背包合并，即 $\text{Meet-In-Middle}$ 。求出 $F(x{)}^2$ 对应的分子 $H(x)$，接下来只需要计算
$$[x^s]H(x{)}^2(1-x{)}^{-4}=\sum _{{\kappa }_i+{\kappa }_j⩽s}{c}_{{\kappa }_i}{c}_{{\kappa }_j}\left(\genfrac{}{}{0ex}{}{s-{\kappa }_i-{\kappa }_j+3}{3}\right)$$

其中 $c_{{\kappa }_i}=[x^{{\kappa }_i}]H(x)$ 。由于指数并不形成区间，故用 ${\kappa }_i$ 表示。

${\kappa }_i$ 是可枚举的，怎样快速求出该式对 $j$ 的求和呢？说白了，就是要把那个组合数写成关于 ${\kappa }_i$ 的式子。可以直接按照定义展开成下降幂；也可以稍微高级一点，用二项式卷积，展开为
$$=\sum _{{\kappa }_i+{\kappa }_j⩽s}{c}_{{\kappa }_i}\sum _{p=0}^3\left(\genfrac{}{}{0ex}{}{-{\kappa }_i}{3-p}\right)\left(\genfrac{}{}{0ex}{}{s-{\kappa }_j+3}{p}\right)c_{{\kappa }_j}$$

于是我们只需要枚举 $i,p$，然后对 ${\kappa }_j$ 的前缀求和。时间复杂度 $\mathcal{O}(n^2)$ 。

类似地，研究 $F(x^2)F(x{)}^2$ 怎么求。还是求出两个分子 $H_1(x),H_2(x)$，则只需求出
$$[x^s]\frac{H_1(x)H_2(x)}{(1-x^2)(1-x{)}^2}$$

这个 $(1-x^2)$ 有点碍事。但是别怕，我们还是按照定义把 $(1-x^2)(1-x{)}^2$ 的倒数求出来：
$$\sum _{i=0}^{+\infty }{x}^i\sum _{j=0}^{\lfloor \frac{i}{2}\rfloor }(i-2j+1)$$

显然右侧求和结果是关于 $\lfloor \frac{i}{2}\rfloor$ 和 $i$ 的多项式。当 $2\mid i$ 时，其等于 $(\frac{i+2}{2}{)}^2$，当 $2\nmid i$ 时，其等于 $\frac{(i+1)(i+3)}{4}$ 。无论如何，是个关于 $i$ 的二次式，记为 $\omega (i)$ 。

按照 $i$ 的奇偶性分开计算，大抵是要计算这样的式子：
$$\sum _{{\kappa }_i+{\nu }_j⩽s}[{\kappa }_i+{\nu }_j\equiv \lambda ]c_{{\kappa }_i}{d}_{{\nu }_j}\cdot \omega (s-{\kappa }_i-{\nu }_j)$$

其中 λ ∈ { 0 , 1 } \lambda\in\{0,1\} λ∈{0,1} 用来固定 $\omega$ 自变量奇偶性，$c_{{\kappa }_i},d_{{\nu }_j}$ 仍表示离散的系数。同理 ${\kappa }_i$ 是可枚举的，同理需要将 $\omega (s-{\kappa }_i-{\nu }_j)$ 写成关于 ${\kappa }_i$ 的式子。这次最好还是直接代入，毕竟 $\omega$ 不过是个二次式。而且要按照 ${\nu }_j$ 的奇偶性分类，分别存储两种 $\omega (i)$ 。

总时间复杂度 $\mathcal{O}(n^2\log n)$ 。五个部分真的是独立的，非常讨厌。

代码

#include <cstdio> // JZM yydJUNK!!!
#include <iostream> // XJX yyds!!!
#include <algorithm> // XYX yydLONELY!!!
#include <cstring> // (the STRONG LONG for loneliness)
#include <cctype> // DDG yydDOG & ZXY yydBUS !!!
#include <vector>
# define rep(i,a,b) for(int i=(a); i<=(b); ++i)
# define rep0(i,a,b) for(int i=(a); i!=(b); ++i)
# define drep(i,a,b) for(int i=(a); i>=(b); --i)
typedef long long llong;
inline int readint(){
	int a = 0, c = getchar(), f = 1;
	for(; !isdigit(c); c=getchar())
		if(c == '-') f = -f;
	for(; isdigit(c); c=getchar())
		a = (a<<3)+(a<<1)+(c^48);
	return a*f;
}

const int MOD = 998244353;
inline void modAddUp(int &x, const int &y){
	if((x += y) >= MOD) x -= MOD;
}
typedef std::pair<int,int> PII;
typedef std::vector<PII> vecset;
void unique(vecset &v){
	if(v.empty()) return ;
	vecset::iterator i = v.begin();
	for(vecset::iterator j=i+1; j!=v.end(); ++j)
		if(j->first == i->first) // merge
			modAddUp(i->second,j->second);
		else ++ i, *i = *j; // new element
	++ i; v.resize(i-v.begin());
}
vecset convolution(const vecset &va, const vecset &vb){
	vecset res; res.clear();
	for(const PII &a : va) for(const PII &b : vb)
		res.push_back(PII{a.first+b.first,int(
			llong(a.second)*b.second%MOD)});
	std::sort(res.begin(),res.end());
	unique(res); return res;
}
/// @brief change indepent variable to @p k times x 
vecset mulx(const vecset &v, const int &k){
	vecset res = v; // copy
	for(PII &a : res) a.first *= k;
	return res;
}

struct function{
	int coe[3]; ///< quadratic
	function operator + (const function &b) const {
		function c; memcpy(c.coe,coe,sizeof(coe));
		rep0(i,0,3) modAddUp(c.coe[i],b.coe[i]);
		return c;
	}
	function operator * (const int &x) const {
		function c;
		rep0(i,0,3) c.coe[i] = int(llong(x)*coe[i]%MOD);
		return c;
	}
	void clear(){ memset(coe,0,sizeof(coe)); }
	int at(const int &x) const {
		return int(((llong(coe[2])*x
			+coe[1])%MOD*x+coe[0])%MOD);
	}
};
const int inv2 = (MOD+1)>>1, inv4 = 748683265;
/// @return the function when the argument is even number ( @p x0 - x )
function omega0(int x0){
	function c; c.coe[2] = inv4, x0 += 2;
	x0 = int(llong(x0)*inv2%MOD), c.coe[1] = MOD-x0;
	*c.coe = int(llong(x0)*x0%MOD); return c;
}
function omega1(int x0){
	function c; c.coe[2] = inv4, x0 += 2;
	*c.coe = int((llong(x0)*x0-1)%MOD*inv4%MOD);
	c.coe[1] = int(llong(MOD-x0)*inv2%MOD); return c;
}

vecset fx; ///< overall OGF
function func[2][2]; ///< parity of exponent, type of function
const int inv3 = (MOD+1)/3;
int main(){
	int n = readint(), s = readint();
	int ans = 0; // calc -6 * F(x^4) by definition
	for(int i=0,l,r,lst=-2; i!=n; ++i){
		l = readint(), r = readint();
		if(l == lst+1) // connect
			fx.back().first = r+1;
		else{ // new range
			fx.push_back(PII{l,1}), lst = r;
			fx.push_back(PII{r+1,MOD-1});
		}
		if(!(s&3) && l <= (s>>2) && (s>>2) <= r)
			ans = MOD-6; // easy calculate
	}
	{ // calc 8 * F(x^3) * F(x)
		vecset &&fx3 = mulx(fx,3);
		vecset &&got = convolution(fx3,fx);
		for(const PII &v : got){
			if(v.first > s) break;
			llong d = (s-v.first)/3+1;
			ans = int((ans+(d*v.second<<3))%MOD);
		}
	}
	vecset &&fx2 = mulx(fx,2); ///< F(x^2)
	if(!(s&1)){ // calc 3 * F(x^2) * F(x^2)
		vecset &&got = convolution(fx2,fx2);
		for(const PII &v : got){
			if(v.first > s) break;
			llong d = ((s-v.first)>>1)+1;
			ans = int((ans+3*d*v.second)%MOD);
		}
	}
	vecset &&fx_2 = convolution(fx,fx); ///< F(x)^2
	{ // calc -6 * F(x^2) * F(x)^2
		rep0(p,0,2) rep0(q,0,2) func[p][q].clear();
		vecset::iterator i = fx2.end(), j = fx_2.begin();
		while(i != fx2.begin()){
			if((--i)->first > s) continue;
			for(; j!=fx_2.end()&&j->first+i->first<=s; ++j)
				rep0(q,0,2){ // type of function
					function &fc = func[j->first&1][q];
					if(q) fc = fc+omega1(s-j->first)*j->second;
					else fc = fc+omega0(s-j->first)*j->second;
				}
			rep0(p,0,2) ans = int((ans+llong(MOD-6)*i->second
				%MOD*func[p][(s^p^i->first)&1].at(i->first))%MOD);
		}
	}
	{ // calc F(x)^2 * F(x)^2
		const int inv[] = {1,inv2,inv3,0};
		int xjx[4] = {0,0,0,0}; // coefficient
		vecset::iterator i = fx_2.end(), j = fx_2.begin();
		while(i != fx_2.begin()){
			if((--i)->first > s) continue;
			for(; j!=fx_2.end()&&j->first+i->first<=s; ++j)
				for(int p=0,v=j->second; p<=3; ++p){
					modAddUp(xjx[p],v), v = int(inv[p]*
						llong(s-j->first+3-p)%MOD*v%MOD);
				}
			for(int p=0,v=i->second; p<=3; ++p){
				ans = int((ans+llong(xjx[3-p])*v)%MOD);
				v = int(llong(MOD-i->first-p)*v%MOD*inv[p]%MOD);
			}
		}
	}
	ans = int(llong(ans)*inv4
		%MOD*inv3%MOD*inv2%MOD);
	printf("%d\n",ans); // div 24
	return 0;
}

后记

当我好不容易做出这道题时，我才发现 $\mathsf{O}\mathsf{U}\mathsf{Y}\mathsf{E}$ 已经把今天的三道题都做完了。

$\mathrm{Q}\mathrm{Q}$ 群里一片摆烂之声。不是别的原因，只是单纯地想知道：我为什么要存在于此。