P3707 [SDOI2017] 相关分析

题目大意

详细题目传送门

\(\verb!1 L R!\)

用直线拟合第 \(L\) 组到第 \(R\) 组观测数据。用 \(\overline{x}\) 表示这些观测数据中 \(x\) 的平均数，用 \(\overline{y}\) 表示这些观测数据中 \(y\) 的平均数，

\[a={\displaystyle\sum_{i=L} ^R (x_i-\overline{x})(y_i-\overline{y}) \over \displaystyle\sum _{i=L} ^R (x_i -\overline{x})^2} \]

你需要帮助 Frank 计算 \(a\)。

\(\verb!2 L R S T!\)

Frank 发现测量数据第 \(L\) 组到第 \(R\) 组数据有误差，对每个 \(i\) 满足 \(L \leq i \leq R\)，\(x_i\) 需要加上 \(S\)，\(y_i\) 需要加上 \(T\)。

\(\verb!3 L R S T!\)

Frank 发现第 \(L\) 组到第 \(R\) 组数据需要修改，对于每个 \(i\) 满足 \(L \leq i \leq R\)，\(x_i\) 需要修改为 \((S+i)\)，\(y_i\) 需要修改为 \((T+i)\)。
\(n,m\leq10^5\)

思路

先对分子分母进行分开计算。

\[\sum_{i=L} ^R (x_i-\overline{x})(y_i-\overline{y}) \]

\[=\sum_{i=L}^R(x_iy_i-x_i\overline{y}-y_i\overline{x}+\overline{x}\overline{y}) \]

\[=\sum_{i=L}^Rx_iy_i-\overline{y}\sum_{i=L}^Rx_i-\overline{x}\sum_{i=L}^Ry_i+\sum_{i=L}^R\overline{x}\overline{y} \]

\[=\sum_{i=L}^Rx_iy_i-\frac{3\sum_{i=L}^Ry_i\sum_{i=L}^Rx_i}{R-L+1} \]

那么分子同理

\[=\sum_{i=L}^Rx_ix_i-\frac{3\sum_{i=L}^Rx_i\sum_{i=L}^Rx_i}{R-L+1} \]

只需要线段树维护即可。

代码

#include<bits/stdc++.h>
#define rep(i,a,b) for(int i=(a);i<=(b);++i)
using namespace std;
typedef long long ll;
const ll MAXN=1e5+5;
ll n,m;
long double x[MAXN],y[MAXN];
struct Query{
	ll op,l,r,s,t;
}a[MAXN];
struct node{
	long double sx,sy,xy,xx,tx,ty,fx,fy;
	#define lc(u) (u<<1)
	#define rc(u) (u<<1|1)
}t[MAXN*4];
void push_up(ll u){
	t[u].sx=t[lc(u)].sx+t[rc(u)].sx;
	t[u].sy=t[lc(u)].sy+t[rc(u)].sy;
	t[u].xx=t[lc(u)].xx+t[rc(u)].xx;
	t[u].xy=t[lc(u)].xy+t[rc(u)].xy;
}
long double pfh(ll x){
	return x*(x+1)*(2*x+1)/6;
}
long double zh(ll l,ll r){
	return (l+r)*(r-l+1)/2;
}
void push_down(ll u,ll l,ll r){
	ll mid=(l+r)>>1;
	if(t[u].fx!=-1e18){
		t[lc(u)].fx=t[u].fx;
		t[lc(u)].fy=t[u].fy;
		t[lc(u)].xy=pfh(mid)-pfh(l-1)+t[u].fx*zh(l,mid)+t[u].fy*zh(l,mid)+(mid-l+1)*t[u].fx*t[u].fy;
		t[lc(u)].xx=pfh(mid)-pfh(l-1)+t[u].fx*zh(l,mid)+t[u].fx*zh(l,mid)+(mid-l+1)*t[u].fx*t[u].fx;
		t[lc(u)].sx=(mid-l+1)*t[u].fx+zh(l,mid);
		t[lc(u)].sy=(mid-l+1)*t[u].fy+zh(l,mid);
		t[rc(u)].fx=t[u].fx;
		t[rc(u)].fy=t[u].fy;
		t[rc(u)].xy=pfh(r)-pfh(mid)+t[u].fx*zh(mid+1,r)+t[u].fy*zh(mid+1,r)+(r-mid)*t[u].fx*t[u].fy;
		t[rc(u)].xx=pfh(r)-pfh(mid)+t[u].fx*zh(mid+1,r)+t[u].fx*zh(mid+1,r)+(r-mid)*t[u].fx*t[u].fx;
		t[rc(u)].sx=(r-mid)*t[u].fx+zh(mid+1,r);
		t[rc(u)].sy=(r-mid)*t[u].fy+zh(mid+1,r);
		t[u].fx=t[u].fy=-1e18;
		t[lc(u)].tx=t[lc(u)].ty=t[rc(u)].tx=t[rc(u)].ty=0;
	}
	t[lc(u)].xx+=2*t[u].tx*t[lc(u)].sx+(mid-l+1)*t[u].tx*t[u].tx;
	t[lc(u)].xy+=t[u].ty*t[lc(u)].sx+t[u].tx*t[lc(u)].sy+(mid-l+1)*t[u].tx*t[u].ty;
	t[lc(u)].sx+=t[u].tx*(mid-l+1);
	t[lc(u)].sy+=t[u].ty*(mid-l+1);
	t[lc(u)].tx+=t[u].tx;
	t[lc(u)].ty+=t[u].ty;
	t[rc(u)].xx+=2*t[u].tx*t[rc(u)].sx+(r-mid)*t[u].tx*t[u].tx;
	t[rc(u)].xy+=t[u].ty*t[rc(u)].sx+t[u].tx*t[rc(u)].sy+(r-mid)*t[u].tx*t[u].ty;
	t[rc(u)].sx+=t[u].tx*(r-mid);
	t[rc(u)].sy+=t[u].ty*(r-mid);
	t[rc(u)].tx+=t[u].tx;
	t[rc(u)].ty+=t[u].ty;
	t[u].tx=t[u].ty=0;
}
void modify(ll u,ll l,ll r,ll ql,ll qr,long double tx,long double ty){
	if(ql<=l&&r<=qr){
		t[u].xx+=2*tx*t[u].sx+(r-l+1)*tx*tx;
		t[u].xy+=ty*t[u].sx+tx*t[u].sy+(r-l+1)*tx*ty;
		t[u].sx+=tx*(r-l+1);
		t[u].sy+=ty*(r-l+1);
		t[u].tx+=tx;
		t[u].ty+=ty;
		return;
	}
	push_down(u,l,r);
	ll mid=(l+r)>>1;
	if(ql<=mid){
		modify(lc(u),l,mid,ql,qr,tx,ty);
	}
	if(mid+1<=qr){
		modify(rc(u),mid+1,r,ql,qr,tx,ty);
	}
	push_up(u);
}
void modify2(ll u,ll l,ll r,ll ql,ll qr,long double fx,long double fy){
	if(ql<=l&&r<=qr){
		t[u].fx=fx;
		t[u].fy=fy;
		t[u].tx=t[u].ty=0;
		t[u].xy=pfh(r)-pfh(l-1)+fx*zh(l,r)+fy*zh(l,r)+(r-l+1)*fx*fy;
		t[u].xx=pfh(r)-pfh(l-1)+fx*zh(l,r)+fx*zh(l,r)+(r-l+1)*fx*fx;
		t[u].sx=(r-l+1)*fx+zh(l,r);
		t[u].sy=(r-l+1)*fy+zh(l,r);
		return;
	}
	push_down(u,l,r);
	ll mid=(l+r)>>1;
	if(ql<=mid){
		modify2(lc(u),l,mid,ql,qr,fx,fy);
	}
	if(mid+1<=qr){
		modify2(rc(u),mid+1,r,ql,qr,fx,fy);
	}
	push_up(u);
}
node query(ll u,ll l,ll r,ll ql,ll qr){
	if(ql<=l&&r<=qr){
		return t[u];
	}
	push_down(u,l,r);
	ll mid=(l+r)>>1;
	node val={0,0,0,0,0,0};
	if(ql<=mid){
		val=query(lc(u),l,mid,ql,qr);
	}
	if(mid+1<=qr){
		node rn=query(rc(u),mid+1,r,ql,qr);
		val.sx+=rn.sx;
		val.sy+=rn.sy;
		val.xx+=rn.xx;
		val.xy+=rn.xy;
	}
	return val;
}
void build(ll u,ll l,ll r){
	t[u].tx=t[u].ty=0;
	t[u].fx=t[u].fy=-1e18;
	if(l==r){
		t[u].sx=x[l];
		t[u].sy=y[l];
		t[u].xx=x[l]*x[l];
		t[u].xy=x[l]*y[l];
		return;
	}
	ll mid=(l+r)>>1;
	build(lc(u),l,mid);
	build(rc(u),mid+1,r);
	push_up(u);
}
namespace Baoli{
	double val(ll l,ll r){
		node vv=query(1,1,n,l,r);
		long double xs=0,ys=0,xy=0,xx=0;
		xs=vv.sx,ys=vv.sy,xy=vv.xy,xx=vv.xx;
		long double fz=xy
				-xs*1.0*ys*1.0/(r-l+1)
				-xs*1.0/(r-l+1)*1.0*ys
				+xs*1.0/(r-l+1)*ys*1.0/(r-l+1)*(r-l+1);
		long double fm=xx
				-xs*2.0*xs*1.0/(r-l+1)
				+xs*1.0/(r-l+1)*xs*1.0/(r-l+1)*(r-l+1);
		return fz*1.0/fm;
	};
	void Do(){
		build(1,1,n);
		rep(i,1,m){
			ll op=a[i].op;
			if(op==1){
				printf("%0.10lf\n",val(a[i].l,a[i].r));
				continue;
			}
			if(op==2){
				modify(1,1,n,a[i].l,a[i].r,a[i].s,a[i].t);
			}
			if(op==3){
				modify2(1,1,n,a[i].l,a[i].r,a[i].s,a[i].t);
			}
		};
	};
}
int main(){
	ios::sync_with_stdio(false); 
	cin.tie(0);cout.tie(0);
	cin>>n>>m;
	rep(i,1,n){
		cin>>x[i];
	};
	rep(i,1,n){
		cin>>y[i];
	};
	rep(i,1,m){
		cin>>a[i].op>>a[i].l>>a[i].r;
		if(a[i].op!=1){
			cin>>a[i].s>>a[i].t;
		}
	};
	Baoli::Do();
	return 0;
}