Scout YYF I POJ - 3744（概率+矩阵快速幂）

题目地址：http://poj.org/problem?id=3744
题意：路上有N个雷（踩雷之后你会牺牲），初始时你在第一单元，每次可以走一个单元或者两个，概率分别为p，1-p；路长最大为100000000，求你安全通过的概率。
思路：
第一步：
设dp[i]是安全通过第i单元的概率，那么初识时ddp[0]=0, p[1]=1，递推关系为：dp[i]=p * dp[i-1]+(i-p) * dp[i-2].当i位置有雷时，让dp[i]=0.
第二步：
完成第一步，我们会发现路长最大为100000000，for循环会超时，而N最大只有10，我们可以用矩阵快速幂来优化
$$\left[\begin{array}{c}dp[i]\\ dp[i-1]\end{array}\right]=\left[\begin{array}{cc}1& 1\\ 1& 0\end{array}\right]\left[\begin{array}{c}dp[i-1]\\ dp[i-2]\end{array}\right]\text{\hspace{1em}(1)}$$
$$\left[\begin{array}{c}dp[x2]\\ dp[x2-1]\end{array}\right]=\left[{\begin{array}{cc}1& 1\\ 1& 0\end{array}}^i\right]\left[\begin{array}{c}dp[x1]\\ dp[x1-1]\end{array}\right]\text{\hspace{1em}(2)}$$
i是两个雷的距离,i=x2-x1，然后让dp[x2]=0，因为x2单元上有雷，然后以x2为起点进行下一个快速幂，直到走最后一个雷Xn，那么安全概率为dp[xn-1] * (i-p)

#include <cstdio>
#include <algorithm>
#include <iostream>
#include <cstring>
using namespace std;
 double p;
struct Matrix
{
    double m[2][2];
}P,I;
Matrix matrixmul(Matrix a,Matrix b)///a *b
{
    int i,j,k;
    Matrix c;
    for(i=0;i<2;i++)
    for(j=0;j<2;j++)
    {
        c.m[i][j]=0;
        for(k=0;k<2;k++)
        c.m[i][j]+=a.m[i][k]*b.m[k][j];
    }
    return c;
}
Matrix quick(int n)
{
    Matrix m=P,b=I;
    while(n>0)
    {
        if(n&1)
         b=matrixmul(b,m);
        n=n>>1;
        m=matrixmul(m,m);
    }
    return b;
}
void init()
{
     for(int i=0;i<2;i++)
     for(int j=0;j<2;j++)
      I.m[i][j]=0;
    I.m[0][0]=I.m[1][1]=1;
    P.m[0][0]=p; P.m[0][1]=1-p;
    P.m[1][0]=1;P.m[1][1]=0;
}
int a[12];
int main()
{
   int n;
   Matrix ans,ans1;
    while(scanf("%d",&n)!=EOF)
    {
        scanf("%lf",&p);
        init();
        for(int i=0;i<n;i++)
        scanf("%d",&a[i]);
        sort(a,a+n);
        int x;
        for(int i=0;i<n;i++)
        {
            if(i==0)
              x=a[i]-1;
            else
            x=a[i]-a[i-1];
             ans=quick(x);
             if(i==0)
             {
                 ans1.m[0][0]=1; ans1.m[0][1]=0;ans1.m[1][0]=0;ans1.m[1][1]=0;
             }
             else
             {
                 ans1.m[0][0]=0;
             }
             ans1=matrixmul(ans,ans1);
        }
        double sum=ans1.m[1][0]*(1-p);
        printf("%.7lf\n",sum);
    }
    return 0;
}