Problem213 ノミのサーカス

問題

30×30の格子状に並んだ正方形の中に900匹のノミがいる.
初期状態では, 1つの正方形に1匹のノミが入っている.
一回ベルが鳴るたびに, 各ノミはランダムに隣の正方形にジャンプする
(殆どの場所では4方向を選べる. ただし, 端や角にいるノミは除く.)

50回ベルが鳴ったときに, 空っぽの正方形の数の期待値はいくつだろうか?
10進小数点以下6桁に四捨五入して答えよ.

ソース

using System;
using System.Collections.Generic;
using System.Linq;

class Program
{
    const int UB = 30 - 1;

    static void Main()
    {
        var sw = System.Diagnostics.Stopwatch.StartNew();

        decimal[,] KitaitiArr = new decimal[UB + 1, UB + 1];
        for (int X = 0; X <= UB; X++) {
            for (int Y = 0; Y <= UB; Y++) {
                KitaitiArr[X, Y] = 1M;
            }
        }

        //ノミの座標ごとのループ
        for (int LoopStaX = 0; LoopStaX <= UB / 2; LoopStaX++) {
            for (int LoopStaY = 0; LoopStaY <= UB / 2; LoopStaY++) {
                decimal[,] KakurituBunpuArr = DeriveKakurituBunpuArr(LoopStaX, LoopStaY);

                //確率分布の値で、確率の乗法定理
                for (int BunpuX = 0; BunpuX <= UB; BunpuX++) {
                    for (int BunpuY = 0; BunpuY <= UB; BunpuY++) {
                        KitaitiArr[BunpuX, BunpuY] *= (1 - KakurituBunpuArr[BunpuX, BunpuY]);

                        //開始のY座標が対称なノミの分
                        KitaitiArr[BunpuX, BunpuY] *= (1 - KakurituBunpuArr[BunpuX, UB - BunpuY]);

                        //開始のX座標が対称なノミの分
                        KitaitiArr[BunpuX, BunpuY] *= (1 - KakurituBunpuArr[UB - BunpuX, BunpuY]);

                        //開始のX座標,Y座標が対称なノミの分
                        KitaitiArr[BunpuX, BunpuY] *= (1 - KakurituBunpuArr[UB - BunpuX, UB - BunpuY]);
                    }
                }
            }
        }
        decimal Answer = KitaitiArr.Cast<decimal>().Sum();
        Answer = Math.Round(Answer, 6, MidpointRounding.AwayFromZero);
        Console.WriteLine("Answer={0}。経過時間={1}", Answer, sw.Elapsed);
    }

    //ノミの開始座標を引数として、50回移動後の確率分布の配列を返す
    static decimal[,] DeriveKakurituBunpuArr(int pStaX, int pStaY)
    {
        decimal[,] PrevDP = new decimal[UB + 1, UB + 1];
        PrevDP[pStaX, pStaY] = 1M;

        for (int I = 1; I <= 50; I++) {
            decimal[,] CurrDP = new decimal[UB + 1, UB + 1];
            for (int X = 0; X <= UB; X++) {
                for (int Y = 0; Y <= UB; Y++) {
                    if (PrevDP[X, Y] == 0M) continue;

                    List<PointDef> wk4KinbouPointList =
                        Derive4KinbouPointList(X, Y);

                    int wkCnt = wk4KinbouPointList.Count;
                    foreach (PointDef EachPoint in wk4KinbouPointList) {
                        CurrDP[EachPoint.X, EachPoint.Y] += PrevDP[X, Y] / wkCnt;
                    }
                }
            }
            PrevDP = CurrDP;
        }
        //DebugPrint(PrevDP);
        return PrevDP;
    }

    struct PointDef
    {
        internal int X;
        internal int Y;
    }

    //座標を引数として、4近傍の座標のListを返す
    static List<PointDef> Derive4KinbouPointList(int pX, int pY)
    {
        var WillReturn = new List<PointDef>();

        Action<int, int> AddAct = (pNewX, pNewY) =>
        {
            if (pNewX < 0 || UB < pNewX) return;
            if (pNewY < 0 || UB < pNewY) return;
            WillReturn.Add(new PointDef() { X = pNewX, Y = pNewY });
        };
        AddAct(pX, pY - 1);
        AddAct(pX, pY + 1);
        AddAct(pX - 1, pY);
        AddAct(pX + 1, pY);

        return WillReturn;
    }

    //確率分布を表示
    static void DebugPrint(decimal[,] pDP)
    {
        var sb = new System.Text.StringBuilder();

        for (int Y = 0; Y <= UB; Y++) {
            for (int X = 0; X <= UB; X++) {
                sb.AppendFormat("{0} ", pDP[X, Y].ToString("0.0000"));
            }
            sb.AppendLine();
        }
        Console.WriteLine(sb.ToString());
    }
}

実行結果

Answer=330.721154。経過時間=00:00:03.2689177

解説

AB
CD

の4マスで、各ノミが2回移動後の、空の正方形の数の期待値を考えると、

Aのマスが空白である確率は、
  (1 - ノミAが2回移動後にAのマスにいる確率)
* (1 - ノミBが2回移動後にBのマスにいる確率)
* (1 - ノミCが2回移動後にCのマスにいる確率)
* (1 - ノミDが2回移動後にDのマスにいる確率)

ですので、
空の正方形の数の期待値 = Aのマスが空白である確率
                        +Bのマスが空白である確率
                        +Cのマスが空白である確率
                        +Dのマスが空白である確率
となります。