Isa@Diary

ソフトウェア開発やってます。プログラミングとか、US生活とかについて書きます。

AOJ 1056 Ben Toh

方針

DP。誤差制限が緩いので
dp[n][k] := n回目にゲットできる確率が2^-kである確率
として

dp[i+1][0] += dp[i][k] * (1 - pow(0.5,k)) //取れなかったら次の確率は1
dp[i+1][k+1] += dp[i][k] * pow(0.5,k);    //取れたら次の確率は更に1/2

とすればよい。

dp[100001][30]*1とするとREになる*2のでdp[2][30]と逐次で計算していき、テーブル更新と同時に期待値も計算していけばいい。
+=を使う場合は足す先を初期化しないといけないことに注意。

以下ソース

#include <iostream>
#include <algorithm>
#include <cstdio>
#include <cmath>

#define FILL(ptr, value) FILL_((ptr), sizeof(ptr)/sizeof(value), (value))
 
template <typename T>
void FILL_(void * ptr, size_t size, T value){
  std::fill((T*)ptr, (T*)ptr+size, value);
}

using namespace std;

int main(void){
	double p[2][100],ret;
	int n;
	while(cin >> n && n){
		int cur=0,nxt=1;
		ret=0;
		FILL(p,0);
		p[cur][0] = 1;
		for(int i=0;i<n;i++){
			for(int j=0;j<100;j++){
				if(p[cur][j] == 0) continue;
				p[nxt][0] += p[cur][j] * (1 - pow(0.5,j));
				p[nxt][j+1] += p[cur][j] * pow(0.5,j);
				ret += p[cur][j] * pow(0.5,j);
			}
			cur = 1 - cur;
			nxt = 1 - nxt;
			for(int j=0;j<100;j++) p[nxt][j] = 0;
		}
		printf("%.8lf\n", ret);
	}
}

*1:2^-30まで計算すればだいたい桁は足りるはず。

*2:10^5*30*8 = 2.4 * 10^7 byteで足りるはずなのですが