Problem: lemonade_life (2022 round1)

There's a famous saying about what to do when life gives you lemons. As a traveling lemonade salesman who's never been given any lemons, you sadly can't relate. You must shamefully concoct your lemonade from store-bought powder with your back turned to the world, lest someone see your dirty secret.

Your sales route can be mapped out on a Cartesian plane with \(N\) houses, the \(i\)th of which is at coordinates \((X_i, Y_i)\). Your journey starts at house \(1\), the leftmost house, and ends at house \(N\), the rightmost house. Along the way, you may stop at zero or more other houses to sell lemonade.

You may only stop at a house \(h\) if:

1) standing at house \(h\), there exists some direction you can face in which all other houses are *strictly more behind you than they are in front of you* (formally, if there exists a [half-plane](https://mathworld.wolfram.com/Half-Plane.html#:~:text=A%20half%2Dplane%20is%20a,called%20an%20open%20half%2Dplane.) containing only house \(h\)), and
2) house \(h\) is at most Euclidean distance \(D\) from the previous house you were at.

Your brand image is hurt if you go too long without selling lemonade. The *brand damage* incurred by traveling from one house to another is the larger of \(K\) and the squared Euclidean distance between them. Formally, if your journey consists of \(M\) \((2 \le M \le N)\) houses with the \(i\)th being house \(H_i\) \((H_1 = 1, H_M = N)\), the total brand damage is:

\[\sum_{i=1}^{M-1} \max(K, (X_{H_i} - X_{H_{i + 1}})^2 + (Y_{H_i} - Y_{H_{i + 1}})^2)\]

Is it possible to make the journey? If so, what is the minimum possible total brand damage to do so? Note that the answer may be large, but will fit in a 64-bit integer.


# Constraints

\(1 \le T \le 90\)
\(2 \le N \le 1{,}000{,}000\)
\(0 \le K, D \le 10^9\)
\(0 \le X_i, Y_i \le 1{,}000{,}000\)
\(X_1\) is strictly less than all other \(X_i\).
\(X_N\) is strictly greater than all other \(X_i\).
All \((X_i, Y_i)\) are distinct within a given test case.

The sum of \(N\) across all test cases is at most \(4{,}000{,}000\).
There are at most \(15\) test cases in which \(N > 5{,}000\).


# Input Format

Input begins with a single integer \(T\), the number of test cases. For each test case, there is first a line containing three space-separated integers \(N\), \(K\), and \(D\). Then, \(N\) lines follow, the \(i\)th of which contains two space-separated integers \(X_i\) and \(Y_i\).


# Output Format

For the \(i\)th test case, print a line containing `"Case #i: "` followed a single integer, the minimum total brand damage that must be incurred to make the journey, or \(-1\) if it's impossible to do so.


# Sample Explanation

The first three sample cases are depicted below, with the optimal paths given in blue.

{{PHOTO_ID:654377492455969|WIDTH:700}}

In the first case, going from one house to another takes at least \(K = 25\) brand damage and must not exceed a distance of \(D = 8\). The total brand damage is \(25+50+40 = 115\). Note that you cannot stop at house \((6, 7)\) because there is no direction you could face from there in which your back is at least slightly facing both house \((1,6)\) and \((11,8)\).

In the second case, you can stop at house \((4, 1)\) because for instance, the line \(y = 0.3x - 0.2\) contains only \((4, 1)\), and the half-plane below it contains no other houses.

In the third case, you cannot stop at house \((4, 1)\). There are no other houses within \(D=7\) units of your starting house, so it is not possible to reach the house at \((8, 2)\).


Sample input:
5
9 25 8
0 5
1 6
6 3
6 7
3 4
9 2
2 1
1 2
11 8
3 100 7
0 0
4 1
7 2
3 100 7
0 0
4 1
8 2
6 0 1000000000
0 10
2 5
1 7
7 4
8 1
10 0
12 1600 2000
0 30
16 48
36 57
951 45
397 63
447 63
185 16
362 10
432 9
507 11
643 16
1000 30


Sample output:
Case #1: 115
Case #2: 200
Case #3: -1
Case #4: 56
Case #5: 184654


Solution: Another way to view the half-plane constraint is that you may only stop at points that are vertices of the convex hull of all \(N\) points, which can be computed in \(\mathcal{O}(N \log N)\) time.

There are many incorrect approaches which are broken by the samples. For instance, it is not always optimal to travel across either the upper or lower hull, nor is it obvious how a two pointer solution would work correctly in the general case with the nonlinear cost function involving \(\max()\).

One potential solution is to construct the full graph of all \(N\) nodes with edge weights set to the brand damage between houses, and run Dijkstra's algorithm. This solution is correct, but at first seems to have a time complexity of \(\mathcal{O}(N^2 \log N)\), much too slow for \(N = 1{,}000{,}000\). However, we are only interested in points on the convex hull border — how many vertices can it have?

Let's treat the convex hull edges as vectors lined up tip-to-tail. It turns out that even if these vectors are chosen optimally, one can only fit a convex hull of \(\mathcal{O}(a^{2/3})\) vectors in an \(a \times a\) grid. To see why, consider the vectors in the lower-right quadrant of the hull, sorted radially. In the worst case, the hull will resemble a circle centered at \((500\,000,\,500\,000)\). What is the maximum number of vectors we can generate such that:
- the sum of the vectors' \(x\)-coordinates is at most \(500{,}000\),
- the sum of the vectors' \(y\)-coordinates is at most \(500{,}000\), and
- no vector is a scalar multiple of another?

There aren't that many unique vectors with both small \(x\)- and \(y\)-coordinates. The total number of vectors we can pick is on the order of \((10^6)^{2/3} = 10^4\). In practice, the worst-case convex hull only consists of about \(35{,}000\) vertices. We've included a script, [lemonade_life_count_hull_vertices.py](https://www.dropbox.com/s/nqolapgsm4y6313/lemonade_life_count_hull_vertices.py?dl=0), in the Dropbox for demonstration. This brings the running time of an \(\mathcal{O}(E \log V)\) Dijkstra's implementation to \(35{,}000^2\log_2(35{,}000) \approx 18\) billion steps, which is almost feasible if not for constant factors and the \(\mathcal{O}(E)\) space complexity.

To get the rest of the way there, we can take advantage of the fact that the graph here is complete. On complete graphs, it is more efficient to run Dijkstra's *without a priority queue*. Instead, we can manually search in \(\mathcal{O}(V)\) time to "dequeue" the next node of minimum tentative distance at each step, ridding ourselves of the \(\log(V)\) factor due to heap pushes. We end up with a running time of \(\mathcal{O}(V^2) \approx 1.2\) billion steps, with a low constant factor and \(\mathcal{O}(V)\) space.

With a pessimistic assumption of \(10^8\) steps per second, one would estimate that each of the \(15\) large cases takes under \(12\) seconds, thus comfortably finishing in the \(6\) minute limit. In practice, our official C++ implementation processes the whole file in under \(20\) seconds on a modest laptop.

[See David Harmeyer's solution video here.](https://www.youtube.com/watch?v=D_ZXdemwpgk)


Code: #include <algorithm>
#include <iostream>
#include <utility>
#include <vector>
using namespace std;

const long long INF = 1LL << 60;

using Pt = pair<long long, long long>;
#define x first
#define y second

long long cross(const Pt &a, const Pt &b, const Pt &o = Pt(0, 0)) {
  return (a.x - o.x) * (b.y - o.y) - (a.y - o.y) * (b.x - o.x);
}

template <class It>
vector<Pt> convex_hull(It lo, It hi) {
  int k = 0;
  if (hi - lo <= 1) {
    return vector<Pt>(lo, hi);
  }
  vector<Pt> res(2 * distance(lo, hi));
  sort(lo, hi);
  for (It it = lo; it != hi; ++it) {
    while (k >= 2 && cross(res[k - 1], *it, res[k - 2]) >= 0) {
      k--;
    }
    res[k++] = *it;
  }
  int t = k + 1;
  for (It it = hi - 2; it != lo - 1; --it) {
    while (k >= t && cross(res[k - 1], *it, res[k - 2]) >= 0) {
      k--;
    }
    res[k++] = *it;
  }
  res.resize(k - 1);
  return res;
}

int N;
long long K, D;

long long solve() {
  cin >> N >> K >> D;
  vector<Pt> P(N);
  for (int i = 0; i < N; i++) {
    cin >> P[i].x >> P[i].y;
  }
  auto hull = convex_hull(P.begin(), P.end());
  int V = hull.size();
  // Find start/end.
  int start = -1, dest = -1;
  for (int i = 0; i < V; i++) {
    if (hull[i] == P[0]) {
      start = i;
    } else if (hull[i] == P[N - 1]) {
      dest = i;
    }
  }
  // Dijkstra in O(V^2) time / O(V) space.
  vector<bool> vis(V);
  vector<long long> dist(V, INF);
  dist[start] = 0;
  while (true) {
    int u = -1;
    long long mindist = INF;
    for (int i = 0; i < V; i++) {
      if (!vis[i] && mindist > dist[i]) {
        u = i;
        mindist = dist[i];
      }
    }
    if (mindist == INF || u == dest) {
      break;
    }
    vis[u] = true;
    for (int v = 0; v < V; v++) {
      long long d = (hull[u].x - hull[v].x) * (hull[u].x - hull[v].x) +
                    (hull[u].y - hull[v].y) * (hull[u].y - hull[v].y);
      if (!vis[v] && d <= D*D && dist[v] > dist[u] + max(K, d)) {
        dist[v] = dist[u] + max(K, d);
      }
    }
  }
  return dist[dest] < INF ? dist[dest] : -1;
}

int main() {
  ios_base::sync_with_stdio(false);
  cin.tie(nullptr);
  int T;
  cin >> T;
  for (int t = 1; t <= T; t++) {
    cout << "Case #" << t << ": " << solve() << endl;
  }
  return 0;
}