【洛谷5294】[HNOI2019] 序列（主席树维护单调栈+二分）

点此看题面

大致题意： 给你一个长度为(n)的序列(A)，每次询问修改一个元素（只对当前询问有效），然后让你找到一个不下降序列(B)，使得这两个序列相应位置之差的平方和最小，并输出这个最小平方和。

如何预处理

首先，仔细观察样例解释，我们可以发现一个有趣的性质：对于(B)序列中相同的一段元素，它们在(A)序列中恰好是这一段区间中所有数的平均数。

因此，我们大胆猜测：我们可以把(A)序列划分成若干段，然后求出每段的平均值，就可以求出最后的(B)序列。

那么，现在我们要做的，就是如何划分的问题。

一个显然的结论，肯定是划分出越多的区间越优。

但我们如何判断是否可以划分一个区间呢？

这就可以通过单调栈维护，来预处理出答案了。

我们考虑枚举(i)，然后对于每个(i)，先将区间([i,i])扔入单调栈中。

然后，我们每次比较栈顶区间的平均值和栈顶下面一个区间的平均值，如果栈顶区间平均值较小，则不满足不下降性质，因此我们将这两个区间合并，然后继续与再下面一个区间比较。（注：其实等于也可以合并，不会影响平均值）

于是新问题来了，如何维护并合并区间信息？

设平均值(d=frac{a_1+a_2+...+a_x}x)，则一个区间内的答案是：

[(a_1-d)^2+(a_2-d)^2+...+(a_x-d)^2 ]

然后公式展开得到：

[a_1^2-2a_1d+d^2+a_2^2-2a_2d+d^2+...+a_x^2-2a_xd+d^2 ]

用加法交换律和加法结合律改变一下顺序得到：

[(a_1^2+a_2^2+...+a_x^2)-d(2a_1+2a_2+...+2a_x-xd) ]

将(d)的值代入得到：

[(a_1^2+a_2^2+...+a_x^2)-frac{a_1+a_2+...+a_x}xcdot(2a_1+2a_2+...+2a_x-a_1-a_2-...-a_x) ]

[(a_1^2+a_2^2+...+a_x^2)-frac{(a_1+a_2+...+a_x)^2}x ]

所以，我们只需要维护区间元素个数、区间元素和、区间元素平方和三个值就能完成求解答案和合并啦。

具体实现详见代码。

如何完成询问

设修改的位置为(x,y)，则可以发现，对于前(x-1)个位置，它的单调栈操作过程必然与初始化时一样。

而加入第(xsim n)个数之后，对于第(x-1)个单调栈，它大致可以表示成这个样子：

那么，我们只要找到(L,R)，并确定第(1sim L-1,R+1sim n)个元素的答案，就可以确定最终答案了。

先考虑我们要用到第(x-1)个单调栈，因此我们就可以用一个主席树来维护单调栈，并记录下每个位置此时的答案以及单调栈的大小，然后就能很方便的求出第(1sim L-1)的答案了。

再考虑第(R+1sim n)个元素的答案就很简单了，直接采取与上面预处理类似的方式，倒着做一遍单调栈，预处理出后缀的全部信息即可。

于是，现在仅剩的问题就是，如何找到(L,R)。

可以发现，这显然具有可二分性。

因此我们先二分一个右端点，然后对于这个右端点，我们再在主席树上二分求出其对应的左端点，并以此来验证这一段的平均值是否小于等于右端点的后一段的平均值即可。

注意对于栈中的元素我们肯定是整段整段取的，因为我们其实模拟了一个将元素加入栈的过程，而这一过程只会将某些段的元素合并，因此只能取整段。

具体实现详见代码。

代码

#include<bits/stdc++.h>
#define Tp template<typename Ty>
#define Ts template<typename Ty,typename... Ar>
#define Reg register
#define RI Reg int
#define RL Reg LL
#define Con const
#define CI Con int&
#define CL Con LL&
#define I inline
#define W while
#define N 100000
#define Log 20
#define LL long long
#define X 998244353
#define Inc(x,y) ((x+=(y))>=X&&(x-=X))
#define Dec(x,y) ((x-=(y))<0&&(x+=X))
using namespace std;
int n,a[N+5],Inv[N+5],rt1[N+5],rt2[N+5],top1[N+5],top2[N+5],ans1[N+5],ans2[N+5];
I int XSum(CI x,CI y) {return x+y>=X?x+y-X:x+y;}
I int XSub(CI x,CI y) {return x<y?x-y+X:x-y;}
struct data
{
	int t,ps;LL s;I data(CI x=0,CL y=0,CI z=0):t(x),s(y),ps(z){}
	I data operator + (Con data& o) {return data(t+o.t,s+o.s,XSum(ps,o.ps));}
	I data operator - (Con data& o) {return data(t-o.t,s-o.s,XSub(ps,o.ps));}
	I bool operator < (Con data& o) Con {return o.t?(t?1.0*s/t<1.0*o.s/o.t:1):0;}
	I bool operator <= (Con data& o) Con {return t?(o.t?1.0*s/t<=1.0*o.s/o.t:0):1;}
	I int GV() {return XSub(ps,1LL*(s%X)*(s%X)%X*Inv[t]%X);}
}s[N+5],S[N+5];
class FastIO
{
	private:
		#define FS 100000
		#define tc() (A==B&&(B=(A=FI)+fread(FI,1,FS,stdin),A==B)?EOF:*A++)
		#define pc(c) (C^FS?FO[C++]=c:(fwrite(FO,1,C,stdout),FO[(C=0)++]=c))
		#define tn (x<<3)+(x<<1)
		#define D isdigit(c=tc())
		int T,C;char c,*A,*B,FI[FS],FO[FS],S[FS];
	public:
		I FastIO() {A=B=FI;}
		Tp I void read(Ty& x) {x=0;W(!D);W(x=tn+(c&15),D);}
		Tp I void write(Ty x) {W(S[++T]=x%10+48,x/=10);W(T) pc(S[T--]);}
		Ts I void read(Ty& x,Ar&... y) {read(x),read(y...);}
		Tp I void writeln(Con Ty& x) {write(x),pc('
');}
		I void clear() {fwrite(FO,1,C,stdout),C=0;}
}F;
class ChairmanTree//主席树维护单调栈
{
	private:
		#define L l,mid,O[rt].S[0]
		#define R mid+1,r,O[rt].S[1]
		#define PushUp(x)
		(
			O[x].l=O[O[x].S[0]].l,O[x].r=O[O[x].S[O[x].S[1]>0]].r,
			O[x].k=O[O[x].S[0]].k
		)//上传信息
		int tot;struct node {int l,r,k,S[2];}O[N*Log<<1];
	public:
		I void Update(CI l,CI r,int& rt,CI p,CI vl,CI vr)//修改单调栈中的值
		{
			if(O[++tot]=O[rt],rt=tot,l==r) return (void)(O[rt].l=vl,O[rt].r=O[rt].k=vr);
			RI mid=l+r>>1;p<=mid?Update(L,p,vl,vr):Update(R,p,vl,vr),PushUp(rt);
		}
		I int Query1(CI l,CI r,CI rt,CI p,data& v)//主席树上二分左端点，用v存储下取到的段的信息
		{
			if(r<=p)
			{
				data dl=s[O[rt].k]-s[O[rt].l-1],dr=s[O[rt].r]-s[O[rt].k];
				if(dr+v<=dl) return v=v+s[O[rt].r]-s[O[rt].l-1],0;
				if(l==r) return O[rt].r;
			}
			RI mid=l+r>>1,res;return p>mid&&(res=Query1(R,p,v))?res:Query1(L,p,v);
		}
		I int Query2(CI l,CI r,CI rt,CI p)//求出第p个位置右端点的值
		{
			if(l==r) return O[rt].r;RI mid=l+r>>1;
			return p<=mid?Query2(L,p):Query2(R,p);
		}
}C;
int main()
{
	RI Qtot,i,x,y,t,T,l,r,mid,p1,p2,ans=0;data v;
	for(F.read(n,Qtot),i=1;i<=n;++i) F.read(a[i]),s[i]=s[i-1]+data(1,a[i],1LL*a[i]*a[i]%X);//预处理信息前缀和
	for(Inv[1]=1,i=2;i<=n;++i) Inv[i]=1LL*Inv[X%i]*(X-X/i)%X;//预处理逆元
	for(T=0,i=1;i<=n;++i)//正着做一遍单调栈
	{
		S[++T]=data(1,a[i],1LL*a[i]*a[i]%X),rt1[i]=rt1[i-1];//加入栈
		W(T^1&&S[T]<=S[T-1]) S[T-1]=S[T-1]+S[T],--T;//对于栈顶不合法元素进行合并
		C.Update(1,n,rt1[i],top1[i]=T,i-S[T].t+1,i),ans1[i]=XSum(ans1[i-S[T].t],S[T].GV());//主席树上修改，记录下每个位置此时的答案以及单调栈的大小
	} 
	for(T=0,i=n;i;--i)//倒着做一遍单调栈，与上面类似
	{
		S[++T]=data(1,a[i],1LL*a[i]*a[i]%X),rt2[i]=rt2[i+1];
		W(T^1&&S[T-1]<=S[T]) S[T-1]=S[T-1]+S[T],--T;//注意比较顺序相反
		C.Update(1,n,rt2[i],top2[i]=T,i,i+S[T].t-1),ans2[i]=XSum(ans2[i+S[T].t],S[T].GV());
	}
	F.writeln(ans1[n]);W(Qtot--)//处理询问
	{
		F.read(x,y),t=a[x],a[x]=y,l=0,r=top2[x+1]-1;//读入，确定二分上下界
		#define Work(pos)
			p2=(pos)?C.Query2(1,n,rt2[x+1],top2[x+1]-(pos)+1):x,
			v=data(1,a[x],1LL*a[x]*a[x]%X)+s[p2]-s[x],
			p1=x^1?C.Query1(1,n,rt1[x-1],top1[x-1],v):x
		W(l<=r) mid=l+r>>1,Work(mid),//二分
			v<s[C.Query2(1,n,rt2[x+1],top2[x+1]-mid)]-s[p2]?r=mid-1:l=mid+1;//验证
		Work(r+1),F.writeln(XSum(v.GV(),XSum(ans1[p1],ans2[p2+1]))),a[x]=t;//求解并输出答案
		//这里p1不用减1是因为它本来就是取左端点前一段的r值，而p2刚好是取这一段的r值因此要加1
	}return F.clear(),0;
}