TopCoder SRM 614, Division 2, Level 3 : TorusSailingEasy

問題文

http://community.topcoder.com/stat?c=problem_statement&pm=13085&rd=15847

問題概要

　幅が N 、高さが M であって、上端と下端・左端と右端がそれぞれ繋がっている二次元のグリッド状の盤面を考える。最初、きつねのしえるが位置 ( 0, 0 ) にいる。以降、しえるは等確率で以下の二つの内いずれかの行動をとる。しえるの現在位置を ( m, n ) として、

• 位置 ( ( n + 1 ) mod N, ( m + 1 ) mod M ) に移動する
• 位置 ( ( n - 1 ) mod N, ( m - 1 ) mod M ) に移動する

位置 ( goalX, goalY ) に移動するまでに要する移動回数の期待値を求めよ。到達不可能な場合は -1 で示せ。

解法

　しえるの移動を全パターン列挙するのは無理なので、もっと効率的な方法を考える必要があります。まず、しえるの現在位置とそれまでの移動回数が等しいような移動の仕方が複数あったとき、それらを区別する必要はありません。従って、そのような状態は一つにまとめることにして、次のような DP を考えることができます。

• dp[ i 解移動した ][ X 座標が j ][ Y 座標が k ] := 該当する状態に至る確率

基底は dp[ 0 ][ 0 ][ 0 ] = 1 で、それ以外は全て 0 です。
　dp[ i ][ j ][ k ] からの遷移は題意より二通りで、

• dp[ i + 1 ][ ( j + 1 ) mod N ][ ( k + 1 ) mod M ] に dp[ i ][ j ][ k ] / 2 を加算
• dp[ i + 1 ][ ( j - 1 ) mod N ][ ( k - 1 ) mod M ] に dp[ i ][ j ][ k ] / 2 を加算

となります。ただし、( goalX, goalY ) に到達すれば手順が終了するので、そこからは更新しないようにします。
　必要な移動回数は無限に大きくなり得ますが、ある状態に至る確率は移動回数が増えるに従って極端に小さくなるので、ある程度の移動回数で打ち切ってしまってよいことになります。（確率が減らないとすれば、それはゴールに辿りつけていないので、いずれにしても無視してよいです）
　十分な移動回数分の DP をして表を埋めたあと、有効な i について i * dp[ i ][ goalX ][ goalY ] の総和を計算すると答えが得られます。この値が 0 であれば、到達不可能なので、代わりに -1 を返します。

コード

template < typename T = int > using VT = vector<T>;
template < typename T = int > using VVT = VT< VT<T> >;

template < typename T > inline T fromString( const string &s ){ T res; istringstream iss( s ); iss >> res; return res; };
template < typename T > inline string toString( const T &a ){ ostringstream oss; oss << a; return oss.str(); };

#define REP( i, m, n ) for ( int i = (int)( m ); i < (int)( n ); ++i )

const int ITERATE_NUM = 100000;

class TorusSailingEasy
{
public:
	double expectedTime( int N, int M, int goalX, int goalY )
	{
		VT<VVT<double>> dp( ITERATE_NUM, VVT<double>( N, VT<double>( M ) ) );
		dp[0][0][0] = 1;

		REP( i, 0, ITERATE_NUM - 1 )
		{
			REP( j, 0, N )
			{
				REP( k, 0, M )
				{
					if ( j == goalX && k == goalY )
					{
						continue;
					}
					dp[ i + 1 ][ modAdd( j, 1, N ) ][ modAdd( k, 1, M ) ] += dp[i][j][k] / 2;
					dp[ i + 1 ][ modAdd( j, -1, N ) ][ modAdd( k, -1, M ) ] += dp[i][j][k] / 2;
				}
			}
		}

		double res = 0;
		REP( i, 0, ITERATE_NUM )
		{
			res += i * dp[i][ goalX ][ goalY ];
		}
		return res <= 1e-8 ? -1 : res;
	}

	int modAdd( const int n, const int d, const int MOD )
	{
		return ( n + d + MOD ) % MOD;
	}
};