SRM 513 500pt: PerfectMemory

問題概要

N*M(N,M<50, N*Mは偶数)のボードがあり、各マスにはカードが裏向きにおいてある。カードは何か印がついていて、完全マッチングになっている。完全な記憶力と戦略を持って神経衰弱するとき、全部回収するのに必用なターン数の期待値を求める問題。

考えたこと

  • これN*M以外の値関係あるのか?長方形に配置とか問題設定上意味ない様な気がするんだけど…。
    • N*Mがあまり大きな素数を含まないとかが効いてくるのかなあ?
  • 細かいことを気にせず普通に解き方考えてみると、状態数少ないし何枚オープンになっているかと残り枚数組にしてDPしたら上手く行くだろう。
    • 神経衰弱で完全な記憶力をもっている、というのは一度オープンにしたカードを再び裏返しにしない、と考えるとイメージしやすい。
  • 何かそれだと簡単過ぎるように見えるけど、でもそれで行けるはずだよなあ。
  • 最適な戦略は、オープンになってないカードをまずひっくり返して、ペアが取れるなら取ればいい、とか。というか当たり前に思える戦略を取れば良さげ。
  • 更新順序が微妙に面倒くさいし、DPじゃなくてメモ化再帰で書くか。
    • 一応メモリはギリギリ足りてるはず。
    • というかこれって状態ループしないよね?一応確認しとこう。
  • 避けようのない場合分けが出たので、コメントつけて注意しておく。
  • 実装したら-nanとかが返される。これは間違いなく0割。
  • メモ化の終了条件とか、遷移確率とかいじったら(訂正したら)サンプル通った。サーバ上で最大ケースが間に合うことを確認してからサブミット。
  • include無しで書けたのでちょっと嬉しい。

本番で通ったコード

bool visited[2509][2509];
double memo[2509][2509];

struct PerfectMemory {

	double getExpectation(int N, int M) {
		return rec(N*M/2,0);
	}

	double rec(int N, int M){ //N: 残りペア、 M: 片方がオープン
		if(visited[N][M]) return memo[N][M];

		visited[N][M] = true;

		if(N==0) return memo[N][M] = 0.0;
		if(M>N) return memo[N][M] = 0.0;
		if(M==0){
			double pp = (1.0/(2*N-1));
			return memo[N][M] = pp*(1 + rec(N-1,0)) + (1-pp)*(1 + rec(N,2));
		}
		if(N==M){
			return memo[N][M] = N;
		}

		double ret = 0, p = (double)M/(2*N-M); //pは1枚目で当たりを引く確率

		//1枚目が既にオープン
		ret += p*(1 + rec(N-1,M-1));

		if(2*N-(M+1)>0){
			//1枚目と2枚目が同じ
			double r = (1.0/(2*N-(M+1)));
			ret += (1-p) * r * (1 + rec(N-1,M));

			double q = (double)M/(2*N-M-1); //qは2枚目で当たりを引く確率
			//1枚目と2枚目が違って、2枚目も外れ
			ret += (1-p) * (1-q-r) * (1 + rec(N, M+2));

			//1枚目と2枚目が違って、2枚目が当たり
			ret += (1-p) * q * (2 + rec(N-1, M));
		}
		else{
			ret += (1-p)*(1 + rec(N,M+1));
		}

		return memo[N][M] = ret;
	}

};

後から思ったこと

  • メモ化再帰でスタックオーバーフローするかどうか検証してなかった。結果的にそんなに深くならないような再帰の書き方をしていたので大丈夫だったけど、再帰関数の呼び出し順によっては深さ1250*1250/2くらいになってしまうこともあった。これは通ったのは偶然と言わざるを得ない(実際はテストで気づくだろうが)。
    • と、思ったら呼び出しの順序を一番深くなるようにしても大丈夫だったっぽい。たまたま運が良かっただけではあるけど。
  • ペアの数をキーにしてたんだからメモリは2500じゃなくて2500/2だけでよかった。
  • 結局N,Mの形式で引数が与えられるのは意味がなかった。TopCoderでこの手の意味のない引数の与え方されるのは珍しい気がする。しかも「実は関係ありません」じゃなくて「どう見ても関係ありません」だし。

バグを減らす工夫

  • 実は本番中に提出した上のコードは、一番最後のelseは通らない。というか再帰呼び出しの意味が我ながらよく分からない。
    • こういうことが起きるのは、目的が正しいコードを書くことではなくとにかくサンプルを通すことになっているせいな気がする。
  • しかも終了条件も何か場合分け多すぎで意味が分かっているのか怪しい。というわけで、終了条件の部分にもコメントをつけておくべきだろう。
  • 状態の1つ目のキーを、「未回収のペア数」にするか「未回収の枚数」のどちらにするかは好み次第かな。
    • 一応、「未回収のペア数」にした方がメモリは節約できる。この問題じゃどちらでも問題ないから分かりやすさ優先で良いけど。
  • 状態の2つ目のキーは、「オープンになっているカードの枚数」にするのが流石に自然な気がする。
    • こうすると、DPで書いたとき舐める順番が降順になるのが気になるかもだけど。
  • 絶対通らない状態に入ったときはassertするのがいいのかなあ?
    • でもそれだと、その状態への遷移確率が0みたいな感じで場合分けを減らせるからやっぱり即座に0を返す方が良いのかもしれない。

以上を考慮した上で書き直したのがこちら。

	double rec(int N, int M){ //N: 残りペア、 M: 片方がオープン
		if(M>N || M<0 || N<0) return 0.0;

		//メモ化処理
		if(visited[N][M]) return memo[N][M];
		visited[N][M] = true;

		//全部とり終わった場合
		if(N==0) return memo[N][M] = 0.0;

		//この場合分けがないと、一番最後の再帰呼び出しでおかしいことが起きてしまう
		if(N==M){
			return memo[N][M] = N;
		}

		double ret = 0, p = (double)M/(2*N-M); //pは1枚目で当たりを引く確率

		//1枚目が既にオープン
		ret += p*(1 + rec(N-1,M-1));

		//1枚目と2枚目が同じ
		double r = (1.0/(2*N-(M+1))); //rは1枚目で引いたもののペアを引く確率
		ret += (1-p) * r * (1 + rec(N-1,M));

		//1枚目と2枚目が違って、2枚目が当たり(次ターンで回収)
		double q = (double)M/(2*N-M-1); //qは2枚目で当たりを引く確率
		ret += (1-p) * q * (2 + rec(N-1, M));

		//1枚目と2枚目が違って、2枚目も外れ
		ret += (1-p) * (1-q-r) * (1 + rec(N, M+2));

		return memo[N][M] = ret;
	}
  • でもこれでもちょっと気持ち悪い。何が気持ち悪いって、if(N==M)の部分がわざわざ特別扱いしているのがとてもイヤだ。
    • 自分の実力を考えると、この場合分けは実行してRE起こすまで絶対に気づけない。
  • こんな事をしなくちゃいけないのは、遷移確率0.0なのに遷移先の状態を計算しに行ってしまうからだ。
  • なので、遷移確率0.0のときは遷移先の状態を調べないことにする。
const double EPS = 1e-9;

	double rec(int N, int M){ //N: 残りペア、 M: 片方がオープン
		if(M>N || M<0 || N<0) return 0.0;

		//メモ化処理
		if(visited[N][M]) return memo[N][M];
		visited[N][M] = true;

		//全部とり終わった場合
		if(N==0) return memo[N][M] = 0.0;

		double ret = 0, p = (double)M/(2*N-M), coef; //pは1枚目で当たりを引く確率

		//1枚目が既にオープン
		if((coef = p) > EPS){
			ret += coef * (1 + rec(N-1,M-1));
		}

		//1枚目と2枚目が同じ
		double r = (1.0/(2*N-(M+1))); //rは1枚目で引いたもののペアを引く確率
		if((coef = (1-p) * r) > EPS){
			ret += coef * (1 + rec(N-1,M));
		}

		//1枚目と2枚目が違って、2枚目が当たり(次ターンで回収)
		double q = (double)M/(2*N-M-1); //qは2枚目で当たりを引く確率
		if((coef = (1-p) * q) > EPS){
			ret += coef * (2 + rec(N-1, M));
		}

		//1枚目と2枚目が違って、2枚目も外れ
		if((coef = (1-p) * (1-q-r)) > EPS){
			ret += coef * (1 + rec(N, M+2));
		}

		return memo[N][M] = ret;
	}
  • こう書くと全ての場合分けがちゃんと意味を持ってる感じで気分が良い。
    • 実際は1行目のif文はなくても通る。でもつけといた方が精神安定上よろしい。
DP解法

まずは普通に。今回は状態の1つ目のキーを「未回収の枚数」にしてみる(特に理由はない)。

	double solve(int num_cards){
		for(int i=2; i<=num_cards; i+=2){
			for(int j=i/2; j>=0; j--){ //回収できずにオープンなカード枚数だけが増えることがあるので、降順にする
				double p = (double)j/(i-j), coef;

				//1枚目が既にオープン
				if((coef = p) > EPS){
					dp[i][j] += coef * (1.0 + dp[i-2][j-1]);
				}

				//1枚目と2枚目が同じ
				double r = 1.0/(i-(j+1)); //rは1枚目で引いたもののペアを引く確率
				if((coef = (1-p) * r) > EPS){
					dp[i][j] += coef * (1.0 + dp[i-2][j]);
				}

				//1枚目と2枚目が違って、2枚目が当たり(次ターンで回収)
				double q = (double)j/(i-(j+1)); //qは2枚目で当たりを引く確率
				if((coef = (1-p) * q) > EPS){
					dp[i][j] += coef * (2.0 + dp[i-2][j]);
				}

				//1枚目と2枚目が違って、2枚目も外れ
				if((coef = (1-p) * (1-q-r)) > EPS){
					dp[i][j] += coef * (1.0 + dp[i][j+2]);
				}
			}
		}
		return dp[num_cards][0];
	}
メモリ節約判

空間計算量O((N*M)^2)からO(N*M)へ。時間計算量はO((N*M)^2)で変わらず。

#include <vector>
using namespace std;

	double solve(int num_cards){
		vector<double> next(num_cards+1, 0.0);
		for(int i=2; i<=num_cards; i+=2){
			vector<double> old = next;
			next.assign(num_cards+1, 0.0);
			for(int j=i/2; j>=0; j--){ //回収できずにオープンなカード枚数だけが増えることがあるので、降順にする
				double p = (double)j/(i-j), coef;

				//1枚目が既にオープン
				if((coef = p) > EPS){
					next[j] += coef * (1.0 + old[j-1]);
				}

				//1枚目と2枚目が同じ
				double r = 1.0/(i-(j+1)); //rは1枚目で引いたもののペアを引く確率
				if((coef = (1-p) * r) > EPS){
					next[j] += coef * (1.0 + old[j]);
				}

				//1枚目と2枚目が違って、2枚目が当たり(次ターンで回収)
				double q = (double)j/(i-(j+1)); //qは2枚目で当たりを引く確率
				if((coef = (1-p) * q) > EPS){
					next[j] += coef * (2.0 + old[j]);
				}

				//1枚目と2枚目が違って、2枚目も外れ
				if((coef = (1-p) * (1-q-r)) > EPS){
					next[j] += coef * (1.0 + next[j+2]);
				}
			}
		}
		return next[0];
	}

おまけ

  • メモ化をDPに書き直したとき、何故自分がDPよりメモ化の方が書きやすいのか分かった気がする。
    • 1つは、ネストが深くならないこと。
    • 2つめは、DPでループ回すときにいつも変数名をi,jとしてしまって、意味のある変数名になっていない。
    • 後者は悪癖なので治したい。
  • 配列使いまわすタイプのDPは、一度普通のDP書いてから後で書き直すとやりやすい。最初から配列2本で書くより安心感があった。