P13851 [CERC 2023] Going to the Moon 题解

imnotcfz · 2025-08-29 13:49:37 · 题解

0x00 省流

将军饮马，但是小河变成了一个圆。

特别地，经过圆的内部也算。

0x01 分析

显然这是数学题而不是单源最短路。

第一种情况

首先考虑最简单的情况：线段 AB 与圆有交点，或者 A，B 都在圆内。

两点之间线段最短，显然直接连接 AB 就能穿过圆，答案即为 AB 的长度。

第二种情况

如果线段 AB 与圆没有交点怎么办？

如图，以第一个样例为例，我们假设 A 经过圆上一点 P 到达 B，我们要求 AP+BP 的最小值。

作射线 CP，不难发现，根据将军饮马问题，当 CP 是 \angle APB 的角平分线时， AP+BP 最小。

:::align{center} :::

然而这玩意一点也不好求（我推了半天也没推出式子），但是别忘了我们的手里有一个十分强大的计算器，而且只要你看一眼数据范围就会发现，T\le 10^3，说明这个问题不一定需要 O(1) 解决。

考虑枚举点 P，由于 P 是圆上一点，我们设一个角 \theta（0\le\theta\le2\pi）：

:::align{center}

:::

运用三角函数可以算出 P 的坐标，这里记作 P_\theta：

设答案为 f(\theta)，则有：

我们只要求出这个函数的最小值即可。

然而题目要求 10^{-6} 的精度，直接枚举显然是不现实的，于是我们可以利用三分法求单峰函数的最值。

但前提是 f 得是个单峰函数，画出函数图象，我们会发现 f 并不满足这一点，怎么办呢？

:::align{center}

:::

可以证明，当线段 AB 与圆不相交时，f 有且仅有一个最大值和一个最小值。

因此，我们可以把它切成单峰函数，只需要平均分成两部分，分别做三分查找就好了。

0x02 实现

Talk is cheap, show me the code!

#include <iostream>
#include <cmath>
using namespace std;

const double eps = 1e-7;
const double PI = acos(-1);  // 注意到 cos(-1) 等于 π
const double TAU = PI * 2;

int T;

// 向量（或点）
// 此 Vector2 非彼 vector
struct Vector2 {
  double x, y;
  // 只用到了这三种运算
  Vector2 operator+ (Vector2 b) { return {x + b.x, y + b.y}; }
  Vector2 operator- (Vector2 b) { return {x - b.x, y - b.y}; }
  Vector2 operator* (double b) { return {x * b, y * b}; }
  // 向量的长度以及长度的平方
  double magnitude() { return sqrt(x * x + y * y); }
  double sqrMagnitude() { return x * x + y * y; }
} A, B;

// 圆
struct Circle {
  Vector2 center;  // 圆心
  double radius;  // 半径
} C;

// 求两点距离
double distance(Vector2 a, Vector2 b) { return (b - a).magnitude(); }

// 线段
struct Segment {
  Vector2 a, b;  // 起讫点（其实没有区别）

  double length() { return distance(a, b); }
};

// 圆与线段是否有交点（其实交点坐标你根本不想知道对吧）
bool intersects(Circle circle, Segment segment) {
  // 提取一些常量，为了使表达式更简洁
  auto O = circle.center;
  auto r = circle.radius;
  auto dir = segment.b - segment.a;
  auto P = segment.a - O;
  // 一元二次方程系数
  double a = dir.sqrMagnitude();
  double b = 2 * (P.x * dir.x + P.y * dir.y);
  double c = P.sqrMagnitude() - r * r;
  // 判别式
  double delta = b * b - 4 * a * c;
  // 无交点
  if (delta < 0) return false;
  // 求解
  delta = sqrt(delta);
  // 注意：是圆与线段，所以还要判断交点在不在线段上
  double t = (-b + delta) / (2 * a);
  return t >= 0 && t <= 1;
}

// 待求解的 f 函数
double f(double theta) {
  Vector2 P = C.center + Vector2{cos(theta), sin(theta)} * C.radius;
  return distance(A, P) + distance(B, P);
}

// 三分法求解 f 函数区间最小值
double ternarySearch(double l, double r) {
  while (r - l > eps) {
    double lmid = l + (r - l) / 3;
    double rmid = r - (r - l) / 3;
    if (f(lmid) < f(rmid)) r = rmid;
    else l = lmid;
  }
  return f(l);
}

int main() {
  cin >> T;
  cout.flags(ios::fixed);
  cout.precision(8);
  while (T--) {
    cin >> A.x >> A.y >> B.x >> B.y >> C.center.x >> C.center.y >> C.radius;
    Segment AB{A, B};
    // 线段 AB 与圆相交或 A、B 的其中一个在圆内部
    if (distance(A, C.center) <= C.radius
      || distance(B, C.center) <= C.radius
      || intersects(C, AB)) {
      cout << AB.length() << endl;
      continue;
    }
    // 分成两段进行搜索
    cout << min(ternarySearch(0, PI), ternarySearch(PI, TAU)) << endl;
  }
  return 0;
}

这个方法速度也是可以接受的，最慢 11ms 出解。