Topcoder Open 2016 Marathon Match(Round 3)の参加日記

制約が多くて面倒な問題だった。そのせいか、レート中位~下位の参加者がずいぶん少なかった。
結果は暫定62位/101人の定位置。

問題概要

  • 長さ1の正方形セルSxS個で構成されるマップがある
  • 各セルはタイプ(0~9)が設定されている
  • マップ上にはN個のアイテムとN個のターゲットがある
  • アイテムの近く(距離10^-3以内)に止まると、該当アイテムを拾う
  • ターゲットの近くで止まると、持っているアイテムを落とす
  • 1つのターゲットに落とすのは1つのアイテム
  • 持てるアイテムの最大値はcapacityで与えられる
  • すべてのアイテムをターゲットに落とさなければならない
  • スタート地点はマップの外周から距離10^-3以内ならどこでもよい
  • このとき、なるべく小さいコストで全アイテムをターゲットに落とすのが目的
  • コストは、セルタイプx移動距離。タイプの違うセルをまたぐときは、さらにタイプ差の二乗が可算される
  • 1回のステップでは、1つ以内のセルしかまたげない
  • セルボーダーから距離10^-3以内には止まってはならない。
  • ステップ間は10^-3以上の距離をおくこと
  • S:10~50、タイプ数2~10、アイテム数:5~S*S/10、capacity:1~10

https://community.topcoder.com/longcontest/?module=ViewProblemStatement&rd=16704&pm=14283


方針

  1. Dijkstraでアイテム・ターゲット・エッジセル間の距離を求める
  2. Dijkstraの結果を使い、アイテムとターゲットでTSP(ビームサーチ)
  3. TSP結果から有効解を作成(経路復元)
  4. 有効解を改善(貪欲+焼きなまし)


上位との比較

  • 全体的な方針は一緒
  • TS部分で大きな差がついている
  • TSは焼きなましの人と、2optなどTSP特有のアルゴリズムの人がいる
  • Dijkstraも、単純なセル間ではなくセルをさらに分割してる人もいた


反省
TSP部分で大きなスコア差が付いている。TSPがイマイチなのは気づいていたが、それより4の改善部分に時間を使ってしまった。
目についた部分を適当に改善するのではなく、どこに手を付けるべきかまず分析するべきだった。

Topcoder Marathon Match Round90 勉強メモ

2位のEvbCFfp1XBさんのコードを読んでみる。
https://community.topcoder.com/longcontest/?module=ViewProblemSolution&pm=14094&rd=16495&cr=23100980&subnum=3

問題文はこちら。
https://community.topcoder.com/longcontest/?module=ViewProblemStatement&rd=16495&pm=14094


方針(Forumより)

  1. ランダムにポイント0ポイント1未満のターゲットを選ぶ ※コードより修正
  2. ターゲットにポイント1が入るか、ループ回数が10になるまで以下を繰り返す
  • ターゲット近くの同色ボールを選ぶ
  • そのボールの近くのボールを1~10個選ぶ ※上で選んだボール含む
  • 選んだボールでビームサーチ(Max Depth=20、幅250)


情報の持ち方

class Board {
    public static final int WALL = 1 << 27;
    public static final int EMPTY = 1 << 28;
    public static final int[] DR = { 0, 1, 0, -1, };
    public static final int[] DC = { -1, 0, 1, 0, };
    public static final int[] DELTA = { -1, 64, 1, -64, };
    public static final int SHIFT = 6;
    public static final int MASK = (1 << SHIFT) - 1;
    public static final int[] ballIndex = new int[64 * 64];
    public static final long[][] hashes = new long[64 * 64][10];
    public int[][] distanceFromTarget;
    public static final int[] targetIndex = new int[64 * 64];
    public List<ColoredPoint> balls = new ArrayList<>();
    public List<ColoredPoint> targets = new ArrayList<>();
    public int R;
    public int C;
    public long hash = 0;
    int countPoint1;
    private ArrayList<State> moveHistory = new ArrayList<>();  // ボールを動かした履歴

    /// 省略 ///

    // ボールをある方向へ動かす
    public void next(int ballIndex, int direction) {
        ColoredPoint ball = balls.get(ballIndex);

        State current = new State(null, ballIndex, ball.z, direction, 0);
        moveHistory.add(current);
        {
            int index1 = targetIndex[ball.z];
            if (index1 != EMPTY && getTargetColor(index1) == ball.color) {
                countPoint1--;
            }
        }
        hash ^= hashes[ball.z][ball.color];

        // ボールを転がす。この更新のやりかたは参考になる  DELTA = { -1, 64, 1, -64, };
        int delta = DELTA[direction];
        int z = ball.z;
        this.ballIndex[z] = EMPTY;
        do {
            z += delta;
        } while (this.ballIndex[z] == EMPTY);
        z -= delta;
        ball.z = z;
        this.ballIndex[z] = ballIndex;

        hash ^= hashes[ball.z][ball.color];
        {
            int index1 = targetIndex[ball.z];
            if (index1 != EMPTY && getTargetColor(index1) == ball.color) {
                countPoint1++;
            }
        }
    }

    /// 省略 ///
}

class State {
    State parent;
    public int index;
    public int z;
    public int direction;
    public double score;
    public int numMoves;

    /// 省略 ///
}

class State2 implements Comparable<State2> {
    State2 parent;
    public int index;
    public int z;
    public int direction;
    public double score;
    public double distance;  // ビームサーチ用。ターゲットまでの(最短)距離合計
    public int numMoves;

    /// 省略 ///
}

巨大なboardクラスを用意している。履歴はStateクラスで持つ。State2はビームサーチ用。


メインルーチン

private void solve() {
    for (; timeManager.getSecond() < 9.5;) {

        if (board.getScorePoint1() > 0.999) {
            break;
        }

        // スコアが1未満のターゲット
        ArrayList<Integer> targetsScore0 = getTargetsScore0();

        // ランダムに1つ選ぶ
        IntArray targetIndexes = new IntArray();
        int targetIndex = targetsScore0.get((int) (targetsScore0.size() * rng.nextDouble()));
        targetIndexes.add(targetIndex);

        // 近隣ボールの選ぶ数。ランダムに試す
        int[] ns = new int[] { 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, };
        Utils.shuffle(ns, rng);
        for (int maxN : ns) {
            if (board.getTargetScore(targetIndex) > 0.999) {
                break;
            }
            if (timeManager.getSecond() > 9.5) {
                break;
            }

            // ターゲットから最寄りの同色ボール
            int nearestSameColorBallIndex = getSameColorScore0BallsDistanceFromTargetAscending(targetIndex).get(0).second;

            ArrayList<Pair<Double, Integer>> balls = getBallsDistanceFromBallAscending(nearestSameColorBallIndex);

            // 同色ボールの近くのボールmaxN個
            IntArray ballIndexes = new IntArray();
            int n = Math.min(maxN, balls.size());
            for (int i = 0; i < n; i++) {
                ballIndexes.add(balls.get(i).second);
            }

            int currentNumMoves = board.numMoves();
            int maxStep = 20;
            int maxBeamWidth = (int) (5 * 1000 / maxStep);
            countBeamsearch++;

            // ビームサーチ
            ArrayList<State> moves = beamsearch(ballIndexes, maxStep, maxBeamWidth, targetIndexes);

            if (moves == null) {
                continue;
            }
            for (State state : moves) {
                assert board.canMove(state.index, state.direction);
                board.next(state.index, state.direction);
            }

            // ベストスコアなら保持
            if (board.getScorePoint1() > bestScore && board.numMoves() < board.maxNumRolls()) {
                bestScore = board.getScorePoint1();
                bestSolution = board.getSolution();
            } else {
                // ベストでなく、ターン数に余裕がないようなら棄却
                if ((double) board.numMoves() / board.maxNumRolls() > 0.99 * timeManager.getSecond() / 9.5) {
                    for (; board.numMoves() > currentNumMoves;) {
                        board.previous();
                    }
                }
            }
        }
    }
}

単純にN個のボールを選ぶのではなく、1個~10個と少ないほうから試している。後の参考になりそう。


ビームサーチの部分

// 動かすボール、最大深さ、ビーム幅、対象ターゲット(1つだけ)
private ArrayList<State> beamsearch(IntArray ballIndexes, int maxStep, int maxBeamWidth, IntArray targetIndexes) {
    ArrayList<State2> currents = new ArrayList<>();
    ArrayList<State2> nexts = new ArrayList<>();
    ArrayList<State2> all = new ArrayList<>();  // 探索にでてきた全ノード。最後に一番いいやつを選ぶ

    int currentNumMoves = board.numMoves();

    used.clear();
    used.add(board.hash);

    // 初手の集合currentsを得る
    for (int i = 0; i < ballIndexes.length; i++) {
        int ballIndex = ballIndexes.values[i];
        ColoredPoint ball = board.balls.get(ballIndex);
        int z = ball.z; // zはボール位置のIntパック
        for (int direction = 0; direction < 4; direction++) {
            if (!board.canMove(ballIndex, direction)) {
                continue;
            }

            board.next(ballIndex, direction);
            if (used.add(board.hash)) {
                double score = getSumScore(ballIndexes);
                double distance = getSumMinDistance(ballIndexes, targetIndexes);
                distance += 1e-2 * rng.nextDouble();
                currents.add(new State2(null, ballIndex, z, direction, score, distance, 1));
            }
            board.previous();
        }
    }


    for (int step = 0; step < maxStep; step++) {
        if (currents.size() == 0) {
            break;
        }
        int beamWidth = Math.min(maxBeamWidth, currents.size());            
        Utils.select(currents, 0, currents.size(), beamWidth);  // currentsをスコア>距離でソートか?最初のビーム幅分のみ抽出?
        for (int beam = 0; beam < beamWidth; beam++) {
            State2 currentState = currents.get(beam);

            all.add(currentState);

            {
                ArrayList<State2> currentStateList = reverse2(toList(currentState));    // 初手からcurrentStateまでの手のリスト
                ArrayList<State> res = new ArrayList<>();
                for (State2 state : currentStateList) {
                    res.add(new State(null, state.index, state.z, state.direction, state.score));
                }

                board.set(res, currentNumMoves);   // boardをcurrentStateにする
            }

            {
                assert getSumScore(ballIndexes) == currentState.score;
                if (currentState.score > ballIndexes.length - 1e-9) {
                    break;
                }
            }
            for (int i = 0; i < ballIndexes.length; i++) {
                int ballIndex = ballIndexes.values[i];
                ColoredPoint ball = board.balls.get(ballIndex);
                int z = ball.z;
                for (int direction = 0; direction < 4; direction++) {
                    if (!board.canMove(ballIndex, direction)) {
                        continue;
                    }

                    int nextZ = board.simulateNext(ball, direction);    // 該当の方向へ動かしてみたときのz

                    long currentHash = board.hash;

                    boolean add = used.add(currentHash ^ (board.hashes[z][ball.color] ^ board.hashes[nextZ][ball.color]));

                    // 使ってないHashになるなら
                    if (add) {
                        // scoreとdistanceを計算
                        double nextScore = currentState.score;
                        {
                            int index1 = board.targetIndex[z];
                            if (index1 != Board.EMPTY && board.getTargetColor(index1) == ball.color) {
                                nextScore--;
                            }
                        }
                        {
                            int index1 = board.targetIndex[nextZ];
                            if (index1 != Board.EMPTY && board.getTargetColor(index1) == ball.color) {
                                nextScore++;
                            }
                        }
                        double score = nextScore;

                        double nextDistance = currentState.distance;
                        int targetIndex = targetIndexes.values[0];
                        if (board.distanceFromTarget[targetIndex][z] == 0 && board.getTargetColor(targetIndex) != ball.color) {
                            nextDistance -= 1e9;
                        } else {
                            nextDistance -= board.distanceFromTarget[targetIndex][z];
                        }
                        if (board.distanceFromTarget[targetIndex][nextZ] == 0 && board.getTargetColor(targetIndex) != ball.color) {
                            nextDistance += 1e9;
                        } else {
                            nextDistance += board.distanceFromTarget[targetIndex][nextZ];
                        }
                        double distance = nextDistance;

                        distance += 1e-2 * rng.nextDouble();    // これは・・・?
                                                                // ちょっとバラけさせると探索が良くなるのか?

                        // 次手の集合nextsに加える
                        nexts.add(new State2(currentState, ballIndex, z, direction, score, distance, currentState.numMoves + 1));
                    }
                }
            }
        }
        {
            ArrayList<State2> swap = currents;
            currents = nexts;
            nexts = swap;
        }
        nexts.clear();

    }
    for (; board.numMoves() > currentNumMoves;) {
        board.previous();
    }

    if (all.size() == 0) {
        return null;
    }

    Collections.sort(all);
    ArrayList<State2> allList = reverse2(toList(all.get(0)));
    ArrayList<State> res = new ArrayList<>();
    for (State2 state : allList) {
        res.add(new State(null, state.index, state.z, state.direction, state.score));
    }
    return res;

}

ループで幅優先ビームサーチをしている。各NodeはState2クラスで保持していて、これが親へのポインタを持っているようだ。


この人のコードはかなり読みやすかった。
Javaだからか、それとも書き方が自分に似ているからだろうか。いずれにしても、とても参考になる。

かなり強い方なので、他のマッチのコードも読んでみても良さそうだ。

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

4位のnikaさんのコードを読んでみる。
https://community.topcoder.com/longcontest/?module=ViewProblemSolution&pm=10729&rd=16702&cr=20315020&subnum=18

Forumのご本人のポストによると、次の3つのフェーズに分かれているらしい。

パート1(30%)
rootが分割されるまでカットを足していく。足すカットは、いくつかの候補から、(コスト)/(新たに分割されるrootの数)がベストのものを選ぶ。

パート2(60%)
カットの増減を繰り返してスコアを改善する。カットをランダムに選んで取り除き、そしてrootがすべて分割されてないときは、パート1と同じ手段でカットを足す。スコアが改善するようならこれを使い、しないなら元に戻す。

パート3(10%)
カット位置を1ポイントずつ動かして改善をはかる。

カットを足す方法について。なるべく良い候補にするため①始点と終点をランダムに選ぶ②垂直・水平のどちらかを選ぶ③どちらかの方向の全パターンを見て、最も良いものを選ぶ、という処理を行っている。ここでコストの再計算を避けるため、スキャンラインテクニック(Scan Line Technique)というのを使っているらしい。

情報の持ち方

int m, n;               // m: Plantの数、n: Endpointの数
int x[N], y[N]          // x[], y[]: Endpointの座標
int parent[N];          // EndpointのPlant番号
int p[N];               // Endpointの親Endpoint番号
double edgeLength[N];   // Endpointとその親の距離
double edgeSum;         // 距離の総合計
vector<int> plants[M];  // Plantに属するEndpoint
vector<int> hull[M];    // Plantに属するEndpoint。外殻のみ。左下から時計回り
                        // CutがPlantを切断するかどうかを判断するのに使う
vector<int> ret;

int separated[M][M];    // plantが分離されてるかどうか(いくつのカットで分離されたか)
int phase;

int minX[N], maxX[N];
int dead[N];            // 根が死んでいるかどうか
int wouldKill[N];
int alive[N];           // 生きている根
int aliven;             // 生きている音根の数
int wk;
double aliveLength[N];  // 生きている根の長さ
multiset<double> cutRatio[N];   // 根が切られた位置。multiset: sorted、重複OKのset

typedef pair<Line, Mask> Cut;   // Lineは始点xy・終点xyなどをもつ構造体
                                // Maskは120ビットのマスク。ビット位置ごとのXor、Orなどができる

double bestSaved = -1e100;      // 切り取られた根の長さ合計ベスト
vector<Cut> bestSolution;
double savedEdgeSum;
vector<Cut> solution;

// グループ: カットされてないPlantの集合
int group_n;            // グループの総数
int group_cnt[M]        // あるグループ番号に属するPlant数
int group[M][M];        // group[i][j]: グループ番号iに属する、j番目のPlant番号
int group_vis[M];       // ワーク

やはりClassを作るのは避けて、なるべく単純な配列等で保持している。
マスクは該当Cutが分割しているPlantの情報。Cut線からLeft集合とRight集合に分けている(ビット位置で保持)。この方法は汎用性がありそうだ。

メインルーチン

vector<int> makeCuts(int NP, vector<int> points, vector<int> roots) {
    Init();

    m = NP;
    n = SZ(points) / 2;

    for (int i = 0; i < n; i++) {
        x[i] = points[2 * i];
        y[i] = points[2 * i + 1];
    }
    edgeSum = 0;

    for (int i = 0; i < n; i++) {
        if (i < m) {
            parent[i] = i;
            p[i] = -1;
        }
        else {
            p[i] = roots[2 * (i - m)];
            parent[i] = parent[p[i]];
            edgeLength[i] = hypot(1.0*x[i] - x[p[i]], 1.0*y[i] - y[p[i]]);
            edgeSum += edgeLength[i];
        }
        plants[parent[i]].push_back(i);
    }

    F0(i, m) BuildHull(i);

    CL(0, separated); F0(i, m) separated[i][i] = 1; group_n = 1; group_cnt[0] = m; F0(i, m) group[0][i] = i;
    FirstPhase();
    SecondPhase();
    FinalPhase();
    if (savedEdgeSum > bestSaved + 1e-4) {
        UpdateBestSolution();
    }
    for (auto v : bestSolution) {
        v.first.AddToSolution();
    }
    fprintf(stderr, "My score: %.2lf\n", bestSaved);

    c_o[9] = SZ(ret) / 4;
    Report();
    return ret;
}

BuildHullまでで必要な情報を取得したあと、フェーズ1~3を実行している。BuildHullのコードは今後の参考になりそうだった。

BuildHull

bool cmp_start(int a, int b) {
    int t = cw(startpoint, a, b);
    if (t != 0) return t>0;
    t = x[a] - x[b];
    if (t != 0) return t<0;
    t = y[a] - y[b];
    if (t != 0) return t<0;
    return 0;
}

int h[N+5];
void BuildHull(int k)
{
    startpoint = plants[k][0];
    for (int i : plants[k]) {
        if (x[i] < x[startpoint] || (x[i] == x[startpoint] && y[i] < y[startpoint])) startpoint = i;
    }
    
    // a[] はこのPlantに属するすべてのEndpoint
    // 左下のStartpointから、CWでソートされている
    vector<int> a;
    a.push_back(startpoint);
    for (int i : plants[k]) {    // C++11の構文らしい
        if (x[i] != x[startpoint] || y[i] != y[startpoint]) a.push_back(i);
    }
    sort(a.begin() + 1, a.end(), cmp_start);

    // 外殻だけを取り出す
    vector<int>& h = hull[k];
    h.push_back(a[0]); h.push_back(a[1]);
    int hn = 2;
    for (int i = 2; i < SZ(a); i++) {
        while (hn >= 2 && cw(h[hn - 2], h[hn - 1], a[i]) <= 0) {
            hn--; h.pop_back();
        }
        h.push_back(a[i]);
        hn++;
    }
    h.push_back(a[0]);
}

(1)左下を起点に、時計回りにソート(2)そのうち、cwが負になる点を除く。

パート1

void FirstPhase() {
    phase = 1;
    savedEdgeSum = edgeSum;
    while (1) {
        // カットされてないPlantの集合がなくなったら、調整してから終了
        if (!group_n) {
            int sg = 1;
            while (sg) {
                sg = 0;
                // 該当カット単独でのPlant分割がなかったら、無駄なので削除
                F0(i,SZ(solution)) if (!BreaksSolution(solution[i].second)) {
                    ApplyMaskToSeparated(solution[i].second, -1);
                    savedEdgeSum += ApplyLine(solution[i].first, -1, true);;  // Lineを削除
                    solution.erase(solution.begin() + i);
                    sg = 1;
                    break;
                }
            }

            if (bestSaved < savedEdgeSum) {
                bestSaved = savedEdgeSum;
                bestSolution = solution;
            }
            return;
        }
        double t = Elapsed();
        bestRatioFound = 1e10;

        int it = 0;
        aliven = 0;
        for (int i = m; i < n; i++) if (!dead[i]) {
            alive[aliven++] = i;
            aliveLength[i] = cutRatio[i].empty() ? edgeLength[i] : *cutRatio[i].begin();
        } 

        // ある程度の時間、繰り返してベストなカット(+マスク)を得る
        int needs = 0; F0(i,group_n) needs += group_cnt[i] - 1;
        while (1) {
            it++;
            FindBestStraightCut();  // まあまあ良いCutを探し出して足す
                                    // ランダムな始点終点を選び、その後xかy方向に1ビットずつスライドさせて
                                    // ベストを選んだもの
            c_o[8]++;
            double portion = (Elapsed() - t) / (LTL - t);
            if (pow(portion, 1) * needs >= 1 || t >= LTL) break;    // なぜこの条件が良く分からない・・・
        }
        solution.push_back(bestCutFound);
        savedEdgeSum += ApplyLine(bestCutFound.first);  // Lineを足す
        ApplyMaskToSeparated(bestCutFound.second);  // separated[][], group_n,
                                                    // group_vis[], group[][], group_cnt[] を更新
    }
}

すべてのPlantが切り離されるようCutを作成。
FindBestStraightCut()だけである程度良いCutになっているはずだが、さらにしつこく時間でループしてベストな候補を選んでいる。
BreakSolution()の実装はこれ。該当カットがソロで分割してるPlantがあればTrueを返している。

bool BreaksSolution(const Mask& mask) {
    vector<int> leftSide, rightSide;
    F0(i, m) if (mask.setAt(i)) leftSide.push_back(i); else rightSide.push_back(i);
    for (int ls : leftSide) for (int rs : rightSide) if (separated[ls][rs] == 1) return true;
    return false;
}


パート2

void SecondPhase() {
    phase = 2;
    
    for (int i = m; i < n; i++) aliveLength[i] = edgeLength[i];

    // このzの使いかたがよく分からない
    int z = min(n, 2000);
    z = max(z, n / 10);
    
    while (Elapsed() < GTL) {
        int sg = 1;

        // 無駄Cutを削除
        while (sg) {
            sg = 0;
            for (int i = SZ(solution) - 1; i >= 0; i--) if (!BreaksSolution(solution[i].second)) {
                ApplyMaskToSeparated(solution[i].second, -1);
                bool needDelete = true;
                F0(j,SZ(bestSolution)) if (bestSolution[j].first == solution[i].first) { needDelete = false; break; }
                savedEdgeSum += ApplyLine(solution[i].first, -1, needDelete);
                solution.erase(solution.begin() + i);
                sg = 1;
                break;
            }
        }

        // カットされてないPlantがあるのはいいが、スコア悪化は嫌う
        if (savedEdgeSum < bestSaved + eps && group_n) {
            RevertSolution();
        }

        if (!group_n) {
            // 全Plantがカット済で、かつベストスコア更新
            if (savedEdgeSum > bestSaved + eps) {
                UpdateBestSolution();
            // 全Plantはカット済だが、スコア悪化
            } else if (savedEdgeSum < bestSaved - eps) {
                RevertSolution();
            } 

            // 適当に1~4本、Cutを削除
            if (SZ(solution) > 0) {
                int k = my.nextInt(SZ(solution));
                auto v = solution[k];
                savedEdgeSum += ApplyLine(v.first, -1);
                ApplyMaskToSeparated(v.second, -1);
                Mask mask = v.second;
                solution.erase(solution.begin() + k);
                
                int sg = my.nextInt(15);
                if (sg == 14) sg = 3;
                else if (sg >= 12) sg = 2;
                else if (sg >= 8) sg = 1;
                else sg = 0;
                if (sg)
                {
                    F0(k,SZ(solution)) {
                        int cnt = 0;
                        F0(i,m) if (mask.setAt(i) != solution[k].second.setAt(i)) cnt++;
                        if (cnt <= sg || cnt >= m - sg) {
                            auto v = solution[k];
                            savedEdgeSum += ApplyLine(v.first, -1);
                            ApplyMaskToSeparated(v.second, -1);
                            solution.erase(solution.begin() + k);
                            break;
                        }
                    }
                }
            }
        }
        bestRatioFound = 1e10;

        int it = 0;
        aliven = 0;

        for (int i = m; i < z; i++) if (!dead[i] && cutRatio[i].empty()) {
            alive[aliven++] = i;
        }
        
        // パート1と同じようにCutを足していく
        int z = 8 + my.nextInt(16);
        while (1) {
            it++;
            FindBestStraightCut();
            c_o[7]++;
            if ((it >= z) || Elapsed() > GTL) break;
        }
        if (bestRatioFound > 1e9) {
            RevertSolution();
            continue;
        }
        solution.push_back(bestCutFound);
        savedEdgeSum += ApplyLine(bestCutFound.first);
        ApplyMaskToSeparated(bestCutFound.second);
    }
}

うまくパート1とコードを共用している。

パート3は影響が小さいので省略。

FindBestStraightCut()とApplyLine()が重要なところだが、難しくて理解できなかった・・・。

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万回ほどだった。

Topcoder Open 2016 Marathon Match(Round 2)の参加日記

62位/183人の定位置だった。そろそろ一皮剥けたいものだ。

問題概要

  • 二次元平面(1024x1024)に星が100~2000個ある
  • 星の場所に、ランダムに船が1~10隻ある
  • 星の場所に、ランダムにUFOが0~(星の数/100)個ある
  • 1ターンごとに船を動かすことができる
  • 1つのターンでいくつの船を、どの星に動かしても良い(同じ星もOK)
  • 1つの星に何度でも到着してよい
  • 星は、どの船でも1回以上訪れていれば、到着済みになる
  • ただし、最初のターンに船がいる星は未到着扱い
  • UFOも1ターンごとにランダムに動く
  • 1ターンごと、UFOの現在地・次の星・次の次の星が与えられる
  • ターン数の上限は星の数x4
  • ターン上限までにすべての星が到着済になるよう、ターンごとに船に指示を与えるのが目的
  • ただし、船の移動コスト合計はなるべく小さくしたい
  • 移動コストは基本的に移動距離と同じ
  • ただし、もし船の航路とUFOの航路が同じなら、移動コストは距離x(1/1000^同航路のUFOの数)になる
  • 星の位置は完全にランダムではない。いくつか(1~16)の銀河が存在し、星は各銀河の中心から正規分布している
  • またUFOごとに飛び方に特徴がある。UFOは星を全体からN個ランダムに選んで、その中の最も近いものに飛ぶが、Nの値はUFO固有(10~10+星の数/10)

https://community.topcoder.com/longcontest/?module=ViewProblemStatement&rd=16703&pm=14268

方針

  1. UFOが近く(距離100くらい)にくるのを待ってから乗る。なかなか来なかったら最も近いものに乗ってしまう
  2. 残りターンが少なくなるまで乗り続ける。同じ星に良さそうな別UFOがいたら乗り換える
  3. UFO未到着の星で巡回セールスマン(Nearest Neighborのあと、貪欲的に改善)

巡回セールスマンは焼きなましで組んでみたが、どうもうまく動かず断念した。

上位との比較

  • 全体的な方針は一緒
  • TS部分で大きな差がついていそう
  • TSは焼きなましの人と、NNのあと2.5optなどTSP特有のアルゴリズムがいる
  • 最上位の方たちは、UFOランデブーもかなり注意深いロジックを組んでいる。UFOの挙動から特徴を推定・シミュレーションして、gain/costを判断してたり

Psyhoさん(暫定1位)のコード http://pastebin.com/GQV4yrDj
(UFOランデブーではモンテカルロシミュレーション、TSでは焼きなましを使ってるらしい)

気づかなかったこと+実装できなかったこと

  • TSパートでもUFOが使えるのに気づかなかった(特に最初に1手)
  • UFOから降りるとき、未到着の星の近くで降りのも良い考えだった
  • 別UFOと重なったときの判断に、UFOごとの統計情報が使えた
  • 最初の1手で船を動かさない場合、その星は到着済扱いだった

反省点
焼きなましを練習しよう。

ベルマンフォード法とダイカストラ法について

メモ。
AOJのこの問題で、入力数が大きいパターンでタイムアウトになっていた。

有向グラフが与えられる。このとき、頂点rから各頂点までの最短経路を求めよ。

単一始点最短経路 | グラフ | Aizu Online Judge

ダイカストラではなくベルマンフォードで解いていたのが原因。この問題では、Vが最大100000、Eが最大500000と大きいので、O(VE)であるベルマンフォードだと制限時間内(3s)に処理が終わらない。 追記:ベルマンフォードを適切に実装したら問題なかったようだ。コードも差し替えておく。
ダイカストラ: O(E logV)
ベルマンフォード: O(VE)

教訓:負の辺がない場合はダイカストラの利用を考える

  • ベルマンフォードでの解答(抜粋)
/// <summary>
/// find the shortest path from start to all destination
/// </summary>
public class BellmanFord
{
    public List<edge> Edge = new List<edge>();
    public int V;   // number of vertex

    public class edge
    {
        public int from, to, cost;
        public edge(int _from, int _to, int _cost)
        {
            from = _from;
            to = _to;
            cost = _cost;
        }
    }

    public BellmanFord(int v)
    {
        V = v;
    }

    /// <summary>
    ///  return shortestPath[V] represents distance from startIndex
    /// </summary>
    public int[] GetShortestPath(int startIndex)
    {
        int[] shortestPath = new int[V];
        var INF = Int32.MaxValue;

        for (int i = 0; i < V; i++) shortestPath[i] = INF;
        shortestPath[startIndex] = 0;
        while (true)
        {
            bool update = false;
            foreach (edge e in Edge)
            {
                if (shortestPath[e.from] != INF && shortestPath[e.to] > shortestPath[e.from] + e.cost)
                {
                    shortestPath[e.to] = shortestPath[e.from] + e.cost;
                    update = true;
                }
            }
            if (!update) break;
        }

        return shortestPath;
    }

    /// <summary>
    ///  return true if it has negative close loop
    /// </summary>
    public bool HasNegativeLoop()
    {
        int[] d = new int[V];
        for (int i = 0; i < V; i++)
        {
            foreach (edge e in Edge)
            {
                if (d[e.to] > d[e.from] + e.cost)
                {
                    d[e.to] = d[e.from] + e.cost;
                    if (i == V - 1) return true;
                }
            }
        }
        return false;
    }
}
  • ダイカストラでの解答(抜粋)
/// <summary>
/// Get min cost between two points
/// </summary>
public class Dijkstra
{
    private class Edge
    {
        public int T { get; set; }
        public int C { get; set; }
    }

    private int maxIndex = -1;
    private const int INF = Int32.MaxValue;

    private Dictionary<int, List<Edge>> _edgeArray = new Dictionary<int, List<Edge>>();

    /// <summary>
    /// Dijkstra, get min cost between two points
    /// should not contain negatie cost path
    /// </summary>
    /// <param name="size">max index of vertice</param>
    public Dijkstra(int size)
    {
        maxIndex = size + 1;
    }

    // Add a path with its cost
    public void AddPath(int s, int t, int cost)
    {
        if (!_edgeArray.ContainsKey(s)) _edgeArray.Add(s, new List<Edge>());
        _edgeArray[s].Add(new Edge() { T = t, C = cost });
    }

    //Get the min cost from s
    //return Int32.MaxValue if no path
    public int[] GetMinCost(int s)
    {
        int[] cost = new int[maxIndex];
        for (int i = 0; i < cost.Length; i++) cost[i] = INF;
        cost[s] = 0;

        var priorityQueue = new PriorityQueue<ComparablePair<int, int>>(_edgeArray.Sum( e => e.Value.Count()) + 1);
        priorityQueue.Push(new ComparablePair<int, int>(0, s));

        while (priorityQueue.Count() > 0)
        {
            var costDestinationPair = priorityQueue.Pop();
            if (cost[costDestinationPair.Item2] < costDestinationPair.Item1) continue;
            if (!_edgeArray.ContainsKey(costDestinationPair.Item2)) continue;

            foreach (var edge in _edgeArray[costDestinationPair.Item2])
            {
                var newCostToi = costDestinationPair.Item1 + edge.C;
                if (newCostToi < cost[edge.T])
                {
                    cost[edge.T] = newCostToi;
                    priorityQueue.Push(new ComparablePair<int, int>(newCostToi, edge.T));
                }
            }
        }

        return cost;
    }
}

PriorityQueueを自作する必要があるのと、Dijkstra内を単純な配列にするとメモリが足りなくなるのが面倒だった。
こういう問題だとC++を使ったほうが楽そうだ。

Topcoder Open 2016 Marathon Match(Round 1)の参加日記

結果は55位/136人と振るわなかったが、今後のために記録を残しておく。

問題概要
https://community.topcoder.com/longcontest/?module=ViewProblemStatement&rd=16702&pm=10729

  • 平面上に、木の幹がランダムに点として存在する
  • それぞれの幹からは、ランダムに根が生えている
  • 根の終端からまた根もランダムに生えている
  • このとき、すべての木の幹が別のグループになるよう、平面を切りたい
  • 切断線は無限大の長さになる(=ある点からある点まで、とは切れない)
  • なるべく切られてしまう根を減らして、切断線を求めよ

実装内容
雑な貪欲。すべての幹の組み合わせを求め、幹間の距離が短い順に

  1. 幹間に線を引き、30等分くらいで切断開始点を決める
  2. それぞれの開始点から、数度刻みでぐるっと一周で切断面を作ってみる
  3. 一番、切ってしまった根の長さ合計が少なかったものを採用

を繰り返しただけ。
あとは、すでに切られている幹間は無視したり、根が多すぎるときは単純化して計算量を減らしたりした。

感想と反省
開始3日くらいで実装・提出してから、終了までほぼ改善できずに終わってしまった。
7日目くらいで焼きなましをやってみようとしたのだが、

  • 体調などのせいで時間が取れない(そして睡眠不足の悪循環)
  • データ構造が良くなくほぼ作り直し

のせいで日の目を見ずにタイムアップ。
まあ3日で何とか黄色を維持できるコードが書けたし、独自Visualizerも悪くない感じにできたので、上達はしているのかもしれない。

教訓

  • 時間がない時こそ、実装を急がず戦略を立てて行動する
  • プロトタイプといえど、データ構造はなるべく柔軟性を持たせる