POJ3744 Scout YYF I ( 矩阵快速幂 + 期望概率DP )

POJ3744 Scout YYF I

原题地址：http://poj.org/problem?id=3744

题意：
YYF在1点，对于之后的每步，YYF将以p的概率走一步，或以1-p的概率跳跃两步。 有n个地雷，给出坐标。请问YYF安全走过地雷区的概率为多少，保留7位小数。
多组数据，处理到EOF为止。

数据范围
1 ≤ N ≤ 10， 0.25 ≤ p ≤ 0.75，地雷坐标∈[1, 100000000].

题解：
矩阵快速幂优化DP。

O(n)递推很好写：
dp[i]=p * dp[i-1]+(1-p) * dp[i-2]
但是，由于地雷坐标∈[1, 100000000]，这样肯定TLE。
考虑每两个地雷之间，用矩阵快速幂优化到 log
可以写出这样的转移方程：
$\bigl( \begin{smallmatrix} dp[i] \\ dp[i-1] \end{smallmatrix} \bigr)$ =$\bigl( \begin{smallmatrix} p& 1-p \\ 1&0\end{smallmatrix} \bigr)$ X $\bigl( \begin{smallmatrix} dp[i-1] \\ dp[i-2] \end{smallmatrix} \bigr)$

其中，初始化 把dp[0]当做有地雷 ，dp[0]=0，dp[1]=1。

坑点：给出雷的坐标不是有序的，需要先排序。

——————————————————
顺便复习下矩阵快速幂。

$\ f(n) = a_1 *f(n-1)+a_2*f(n-2)+...+a_k *f(n-k)+c$

$$\begin{bmatrix} f(n) \\ f(n-1) \\ \vdots \\ f(n-k+1) \\ 1 \\ \end{bmatrix} \tag{1}$$

$$\left[ \begin{matrix} a_1 & a_2 & \cdots &a_{k-1}& a_k&c \\ 1 & 0 & \cdots & 0&0&0 \\ 0 &1 \\ \vdots & \vdots & \ddots & \vdots& \vdots& \vdots \\ 0 & 0 & \cdots & 1 &0&0 \\ 0 &0 & \cdots & 0 &0&1 \\ \end{matrix} \right]\tag{2}$$
$$\begin{bmatrix} f(n-1) \\ f(n-2) \\ \vdots \\ f(n-k) \\ 1 \\ \end{bmatrix} \tag{3}$$
$$(1)=(2) * (3)$$
代码：

#include<cstdio>
#include<iostream>
#include<algorithm>
#include<cstring>
using namespace std;
const int N=15;
int n,pos[N];
double p;
struct node
{
    double mat[2][2];
    void clean()
    {
        mat[0][0]=1; mat[0][1]=0;
        mat[1][0]=0; mat[1][1]=1;
    }
    void init()
    {
        mat[0][0]=p; mat[0][1]=1-p;
        mat[1][0]=1; mat[1][1]=0;
    }
    void mem()
    {
        memset(mat,0,sizeof(mat));
    }
};
node mul(node A,node B,int r,int c)
{
    node C; C.mem();
    for(int i=0;i<=r;i++)
    for(int j=0;j<=c;j++)
    for(int k=0;k<=1;k++)
    C.mat[i][j]+=A.mat[i][k]*B.mat[k][j];   
    return C;   
}
node mulpow(node A,int B)
{
    node ans; ans.clean();
    node base; base.init();
    for(;B;B>>=1)
    {
        if(B&1) 
        ans=mul(ans,base,1,1);
        base=mul(base,base,1,1);
    }
    return mul(ans,A,1,0);
}
double solve()
{   
    for(int i=2;i<=n;i++)
    if(pos[i]==pos[i-1]+1) return 0;

    double ans=0;
    if(pos[1]==1) return 0;
    int i=1; int tmp=1;
    node A; 
    double last[2]={1.0,0.0};
    while(1)
    {
        A.mat[0][0]=last[0]; A.mat[1][0]=last[1];
        int nxt=pos[tmp]-1;
        A=mulpow(A,nxt-i);
        last[1]=A.mat[0][0];
        last[0]=0;
        if(tmp==n) break;
        i=pos[tmp],tmp++; 
    }
    A.mat[0][0]=last[0]; A.mat[1][0]=last[1];
    A=mulpow(A,1);
    return A.mat[0][0];
}
int main()
{   
    while(scanf("%d",&n)!=EOF)
    {
        scanf("%lf",&p);
        for(int i=1;i<=n;i++)
        scanf("%d",&pos[i]);
        sort(pos+1,pos+n+1);
        printf("%0.7lf\n",solve());
    }   
    return 0;
}