概要

母船・デブリ・宇宙ステーションを頂点とし、距離 \(D\) 以下の場所へ辺を張ったグラフで、最短経路を求める問題です。

ただし、すべての辺を調べると間に合わないため、格子による空間分割と円弧の上側包絡線を使い、幅優先探索（BFS）を高速化します。

考察

1. 最小ジャンプ回数は BFS で求められる

どのジャンプもコストは \(1\) なので、BFS によって最小ジャンプ回数を求められます。

ただし、母船からステーションへの直接移動は禁止されています。そのため、最初に母船から届くデブリを距離 \(1\) の頂点として登録し、以後はデブリを経由して探索します。

問題は、各デブリから移動できる相手をどう探すかです。すべてのデブリの組について距離を計算すると \(O(N^2)\) となり、\(N \leq 10^5\) では間に合いません。

2. 一辺 \(D/2\) の格子に分ける

平面を一辺 \(D/2\) の正方形に分け、各デブリを所属するセルにまとめます。

点 \((x,y)\) のセル番号は、次のように求めます。

\((g_x,g_y)=\left(\left\lfloor \frac{2x}{D}\right\rfloor,\left\lfloor \frac{2y}{D}\right\rfloor\right)\)

この分割には、重要な性質が二つあります。

同じセル内では、どのデブリ同士も移動できる

セルの対角線の長さは、

\(\sqrt{(D/2)^2+(D/2)^2}=\frac{D}{\sqrt{2}}<D\)

です。したがって、同じセル内のデブリ同士は必ず \(1\) 回で移動できます。

あるセルで初めてデブリに到達したとき、その最短距離を \(k\) とします。そのセルの残りのデブリにも、遅くとも \(k+1\) 回で到達できます。

よって、一つのセルに含まれるデブリの最短距離は、\(k\) と \(k+1\) の高々二種類です。これにより、一つのセルが BFS の探索元になるのも高々 \(2\) 層に限られます。

調べる必要があるのは周囲の \(24\) セルだけ

セル番号がいずれかの軸で \(3\) 以上離れていれば、点同士の距離は \(D\) より大きくなります。

したがって、別セルへの移動では、

\(-2 \leq \Delta g_x \leq 2,\qquad -2 \leq \Delta g_y \leq 2\)

を満たす周囲の \(24\) セルだけを調べれば十分です。

3. 格子に分けるだけでは、まだ不十分

一つのセルに多数のデブリが集まる場合があります。

探索元セルに \(a\) 個、移動先セルに \(b\) 個のデブリがあるとき、すべての組を確認してしまうと \(O(ab)\) です。格子に分けただけでは、最悪 \(O(N^2)\) を避けられません。

そこで、次の問い合わせを高速に処理します。

現在の BFS 層にある探索元の集合のうち、移動先の点に届くデブリが一つでもあるか？

この判定に使うのが、円弧の上側包絡線です。

アルゴリズム

1. 円の右端を使って到達可能性を判定する

まず、移動先セルが探索元セルより右側にある場合を考えます。

探索元のデブリを \(p=(x_p,y_p)\) とします。そこから到達できる範囲は、中心 \(p\)、半径 \(D\) の円の内部です。

高さ \(t\) における、この円の右端の座標は、

\(f_p(t)=x_p+\sqrt{D^2-(t-y_p)^2}\)

です。ただし、\(|t-y_p|\leq D\) の範囲でのみ定義されます。

移動先の点 \(q=(x_q,y_q)\) は、すべての探索元より右側にあります。そのため、ある探索元から到達できる条件は、

\(x_q\leq \max_p f_p(y_q)\)

となります。

つまり、高さ \(y_q\) で最も右まで届く円を一つ選べれば、その円の中心との距離を調べるだけで判定できます。

2. 円弧の上側包絡線を構築する

複数の関数 \(f_p(t)\) に対し、その最大値を集めたものを上側包絡線と呼びます。

コードの Envelope は、次の情報を保持しています。

• 各範囲で最も右まで届く円弧
• その円弧が最良になる最初の整数座標 start

半径が等しい円の右半分の円弧を考えると、中心の高さが異なる二つの円弧の優劣は、共通の定義域で高々一度しか切り替わりません。高さを増やしたとき、中心が低い円弧から高い円弧へと切り替わります。

この性質を利用し、次の手順で構築します。

1. 探索元を中心の高さ順にソートする。
2. 高さが同じなら、最も右にある中心だけを残す。
3. 円弧を順に追加する。
4. 末尾の円弧が不要になったら取り除く。
5. 新しい円弧に切り替わる最初の整数座標を二分探索する。

問い合わせでは、start の列を二分探索することで、その高さで最良の円弧を求められます。

得られた中心と移動先の点について、最後に通常の距離判定を行います。どの円もその高さに届かない場合も、この距離判定によって到達不能と判定できます。

3. 左・上・下方向も座標変換で処理する

右方向以外についても、座標を反転・交換すれば同じ処理を使えます。

移動方向	変換後の座標
右	\((x,y)\)
左	\((-x,y)\)
上	\((y,x)\)
下	\((-y,x)\)

移動先セルとの横方向のセル番号が異なる場合は、右または左として処理します。同じ場合は、上または下として処理します。

各探索元セルについて、必要な方向の包絡線だけを構築し、同じ方向の移動先セルで使い回します。

4. セル単位で BFS を進める

各セルについて、デブリを次の三つに分けて管理します。

• remaining：まだ到達していないデブリ
• frontier：現在の BFS 層で到達したデブリ
• next：次の BFS 層で到達するデブリ

最初に、母船から距離 \(D\) 以下のデブリを frontier に入れ、現在のジャンプ回数を depth = 1 とします。

その後、frontier が空でないセルごとに、次を行います。

1. 同じセルの未到達デブリをすべて発見する。
   同じセル内なら必ず移動できるため、すべて next に入れます。

2. 周囲の \(24\) セルを調べる。
   現在のセルの frontier から包絡線を構築し、移動先セルの未到達デブリごとに到達可能性を判定します。

3. 発見したデブリを未到達集合から取り除く。
   BFS なので、初めて発見したときのジャンプ回数が最短です。

現在の層から発見したデブリには depth + 1 回で到達できます。そこからステーションに届くなら、答えは depth + 2 です。

現在の層の処理がすべて終わったら、next を新しい frontier にして次の層へ進みます。探索が尽きてもステーションに届かなければ -1 を出力します。

5. この方法で最短距離を求められる理由

• 同じセル内への移動は、すべて正しく列挙しています。
• 別セルへの移動は、届く可能性がある周囲の \(24\) セルをすべて調べています。
• 各移動先について、包絡線によって「現在の層のどれかのデブリから届くか」を正確に判定しています。
• 新しく発見したデブリは、必ず次の層で探索元にします。

したがって、辺を明示的に作った通常の BFS と同じ順序で到達可能なデブリを発見でき、最小ジャンプ回数が得られます。

計算量

座標の二分探索範囲の大きさを \(C\) とします。コードでは \(C=2\times 10^9+1\) です。

• 時間計算量: 期待 \(O(N\log N+N\log C)\)
• 空間計算量: \(O(N)\)

各セルが探索元になるのは高々 \(2\) 層で、隣接セル数も \(24\) と定数です。したがって、各未到達デブリが問い合わせ対象になる回数は定数回に抑えられます。

また、各デブリが frontier に入るのは一度だけです。包絡線のソート・構築・問い合わせを全体で合計すると、上記の計算量になります。期待計算量としているのは、セルの検索にハッシュ表を使うためです。

実装のポイント

負の座標でも、セル番号は切り下げで求める

C++ の整数除算は負の値を \(0\) に近づける方向へ丸めますが、必要なのは床関数です。

例えば \(D=4,\ x=-1\) なら、

\(\left\lfloor \frac{2x}{D}\right\rfloor=\left\lfloor-\frac12\right\rfloor=-1\)

です。コードの gridCoordinate では、負の値を別処理して正しく計算しています。

距離判定に平方根を使わない

距離が \(D\) 以下かどうかは、

\((x_1-x_2)^2+(y_1-y_2)^2\leq D^2\)

で判定します。整数のまま計算することで、境界上の点も正確に扱えます。

円弧の比較も整数演算で行う

包絡線の構築では、二つの平方根を含む式を比較します。浮動小数点数で比較すると、切り替わる座標がずれるおそれがあります。

better では、符号を確認しながら移項・平方を行い、整数演算だけで比較しています。その途中には 64 ビット整数を超える積が現れるため、__int128_t を使っています。

現在の層と次の層を混ぜない

発見したデブリをその場で frontier に追加すると、同じ BFS 層で複数回ジャンプしたことになってしまいます。必ず next に入れ、現在の層の処理が終わってから入れ替えます。

発見済みのデブリは定数時間で削除する

remaining から要素を削除するときは、末尾の要素を削除位置に移して pop_back() します。順序は不要なので、この方法で一回の削除を \(O(1)\) にできます。

ソースコード

#include <bits/stdc++.h>
using namespace std;

using ll = long long;
using i128 = __int128_t;
using ull = unsigned long long;

struct Hash {
    static ull splitmix64(ull x) {
        x += 0x9e3779b97f4a7c15ULL;
        x = (x ^ (x >> 30)) * 0xbf58476d1ce4e5b9ULL;
        x = (x ^ (x >> 27)) * 0x94d049bb133111ebULL;
        return x ^ (x >> 31);
    }
    size_t operator()(ull x) const {
        static const ull seed =
            chrono::steady_clock::now().time_since_epoch().count();
        return splitmix64(x + seed);
    }
};

struct Point {
    ll x, y;
};

struct Cell {
    int gx, gy;
    vector<int> remaining, frontier, next;
    array<int, 24> neighbor;
};

struct Arc {
    ll x, y;
    int id;
    int start;
};

class Envelope {
    static constexpr int MIN_T = -1000000000;
    static constexpr int MAX_T = 1000000000;

    ll d, d2;
    vector<Arc> hull;

    bool better(const Arc& a, const Arc& b, ll t) const {
        if (t < a.y - d) return false;
        if (t > b.y + d) return true;

        ll da = t - a.y;
        ll db = t - b.y;
        ll va = d2 - da * da;
        ll vb = d2 - db * db;
        ll delta = a.x - b.x;
        ll delta2 = delta * delta;

        if (delta >= 0) {
            if (va >= vb) return true;
            ll z = vb - va - delta2;
            if (z <= 0) return true;
            return (i128)4 * delta2 * va >= (i128)z * z;
        } else {
            if (va < vb) return false;
            ll z = va - vb - delta2;
            if (z < 0) return false;
            return (i128)z * z >= (i128)4 * delta2 * vb;
        }
    }

public:
    Envelope(ll distance) : d(distance), d2(distance * distance) {}

    void build(const vector<int>& ids, const vector<Point>& points, int direction) {
        vector<Arc> arcs;
        arcs.reserve(ids.size());

        for (int id : ids) {
            const auto& p = points[id];
            ll x, y;
            if (direction == 0) {
                x = p.x;
                y = p.y;
            } else if (direction == 1) {
                x = -p.x;
                y = p.y;
            } else if (direction == 2) {
                x = p.y;
                y = p.x;
            } else {
                x = -p.y;
                y = p.x;
            }
            arcs.push_back({x, y, id, MIN_T});
        }

        sort(arcs.begin(), arcs.end(), [](const Arc& a, const Arc& b) {
            if (a.y != b.y) return a.y < b.y;
            return a.x > b.x;
        });

        hull.reserve(arcs.size());

        for (size_t i = 0; i < arcs.size();) {
            Arc a = arcs[i];
            size_t j = i + 1;
            while (j < arcs.size() && arcs[j].y == a.y) ++j;
            i = j;

            while (!hull.empty() && better(a, hull.back(), hull.back().start)) {
                hull.pop_back();
            }

            if (hull.empty()) {
                a.start = MIN_T;
                hull.push_back(a);
                continue;
            }

            if (!better(a, hull.back(), MAX_T)) continue;

            int lo = hull.back().start + 1;
            int hi = MAX_T;
            while (lo < hi) {
                int mid = lo + (hi - lo) / 2;
                if (better(a, hull.back(), mid)) hi = mid;
                else lo = mid + 1;
            }

            a.start = lo;
            hull.push_back(a);
        }
    }

    int query(ll t) const {
        int lo = 0, hi = (int)hull.size();
        while (lo + 1 < hi) {
            int mid = (lo + hi) / 2;
            if (hull[mid].start <= t) lo = mid;
            else hi = mid;
        }
        return hull[lo].id;
    }
};

int main() {
    ios::sync_with_stdio(false);
    cin.tie(nullptr);

    int N;
    ll W, D;
    cin >> N >> W >> D;
    const ll D2 = D * D;

    vector<Point> points(N);
    vector<Cell> cells;
    cells.reserve(N);

    unordered_map<ull, int, Hash> index;
    index.max_load_factor(0.7f);
    index.reserve(N);

    auto key = [](int x, int y) -> ull {
        return (ull)(uint32_t)x << 32 | (uint32_t)y;
    };

    auto gridCoordinate = [&](ll x) -> int {
        ll v = 2 * x;
        if (v >= 0) return (int)(v / D);
        return (int)(-((-v + D - 1) / D));
    };

    for (int i = 0; i < N; ++i) {
        cin >> points[i].x >> points[i].y;
        int gx = gridCoordinate(points[i].x);
        int gy = gridCoordinate(points[i].y);
        ull k = key(gx, gy);

        auto [it, inserted] = index.emplace(k, (int)cells.size());
        if (inserted) {
            cells.emplace_back();
            cells.back().gx = gx;
            cells.back().gy = gy;
        }
        cells[it->second].remaining.push_back(i);
    }

    auto reachesGoal = [&](int id) {
        ll x = points[id].x;
        ll y = W - points[id].y;
        return x * x + y * y <= D2;
    };

    vector<int> active, nextActive;
    active.reserve(cells.size());
    nextActive.reserve(cells.size());

    for (int c = 0; c < (int)cells.size(); ++c) {
        auto& cell = cells[c];
        auto& rem = cell.remaining;
        size_t i = 0;
        while (i < rem.size()) {
            int id = rem[i];
            ll x = points[id].x, y = points[id].y;
            if (x * x + y * y <= D2) {
                if (reachesGoal(id)) {
                    cout << 2 << '\n';
                    return 0;
                }
                cell.frontier.push_back(id);
                rem[i] = rem.back();
                rem.pop_back();
            } else {
                ++i;
            }
        }
        if (!cell.frontier.empty()) active.push_back(c);
    }

    if (active.empty()) {
        cout << -1 << '\n';
        return 0;
    }

    array<int, 24> dxs, dys, directions;
    int count = 0;
    for (int dx = -2; dx <= 2; ++dx) {
        for (int dy = -2; dy <= 2; ++dy) {
            if (dx == 0 && dy == 0) continue;
            dxs[count] = dx;
            dys[count] = dy;
            directions[count] = dx > 0 ? 0 : dx < 0 ? 1 : dy > 0 ? 2 : 3;
            ++count;
        }
    }

    for (auto& cell : cells) {
        for (int j = 0; j < 24; ++j) {
            auto it = index.find(key(cell.gx + dxs[j], cell.gy + dys[j]));
            cell.neighbor[j] = it == index.end() ? -1 : it->second;
        }
    }

    int depth = 1;

    auto enqueue = [&](int c, int id) {
        if (cells[c].next.empty()) nextActive.push_back(c);
        cells[c].next.push_back(id);
    };

    while (!active.empty()) {
        for (int c : active) {
            auto& cell = cells[c];

            for (int id : cell.remaining) {
                if (reachesGoal(id)) {
                    cout << depth + 2 << '\n';
                    return 0;
                }
                enqueue(c, id);
            }
            cell.remaining.clear();

            array<unique_ptr<Envelope>, 4> envelopes;

            for (int j = 0; j < 24; ++j) {
                int target = cell.neighbor[j];
                if (target == -1 || cells[target].remaining.empty()) continue;

                int direction = directions[j];
                if (!envelopes[direction]) {
                    envelopes[direction] = make_unique<Envelope>(D);
                    envelopes[direction]->build(cell.frontier, points, direction);
                }

                auto& envelope = *envelopes[direction];
                auto& rem = cells[target].remaining;

                size_t k = 0;
                while (k < rem.size()) {
                    int id = rem[k];
                    ll t = direction < 2 ? points[id].y : points[id].x;
                    int source = envelope.query(t);

                    ll dx = points[id].x - points[source].x;
                    ll dy = points[id].y - points[source].y;

                    if (dx * dx + dy * dy <= D2) {
                        if (reachesGoal(id)) {
                            cout << depth + 2 << '\n';
                            return 0;
                        }
                        enqueue(target, id);
                        rem[k] = rem.back();
                        rem.pop_back();
                    } else {
                        ++k;
                    }
                }
            }

            cell.frontier.clear();
        }

        for (int c : nextActive) {
            cells[c].frontier.swap(cells[c].next);
        }
        active.swap(nextActive);
        nextActive.clear();
        ++depth;
    }

    cout << -1 << '\n';
    return 0;
}

---

この解説は gpt-6-astra-medium によって生成されました。