AOJ 1334 Cubic Colonies

これはAOJ-ICPC Advent Calendar 2015の2日目の記事です。

AOJ自体長らく触っていませんが、自分が解いたことのある問題の中で「AOJ (問題番号)」でぐぐっても解説記事が出てこない問題がないか探して、この問題について書くことにしました。

問題: Cubic Colonies | Aizu Online Judge

2012年のICPC地区予選の最終防衛問題で、本番中に解かれた問題です。
複数の立方体からなる物体があり、物体の表面にはケーブルを敷くことができます。物体の表面上の2点が与えられるので、その間を結ぶ長さ最小のケーブルを求めて下さい。
内容把握したと思いつつ問題文の下にある図が目に入ると「そんなところも通れるんかーい!」とツッコまずにはいられません。
物理的に可能かどうかはさておき、表面であれば自由にケーブルを敷けるということで、表面に頂点と辺を追加して最短経路問題にして解けます。

解法

どこに面があるか調べる

隣り合ったマスの片方のみに物体があれば面ができます。

面に頂点と辺を加える

ケーブルがどのように敷かれるか展開図で考えてみると、2つの格子点間を結ぶ線分の組み合わせになることが分かります。サイズも 3\times 3\times 3と小さいので、一つの線分もせいぜい 10 \times 10程度に収まりそうです。で、ひとつの面だけ取り出してみると、長さ1の線分を[0,1]とみなして0,1,{n/m; 0 < n < m <= D}に頂点を置く、というのを全ての辺についてやって、異なる辺にある2つの頂点間に辺を張ればよいです。分母の最大値をDと置いてますが、そんなに大きくなりそうにないのでコード内ではとりあえず9にしています。

単一始点最短路を解く

スタートに対応する頂点からゴールに対応する頂点までの最短路を計算します。
面の数は108、頂点の数は粗く見積もって 12\times 4\times 3 \times D^2個未満、辺の数は粗く見積もって 108\times 3 \times D^4個未満、ダイクストラ法で間に合います。(一つのインプットが複数データセットからなるなら、ちゃんとその数も書いてあるといいのになぁ)

2年前にAOJに提出したコードではダイクストラをしてました。(295ms)

今年向けに書き直すついでに高速化しようと思ったのですが、バグが取れませんでした。日付が変わる前に一旦投稿します・・・

感想

やることを挙げると大した事はないですが、実装量は意外とあって、その割にはライブラリで殴る問題でもない不思議な問題でした。デバッグがつらい。

コード

2年前のクソコード

#include <bits/stdc++.h>
#define iter(i,s,n) for(int i=int(s);i<int(n);i++)
#define rep(i,n) iter(i,0,n)
#define all(c) (c).begin(), (c).end()
using namespace std;
typedef long long ll;
template<class T>void chmin(T &t, T f){if(t > f)t = f;}
int gcd(int a,int b){return b?gcd(b,a%b):a;}

vector<double> frac;

struct edge{int to;double cost;};
struct edgeid{int x,y,z,d;};
vector<edge> G[8000];

int edge_id(int x,int y,int z,int d){return ((x*4+y)*4+z)*3+d;}
int edge_ide(edgeid e){return edge_id(e.x,e.y,e.z,e.d);}
int node_id(int x,int y,int z, int d,int f){return edge_id(x,y,z,d)*frac.size()+f;}
int node_ide(edgeid e, int f){return node_id(e.x,e.y,e.z,e.d,f);}

struct node{double x,y,z;};
node no[8000];
void print_node(node a){printf("%f %f %f\n",a.x,a.y,a.z);}
void print(int id){print_node(no[id]);}

double d[8000];

void connect(int x,int y,int z,int d){
    int S = frac.size();
    int added = 0;
    //printf("<%d %d %d %c>\n",x,y,z,d+'x');
    //x
    //e[x,y,z,2], e[x,y,z,3], e[x,y+1,z,2], e[x,y,z+1,3]
    vector<edgeid> f;
    if(d==0){
        //printf("<%d %d %d x>\n",x,y,z);
        f.push_back({x,y,z,1});
        f.push_back({x,y,z,2});
        f.push_back({x,y,z+1,1});
        f.push_back({x,y+1,z,2});
    }
    //y
    if(d==1){
        //printf("<%d %d %d y>\n",x,y,z);
        f.push_back({x,y,z,0});
        f.push_back({x,y,z,2});
        f.push_back({x,y,z+1,0});
        f.push_back({x+1,y,z,2});
    }
    //z
    if(d==2){
        //printf("<%d %d %d z>\n",x,y,z);
        f.push_back({x,y,z,0});
        f.push_back({x,y,z,1});
        f.push_back({x,y+1,z,0});
        f.push_back({x+1,y,z,1});
    }

    //mukaiau
    rep(i1,2){
        //printf("[[%s++]]\n",d==1?i1?"z":"y":d==2?i1?"z":"x":i1?"y":"x");
        rep(j1,S){
            int s = node_ide(f[i1],j1);
            rep(j2,S){
                int t = node_ide(f[i1+2],j2);
                double d = hypot(1.0,frac[j1]-frac[j2]);
                G[s].push_back({t,d}); G[t].push_back({s,d});
                //print(s); print(t); printf("%f\n",d);
                added++;
            }
        }
    }

    //tonariau
    rep(j1,S){
        int s = node_ide(f[0],j1);
        rep(j2,S){
            int t = node_ide(f[1],j2);
            double d = hypot(frac[j1],frac[j2]);
            G[s].push_back({t,d}); G[t].push_back({s,d});
            added++;
        }
    }

    rep(j1,S){
        int s = node_ide(f[0],j1);
        rep(j2,S){
            int t = node_ide(f[3],j2);
            double d = hypot(1-frac[j1],frac[j2]);
            G[s].push_back({t,d}); G[t].push_back({s,d});
            added++;
        }
    }

    rep(j1,S){
        int s = node_ide(f[1],j1);
        rep(j2,S){
            int t = node_ide(f[2],j2);
            double d = hypot(1-frac[j1],frac[j2]);
            G[s].push_back({t,d}); G[t].push_back({s,d});
        }
    }

    rep(j1,S){
        int s = node_ide(f[2],j1);
        rep(j2,S){
            int t = node_ide(f[3],j2);
            double d = hypot(1-frac[j1],1-frac[j2]);
            G[s].push_back({t,d}); G[t].push_back({s,d});
        }
    }

    //cout << added << endl;
}

int main(){
    //init
    for(int i=2; i<=9; i++){
        for(int j=1; j<i; j++){
            if(gcd(i,j) > 1) continue;
            frac.push_back((double)j/i);
        }
    }
    frac.push_back(0); frac.push_back(1);
    sort(all(frac));

    //cout << frac.size() << endl;

    rep(x,4) rep(y,4) rep(z,4) rep(f,frac.size()){
        //printf("%d %d %d\n",node_id(x,y,z,0,f),node_id(x,y,z,1,f),node_id(x,y,z,2,f));
        no[node_id(x,y,z,0,f)] = {x+frac[f],(double)y,(double)z};
        no[node_id(x,y,z,1,f)] = {(double)x,y+frac[f],(double)z};
        no[node_id(x,y,z,2,f)] = {(double)x,(double)y,z+frac[f]};
    }
    //repl
    int sx,sy,sz,tx,ty,tz;
    while(cin >> sx >> sy >> sz >> tx >> ty >> tz, sx|sy|sz|tx|ty|tz){
        char bl[5][5][5];
        memset(bl,'.',sizeof(bl));
        iter(z,1,4) iter(y,1,4) iter(x,1,4){
            int c = getchar();
            while(c == 32 || c == 10) c = getchar();
            bl[x][y][z] = c;
        }

        //rep(i,5){
        //    rep(j,5){
        //        rep(k,5) putchar(bl[i][j][k]);
        //        puts("");
        //    }
        //    puts("------------");
        //}

        //edge_init
        rep(i,8000) G[i].clear();

        rep(x,4) rep(y,4) rep(z,4){
            vector<int> same;
            same.push_back(node_id(x,y,z,0,0));
            same.push_back(node_id(x,y,z,1,0));
            same.push_back(node_id(x,y,z,2,0));
            if(x>0) same.push_back(node_id(x-1,y,z,0,frac.size()-1));
            if(y>0) same.push_back(node_id(x,y-1,z,1,frac.size()-1));
            if(z>0) same.push_back(node_id(x,y,z-1,2,frac.size()-1));
            rep(i,same.size()) rep(j,same.size()){
                if(i==j) continue;
                G[same[i]].push_back({same[j],0.0});
            }
            //G[node_id(x,y,z,0,0)].push_back({node_id(x,y,z,0,frac.size()-1),1.0});
            //G[node_id(x,y,z,1,0)].push_back({node_id(x,y,z,1,frac.size()-1),1.0});
            //G[node_id(x,y,z,2,0)].push_back({node_id(x,y,z,2,frac.size()-1),1.0});
        }

        rep(x,4) rep(y,4) rep(z,4){
            if(bl[x+1][y][z] != bl[x][y][z]) connect(x,y-1,z-1,0);
            if(bl[x][y+1][z] != bl[x][y][z]) connect(x-1,y,z-1,1);
            if(bl[x][y][z+1] != bl[x][y][z]) connect(x-1,y-1,z,2);
        }

        //dist init
        rep(i,8000) d[i] = 1e9;
        d[node_id(sx,sy,sz,0,0)] = 0;
        //G[node_id(sx+1,sy+1,sz+1,0,0)].push_back({node_id(sx+1,sy+1,sz+1,1,0),0.0});
        priority_queue<pair<double,int>, vector<pair<double,int>>,greater<pair<double,int>>> pq;
        pq.push({0.0, node_id(sx,sy,sz,0,0)});
        int goal = node_id(tx,ty,tz,0,0);
        int zz = 0;
        while(!pq.empty()){
            zz++;
            //printf("[%d]",zz);
            auto p = pq.top(); pq.pop();
            int v = p.second;
            if(d[v] < p.first) continue;

            if(v == goal){
                cout << fixed << setprecision(12) << d[goal] << endl;
                break;
            }

            //printf("((%d))",G[v].size());
            rep(i,G[v].size()){
                edge e = G[v][i];
                if(d[e.to] > d[v] + e.cost){
                    d[e.to] = d[v] + e.cost;
                    pq.emplace(d[e.to],e.to);
                }
            }
        }

        //cout << fixed << setprecision(12) << d[goal] << endl;
    }
}