Topcoder Open 2016 Marathon Match(Round 2)勉強メモ

前回の続き。上位者のコードで勉強したメモ。

  • Psyhoさんのコード

http://pastebin.com/GQV4yrDj
UFOモードではモンテカルロシミュレーション、TSPモードではGreedyに経路を決めたあとに焼きなまししているようだ。TSP部分だけ読んでみる。→Greedy部分は(実行はしてるけど)使ってない。ローカル用か?

情報の持ち方(TSP関連のみ)

int PNO;                                // Starの数
PII PP[MAX_STARS];                      // Star
double PD[MAX_STARS][MAX_STARS+4];      // Star i と j の距離。UFOが通るところは0にする
int PN[MAX_STARS][MAX_STARS];           // Star i から j 番目に近い星
int PO[MAX_STARS][MAX_STARS];           // Star i からみて j が 何番目に近いか

double OPD[MAX_STARS][MAX_STARS+4];     // PDと同じだがUFO経路は無視
int PC[MAX_STARS][MAX_NEIGHBOURS];      // Star i から近い順の星、未到着のもののみ
PII PPOS[MAX_STARS];                    // 各星について、[船ID、航路Idx]のペア

int SNO;                                // Shipの数
int UNO;                                // UFOの数

int vs[MAX_STARS];                      // すでに訪れた星のフラグ
int vscopy[MAX_STARS];

int TPATH[MAX_STARS+10];                // ワーク

int PATH[MAX_SHIPS][MAX_STARS+10];      // 船ごとのPath
int PATHLEN[MAX_SHIPS];                 // 船ごとのPathの長さ

int BPATH[MAX_SHIPS][MAX_STARS+10];     // PATHと同じ。焼きなましのベストのやつ
int BPATHLEN[MAX_SHIPS];                // PATHLENと同じ。焼きなましのベスト

StarやUFOのクラスを作ったりせず、単純に配列の組み合わせで保存している。また、距離情報+近い星の情報は初期化時にすべて作ってしまっているようだ。

Greedyの部分

// pos: 船の位置
// planets: 未到着の星、シャッフル済
// 返り値は距離の合計
double greedy(VI pos, VI planets) {
    double rv = 0;
    REP(s, SNO) {
        PATHLEN[s] = 1;
        PATH[s][0] = pos[0];
    }
    
    while (planets.S) {
        double bv = 1e30;   // ベストな増加距離
        int bs = -1;    // 更新する船
        int bp = -1;    // 挿入するPathの位置
        int bi = -1;    // 挿入する星
        REP(s, SNO) FOR(p, 1, PATHLEN[s]+1) REP(i, planets.S) {
            double av = 0;
            av += PD[PATH[s][p-1]][planets[i]];
            if (p < PATHLEN[s]) av += PD[PATH[s][p]][planets[i]];
            if (av < bv) {
                bv = av;
                bs = s;
                bp = p;
                bi = i;
            }
        }
        PATHLEN[bs]++;
        for (int i = PATHLEN[bs]; i > bp; i--) PATH[bs][i] = PATH[bs][i-1];
        PATH[bs][bp] = planets[bi];
        rv += bv;
        planets[bi] = planets.back();
        planets.pop_back();
    }
    return rv;
}

距離増加が最小になるように、星を船の航路に挿入している。

焼きなましの部分

while (true) {
    if ((steps & 1023) == 0) {
        if (notVisited.S == 1) break;
        timePassed = (getTime() - startTime) / (TIME_LIMIT - totalTime);
        temp = (1.0 - pow(timePassed, 0.15)) * 160;
        if (timePassed >= 1.0) break;
    }
    steps++;
    
    // s0: ランダムに選んだ船、p0: s0のPathのIdx(ランダム)
    int s0 = rng.next(SNO);
    if (PATHLEN[s0] == 1) continue;
    int p0 = rng.next(1, PATHLEN[s0]);

    // s1: 船、p1: s1のPathのIdx
    int p1, s1;
    int mode = rng.rand() & 1;
    if (mode == 0) {    // 半分の確率でs0/p0の近くの星から選ぶ
        PII p = PPOS[PC[PATH[s0][p0]][rng.next(NEIGHBOURS)]];
        s1 = p.X;
        p1 = p.Y;
    } else {
        p1 = 0;         // 半分の確率で適当な船を選ぶ。p1はあとで
        s1 = rng.next(SNO);
    }
    int p0b, p1b;
    int type;
    
    //ちなみにp0、p1はかならず1以上。Path[0]は現在地なので無視。->Path操作がラクになる

    double diff = 0;
    if (s0 == s1) {
        if (p1 == 0) p1 = rng.next(1, PATHLEN[s0]);
        type = rng.next(2) == 0;
        if (type == 0) {    // p0とp1を交換
            if (p0 == p1) continue;
            if (p0 > p1) swap(p0, p1);
            diff -= pathDist(s0, p0-1, PATH[s0][p0]);
            diff -= pathDist(s0, p1+1, PATH[s0][p1]);
            diff += pathDist(s0, p0-1, PATH[s0][p1]);
            diff += pathDist(s0, p1+1, PATH[s0][p0]);
        } else {            // p0をp1の前に挿入
            if (abs(p0-p1)<=2) continue;
            diff -= pathDist(s0, p0-1, PATH[s0][p0]);
            diff -= pathDist(s0, p0+1, PATH[s0][p0]);
            diff += pathDist(s0, p0+1, PATH[s0][p0-1]);
            diff -= pathDist(s0, p1, PATH[s0][p1-1]);
            diff += pathDist(s0, p1-1, PATH[s0][p0]);
            diff += pathDist(s0, p1, PATH[s0][p0]);
        }
    } else {
        type = rng.next(2) == 0;
        if (type == 0) {    // 半分の確率で(p0b, p0) と (p1b, p1) を入れ替える
            // { p0b -> a -> b -> p0 }, { p1b -> x -> y -> p1 }  これを
            // { p0b -> x -> y -> p0 }, { p1b -> a -> b -> p1 }  こうする
            if (PATHLEN[s1] == 1) continue;
            if (p1 == 0) p1 = rng.next(1, PATHLEN[s1]);
            p0b = rng.next(1, PATHLEN[s0]);
            p1b = rng.next(1, PATHLEN[s1]);
            if (p0 > p0b) swap(p0, p0b);
            if (p1 > p1b) swap(p1, p1b);
            diff -= pathDist(s0, p0-1, PATH[s0][p0]);
            diff -= pathDist(s0, p0b+1, PATH[s0][p0b]);
            diff += pathDist(s0, p0-1, PATH[s1][p1]);
            diff += pathDist(s0, p0b+1, PATH[s1][p1b]);
            diff -= pathDist(s1, p1-1, PATH[s1][p1]);
            diff -= pathDist(s1, p1b+1, PATH[s1][p1b]);
            diff += pathDist(s1, p1-1, PATH[s0][p0]);
            diff += pathDist(s1, p1b+1, PATH[s0][p0b]);
        } else {    // 半分の確率で、2つの航路を途中で切り替えている
            p1 = p1 == 0 ? rng.next(1, PATHLEN[s1]+1) : p1 + 1;
            if (rng.next(2)) {
                // ここはよくわからない・・・未使用の星ができてしまうような・・・
                p0b = rng.next(1, PATHLEN[s0]);
                if (p0 > p0b) swap(p0, p0b);
            } else {
                // { a -> b -> p0 -> c -> d }, { p -> q -> p1 -> r -> s }  これを
                // { a -> b -> p0 -> r -> s }, { p -> q -> p1 -> c -> d }  こうする
                p0b = p0;
            }
            p1b = p1-1;
            
            diff -= pathDist(s0, p0-1, PATH[s0][p0]);
            diff -= pathDist(s0, p0b+1, PATH[s0][p0b]);
            diff += pathDist(s0, p0b+1, PATH[s0][p0-1]);
            diff -= pathDist(s1, p1, PATH[s1][p1-1]);
            diff += pathDist(s1, p1-1, PATH[s0][p0]);
            diff += pathDist(s1, p1, PATH[s0][p0b]);
        }
    }
    
    if (diff <= rng.nextDouble() * temp) {
        if (s0 == s1) {
            if (type == 0) {    // p0とp1を交換ということは、[p0, p1]を反転させるということ
                reverse(PATH[s0]+p0, PATH[s0]+p1+1);
                updatePPOS(s0, p0, p1);
            } else {
                int tmp = PATH[s0][p0];
                if (p0 < p1) {
                    FOR(i, p0, p1-1) PATH[s0][i] = PATH[s0][i+1];
                    PATH[s0][p1-1] = tmp;
                } else {
                    for (int i = p0; i > p1; i--) PATH[s0][i] = PATH[s0][i-1];
                    PATH[s0][p1] = tmp;
                }
                updatePPOS(s0, p0, p1);
            }
        } else {
            if (p0b-p0>p1b-p1) {
                swap(s0, s1);
                swap(p0, p1);
                swap(p0b, p1b);
            }
            int move = (p1b-p1)-(p0b-p0);
            assert(move >= 0);
            FOR(i, p0, p0b+1) TPATH[i] = PATH[s0][i];
            PATHLEN[s0] += move;
            PATHLEN[s1] -= move;
            assert(PATHLEN[s0] <= turnsLeft+1);
            assert(PATHLEN[s1] <= turnsLeft+1);
            for (int i = PATHLEN[s0]-1; i >= p0b+1+move; i--) PATH[s0][i] = PATH[s0][i-move]; //TODO: fix
            FOR(i, p1, p1b+1) PATH[s0][i-p1+p0] = PATH[s1][i];
            FOR(i, p0, p0b+1) PATH[s1][i-p0+p1] = TPATH[i];
            FOR(i, p1b+1-move, PATHLEN[s1]) PATH[s1][i] = PATH[s1][i+move];
            updatePPOS(s0, p0, PATHLEN[s0]-1);
            updatePPOS(s1, 1, PATHLEN[s1]-1);
        }
        PATH[s0][PATHLEN[s0]] = PNO;
        PATH[s1][PATHLEN[s1]] = PNO;
        bv += diff;
        acc++;
        // ベスト経路を取っておく
        if (bv < xv) {
            xv = bv;
            bestacc++;
            REP(s, SNO) {
                BPATHLEN[s] = PATHLEN[s];
                REP(i, PATHLEN[s]) BPATH[s][i] = PATH[s][i];
            }
        }
    }
}

コスト差分は予め計算しておいた距離を使ってやっている。経路は静的配列で持ち、ループで更新している。

Seed=2(船が1槽でUFOなし)で挙動を計測してみた。

f:id:yambe2002:20160606124734p:plain
温度と更新頻度のプロット。後半はほとんど更新されなくなっている。

f:id:yambe2002:20160606124741p:plain
およびスコアのプロット。前半に一気に改善して、だんだんなだらかになっている。スコアが悪化する場合もあるはずだが、プロット数をかなり減らしているので見えていない(データを取ってるのが0xFFF回に一回だけなので)。


航路変更によるスコア(距離)差分が、平均190(最小-130、最大1100、標準偏差180)くらい。

if (diff <= rng.nextDouble() * temp) {

遷移条件がこれ。コスト差分量によって遷移確率が異なる。
ループ総数3800万回のうち、遷移したのが130万回ほど。悪化したのに遷移したのは70万回ほどだった。