[ARC124E]Pass to Next

最新推荐文章于 2024-09-13 16:27:15 发布

OneInDark

最新推荐文章于 2024-09-13 16:27:15 发布

阅读量303

点赞数

分类专栏： C++ 动态规划文章标签：算法动态规划线性代数

本文链接：https://blog.csdn.net/qq_42101694/article/details/120817682

版权

动态规划线性变换矩阵乘法编程竞赛复杂度优化

关键词由CSDN通过智能技术生成

C++ 同时被 2 个专栏收录

547 篇文章 3 订阅

订阅专栏

动态规划

220 篇文章 2 订阅

订阅专栏

题目

传送门 to AtCoder

思路

说个笑话：我一上来就是 传递糖果，于是一个人可以传递 $\sum a$ 个糖果……

事实上普通 $\tt dp$ 应该很简单。设 $dp_{i,j}$ 为第 $i$ 个人给出 $j$ 个糖果，有转移
$dp_{i,j}=\sum_{k\leqslant a_{i-1}}(a_i+k-j)\cdot dp_{i-1,k}$

由于 $j, k$ 实际上是独立的，我们很有可能直接从 $\sum_{k\leqslant a_{i-1}}dp_{i-1,k}$ 和 $\sum_{k\leqslant a_{i-1}}dp_{i-1,k}\cdot k$ 得到 $\sum_{j\leqslant a_i}dp_{i,j}$ 。
$\sum_{j\leqslant a_i}dp_{i,j}=s^{[0]}(a_i)\sum_{k\leqslant a_{i-1}}(a_i+k)\cdot dp_{i-1,k}-s^{[1]}(a_i)\sum_{k\leqslant a_{i-1}}dp_{i-1,k}\\ \sum_{j\leqslant a_i}dp_{i,j}\cdot j=s^{[1]}(a_i)\sum_{k\leqslant a_{i-1}}(a_i+k)\cdot dp_{i-1,k}-s^{[2]}(a_i)\sum_{k\leqslant a_{i-1}}dp_{i-1,k}$

其中 $s^{[p]}(v)=\sum_{i=0}^{v}i^p$ 。这很容易理解。

但是这个会算重，因为不同的给出量不代表不同的剩余量。事实上，决定剩余量的是差分值，所以我们强迫一个值为 $0$ 即可。枚举它，时间复杂度 $\mathcal O(n^2)$ 。

然而这是线性变化，可以写出转移矩阵
$\begin{bmatrix} s^{[0]}\cdot a_i-s^{[1]} & s^{[1]}\cdot a_i-s^{[2]} \\ s^{[0]} & s^{[1]} \end{bmatrix}$

注：状态矩阵是行向量。那么维护矩阵的前后缀乘积即可 $\mathcal O(1)$ 计算。总复杂度 $\mathcal O(n)$ 。

代码

#include <cstdio>
#include <iostream>
#include <vector>
#include <cstring>
#include <algorithm>
using namespace std;
# define rep(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 int_;
inline int readint(){
	int a = 0, c = getchar(), f = 1;
	for(; '0'>c||c>'9'; c=getchar())
		if(c == '-') f = -f;
	for(; '0'<=c&&c<='9'; c=getchar())
		a = (a<<3)+(a<<1)+(c^48);
	return a*f;
}
void writeUnsigned(const unsigned &x){
	if(x > 9) writeUnsigned(x/10);
	putchar((x-x/10*10)^48);
}
inline void writeint(const int &x){
	if(x < 0){
		putchar('-');
		writeUnsigned(-x);
	}
	else writeUnsigned(x);
}

const int Mod = 998244353;
const int inv3 = (Mod+1)/3;
struct Matrix{
	int a[2][2];
	void clear(){ memset(a,0,2*2<<2); }
	Matrix operator * (const Matrix &b) const {
		Matrix c; c.clear();
		rep(i,0,1) rep(j,0,1) rep(k,0,1)
			c.a[i][k] = static_cast<int>(
				(c.a[i][k]+int_(a[i][j])*b.a[j][k])%Mod);
		return c;
	}
	static Matrix I(){
		Matrix c;
		c.a[0][1] = c.a[1][0] = 0;
		c.a[0][0] = c.a[1][1] = 1;
		return c;
	}
};

const int MAXN = 1000005;
/// @param f if x can be 0
Matrix fuck(int x,bool f=true){
	Matrix res;
	int &s0 = res.a[1][0] = x+f;
	int &s1 = res.a[1][1] = int((int_(x+1)*x>>1)%Mod);
	int s2 = int(s1*int_(x<<1|1)%Mod*inv3%Mod);
	res.a[0][0] = int((int_(s0)*x+Mod-s1)%Mod);
	res.a[0][1] = int((int_(s1)*x+Mod-s2)%Mod);
	return res;
}
Matrix pre[MAXN], suf[MAXN];
int a[MAXN];
int main(){
	int n = readint(), ans = 0;
	pre[0] = suf[n+1] = Matrix::I();
	rep(i,1,n){
		a[i] = readint();
		pre[i] = pre[i-1]*fuck(a[i],false);
	}
	drep(i,n,1){
		suf[i] = fuck(a[i],true)*suf[i+1];
		Matrix mat; mat.clear();
		mat.a[0][0] = 1; // f(i,0) = 1
		mat = mat*suf[i+1]*pre[i-1];
		ans = int((ans+int_(a[i])*mat.a[0][0])%Mod);
		ans = (ans+mat.a[0][1])%Mod;
	}
	printf("%d\n",ans);
	return 0;
}