Hatena::Grouptopcoder

kuuso1の日記

2017-07-21

TCO17 Round 3 PoisonedWine の少数毒入りワインケースのソルバーについて

| 05:48 |  TCO17 Round 3 PoisonedWine の少数毒入りワインケースのソルバーについて - kuuso1の日記 を含むブックマーク はてなブックマーク -  TCO17 Round 3 PoisonedWine の少数毒入りワインケースのソルバーについて - kuuso1の日記

TCO2017のMarathon Match に参加していました.

https://community.topcoder.com/longcontest/?module=ViewProblemStatement&rd=16944&pm=14636

・N本のワインにB本の毒入りワインが混入している.

・S枚のstrip(試験紙)とT回のテストができる.

・1回のテストは,各stripに確認したいワインを浸して行う.

・stripには複数のワインを同時に浸しても良いが,1つでも毒入りワインが混入しているとstripが反応して以後使用不可になる.

・反応しなかったstripは次以降のテストで再び使用できる.

毒入りと疑わしいワインをすべて廃棄するので,被疑ワインの数を最小化せよ.

(制約) 50 ≤ N ≤ 10000,5 ≤ S ≤ 20,1 ≤ T ≤ 10,1 ≤ B ≤ (N / 50), TLE: 60sec

(スコア) (残ったワイン数 / 正常ワイン数) ^ 2,残ったワインの中に毒入りワインがある場合はスコア0,5000ケースの平均値を最終スコアとする.

方針が立たず,T=1のケースはモンテカルロで各stripに割り当てるのに一番よさそうなサイズを決める期待値ギャンブルに.

またそのほかのケースは,確率pで毒が混入しているとみなして,s枚のstripにそれぞれlen個のワインを割り当てると,

大体生き残るワインの個数の期待値が s * len * (1-p)^len くらいになりそうだというところで,各test前にいくつ割り当てるのがいいかを

ビームサーチで適当に探すという方針で実装.provisional 66位で終了という今一つなスコアで終了.

で,全体的には確率ゲーをどう最大化していくかという問題なのですが,Bが非常に小さいケースだけはexplicitに追い込んでいく方法があるので,この日記ではそれについて書きます.

B=1のケース

ググると出てきたワイン問題詳しく に書かれているそのままですが,2の冪に割り当てていくことでかなり絞り込むことができます.

f:id:kuuso1:20170722045133p:image

2回目のテストができる場合,2回目に残るstripの数に合わせてあらかじめ重みを付けてグルーピングすることでより効率的に追い込んでいけます.

f:id:kuuso1:20170722045241p:image

一般にS個stripがあってT回できる場合,

f:id:kuuso1:20170722045242p:image

出てくる係数が綺麗すぎてびっくりする.二項定理万歳.

Bが複数のケース

B=1のソルバーを使うためには,毒ワインを1つしか含まないクラスタに分ける必要がありますが,これはいい方法が思いつかず.

ただ,例えばB=2の場合だと,全体をS個に分けてテストをすると,高々2枚のstripしか死なないことが分かります.

1枚だけ死んだ場合は,その1枚の中にB=2が2つとも含まれているのですが,その時点で(S-1)枚分はシロだとわかるので,

何回か繰り返すとちょうどB枚死ぬケースを引くケースが確率的に出てきます.その時にパラレルにB=1のソルバーを走るように分岐する方針にしました.

S>Bのケースでしか意味がないうえに,例えばS=19, B=16みたいなケースで 16枚のstripが死ぬ確率は相当低いので,事実上はB<5くらいまでしか機能しません...

で,そのパラレルケースですが,stripの割り振り方に少し改善ポイントがありました.

f:id:kuuso1:20170722045243p:image

実装

B=1ソルバー

long[][] Resolution;  // Resolution[t][s] : Pow((1+t), s) 前計算済み

int[] OneSolver(){
    // 使うstrip数(最大分解できる最小strip数)を決めて UnitSolverに投げる
    List<int> candi = new List<int>();
    for(int i=0;i<N;i++) if(WS[i] != TriState.True) candi.Add(i);
    if(candi.Count == 1){
        T = 0; S = 0; return candi.ToArray();
    }
    int testTurn = T;
    int useStrip = S;
    if(Resolution[testTurn].Length < useStrip) useStrip = Resolution[testTurn].Length - 1;
    while(Resolution[testTurn][useStrip - 1] >= candi.Count) useStrip--;
    var bad = UnitSolver(testTurn,useStrip,candi);
    return bad.ToArray();
}

List<int> UnitSolver(int turn, int strip, List<int> candi){
    // turn/strip での分解能 のグループに candi(被疑メンバのリスト)のメンバをグルーピング
    if(strip == 0 || turn == 0 || candi.Count == 1) return candi;
    long reso = Resolution[turn][strip];
    List<int>[] gp = new List<int>[reso];
    for(int i=0;i<reso;i++) gp[i] = new List<int>();
    for(int i=0;i<candi.Count;i++){
        gp[i % reso].Add(candi[i]);
    }
    
    // 2^strip に分ける (死ぬstrip数はbitCount)
    List<int>[] cruster = new List<int>[1<<strip];
    for(int i=0;i<(1<<strip);i++) cruster[i] = new List<int>();
    int assigned = 0;
    for(int i=0;i<(1<<strip);i++) {
        int bcnt = PopCnt(i);
        int mem = (int) Resolution[turn - 1][strip - bcnt];
        for(int j=0;j<mem;j++){
            cruster[i].AddRange(gp[assigned]);
            assigned++;
        } 
    }

    // 各クラスタをstripにのせる
    List<int>[] stripEntry = new List<int>[strip];
    for(int i=0;i<strip;i++) stripEntry[i] = new List<int>();
    for(int i=0;i<(1<<strip);i++){
        for(int j=0;j<strip;j++){
            if(((i>>j) & 1) == 1){
                stripEntry[j].AddRange(cruster[i]);
            }
        }
    }
    
	// empty queryが発生しないようにひと手間
    List<String> ss = new List<String>();
    int[] Idx = new int[strip];
    for(int i=0;i<strip;i++){
        if(stripEntry[i].Count == 0){
            Idx[i] = -1;
        } else {
            var csv = String.Join(",",stripEntry[i]);
            ss.Add(csv);
            Idx[i] = ss.Count - 1;
        }
    }
    
	// テストクエリ 
    var ret = PoisonTest.useTestStrips(ss.ToArray());
    T--;
    int[] ret0 = new int[strip];
    for(int i=0;i<strip;i++){
        if(Idx[i] == -1){
            ret0[i] = 0;
        } else {
            ret0[i] = ret[Idx[i]];
        }
    }

    // テスト結果から被疑クラスタを決定
    int crusterIdx = 0;
    int deathStrip = 0;
    for(int i=0;i<strip;i++){
        if(ret0[i] == 1){
            crusterIdx |= (1 << i);
            deathStrip++;
        }
    }

    for(int i=0;i<(1<<strip);i++){
        if(i == crusterIdx) continue;
        foreach(var m in cruster[i]) WS[m] = TriState.True;
    }
    S -= deathStrip;
    return UnitSolver(turn - 1,strip - deathStrip, cruster[crusterIdx]);
}

int PopCnt(int bits){
    bits = (bits & 0x55555555) + (bits >> 1 & 0x55555555);
    bits = (bits & 0x33333333) + (bits >> 2 & 0x33333333);
    bits = (bits & 0x0f0f0f0f) + (bits >> 4 & 0x0f0f0f0f);
    bits = (bits & 0x00ff00ff) + (bits >> 8 & 0x00ff00ff);
    return (bits & 0x0000ffff) + (bits >>16 & 0x0000ffff);
}

パラレルソルバー

UnitSolverを使いまわせるかと思いきや,テストクエリは全体をまとめて一度に行わないといけなくて発狂した.

int S;  // 残りstrip数
int T;  // 残りtest数
int B:  // 毒ワイン数
TriState[] WS; // 各ワインの状態 True:シロ / False:クロ / Undef:未定

int[] KSolver() {
    // 全体を S 個に分割, B個死ぬのを引くまで繰り返す
    var rnd = new Xor128(1234567);
    while(T > 0 && S > 0 && S > B){
        List<int> L = new List<int>();
        for(int i=0;i<N;i++) if(WS[i] != TriState.True) L.Add(i);
        for(int i=L.Count-1;i>=0;i--){
            int j = rnd.Next(i);
            var tmp = L[i]; L[i] = L[j]; L[j] = tmp;
        }
        List<int>[] test = new List<int>[S];
        int NN = L.Count;
        for(int i=0;i<S;i++) test[i] = new List<int>();
        for(int i=0;i<NN;i++){
            test[i % S].Add(L[i]);
        }
        List<String> ss = new List<String>();
        for(int i=0;i<S;i++){
            if(test[i].Count != 0) ss.Add(String.Join(",",test[i]));
        } 
        
        var ret = PoisonTest.useTestStrips(ss.ToArray());
        T--;
        int dead = 0;
        for(int i=0;i<ss.Count;i++){
            if(ret[i] == 0){
                foreach(var n in test[i]) WS[n] = TriState.True;
            } else {
                dead++;
            }
        }
        S -= dead;
        if(dead == B){
        // B個死んだ.やったぜ!
Console.Error.WriteLine("Yattaze!! KS rest:T:{0},S:{1}",T,S);
            // 生きてるstripのワインはすべてシロ. 
            var hs = new HashSet<int>();
            List<List<int>> candis = new List<List<int>>();
            for(int i=0;i<ss.Count;i++){
                if(ret[i] == 1){
                    foreach(var n in test[i]) hs.Add(n);
                    candis.Add(new List<int>(test[i]));
                }
            }
            for(int i=0;i<N;i++) if(!hs.Contains(i)) WS[i] = TriState.True;
            
            // 最適なstripの残り割り振りを AssignStripsで計算してParaSolverに投げる.
            int[] stp = AssignStrips(S,candis.Select(l => l.Count).ToArray());
            
            var bads = ParaSolver(T,stp,candis);
        }
    }
    
    List<int> bad = new List<int>();
    for(int i=0;i<N;i++) if(WS[i] != TriState.True) bad.Add(i);
    return bad.ToArray();

}

List<List<int>> ParaSolver(int turn, int[] stp, List<List<int>> candis){
    
    if(turn == 0){
        return candis;
    }
    int NP = candis.Count;
    
    // 基本的にUnitSoverのコピペ+ループ
    long[] reso = new long[NP];
    List<int>[][] gp = new List<int>[NP][];
    List<int>[][] cruster = new List<int>[NP][];
    List<int>[][] stripEntry = new List<int>[NP][];
    for(int k=0;k<NP;k++){
        reso[k] = Resolution[turn][stp[k]];
        gp[k] = new List<int>[reso[k]];
        for(int i=0;i<reso[k];i++) gp[k][i] = new List<int>();
        for(int i=0;i<candis[k].Count;i++){
            gp[k][i % reso[k]].Add(candis[k][i]);
        }
        
        cruster[k] = new List<int>[1<<stp[k]];
        for(int i=0;i<(1<<stp[k]);i++) cruster[k][i] = new List<int>();

        int assigned = 0;
        for(int i=0;i<(1<<stp[k]);i++) {
            int bcnt = PopCnt(i);
            int mem = (int) Resolution[turn - 1][stp[k] - bcnt];
            for(int j=0;j<mem;j++){
                cruster[k][i].AddRange(gp[k][assigned]);
                assigned++;
            } 
        }

        stripEntry[k] = new List<int>[stp[k]];
        for(int i=0;i<stp[k];i++) stripEntry[k][i] = new List<int>();
        for(int i=0;i<(1<<stp[k]);i++){
            for(int j=0;j<stp[k];j++){
                if(((i>>j) & 1) == 1){
                    stripEntry[k][j].AddRange(cruster[k][i]);
                }
            }
        }
    }

    int[][] Idx = new int[NP][];
    List<String> ss = new List<String>();
    for(int k=0;k<NP;k++){
        Idx[k] = new int[stp[k]];
        for(int i=0;i<stp[k];i++){
            if(stripEntry[k][i].Count == 0){
                Idx[k][i] = -1;
            } else {
                var csv = String.Join(",",stripEntry[k][i]);
                ss.Add(csv);
                Idx[k][i] = ss.Count - 1;
            }
        }
    }
    
    // テストクエリはまとめて実施しないとだめ
    var ret = PoisonTest.useTestStrips(ss.ToArray());
    T--;

    int[][] ret0 = new int[NP][];
    int[] crusterIdx = new int[NP];
    int[] deathStrip = new int[NP];
    for(int k=0;k<NP;k++){
        ret0[k] = new int[stp[k]];
        for(int i=0;i<stp[k];i++){
            if(Idx[k][i] == -1){
                ret0[k][i] = 0;
            } else {
                ret0[k][i] = ret[Idx[k][i]];
            }
        }

        crusterIdx[k] = 0;
        deathStrip[k] = 0;
        for(int i=0;i<stp[k];i++){
            if(ret0[k][i] == 1){
                crusterIdx[k] |= (1 << i);
                deathStrip[k]++;
            }
        }

        for(int i=0;i<(1<<stp[k]);i++){
            if(i == crusterIdx[k]) continue;
            foreach(var m in cruster[k][i]) WS[m] = TriState.True;
        }
        S -= deathStrip[k];
    }
    
    List<List<int>> ncandis = new List<List<int>>();
    for(int i=0;i<NP;i++) ncandis.Add(cruster[i][crusterIdx[i]]);
    
    // 再度,最適なstripの残り割り振りを AssignStripsで計算してParaSolverに投げる.
    int[] nstp = AssignStrips(S,ncandis.Select(l => l.Count).ToArray());

    return ParaSolver(turn - 1,nstp, ncandis);
}

int[] AssignStrips(int strip, int[] mcnt){
    // strip:合計strip数
    // mcnt: 各クラスタのメンバーの数.
    
    int n = mcnt.Length;
    int s = strip;
    if(strip == 0){
        return new int[n];
    }
    long XX = (long) 1e18; 
    long[] dp = new long[s+1];
    for(int i=0;i<=s;i++) dp[i] = XX;
    dp[0] = 0;

    List<int>[] dp2 = new List<int>[s+1];
    dp2[0] = new List<int>();
    
    // 各クラスタにいくつのstripを割り振ればよいかのdp(flying Table)
    // 経路復元が面倒なのでListで経路を持って同時に更新(dp2)
    for(int t=0;t<n;t++){
        long[] ndp = new long[s+1];
        for(int i=0;i<=s;i++) ndp[i] = XX;
        List<int>[] ndp2 = new List<int>[s+1];
        for(int i=0;i<=s;i++){
            if(dp[i] == XX) continue;
            for(int j=0;j<Resolution[T].Length;j++){
                if(i + j > s) continue;
                long inc = Ceil(mcnt[t],Resolution[T][j]);
                if(ndp[i+j] > dp[i] + inc){
                    ndp[i+j] = dp[i] + inc;
                    ndp2[i+j] = new List<int>(dp2[i]);
                    ndp2[i+j].Add(j);
                }
            }
        }
        dp = ndp;
        dp2 = ndp2;
    }

    long min = XX;
    int idx = -1;
    for(int i=0;i<=s;i++){
        if(dp[i] < min){
            idx = i;
        }
    }
    return dp2[idx].ToArray();
}

long Ceil (long x, long d){
    return x / d + (x % d == 0 ? 0 : 1);
}

少数ケースだけ結構highest取れてるんじゃないかという気もするのでそこだけは希望をもって最終結果たのしみにしてます.

まぁ,全体の10%に満たないボリュームの重箱なんですが.