龙与地下城
Description
小Q同学是一个热爱学习的人,但是他最近沉迷于各种游戏,龙与地下城就是其中之一。
在这个游戏中,很多场合需要通过掷骰子来产生随机数,并由此决定角色未来的命运,因此骰子堪称该游戏的标志性道具。
骰子也分为许多种类,比如4面骰、6面骰、8面骰、12面骰、20面骰,其中20面骰用到的机会非常多。当然,现在科技发达,可以用一个随机数生成器来取代真实的骰子,所以这里认为骰子就是一个随机数生成器。
在战斗中,骰子主要用来决定角色的攻击是否命中,以及命中后造成的伤害值。举个例子,假设现在已经确定能够命中敌人,那么 YdX (也就是掷出Y个X面骰子之后所有骰子显示的数字之和)就是对敌人的基础伤害。在敌人没有防御的情况下,这个基础伤害就是真实伤害。
众所周知,骰子显示每个数的概率应该是相等的,也就是说,对于一个XX面骰子,显示 0,1,2,…,X−1 中每一个数字的概率都是 1x
更形式地说,这个骰子显示的数W满足离散的均匀分布,其分布列为
除此之外还有一些性质
W的一阶原点矩(期望)为
v1(W)=E(W)=∑X−1i=0iP(W=i)=X−12
W的二阶中心矩(方差)为
μ2(W)=E((W−E(W))2)=∑X−1i=0(i−E(W))2P(W=i)=X2−112
言归正传,现在小Q同学面对着一个生命值为A的没有防御的敌人,能够发动一次必中的
YdX
攻击,显然只有造成的伤害不少于敌人的生命值才能打倒敌人。但是另一方面,小Q同学作为强迫症患者,不希望出现overkill,也就是造成的伤害大于B的情况,因此只有在打倒敌人并且不发生overkill的情况下小Q同学才会认为取得了属于他的胜利。
因为小Q同学非常谨慎,他会进行10次模拟战,每次给出敌人的生命值A以及overkill的标准B,他想知道此时取得属于他的胜利的概率是多少,你能帮帮他吗?
Input
第一行是一个正整数T,表示测试数据的组数,
对于每组测试数据:
第一行是两个整数X , Y,分别表示骰子的面数以及骰子的个数;
接下来10行,每行包含两个整数A , B,分别表示敌人的生命值A以及overkill的标准B。
Output
对于每组测试数据,输出10行,对每个询问输出一个实数,要求绝对误差不超过 0.013579 ,也就是说,记输出为aa,答案为bb,若满足 |a−b|≤0.013579 ,则认为输出是正确的。
Sample Input
1
2 19
0 0
0 1
0 2
0 3
0 4
0 5
0 6
0 7
0 8
0 9
Sample Output
0.000002
0.000038
0.000364
0.002213
0.009605
0.031784
0.083534
0.179642
0.323803
0.500000
Hint
对于 100% 的数据, T≤10 , 2≤X≤20 , 1≤Y≤200000 , 0≤A≤B≤(X−1)Y ,保证满足 Y>800 的数据不超过22组。
被题目名字吸引进来了……
居然没有找到题解……
然后就被坑了好久
思路:
很显然可以有一个DP。
令
f[i][j]
表示丢了
i
个骰子,点数为
那么显然很好转移:
f[i][j]=∑k≤j,k≤x−1f[i−1][j−k]∗1n
然后可以发现这是个卷积,那么FFT成点值表达式后快速幂直接上即可~
才怪。这只有70分。
正当咱完全没有思路的时候,知乎上此题的出题人说,这题的idea 基于“中心极限定理”。
学习一下这个东西可以发现,其实出题人想要让咱们用正态分布:
所谓正态分布,就是满足位置参数为
μ
、尺度参数为
σ
的概率分布的随机变量的分布曲线。
其中,
μ
为期望,
σ
为方差的平方根。
具体详见百度百科,本蒟蒻水平有限~
可以发现出题人已经在题面中将期望和方差给出了。然而并没有什么人知道这两句话有用
同样可以发现然而咱还是没发现,这题的“掷骰子”事件满足正态分布。
那么咱们可以根据正态分布的曲线来估计答案!
不加证明的(不会证明)给出正态分布的解析式:
f(x)=12π√σexp(−(x−μ)22σ2)
同时给出正态分布的图像:
然后现在目标转换为,求图上两点间函数图像的面积。
考虑到求函数的面积,还是这么奇怪的面积,显然需要估算了~
听说考场上有人用拆成小梯形的方法过了,然而明显这样做精度比较令人着急……
那么考虑更靠谱一点的方法:自适应辛普森积分。
辛普森积分嘛,就是用一根迷之二次函数曲线,试图凑出某个奇怪的函数在某一小部分内围成的面积的方法~原理自行百度~
自适应辛普森积分就是,不停二分区间,如果到某一时刻,用辛普森积分对当前的
(l,r)
积分得到的结果
f(l,r)
,与
f(l,mid)+f(mid,r)
之间的差很小,咱就认为它足够精确了。
这算法的原理有种自欺欺人的感觉
然而这样积分出面积,只有80pts。
因为刚才咱们使用的是普通的正态分布,有很多除法,不够精确。
需要化成标准正态分布以提高精度。
具体做法是,对于传入的参数
x
,进行如下处理:
把
y
作为以前的
同时将正态分布的公式换成标准正态分布:
f(x)=12π√σexp(−x22)
这样就少了很多除法,可过100pts!
注意,小数据由于精度误差,仍需要使用dp(咱用的FFT)
#include<bits/stdc++.h>
using namespace std;
inline int read()
{
int x=0;char ch=getchar();
while(ch<'0' || '9'<ch)ch=getchar();
while('0'<=ch && ch<='9')x=x*10+(ch^48),ch=getchar();
return x;
}
typedef double db;
db x,y;
namespace Fast_Fouriet_Transform
{
typedef complex<db> cp;
const int M=800009;
const db pi=acos(-1);
cp a[M],b[M];
db per,ans[M];
int m,l,rev[M];
inline void FFT(cp *a,int n,int f)
{
for(int i=0;i<n;i++)
if(i<rev[i])
swap(a[i],a[rev[i]]);
for(int h=2;h<=n;h<<=1)
{
cp w(cos(2*pi/h),f*sin(2*pi/h));
for(int j=0;j<n;j+=h)
{
cp wn(1,0),x,y;
for(int k=j;k<j+(h>>1);k++)
{
x=a[k],y=wn*a[k+(h>>1)];
a[k]=x+y;
a[k+(h>>1)]=x-y;
wn*=w;
}
}
}
if(f==-1)
for(int i=0;i<n;i++)
a[i].real()/=(db)n;
}
int mian()
{
per=1.0/(db)x;
for(m=1,l=0;m<(y*x);m<<=1)l++;
for(int i=0;i<m;i++)
rev[i]=(rev[i>>1]>>1)|((i&1)<<(l-1));
memset(a,0,sizeof(a));
memset(b,0,sizeof(b));
a[0].real()=1.0;
for(int i=0;i<x;i++)
b[i].real()=per;
FFT(a,m,1);FFT(b,m,1);
for(int i=0;i<m;i++)
{
int cnt=y;
while(cnt)
{
if(cnt&1)
a[i]=a[i]*b[i];
b[i]=b[i]*b[i];
cnt>>=1;
}
}
FFT(a,m,-1);
ans[0]=a[0].real();
for(int i=1;i<m;i++)
ans[i]=a[i].real()+ans[i-1];
for(int i=1;i<=10;i++)
{
int a=read(),b=read();
if(a==0)
printf("%.8lf\n",ans[b]);
else
printf("%.8lf\n",ans[b]-ans[a-1]);
}
return 0;
}
}
namespace simpsons
{
const db pi2=sqrt(acos(-1)*2.0);
const db eps=1e-12;
db sigma,mu;
inline db sqr(db x){return x*x;}
inline db f(db x){return exp(-sqr(x)/2)/pi2;}
inline bool dcmp(db x){return fabs(x)<15*eps;}
inline void init()
{
sigma=sqrt((sqr(x)-1.0)/12.0);
mu=(x-1.0)/2.0;
}
inline db simpson(db a,db b,db fl,db fmid,db fr)
{
return (b-a)*(fl+4.0*fmid+fr)/6.0;
}
inline db solve(db l,db r)
{
db mid=(l+r)/2.0,lm=(l+mid)/2.0,mr=(mid+r)/2.0;
db fl=f(l),fmid=f(mid),fr=f(r),flm=f(lm),fmr=f(mr);
db sipl=simpson(l,mid,fl,flm,fmid);
db sipr=simpson(mid,r,fmid,fmr,fr);
db sipm=simpson(l,r,fl,fmid,fr);
if(dcmp(sipl+sipr-sipm))
return sipl+sipr+(sipl+sipr-sipm)/15;
else return solve(l,mid)+solve(mid,r);
}
int mian()
{
init();
db base=solve(0,(x-1)*y);
for(int i=1;i<=10;i++)
{
db a=(((db)read()-0.5)/y-mu)*sqrt(y)/sigma;
db b=(((db)read()+0.5)/y-mu)*sqrt(y)/sigma;
printf("%.8lf\n",solve(0,b)-solve(0,a));
}
}
}
int main()
{
int T=read();
while(T--)
{
scanf("%lf%lf",&x,&y);
if(x*y<=400000)
Fast_Fouriet_Transform::mian();
else
simpsons::mian();
}
return 0;
}