d p [ i ] = d p [ i − 1 ] ∗ p + d p [ i − 2 ] ∗ ( 1 − p ) dp[i]=dp[i-1]*p+dp[i-2]*(1-p) dp[i]=dp[i−1]∗p+dp[i−2]∗(1−p)
但是由于 i ∈ [ 1 , 100000000 ] i\in[1,100000000] i∈[1,100000000],可以使用矩阵快速幂来加速
比如对于两个地雷 i , j i,j i,j
我们计算出 d p [ i + 1 ] dp[i+1] dp[i+1]的概率,用矩阵加速算到 d p [ j − 1 ] dp[j-1] dp[j−1]
d p [ j + 1 ] = d p [ j − 1 ] ∗ ( 1 − p ) dp[j+1]=dp[j-1]*(1-p) dp[j+1]=dp[j−1]∗(1−p)
然后就可以一直使用这种加速方式了
构造矩阵是一个(1,2)的矩阵
[ d p [ i − 1 ] , d p [ i ] ] \begin{bmatrix} dp[i-1],dp[i] \end{bmatrix} [dp[i−1],dp[i]]
转移矩阵是一个(2,2)的矩阵
[ 0 , 1 − p 1 , p ] \begin{bmatrix} 0,1-p\\ 1,p\\ \end{bmatrix} [0,1−p1,p]
#include <iostream>
#include <algorithm>
using namespace std;
const int maxn=2009;
int a[maxn],n;
double p;
struct rce{
double m[3][3];
rce(){ for(int i=1;i<=2;i++)for(int j=1;j<=2;j++) m[i][j]=0; }
};
void print(rce a)
{
for(int i=1;i<=2;i++)
for(int j=1;j<=2;j++)
{
printf("%.2lf ",a.m[i][j] );
if( j==2 ) cout << endl;
}
}
rce init()
{
rce w;
for(int i=1;i<=2;i++) w.m[i][i]=1.0;
return w;
}
rce operator *(rce a,rce b)
{
rce ans;
for(int i=1;i<=2;i++)
for(int j=1;j<=2;j++)
{
ans.m[i][j]=0;
for(int k=1;k<=2;k++)
ans.m[i][j]+=a.m[i][k]*b.m[k][j];
}
return ans;
}
rce quick_pow(rce a,int n)
{
rce x=init();
while( n )
{
if( n&1 ) x=x*a;
a=a*a;
n>>=1;
}
return x;
}
int main()
{
while( cin >> n >> p )
{
for(int i=1;i<=n;i++) cin >> a[i];
sort(a+1,a+1+n);
if( a[1]==1 )
{
printf("0.0000000\n");
continue;
}
double last=1.0,now=p;
int x=2;//(1,2)
for(int i=1;i<=n;i++)
{
rce chu,gou;
chu.m[1][1]=last,chu.m[1][2]=now;
gou.m[1][1]=0,gou.m[1][2]=1-p;
gou.m[2][1]=1,gou.m[2][2]=p;
gou = quick_pow(gou,a[i]-x );//(a[i-1],a[i])
chu = chu*gou;
last = chu.m[1][1],now = 0;
x=a[i];
}
printf("%.7f\n",last*(1-p));
}
}
