51nod1687 方差求和

题目

题目链接

解题思路

简化公式
∑ i = 1 n ( a i − a ˉ ) 2 = ∑ i = 1 n ( a i 2 − 2 a i a ˉ + a ˉ 2 ) = ∑ i = 1 n a i 2 − 2 a ˉ ∑ i = 1 n a i + n a ˉ 2 = ∑ i = 1 n a i 2 − 2 n a ˉ 2 + n a ˉ 2 = ∑ i = 1 n a i 2 − n a ˉ 2 \sum_{i=1}^n(a_i-\bar a)^2\\=\sum_{i=1}^n(a_i^2-2a_i\bar a+\bar a^2)\\=\sum_{i=1}^na_i^2-2\bar a\sum_{i=1}^na_i+n\bar a^2\\=\sum_{i=1}^na_i^2-2n\bar a^2+n\bar a^2\\=\sum_{i=1}^na_i^2-n\bar a^2 ∑i=1n​(ai​−aˉ)2=∑i=1n​(ai2​−2ai​aˉ+aˉ2)=∑i=1n​ai2​−2aˉ∑i=1n​ai​+naˉ2=∑i=1n​ai2​−2naˉ2+naˉ2=∑i=1n​ai2​−naˉ2
考虑增量法
n a ˉ 2 = ( ∑ i = 1 n a i ) 2 n n\bar a^2=\frac{(\sum_{i=1}^na_i)^2}{n} naˉ2=n(∑i=1n​ai​)2​
( x + a ) 2 n + 1 − x 2 n = n x 2 + 2 n x a + n a 2 − ( n + 1 ) x 2 n ( n + 1 ) = 2 n a x + n a 2 − x 2 n ( n + 1 ) = 2 a x + a 2 n + 1 − x 2 n ( n + 1 ) \frac{(x + a)^2}{n+1}-\frac{x^2}{n}\\=\frac{nx^2+2nxa+na^2-(n+1)x^2}{n(n+1)}\\=\frac{2nax+na^2-x^2}{n(n+1)}\\=\frac{2ax+a^2}{n+1}-\frac{x^2}{n(n+1)} n+1(x+a)2​−nx2​=n(n+1)nx2+2nxa+na2−(n+1)x2​=n(n+1)2nax+na2−x2​=n+12ax+a2​−n(n+1)x2​
无法摆脱与长度的关系，不可行。
考虑枚举长度，同时计算长度相同的区间
我们知道区间长度和左端点L就可以计算右端点R了。

a n s = ∑ l e n = 1 n ∑ L = 1 n − l e n + 1 ( ∑ i = L R a i 2 − l e n a ˉ 2 ) ans=\sum_{len=1}^n\sum_{L=1}^{n-len+1}(\sum_{i=L}^Ra_i^2-len\bar a^2) ans=∑len=1n​∑L=1n−len+1​(∑i=LR​ai2​−lenaˉ2)

其中 ∑ l e n = 1 n ∑ L = 1 n − l e n + 1 ∑ i = L R a i 2 = ∑ i = 1 n a i 2 ∗ i ∗ ( n − i + 1 ) \sum_{len=1}^n\sum_{L=1}^{n-len+1}\sum_{i=L}^Ra_i^2=\sum_{i=1}^na_i^2*i*(n-i+1) ∑len=1n​∑L=1n−len+1​∑i=LR​ai2​=∑i=1n​ai2​∗i∗(n−i+1)可以简单求出。

∑ l e n = 1 n ∑ L = 1 n − l e n + 1 l e n a ˉ 2 \sum_{len=1}^n\sum_{L=1}^{n-len+1}len\bar a^2 ∑len=1n​∑L=1n−len+1​lenaˉ2

= ∑ l e n = 1 n 1 l e n ∑ L = 1 n − l e n + 1 ( ∑ i = L R a i ) 2 =\sum_{len=1}^n\frac{1}{len}\sum_{L=1}^{n-len+1}(\sum_{i=L}^Ra_i)^2 =∑len=1n​len1​∑L=1n−len+1​(∑i=LR​ai​)2

= ∑ l e n = 1 n 1 l e n ∑ L = 1 n − l e n + 1 ( ∑ i = L R a i 2 + 2 ∑ i = L R ∑ j = i + 1 R a i a j ) =\sum_{len=1}^n\frac{1}{len}\sum_{L=1}^{n-len+1}(\sum_{i=L}^Ra_i^2+2\sum_{i=L}^R\sum_{j=i+1}^Ra_ia_j) =∑len=1n​len1​∑L=1n−len+1​(∑i=LR​ai2​+2∑i=LR​∑j=i+1R​ai​aj​)

其中 ∑ l e n = 1 n 1 l e n ∑ L = 1 n − l e n + 1 ∑ i = L R a i 2 \sum_{len=1}^n\frac{1}{len}\sum_{L=1}^{n-len+1}\sum_{i=L}^Ra_i^2 ∑len=1n​len1​∑L=1n−len+1​∑i=LR​ai2​也能够简单维护。

令 f ( l e n ) = ∑ L = 1 n − l e n + 1 ∑ i = L R ∑ j = i + 1 R a i a j f(len)=\sum_{L=1}^{n-len+1}\sum_{i=L}^R\sum_{j=i+1}^Ra_ia_j f(len)=∑L=1n−len+1​∑i=LR​∑j=i+1R​ai​aj​

我们发现$f(len)$的二阶差分可以分为两部分，一部分可以用NTT维护，一部分可以简单维护。
整体复杂度 O ( T n l o g 2 n ) O(Tnlog_2n) O(Tnlog2​n)。

代码

#include 
#include 
#include 
using namespace std;typedef long long ll;void read(int &x) {x = 0; char c = getchar();while (c < '0' || c > '9') c = getchar();while (c >= '0' && c <= '9') x = x * 10 + c - '0', c = getchar();
}const int N = 3e5 + 100;namespace NTT {
#define pw(n) (1<const int N = 262144, P = 998244353, g = 3;//或P=1004535809int n, m, bit, bitnum = 0, a[N + 5], b[N + 5], rev[N + 5];void getrev(int l) {for (int i = 0; i < pw(l); i++) {rev[i] = (rev[i >> 1] >> 1) | ((i & 1) << (l - 1));}}int fastpow(int a, int b) {int ans = 1;for (; b; b >>= 1, a = 1LL * a*a%P) {if (b & 1)ans = 1LL * ans*a%P;}return ans;}void NTT(int *s, int op) {for (int i = 0; i < bit; i++)if (i < rev[i])swap(s[i], s[rev[i]]);for (int i = 1; i < bit; i <<= 1) {int w = fastpow(g, (P - 1) / (i << 1));for (int p = i << 1, j = 0; j < bit; j += p) {int wk = 1;for (int k = j; k < i + j; k++, wk = 1LL * wk*w%P) {int x = s[k], y = 1LL * s[k + i] * wk%P;s[k] = (x + y) % P;s[k + i] = (x - y + P) % P;}}}if (op == -1) {reverse(s + 1, s + bit);int inv = fastpow(bit, P - 2);for (int i = 0; i < bit; i++)s[i] = 1LL * s[i] * inv%P;}}int solve(int *aa, int nn, int *bb, int mm, int *c) {n = nn; m = mm;bit = bitnum = 0;for (int i = 0; i <= n; i++) a[i] = aa[i];for (int i = 0; i <= m; i++) b[i] = bb[i];m += n;for (bit = 1; bit <= m; bit <<= 1)bitnum++;getrev(bitnum);NTT(a, 1);NTT(b, 1);for (int i = 0; i < bit; i++) a[i] = 1LL * a[i] * b[i] % P;NTT(a, -1);for (int i = 0; i < bit; i++) c[i] = a[i];for (int i = 0; i < bit; i++) a[i] = b[i] = 0;return bit;}
}int n;
int sa[N], sum[N], pre[N], suf[N], aa[N], bb[N];int main() {//freopen("0.txt", "r", stdin);int T; read(T);while (T--) {read(n);for (int i = 1; i <= n; i++) read(sa[i]);double ans = 0;for (int i = 1; i <= n; i++) {ans += 1LL * i * (n - i + 1) * sa[i] * sa[i];sum[i] = sum[i - 1] + sa[i] * sa[i];pre[i] = pre[i - 1] + sa[i];}for (int i = n; i >= 1; i--) suf[i] = suf[i + 1] + sa[i];ll r = sum[n];for (int len = 1; len <= n; len++) {ans -= r / double(len);int L = len + 1, R = n - len;if (L <= r) r += sum[R] - sum[L - 1];else r -= sum[L] - sum[R - 1];}for (int i = 1; i <= n; i++) {aa[i] = sa[i];bb[i] = sa[n - i + 1];}NTT::solve(aa, n + 1, bb, n + 1, aa);ll r1 = 0, r2 = 0;for (int len = 2; len <= n; len++) {r1 += aa[n - len + 2] - sa[n - len + 2] * suf[n - len + 3] - sa[len - 1] * pre[len - 2];r2 += r1;ans -= 2 * r2 / double(len);}printf("%.2lf\n", ans);for (int i = 0; i <= n + 1; i++) sum[i] = pre[i] = suf[i] = 0;}return 0;
}

本文来自互联网用户投稿，文章观点仅代表作者本人，不代表本站立场，不承担相关法律责任。如若转载，请注明出处。 如若内容造成侵权/违法违规/事实不符，请点击【内容举报】进行投诉反馈！