相关分析
Time Limit: 10 Sec Memory Limit: 128 MB[Submit][Status][Discuss]
Description
Frank对天文学非常感兴趣，他经常用望远镜看星星，同时记录下它们的信息，比如亮度、颜色等等，进而估算出星星的距离，半径等等。
Frank不仅喜欢观测，还喜欢分析观测到的数据。他经常分析两个参数之间（比如亮度和半径）是否存在某种关系。
现在Frank要分析参数X与Y之间的关系。
他有n组观测数据，第i组观测数据记录了x_i和y_i。
他需要一下几种操作
1 L,R：用直线拟合第L组到底R组观测数据。
用xx表示这些观测数据中x的平均数，用yy表示这些观测数据中y的平均数，即
xx=Σx_i/(R-L+1)(L<=i<=R)
yy=Σy_i/(R-L+1)(L<=i<=R)
如果直线方程是y=ax+b，那么a应当这样计算：
a=(Σ(x_i-xx)(y_i-yy))/(Σ(x_i-xx)(x_i-xx)) (L<=i<=R)
你需要帮助Frank计算a。
2 L,R,S,T：
Frank发现测量数据第L组到底R组数据有误差，对每个i满足L <= i <= R，x_i需要加上S，y_i需要加上T。
3 L,R,S,T：
Frank发现第L组到第R组数据需要修改，对于每个i满足L <= i <= R，x_i需要修改为(S+i)，y_i需要修改为(T+i)。
Input
第一行两个数n,m，表示观测数据组数和操作次数。
接下来一行n个数，第i个数是x_i。
接下来一行n个数，第i个数是y_i。
接下来m行，表示操作，格式见题目描述。
保证1操作不会出现分母为0的情况。
Output
对于每个1操作，输出一行，表示直线斜率a。
选手输出与标准输出的绝对误差不超过10^-5即为正确。
Sample Input
3 5
1 2 3
1 2 3
1 1 3
2 2 3 -3 2
1 1 2
3 1 2 2 1
1 1 3
1 2 3
1 2 3
1 1 3
2 2 3 -3 2
1 1 2
3 1 2 2 1
1 1 3
Sample Output
1.0000000000
-1.5000000000
-0.6153846154
-1.5000000000
-0.6153846154
HINT
1<=n,m<=10^5,0<=|S|,|T|,|x_i|,|y_i|<=10^5
Main idea
维护一个线性回归方程，需要支持区间加，区间覆盖等差数列。
Solution
我们先化一个式子：
然后就只要运用线段树维护 Σx Σy Σxx Σxy 就可以了。
每一个具体怎么维护的话，就是把式子列出来，暴力展开一下看一下其中的关联即可，并不难（BearChild懒得写啦！）。
Code
1 #include<iostream>
2 #include<string>
3 #include<algorithm>
4 #include<cstdio>
5 #include<cstring>
6 #include<cstdlib>
7 #include<cmath>
8 using namespace std;
9 typedef long long s64;
10
11 const int ONE = 200005;
12
13 int n,m,P;
14 int L,R,S,T;
15 double Sumsq[ONE];
16 double Vx[ONE],Vy[ONE];
17
18 struct power
19 {
20 double sumx,sumy,sumxx,sumxy;
21 double addx,addy;
22 double covx,covy;
23 bool cov;
24 }Node[ONE*4];
25
26 struct ans
27 {
28 double x,y,xx,xy;
29 }res;
30
31 inline int get()
32 {
33 int res=1,Q=1; char c;
34 while( (c=getchar())<48 || c>57)
35 if(c=='-')Q=-1;
36 if(Q) res=c-48;
37 while((c=getchar())>=48 && c<=57)
38 res=res*10+c-48;
39 return res*Q;
40 }
41
42 void Covers(int i,int l,int r,double S,double T)
43 {
44 if(l > r) return;
45 double len = r-l+1; double sum = (l+r)*len/2;
46 Node[i].addx = Node[i].addy = 0;
47 Node[i].covx = S; Node[i].covy = T;
48 Node[i].cov = 1;
49 Node[i].sumxx = S*S*len + sum*S + sum*S + Sumsq[r] - Sumsq[l-1];
50 Node[i].sumxy = S*T*len + sum*S + sum*T + Sumsq[r] - Sumsq[l-1];
51 Node[i].sumx = (S+l + S+r)*len / 2;
52 Node[i].sumy = (T+l + T+r)*len / 2;
53 }
54
55 void PC(int i,int l,int r)
56 {
57 if(Node[i].cov)
58 {
59 int mid = l+r>>1;
60 Covers(i<<1,l,mid, Node[i].covx,Node[i].covy);
61 Covers(i<<1|1,mid+1,r, Node[i].covx,Node[i].covy);
62 Node[i].cov = 0;
63 }
64 }
65
66 void Update(int i,int l,int r,double S,double T)
67 {
68 if(l > r) return;
69 PC(i,l,r);
70 double len = r-l+1;
71 Node[i].addx += S; Node[i].addy += T;
72 Node[i].sumxx += 2*S*Node[i].sumx + S*S*len;
73 Node[i].sumxy += S*Node[i].sumy + T*Node[i].sumx + S*T*len;
74 Node[i].sumx += S*len; Node[i].sumy += T*len;
75 }
76
77 void PU(int i,int l,int r)
78 {
79 if(Node[i].addx || Node[i].addy)
80 {
81 int mid = l+r>>1;
82 Update(i<<1,l,mid, Node[i].addx,Node[i].addy);
83 Update(i<<1|1,mid+1,r, Node[i].addx,Node[i].addy);
84 Node[i].addx = Node[i].addy = 0;
85 }
86 }
87
88 void pushdown(int i,int l,int r)
89 {
90 PU(i,l,r); PC(i,l,r);
91 }
92
93 void Renew(int i)
94 {
95 int a = i<<1, b = i<<1|1;
96 Node[i].sumx = Node[a].sumx + Node[b].sumx;
97 Node[i].sumy = Node[a].sumy + Node[b].sumy;
98 Node[i].sumxx = Node[a].sumxx + Node[b].sumxx;
99 Node[i].sumxy = Node[a].sumxy + Node[b].sumxy;
100 }
101
102 void Build(int i,int l,int r)
103 {
104 if(l==r)
105 {
106 Node[i].sumx = Vx[l];
107 Node[i].sumy = Vy[l];
108 Node[i].sumxx = (double)Vx[l] * Vx[l];
109 Node[i].sumxy = (double)Vx[l] * Vy[l];
110 return;
111 }
112 int mid = l+r>>1;
113 Build(i<<1,l,mid); Build(i<<1|1,mid+1,r);
114 Renew(i);
115 }
116
117 void Cov(int i,int l,int r,int L,int R,double S,double T)
118 {
119 if(L<=l && r<=R)
120 {
121 Covers(i,l,r,S,T);
122 return;
123 }
124
125 pushdown(i,l,r);
126 int mid = l+r>>1;
127 if(L<=mid) Cov(i<<1,l,mid,L,R,S,T);
128 if(mid+1<=R) Cov(i<<1|1,mid+1,r,L,R,S,T);
129 Renew(i);
130 }
131
132 void Add(int i,int l,int r,int L,int R,double S,double T)
133 {
134 if(L<=l && r<=R)
135 {
136 Update(i,l,r,S,T);
137 return;
138 }
139
140 pushdown(i,l,r);
141 int mid = l+r>>1;
142 if(L<=mid) Add(i<<1,l,mid,L,R,S,T);
143 if(mid+1<=R) Add(i<<1|1,mid+1,r,L,R,S,T);
144 Renew(i);
145 }
146
147 void Query(int i,int l,int r,int L,int R)
148 {
149 if(L<=l && r<=R)
150 {
151 res.x += Node[i].sumx; res.y += Node[i].sumy;
152 res.xx += Node[i].sumxx; res.xy += Node[i].sumxy;
153 return;
154 }
155
156 pushdown(i,l,r);
157 int mid = l+r>>1;
158 if(L<=mid) Query(i<<1,l,mid,L,R);
159 if(mid+1<=R) Query(i<<1|1,mid+1,r,L,R);
160 }
161
162 int main()
163 {
164 for(int i=1;i<=ONE-1;i++) Sumsq[i] = Sumsq[i-1] + (double)i*i;
165
166 n=get(); m=get();
167 for(int i=1;i<=n;i++) Vx[i]=get();
168 for(int i=1;i<=n;i++) Vy[i]=get();
169 Build(1,1,n);
170
171 while(m--)
172 {
173 P = get(); L = get(); R = get();
174 if(P == 1)
175 {
176 res.x = res.y = res.xx = res.xy = 0;
177 Query(1,1,n,L,R);
178 double len = R-L+1;
179 double Avex = res.x / len;
180 double Avey = res.y / len;
181 printf("%.6lf
", (res.xy - len * Avex * Avey) / (res.xx - len*Avex*Avex));
182 }
183 else
184 {
185 S = get(); T = get();
186 if(P == 2) Add(1,1,n, L,R,S,T);
187 else Cov(1,1,n, L,R,S,T);
188 }
189 }
190 }