「集训队作业2020」Old Problem
「集训队作业2020」Old Problem

2020 年 12 月 17 日

给一个长度为 nn 的序列 aia_i,和 qq 组询问 (l,r,x)(l,r,x),表示求 i=lr(1aix)\displaystyle\prod_{i=l}^r\left(1-\frac{a_i}{x}\right) 的值。实数输出,精度要求 10610^{-6}

n,q6×105, 1ai<x109n,q\le6\times10^5,\ 1\leq a_i < x\leq 10^9

题解

按照 EI 的话说,这是一个误差分析题。

首先需要注意到这个式子可以泰勒展开:

exp(ln(i=lr(1aix)))=exp(i=lrln(1aix))=exp(i=lrk=1Laikkxk)\begin{aligned} &\exp\left(\ln\left(\prod_{i=l}^r\left(1-\frac{a_i}{x}\right)\right)\right) \\ =&\exp\left(\sum_{i=l}^r\ln\left(1-\frac{a_i}{x}\right)\right) \\ =&\exp\left(-\sum_{i=l}^r\sum_{k=1}^{L}\frac{a_i^k}{k\cdot x^k}\right) \\ \end{aligned}

然而,如果直接泰勒展开,需要的 LL 是数十万级别的,无益于我们解决问题。

考虑到导致泰勒展开精度损失的主要原因,是因为当 ai/xa_i/x 较大时,我们对 ln(1ai/x)\ln(1-a_i/x) 的精度要求很高。然而,ai/xa_i/x 较大时,很容易导致答案小于我们要求的精度范围。

故我们不妨设定一个阈值 R=0.5R=0.5,当 ai/xRa_i/x\le R 时,考虑线段树维护泰勒展开;否则,当 ai/x>Ra_i/x>R 时,优先找出这些位置并暴力计算。当答案小于精度要求时就退出。暴力计算的次数显然不会超过 log2106\log_2 10^6 次。

这样就得到了一个 O(nlog2L)O(n\log^2L) 的做法,实践得 LL2020 左右即可。

可以通过本题,但时间较大。实际上,线段树的部分可以直接换为前缀和。为什么精度还在接受范围内呢?注意到 x>aix>a_i 对于任意 xx 和任意 ii 都成立。如果前缀和的部分因为 aia_i 太小被省略,他本身对泰勒展开的影响也是被省略的级别。换句话说,对于泰勒展开的值,前缀和能保证的精度范围,恰为 double 本身的精度范围。

所以,直接将上述做法中的线段树替换为前缀和就能在 O(nlognL)O(n\log nL) 的时间复杂度内解决本题。三个 log 年轻人不讲武德。

代码

#include <bits/stdc++.h>
using namespace std;
const int N = 6e5 + 10, K = 20;
const double eps = 1e-10;
int n, m, x, a[N], lg[N], st[20][N];
double sum, ans;
mt19937 rng(20040725);
inline int better(int i, int j) { return a[i] < a[j] ? j : i; }
inline int query(int l, int r) {
  if (l == r) return l;
  int k = lg[r - l];
  return better(st[k][l], st[k][r - (1 << k)]);
}
struct segment {
  int l, r, mid;
  double s[K];
} p[N << 2];
void build(int u, int l, int r) {
  p[u].l = l, p[u].r = r, p[u].mid = (l + r) >> 1;
  if (l == r) {
    p[u].s[0] = a[l];
    for (int i = 1; i < K; i++) p[u].s[i] = p[u].s[i - 1] * a[l];
    return;
  }
  build(u << 1, l, p[u].mid);
  build(u << 1 | 1, p[u].mid + 1, r);
  for (int i = 0; i < K; i++) {
    p[u].s[i] = p[u << 1].s[i] + p[u << 1 | 1].s[i];
  }
}
double query(int u, int l, int r) {
  if (p[u].l == l && p[u].r == r) {
    double sum = 0;
    for (int i = K - 1; i >= 0; i--) {
      sum = (sum + p[u].s[i] / (i + 1)) / x;
    }
    return sum;
  }
  if (r <= p[u].mid) return query(u << 1, l, r);
  if (l > p[u].mid) return query(u << 1 | 1, l, r);
  return query(u << 1, l, p[u].mid) + query(u << 1 | 1, p[u].mid + 1, r);
}
void solve(int l, int r) {
  if (l > r || ans < eps) {
    return;
  }
  int pos = query(l, r);
  if (a[pos] < (x >> 1)) {
    ans *= exp(-query(1, l, r));
    return;
  }
  ans *= 1 - a[pos] / (double)x;
  if (rng() & 1) {
    solve(l, pos - 1);
    solve(pos + 1, r);
  } else {
    solve(pos + 1, r);
    solve(l, pos - 1);
  }
}
int main() {
#ifdef memset0
  freopen("1.in", "r", stdin);
#endif
  cin.tie(0)->sync_with_stdio(0);
  lg[0] = -1;
  for (int i = 1; i < N; i++) lg[i] = lg[i >> 1] + 1;
  cin >> n >> m;
  for (int i = 1; i <= n; i++) cin >> a[i];
  for (int i = 1; i < n; i++) {
    st[0][i] = better(i, i + 1);
  }
  for (int i = 1; i < 20; i++)
    for (int j = 1; j + (1 << i) <= n; j++) {
      st[i][j] = better(st[i - 1][j], st[i - 1][j + (1 << (i - 1))]);
    }
  build(1, 1, n);
  cout << fixed << setprecision(12);
  for (int l, r, i = 1; i <= m; i++) {
    cin >> l >> r >> x;
    ans = 1;
    solve(l, r);
    cout << 1 - ans << endl;
  }
}

评论

TABLE OF CONTENTS

题解