gmoj 6087. 【GDOI2019模拟2019.3.26】获取名额 题解

题目

https://gmoj.net/senior/#main/show/6087

题解

发现
$$ans=1-\prod _{i=l}^r\left(1-\frac{a_i}{x}\right)$$
$\prod$ 比较难搞，就把它变成 $\sum$ 吧！
$$ans=1-e^{{\sum }_{i=l}^r\ln \left(1-\frac{a_i}{x}\right)}$$
$\ln$ 比较难搞，就把它变成多项式函数吧！

为了方便处理，先弄出 $\ln (x+1)$ 的情况，再把 $-\frac{a_i}{x}$ 代进去。

令 $f(x)=\ln (x+1)$ 。

由 泰勒公式 ,得
$$f(x)=\sum _{i=0}^{\infty }\frac{f^{(i)}(x_0)}{i!}(x-x_0{)}^i$$
其中 $f^{(i)}(x)$ 表示 $f(x)$ 的 $i$ 阶导数。

先求出 $(\ln (x){)}^{\prime }$ 。

设 $y=\ln x$ ，那么
$$(\ln x{)}^{\prime }=\frac{1}{(e^y{)}^{\prime }}=\frac{1}{e^y}=\frac{1}{x}$$
所以
$$\begin{array}{rl}{f}^{\prime }(x)& =[\ln (x+1){]}^{\prime }(x+1{)}^{\prime }=\frac{1}{x+1}=(x+1{)}^{-1}\\ f^{\prime \prime }(x)& =-(x+1{)}^{-2}\\ f^{\prime \prime \prime }(x)& =2(x+1{)}^{-3}\\ & \cdots \\ f^{(n)}(x)& =(-1{)}^{n+1}(n-1)!(x+1{)}^{-n}\end{array}$$

那么
$$\begin{array}{rlr}f(x)& =\sum _{i=0}^{\infty }\frac{f^{(i)}(x_0)}{i!}(x-x_0{)}^i& \\ & =\sum _{i=0}^{\infty }\frac{f^{(i)}(0)}{i!}{x}^i& \text{令 x0=0 }\\ & =\sum _{i=1}^{\infty }\frac{(-1{)}^{i-1}(i-1)!}{i!}{x}^i& \text{ln1=0}\\ & =\sum _{i=1}^{\infty }\frac{(-1{)}^{i-1}{x}^i}{i}& \end{array}$$
将 $-\frac{a_i}{x}$ 代入，得
$$\begin{array}{rl}\sum _{i=l}^r\ln \left(1-\frac{a_i}{x}\right)& =\sum _{i=l}^{r}f\left(\frac{a_i}{x}\right)\\ & =-\sum _{i=l}^r\sum _{j=1}^{\infty }\frac{a_i^j}{j\cdot x^j}\end{array}$$
在实际上做的时候，第二重循环只用枚举到 $20$ 精度就足够了，因此用 前缀和 维护一个 ${\sum }_{i=l}^r{\sum }_{j=1}^{20}\frac{a_i^j}{j}$ 就可以了。

为了防止 $a_i^j$ 爆 long double ，可以把每一个 $a_i$ 和 $x$ 都除以一个 ${\max }_{i=1}^n{a}_i$ 。

还可以发现当 $\frac{a_i}{x}$ 比较大的时候，由于它离 $x_0$ 远，会有较大误差。这时只要暴力地把这些 $\frac{a_i}{x}$ 乘在一起而不把它们放到 $\sum$ 中计算即可，如果这些 $\frac{a_i}{x}$ 较多答案就是 $0$ 了。可以把阈值设置为 $0.5$ 。

---

代码

#include<cstdio>
#include<cmath>
using namespace std;
#define fo(i,l,r) for(i=l;i<=r;++i)
typedef long double LF;
char p[15],buf[100005],*l=buf,*r=buf;
inline char gc()
{return l==r&&(r=(l=buf)+fread(buf,1,100005,stdin),l==r)?EOF:*l++;}
inline void read(int &k)
{
	char ch;while(ch=gc(),ch<'0'||ch>'9');k=ch-48;
	while(ch=gc(),ch>='0'&&ch<='9') k=k*10+ch-48;
}
#define N 600005
const int lim=20;
LF a[N],s[N][lim+1],x,ans1,ans2;
int A[N],f[N][20],highbit[N];
inline int mymax(int x,int y){return a[x]>a[y]?x:y;}
void solve(int l,int r)
{
	if(l>r) return;
	int mid=mymax(f[r-(1<<highbit[r-l+1])+1][highbit[r-l+1]],f[l][highbit[r-l+1]]);
	if(a[mid]/x>0.5)
	{
		ans2*=1-a[mid]/x;
		if(ans2<1e-7) return;
		solve(l,mid-1),solve(mid+1,r);
	}
	else
	{
		LF tmp=x;
		for(register int i=1;i<=lim;++i,tmp*=x)
			ans1-=(s[r][i]-s[l-1][i])/tmp;
	}
}
int main()
{
	freopen("orz.in","r",stdin);
	freopen("orz.out","w",stdout);
	register int i;
	int j,k,n,q,l,r,y,maxx=0;
	LF max,tmp;
	read(n),read(q);
	fo(i,1,n)
	{
		read(A[i]);
		if(A[i]>maxx) maxx=A[i];
	}
	max=maxx;
	fo(i,1,n) a[i]=A[i]/max;
	j=0;
	fo(i,1,n)
	{
		if(i==1<<j+1) ++j;
		highbit[i]=j;
	}
	fo(i,1,n) f[i][0]=i;
	fo(j,1,19)
	{
		k=n-(1<<j)+1;
		fo(i,1,k)
			f[i][j]=mymax(f[i][j-1],f[i+(1<<j-1)][j-1]);
	}
	fo(i,1,n)
	{
		tmp=a[i];
		fo(j,1,lim)
		{
			s[i][j]=s[i-1][j]+tmp/j,
			tmp*=a[i];
		}
	}
	while(q--)
	{
		read(l),read(r),read(y);
		x=y/max,ans1=0,ans2=1;
		solve(l,r);
		printf("%.7Lf\n",1-exp(ans1)*ans2);
	}
	return 0;
}