洛谷P3779 [SDOI2017]龙与地下城（概率论+Simpson+FFT）

心已赠人 2022-01-12 09:25 206阅读 0赞

## 题面 ##

[传送门][Link 1]

## 题解 ##

orz [shadowice][]

### 正态分布 ###

正态分布是随机变量\\(X\\)的一种概率分布形式。它用一个期望\\(\\mu\\)和方差\\(\\sigma^2\\)就可以描述，记为\\(N(\\mu,\\sigma^2)\\)。  
若随机变量\\(X\\)服从一个数学期望为\\(\\mu\\)、方差为\\(\\sigma^2\\)的正态分布，记作\\(X\\)∼\\(N(\\mu,\\sigma^2)\\)，读作\\(X\\)服从\\(N(\\mu,\\sigma^2)\\)。  
当\\(\\mu=0,\\sigma=1\\)时的正态分布称为标准正态分布。

### 概率密度函数 ###

用来描述连续型随机变量的分布情况。随机变量的取值落在某个区域内的概率，为概率密度函数在该区域的积分。（或者就是\\(f(x)\\)在该区域内与\\(x\\)轴围成的图形面积）  
若随机变量\\(X\\)∼\\(N(\\mu,\\sigma^2)\\)，则其概率密度函数为

\\\[f(x)=\\frac\{1\}\{\\sqrt\{2\\pi\}\\sigma\}e^\{-\\frac\{(x-\\mu)^\{2\}\}\{2\\sigma^\{2\}\}\}\\\]

只要证明了一个变量服从正态分布，就可以直接对概率密度函数的这一区间进行积分了

### 中心极限定理 ###

当样本量n逐渐趋于无穷大时，n个抽样样本的均值的频数逐渐趋于正态分布（无论总体是什么分布）。  
该定理说明，设随机变量\\(X\_1,X\_2,…,X\_n\\)独立同分布，它们的期望为\\(\\mu\\)、方差为\\(\\sigma^2\\)，当\\(n\\)足够大时（满足精确度需求时），随机变量

\\\[Y\_\{n\}=\\frac\{\\sum\_\{i=1\}^\{n\}x\_\{i\}-n\\mu\}\{\\sqrt\{n\}\\sigma\}\\\]

\\(Y\_n\\)服从正态分布，求出其范围后就可以直接对正态分布的概率密度函数求积分了。

### 本题 ###

然后关于本题有

\\\[\\mu=\\frac\{n-1\}\{2\}，\\sigma^2=\\frac\{\\sum\_\{i=1\}^n(i-\\mu)^2\}\{n\}=\\frac\{n^2-1\}\{12\}\\\\\\sum\_\{i=1\}^nX\_i\\in\[A,B\]\\\\Y\_n\\in\[\\frac\{A-n\\mu\}\{\\sqrt n\\sigma\},\\frac\{B-n\\mu\}\{\\sqrt n\\sigma\}\]\\\]

那么我们对于\\(Y\_n\\)的值域直接上自适应辛普森法就行了

然而还有一个问题，就是\\(n\\)等于\\(1\\)之类的情况我们显然不能认为它满足精度限制

那么我们构造出概率的生成函数，\\(F(z)=\\sum\_\{i=0\}^\{x-1\}\{1\\over x\}z^i\\)，只要求出\\(F^y(z)\\)的第\\(l\\)到\\(r\\)项系数就行了。而且这里多项式快速幂可以直接把点值给快速幂

听说余奶奶直接用生成函数艹过去了……完全没用到任何概率论的芝士……orz

//minamoto #include<bits/stdc++.h> #define R register #define fp(i,a,b) for(R int i=(a),I=(b)+1;i<I;++i) #define fd(i,a,b) for(R int i=(a),I=(b)-1;i>I;--i) #define go(u) for(int i=head[u],v=e[i].v;i;i=e[i].nx,v=e[i].v) using namespace std; const int N=(1<<19)+5;const double Pi=acos(-1.0),K=1.0/sqrt(2*Pi),eps=1e-7; int read(){ R int res,f=1;R char ch; while((ch=getchar())>'9'||ch<'0')(ch=='-')&&(f=-1); for(res=ch-'0';(ch=getchar())>='0'&&ch<='9';res=res*10+ch-'0'); return res*f; } struct cp{ double x,y; cp(){} cp(R double xx,R double yy):x(xx),y(yy){} inline cp operator +(const cp &b)const{return cp(x+b.x,y+b.y);} inline cp operator -(const cp &b)const{return cp(x-b.x,y-b.y);} inline cp operator *(const cp &b)const{return cp(x*b.x-y*b.y,x*b.y+y*b.x);} inline cp operator *(const double &b)const{return cp(x*b,y*b);} }A[N],O[N]; int r[N],lim,l,X,Y,L;double res; inline cp ksm(R cp x,R int y){ R cp res(1,0); for(;y;y>>=1,x=x*x)if(y&1)res=res*x; return res; } void FFT(cp *A,int ty){ fp(i,0,lim-1)if(i<r[i])swap(A[i],A[r[i]]); for(R int mid=1;mid<lim;mid<<=1) for(R int j=0;j<lim;j+=(mid<<1)) fp(k,0,mid-1){ R cp x=A[j+k],y=A[j+k+mid]*O[mid+k]; A[j+k]=x+y,A[j+k+mid]=x-y; } if(ty==-1){ reverse(A+1,A+lim); R double k=1.0/lim;fp(i,0,lim-1)A[i]=A[i]*k; } } inline double F(R double x){return K*exp(-x*x*0.5);} inline double simpson(R double l,R double r){return (F(l)+F(r)+4*F((l+r)/2))*(r-l)/6;} double calc(double l,double r,double eps,double res){ double mid=(l+r)/2,ql=simpson(l,mid),qr=simpson(mid,r); if(fabs(ql+qr-res)<=15*eps)return ql+qr+(ql+qr-res)/15; return calc(l,mid,eps/2,ql)+calc(mid,r,eps/2,qr); } void solve1(){ fp(i,0,lim-1)A[i]=cp(0,0),r[i]=(r[i>>1]>>1)|((i&1)<<(l-1)); for(R int i=1;i<lim;i<<=1)fp(k,0,i-1)O[i+k]=cp(cos(Pi*k/i),sin(Pi*k/i)); R double xx=1.0/X; fp(i,0,X-1)A[i].x=xx; FFT(A,1); fp(i,0,lim-1)A[i]=ksm(A[i],Y); FFT(A,-1); fp(i,1,lim-1)A[i].x+=A[i-1].x; for(R int i=1,l,r;i<=10;++i)l=read(),r=read(),printf("%.7lf\n",l?A[r].x-A[l-1].x:A[r].x); } void solve2(){ double l,r,mu=1.0*(X-1)/2,sigma=1.0*(X*X-1)/12,a=mu*Y,b=sqrt(sigma*Y); fp(i,1,10){ l=1.0*(read()-a)/b,r=1.0*(read()-a)/b; printf("%.7lf\n",calc(0,r,eps,simpson(0,r))-calc(0,l,eps,simpson(0,l))); } } int main(){ // freopen("testdata.in","r",stdin); for(int T=read();T;--T){ X=read(),Y=read(),L=X*Y,lim=1,l=0; while(lim<=L)lim<<=1,++l; lim<N?solve1():solve2(); } return 0; }

转载于:https://www.cnblogs.com/bztMinamoto/p/10529345.html

[Link 1]: https://www.luogu.org/problemnew/show/P3779
[shadowice]: https://www.luogu.org/blog/ShadowassIIXVIIIIV/solution-p3779#