問題概要
$N$ マスのすごろくを2人で遊んでいる。サイコロは等確率で $1$ から $d$ までの目を出し、出た目の数だけマスを進める。ただし、ゴールを超える場合はその分だけ折り返す。
先手と後手の現在位置はそれぞれ $x,y$ で、今先手の手番である。
先手が後手より先にゴールに到達する確率を求めよ。
制約
$10 \leqq N \leqq 5000$
$1 \leqq d,x,y \leqq N-1$
解法
まずこのルールを観察すると、$N-d, N-d+1, \ldots , N-1$ のマスは同等なものとして扱っていいことが分かる。なぜなら、どのマスでも等確率 $p$ でゴールに到達し、等確率 $q$ でこれらのマスのいずれかに移動するからである。
先手と後手が共に上の位置にいる時の先手の勝利確率を考えていく。
これは、初項 $p$ 公差 $q^{2}$ の無限等比級数の和として表すことができる。
$$p + p \cdot q^{2} + p \cdot q^{4} + \ldots = \frac{p}{1 - q^{2}}$$
$p=\frac{1}{d}, q=\frac{d-1}{d}$ であるから、この値は $\frac{d}{2d-1}$ に等しくなる。
各プレイヤーはサイコロを振るたびに少なくとも1マスずつ進むことが出来るので、高々 $N-d$ 回未満の操作で上で述べた状態、もしくは試合が既に終了している、のいずれかに遷移することが言える。
従って、 $N-d$ 回未満の操作で先手が勝利する確率 $P_{1}$ と、試合の決着が付かない確率 $P _ {2}$ を用いることで、先手の勝利確率は $P _ {1}+P _ {2}\cdot\frac{d}{2d-1}$ と表すことができる。
これで無限回の操作から、有限回の操作の結果における確率が求める対象となり、一気に見通しが良くなった。
ここで、先手と後手の操作を独立して考えてみる。つまり対戦相手のいないすごろくを仮定する。
丁度 $i$ 回の操作でゴールに到達する確率を先手後手それぞれに対して $A _ {i} , B _ {i}$ と定めると、今求めたい $P _ {1} , P _ {2}$ は次のように書ける。
$$P _ {1} = A _ {1} + (1 - B _ {1})A _ {2} + \ldots + (1 - B _ {1} - \ldots - B _ {N-d-2})A _ {N-d-1}$$
$$P _ {2} = (1 - A _ {1} - A _ {2} - \ldots - A _ {N-d-1})(1 - B _ {1} - B _ {2} - \ldots - B _ {N-d-1})$$
各 $i$ に対する $A _ {i} , B _ {i}$ が求まればこれは $O(N)$ で計算可能。
残るステップは $A _ {i} , B _ {i}$ だが、これらも $O(N^{2})$ のDPで求めることができる。
実装
class BouncingDiceGame {
public:
double winProbability(int n, int d, int x, int y) {
x--; y--;
const int N=n-d;
const double p=1.0/(double)d;
vector<double> A(N),B(N);
{// A_{i}を求める
vector<double> dp(n);
vector<double> ndp(n);
dp[x]=1.0;
for(int i=1;i<N;i++){
double sum1=0.0;
double sum2=0.0;
for(int i=0;i<n;i++){
if(i+d>n-1)sum2+=dp[n-1-(i+d-n+1)];
if(i+1==n)sum2=0;
ndp[i]=sum1*p+sum2*p;
sum1+=dp[i];
if(i-d>=0)sum1-=dp[i-d];
}
A[i]=ndp[n-1];
swap(dp,ndp);
}
}
{// B_{i}を求める
vector<double> dp(n);
vector<double> ndp(n);
dp[y]=1.0;
for(int i=1;i<N;i++){
double sum1=0.0;
double sum2=0.0;
for(int i=0;i<n;i++){
if(i+d>n-1)sum2+=dp[n-1-(i+d-n+1)];
if(i+1==n)sum2=0;
ndp[i]=sum1*p+sum2*p;
sum1+=dp[i];
if(i-d>=0)sum1-=dp[i-d];
}
B[i]=ndp[n-1];
swap(dp,ndp);
}
}
double P1=0.0,P2=0.0;
double sumA=0.0,sumB=0.0;
for(int i=1;i<N;i++){
P1+=(1.0-sumB)*A[i];
sumA+=A[i],sumB+=B[i];
}
P2=(1.0-sumA)*(1.0-sumB);
double res=P1+P2*d/(double)(2*d-1);
return res;
}
};
