LOJ3059. 「HNOI2019」序列

给出序列(a_i),你要决定一个实数不下降序列(b_i),最小化(sum(a_i-b_i)^2)

有若干个询问,表示修改一个(a_i)之后答案。

(nle 10^5)


参考题解:https://www.luogu.com.cn/blog/zhongyuwei/solution-p5294

考察对保序回归的进一步理解。

关于解的唯一性的证明:可以把({a_i},{b_i})分别看做高维向量,要最小化它们之间的距离。(b_ile b_{i+1})得到的空间是凸的,所以解唯一。

假设当前位置是(x)。考虑求出在搞完之后,单调栈中(x)所在区间的左右端点。求出区间([1,x-1])([x+1,n])(从后往前)的单调栈。现在要把这两个单调栈和(a_x)合并。

一个算法是:枚举(R_0)表示后面单调栈的第(R_0)个,把(x)([R_0+1,tp_r])先合并,合并之后尝试去和前面的栈合并。如果合并出来了之后小于(v(R_0)),则合法。找到最大的(R_0),并且二分出最大的(L_0)(表示([L_0+1,tp_l],x,[R_0+1,tp_r])合并起来大于(v(L_0))),如果小于(v(R_0)),则得到了我们想要的东西,即([L_0+1,tp_l],x,[R_0+1,tp_r])是最终(x)所在单调栈的区间。

接下来证明(R_0)的可二分性:也就是说如果(R_0)合法,(R_0-1)也合法。设两者求出的(L_0)分别为(L,L')

(v(l,r))表示([l+1,tp_l],x,[r+1,tp_r])合并之后的值。(L_0)是最大的(L_0)满足(v(L_0,R_0)>v(L_0)),等价于使(v(L_0,R_0))最大的(L_0)

因为(v(L',R_0)<v(L,R_0)<v(R_0)<v(R_0-1))(v(L',R_0-1))(v(L',R_0),v(R_0))的带权平均数,所以(v(L',R_0-1)<v(R_0-1))

于是可以二分套二分,先二分(R_0)再二分(L_0)。时间(O(nlg^2 n))


using namespace std;
#include <bits/stdc++.h>
#define N 100005
#define ll long long
#define ull unsigned long long
#define mo 998244353
ll sqr(ll x){return x*x%mo;}
ll qpow(ll x,ll y=mo-2){
	ll r=1;
	for (;y;y>>=1,x=x*x%mo)
		if (y&1)
			r=r*x%mo;
	return r;
}
ll inv[N];
int n,m;
int a[N];
struct Query{int x,c,id;} q[N];
bool cmpq(Query a,Query b){return a.x<b.x;}
bool cmpf(ll s0,ll w0,ll s1,ll w1){return (ull)s0*(ull)w1<=(ull)s1*(ull)w0;}//用<=而不是<(可能会出现分子分母为0,此时直接返回true) 
struct info{
	ll w,v,v2,s,sw,sv,sv2;
} sl[N],sr[N];
int tpl,tpr;
info merge(info a,info b){
	a.w+=b.w;
	a.v+=b.v;
	(a.v2+=b.v2)%=mo;
	return a;
}
struct oper{
	int t,ty;
	info c;
} o[N*2];
int cnt;
int calcl(int R0,int c){
	int l=0,r=tpl,_L0=0;
	while (l<=r){
		int L0=l+r>>1;
		ll tmpv=c+sl[tpl].sv-sl[L0].sv+sr[tpr].sv-sr[R0].sv,tmpw=1+sl[tpl].sw-sl[L0].sw+sr[tpr].sw-sr[R0].sw;
		if (cmpf(sl[L0].v,sl[L0].w,tmpv,tmpw)){
			_L0=L0;
			l=L0+1;
		}
		else
			r=L0-1;
	}
	return _L0;
}
ll calcr(int i,int c){
	int l=0,r=tpr,_R0=0,_L0=0;
	while (l<=r){
		int R0=l+r>>1,L0=calcl(R0,c);
		ll tmpv=c+sl[tpl].sv-sl[L0].sv+sr[tpr].sv-sr[R0].sv,tmpw=1+sl[tpl].sw-sl[L0].sw+sr[tpr].sw-sr[R0].sw;
		if (cmpf(tmpv,tmpw,sr[R0].v,sr[R0].w)){
			_R0=R0,_L0=L0;
			l=R0+1;
		}
		else
			r=R0-1;
	}
	ll tmpv=c+sl[tpl].sv-sl[_L0].sv+sr[tpr].sv-sr[_R0].sv,tmpw=1+sl[tpl].sw-sl[_L0].sw+sr[tpr].sw-sr[_R0].sw;
	ll tmpv2=sqr(c)+sl[tpl].sv2-sl[_L0].sv2+sr[tpr].sv2-sr[_R0].sv2;
	ll res=(sl[_L0].s+sr[_R0].s-sqr(tmpv%mo)*inv[tmpw]%mo+tmpv2)%mo;
	return res;
}
ll ans[N];
int main(){
	freopen("in.txt","r",stdin);
	scanf("%d%d",&n,&m);
	for (int i=1;i<=n;++i)
		scanf("%d",&a[i]);
	for (int i=1;i<=m;++i)
		scanf("%d%d",&q[i].x,&q[i].c),q[i].id=i;
	inv[1]=1;
	for (int i=2;i<=n;++i)
		inv[i]=(mo-mo/i)*inv[mo%i]%mo;
	for (int i=n;i>=1;--i){
		info tmp={1,a[i],sqr(a[i])};
		for (;tpr && cmpf(sr[tpr].v,sr[tpr].w,tmp.v,tmp.w);--tpr){
			tmp=merge(tmp,sr[tpr]);
			o[++cnt]={i,-1,sr[tpr]};
		}
		tmp.s=(sr[tpr].s-sqr(tmp.v%mo)*inv[tmp.w]+tmp.v2)%mo;
		tmp.sw=sr[tpr].sw+tmp.w;
		tmp.sv=sr[tpr].sv+tmp.v;
		tmp.sv2=(sr[tpr].sv2+tmp.v2)%mo;
		sr[++tpr]=tmp;
		o[++cnt]={i,1};
	}
	ans[0]=sr[tpr].s;
	sort(q+1,q+m+1,cmpq);
	for (int i=1,j=1;i<=n;++i){
		for (;cnt && o[cnt].t<=i;--cnt)
			if (o[cnt].ty==-1)
				sr[++tpr]=o[cnt].c;
			else
				--tpr;
		for (;j<=m && q[j].x==i;++j)
			ans[q[j].id]=calcr(i,q[j].c);
		info tmp={1,a[i],sqr(a[i])};
		for (;tpl && cmpf(tmp.v,tmp.w,sl[tpl].v,sl[tpl].w);--tpl)
			tmp=merge(tmp,sl[tpl]);
		tmp.s=(sl[tpl].s-sqr(tmp.v%mo)*inv[tmp.w]+tmp.v2)%mo;
		tmp.sw=sl[tpl].sw+tmp.w;
		tmp.sv=sl[tpl].sv+tmp.v;
		tmp.sv2=(sl[tpl].sv2+tmp.v2)%mo;
		sl[++tpl]=tmp;
	}
	for (int i=0;i<=m;++i)
		printf("%lld
",(ans[i]+mo)%mo);
	return 0;	
}
原文地址:https://www.cnblogs.com/jz-597/p/14622665.html