% Generated by scripts/build_code_appendix.py from mcm26b-核心代码.ipynb. % Do not edit by hand: rerun the generator after the notebook changes. \begingroup\linespread{1.0}\selectfont \subsection*{动作计时与开放路径排序} 三段耗时按整数微秒累计；后备清除的走格顺序用最近邻打底并做 2-opt 改进。 \noindent\small \texttt{CostBreakdown}：一次合法动作的三段耗时：移动 + 换频道（仅检测）+ 动作本身。 \begin{lstlisting}[style=corecode] class CostBreakdown: """一次合法动作的三段耗时 [题面]：移动 + 换频道（仅检测）+ 动作本身。""" move_s: float switch_s: float action_s: float @property def total_s(self) -> float: return self.total_us / 1_000_000 @property def total_us(self) -> int: """逐段量化后按整数微秒累计，与官方虚拟时钟一致。""" return sum(round(value * 1_000_000) for value in (self.move_s, self.switch_s, self.action_s)) \end{lstlisting} \noindent\small \texttt{measure\_cost}：合法 /measure 的耗时；移动按直线距离，换频道只看前后频道是否不同。 \begin{lstlisting}[style=corecode] def measure_cost(prev_pos: Vec2, pos: Vec2, prev_channel: int, channel: int, cost: CostConfig = CostConfig()) -> CostBreakdown: """合法 /measure 的耗时。移动按直线距离；换频道只看前后频道是否不同。""" return CostBreakdown( move_s=dist(prev_pos, pos) / cost.speed_mps, switch_s=0.0 if channel == prev_channel else cost.switch_s, action_s=cost.measure_s, ) \end{lstlisting} \noindent\small \texttt{clear\_cost}：合法 /clear 的耗时；清除不切换测向机频道，故不含换频项。 \begin{lstlisting}[style=corecode] def clear_cost(prev_pos: Vec2, pos: Vec2, success: bool, cost: CostConfig = CostConfig()) -> CostBreakdown: """合法 /clear 的耗时。""" return CostBreakdown( move_s=dist(prev_pos, pos) / cost.speed_mps, switch_s=0.0, action_s=cost.clear_success_s if success else cost.clear_fail_s, ) \end{lstlisting} \noindent\small \texttt{move\_time}：任意两点间的移动耗时，供规划使用。 \begin{lstlisting}[style=corecode] def move_time(a: Vec2, b: Vec2, cost: CostConfig = CostConfig()) -> float: """任意两点间的移动耗时（规划用）。""" return dist(a, b) / cost.speed_mps \end{lstlisting} \noindent\small \texttt{nearest\_neighbor\_order}：最近邻开放路径，作为路径排序的启发式上界。 \begin{lstlisting}[style=corecode] def nearest_neighbor_order(points: list[Vec2], start: Vec2 = (0.0, 0.0)) -> list[Vec2]: """最近邻开放路径（O(n²)，启发式上界）。""" remaining = list(points) order: list[Vec2] = [] position = start while remaining: nxt = min(remaining, key=lambda q: math.dist(position, q)) remaining.remove(nxt) order.append(nxt) position = nxt return order \end{lstlisting} \noindent\small \texttt{two\_opt\_open}：开放路径 2-opt 改进；终点自由，只反转内部片段。 \begin{lstlisting}[style=corecode] def two_opt_open(order: list[Vec2], start: Vec2 = (0.0, 0.0), max_passes: int = 8) -> list[Vec2]: """开放路径 2-opt（终点自由：只反转内部片段，不做首尾相接）。""" best = list(order) n = len(best) if n < 3: return best for _ in range(max_passes): improved = False for i in range(n - 1): prev = best[i - 1] if i > 0 else start pi = best[i] d_prev_i = math.dist(prev, pi) for j in range(i + 1, n): pj = best[j] nxt = best[j + 1] if j + 1 < n else None d_prev_j = math.dist(prev, pj) removed = d_prev_i added = d_prev_j if nxt is not None: removed += math.dist(pj, nxt) added += math.dist(pi, nxt) if added < removed - 1e-9: best[i:j + 1] = best[i:j + 1][::-1] improved = True pi = best[i] d_prev_i = math.dist(prev, pi) if not improved: break return best \end{lstlisting} \subsection*{问题一：角域交会、直径与最小包围圆} 半平面交求定位区域，顶点对枚举直径，Welzl 算法求最小包围圆，并输出直径圆覆盖判定与清除半径判据。 \noindent\small \texttt{Measurement}：一次示向度观测：测点坐标、示向度与误差半角（题设为 $1^\circ$）。 \begin{lstlisting}[style=corecode] class Measurement: x: float y: float bearing_deg: float half_angle_deg: float = 1.0 \end{lstlisting} \noindent\small \texttt{halfplanes}：一次观测的三个半平面约束：角域两侧边界与前向约束，法向量已单位化。 \begin{lstlisting}[style=corecode] def halfplanes(m: Measurement) -> list[tuple[float, float, float]]: """Return a*x+b*y<=c with unit normals, including the forward…""" theta = m.bearing_deg % 360 lower, upper = radians(theta - m.half_angle_deg), radians(theta + m.half_angle_deg) middle = radians(theta) normals = [(sin(lower), -cos(lower)), (-sin(upper), cos(upper)), (-cos(middle), -sin(middle))] return [(a, b, a * m.x + b * m.y) for a, b in normals] \end{lstlisting} \noindent\small \texttt{convex\_hull}：单调链凸包（逆时针、无重复末点），输入坐标先平移缩放到同一量级。 \begin{lstlisting}[style=corecode] def convex_hull(points: list[Point], eps: float = EPS) -> list[Point]: """Monotone chain, CCW, no repeated final vertex; retain segm…""" ordered = sorted(points) unique = [] for p in ordered: if not unique or distance(p, unique[-1]) > eps: unique.append(p) if len(unique) <= 2: return unique def chain(seq): out = [] for p in seq: while len(out) >= 2: a = (out[-1][0] - out[-2][0], out[-1][1] - out[-2][1]) b = (p[0] - out[-1][0], p[1] - out[-1][1]) # Distance-to-line tolerance, not an unscaled area tolerance. if cross(a, b) > eps * max(hypot(*a), hypot(*b), eps): break out.pop() out.append(p) return out return chain(unique)[:-1] + chain(reversed(unique))[:-1] \end{lstlisting} \noindent\small \texttt{diameter}：顶点对枚举区域直径；有界凸集的最远点必为顶点。 \begin{lstlisting}[style=corecode] def diameter(points: list[Point]): if not points: raise ValueError("空点集没有本程序定义下的直径") a = b = points[0] value = 0.0 for p, q in combinations(points, 2): d = distance(p, q) if d > value: value, a, b = d, p, q return value, (a, b) \end{lstlisting} \noindent\small \texttt{smallest\_enclosing\_circle}：Welzl 随机增量最小包围圆；固定随机种子保证结果可复现。 \begin{lstlisting}[style=corecode] def smallest_enclosing_circle(points: list[Point], eps: float = EPS): """Welzl's randomized incremental algorithm, using an explici…""" if not points: raise ValueError("空集不计算最小包围圆") if not isfinite(eps) or eps <= 0: raise ValueError("eps必须为有限正数") if any(len(p) != 2 or any(isinstance(x, bool) or not isinstance(x, (int, float)) or not isfinite(x) for x in p) for p in points): raise ValueError("点集必须由有限二维坐标组成") origin = points[0] scale = max(distance(origin, p) for p in points) if scale == 0: return tuple(origin), 0.0 local = [((p[0]-origin[0])/scale, (p[1]-origin[1])/scale) for p in points] tolerance = min(eps/scale, EPS) Random(20260910).shuffle(local) # Frame = (remaining prefix length, forced boundary, returning from child). stack = [(len(local), (), False)] circle = None while stack: n, boundary, returning = stack.pop() if n == 0 or len(boundary) == 3: # Forced support points must lie ON the boundary. For three points # use their circumcircle, even when their unconstrained MEC is smaller. circle = _circumcircle(*boundary) if len(boundary) == 3 else None if circle is None and boundary: circle = _enclosing_circle_reference(list(boundary), tolerance) elif not returning: stack.append((n, boundary, True)) stack.append((n-1, boundary, False)) elif not _circle_covers(circle, local[n-1], tolerance): stack.append((n-1, boundary + (local[n-1],), False)) assert circle is not None center = (origin[0]+circle[0][0]*scale, origin[1]+circle[0][1]*scale) # Recompute after transforming back: returned circle contains all input points # to distance() precision, including any roundoff from the large translation. return center, max(distance(center, p) for p in points) \end{lstlisting} \noindent\small \texttt{solve}：问题一主入口：半平面交求定位区域，区分空集、无界与退化情形，再给出直径、直径圆覆盖判定与最小包围圆。 \begin{lstlisting}[style=corecode] def solve(measurements: list[Measurement], eps: float = EPS) -> dict: """Solve Q1's pure angular intersection (no disk/range prior…""" if not isfinite(eps) or eps <= 0: raise ValueError("eps必须为有限正数") if not measurements: return {"status": "unbounded", "vertices": [], "diameter_m": None, "diameter_infinite": True, "recession_direction": [1.0, 0.0], "message": "没有观测，整个平面均可行"} if len(measurements) > 100: raise ValueError("参考实现最多接收100个检测点；更大规模应使用半平面交算法") ox = sum(m.x for m in measurements) / len(measurements) oy = sum(m.y for m in measurements) / len(measurements) scale = max(1.0, max(hypot(m.x - ox, m.y - oy) for m in measurements)) local = [Measurement((m.x - ox) / scale, (m.y - oy) / scale, m.bearing_deg, m.half_angle_deg) for m in measurements] planes = [h for m in local for h in halfplanes(m)] def inside(p): tol = eps * max(1.0, abs(p[0]), abs(p[1])) return all(a * p[0] + b * p[1] <= c + tol for a, b, c in planes) feasible_vertices = [] for (a, b, c), (d, e, f) in combinations(planes, 2): det = a * e - b * d if det == 0.0: continue # 精确平行：无交点（近平行对的交点是合法远端顶点，保留） p = ((c * e - b * f) / det, (a * f - c * d) / det) if isfinite(p[0]) and isfinite(p[1]) and inside(p): feasible_vertices.append(p) if not feasible_vertices: return {"status": "empty", "vertices": [], "diameter_m": None, "diameter_infinite": False, "tolerance_m": eps * scale, "message": "测向角域没有公共点；检查频道、数据及误差界"} # A nontrivial 2D recession cone has a ray tangent to an input halfplane. # Use a much tighter angular threshold than the positional feasibility epsilon. ray = None for a, b, _ in planes: for direction in ((b, -a), (-b, a)): if all(d * direction[0] + e * direction[1] <= 2e-14 for d, e, _ in planes): ray = direction break if ray is not None: break hull = convex_hull(feasible_vertices, eps) def world(p): return [ox + p[0] * scale, oy + p[1] * scale] result = {"status": "unbounded" if ray else ("point" if len(hull) == 1 else "segment" if len(hull) == 2 else "polygon"), "vertices": [world(p) for p in hull], "tolerance_m": eps * scale, "diameter_infinite": ray is not None} if ray is not None: result.update(diameter_m=None, recession_direction=list(ray), message="存在非零衰退方向，区域无界，直径为无穷大") return result d, pair = diameter(hull) midpoint = ((pair[0][0] + pair[1][0]) / 2, (pair[0][1] + pair[1][1]) / 2) distances = [distance(midpoint, p) for p in hull] witness_index = max(range(len(hull)), key=distances.__getitem__) max_distance = distances[witness_index] center, radius = smallest_enclosing_circle(hull, eps) area = abs(sum(cross(hull[i], hull[(i + 1) % len(hull)]) for i in range(len(hull)))) / 2 result.update(diameter_m=d * scale, diameter_pair=[world(p) for p in pair], area_m2=area * scale ** 2, diameter_circle={"center": world(midpoint), "radius_m": d * scale / 2, "max_vertex_distance_m": max_distance * scale, "covers": max_distance <= d / 2 + eps * 4, "coverage_tolerance_m": 4 * eps * scale, "max_excess_m": (max_distance-d/2) * scale, "outside_vertex_indices": [i for i, value in enumerate(distances) if value > d/2 + 4*eps], "witness_vertex": world(hull[witness_index]) if max_distance > d/2 + 4*eps else None}, minimum_enclosing_circle={"center": world(center), "radius_m": radius * scale}, guaranteed_clearable_at_mec_center=radius * scale <= 20 - eps * scale * 4, clearance_certificate={"model": "angular_intersection_only", "clear_radius_m": 20.0, "numerical_margin_m": 4*eps*scale, "radius_with_margin_m": (radius+4*eps)*scale, "remaining_margin_m": 20-(radius+4*eps)*scale, "center": world(center)}) return result \end{lstlisting} \subsection*{问题二：候选集与保证接收的第二检测点} 候选集用圆盘外切多边形与角域半平面裁剪表示；安全盘半径取$1000-\rho$，盘内候选按最坏后验包围半径、同分按移动代价取优。 \noindent\small \texttt{CandidateRegion}：某频道候选集合的保守凸外包多边形（顶点按 CCW；允许退化为 1-2 点）。 \begin{lstlisting}[style=corecode] class CandidateRegion: """某频道候选集合的保守凸外包多边形（顶点按 CCW；允许退化为 1-2 点）。""" @classmethod def disk(cls, center: Vec2, radius_m: float, sides: int, slack_m: float | None = None) -> "CandidateRegion": if slack_m is None: slack_m = radius_m * (1 - math.cos(math.pi / sides)) return cls(circumscribed_polygon(center, radius_m, sides), slack_m) @classmethod def initial(cls, station: Vec2, config: ProblemConfig) -> "CandidateRegion": """首个 direction 后的基础集：Ω ∩ B(station, r_max) 的保守外包。""" field_cfg, sides = config.field, config.search.polygon_sides slack = field_cfg.omega_radius_m * (1 - math.cos(math.pi / sides)) region = cls.disk((0.0, 0.0), field_cfg.omega_radius_m, sides, slack) region.clip_disk(station, field_cfg.r_eff_max_m, sides) return region def clip_disk(self, center: Vec2, radius_m: float, sides: int) -> None: """与圆盘 B(center, radius) 的外切多边形求交（近似误差入账）。""" for a, b, c in polygon_halfplanes(circumscribed_polygon(center, radius_m, sides)): self._clip(a, b, c) self.slack_m += radius_m * (1 - math.cos(math.pi / sides)) def clip_wedge(self, station: Vec2, bearing_deg: float, half_angle_deg: float) -> None: """与测向角域 W(station, bearing±half_angle) 求交（半角须 <90°）。 两个边界半平面与 localization.halfplanes 同口径，0°/360° 跨界由 sin/cos 直接处理，无需分支。 """ if not 0 < half_angle_deg < 90: raise ValueError("角域半角须在 (0, 90) 度") lower = math.radians(bearing_deg - half_angle_deg) upper = math.radians(bearing_deg + half_angle_deg) sx, sy = float(station[0]), float(station[1]) self._clip(math.sin(lower), -math.cos(lower), math.sin(lower) * sx - math.cos(lower) * sy) self._clip(-math.sin(upper), math.cos(upper), -math.sin(upper) * sx + math.cos(upper) * sy) def mec(self) -> tuple[Vec2, float]: """最小包围圆（圆心, 半径）。半径是真实候选集 ρ 的上界。""" if self._mec is None: if not self.polygon: raise InconsistentBelief("候选区域为空") self._mec = smallest_enclosing_circle(self.polygon) return self._mec def witnesses(self, spacing_m: float) -> list[Vec2]: """顶点 + 边界按间距细分的采样点（最坏后验评估用）。""" if len(self.polygon) < 3: return list(self.polygon) points = list(self.polygon) count = len(self.polygon) for index in range(count): p, q = self.polygon[index], self.polygon[(index + 1) % count] length = math.dist(p, q) steps = int(length // spacing_m) for k in range(1, steps + 1): ratio = k / (steps + 1) points.append((p[0] + ratio * (q[0] - p[0]), p[1] + ratio * (q[1] - p[1]))) return _dedupe(points, tol=1e-9) def contains(self, point: Vec2, tol_m: float = 0.05) -> bool: """点是否在候选集内（带容差；用于真源保留性检验，不做删点依据）。""" if not self.polygon: return False if len(self.polygon) == 1: return math.dist(point, self.polygon[0]) <= tol_m + self.slack_m if len(self.polygon) == 2: p, q = self.polygon return _point_segment_distance(point, p, q) <= tol_m + self.slack_m for index, p in enumerate(self.polygon): q = self.polygon[(index + 1) % len(self.polygon)] if (q[0] - p[0]) * (point[1] - p[1]) - (q[1] - p[1]) * (point[0] - p[0]) \ < -tol_m - self.slack_m: return False return True \end{lstlisting} \noindent\small \texttt{ClearCertificate}：清除充分条件：候选集外包半径与近似误差之和不超过清除半径。 \begin{lstlisting}[style=corecode] class ClearCertificate: """清除证书：候选集外包半径（含近似误差）不超过 20m 阈值。""" channel: int center: Vec2 radius_m: float limit_m: float certified: bool \end{lstlisting} \noindent\small \texttt{clear\_certificate}：由候选集最小包围圆判定是否可以就地清除。 \begin{lstlisting}[style=corecode] def clear_certificate(region: CandidateRegion, channel: int, config: ProblemConfig = ProblemConfig()) -> ClearCertificate: """ρ_poly <= 20 - 外包近似误差 - 工程余量 才允许发起清除。""" center, radius = region.mec() limit = (config.field.clear_radius_m - region.slack_m - config.search.clear_margin_m) return ClearCertificate(channel, center, radius, limit, radius <= limit) \end{lstlisting} \noindent\small \texttt{candidate\_points}：保证接收盘内的候选点：盘心加各比例半径上的均匀方向点。 \begin{lstlisting}[style=corecode] def candidate_points(center: Vec2, safe_radius_m: float, rings, angles: int) -> list[Vec2]: """保证接收盘内的候选点：盘心 + 各比例半径 × 均匀方向。""" points = [tuple(center)] for fraction in rings: radius = fraction * safe_radius_m if radius <= 1e-9: continue step = 2 * math.pi / angles for k in range(angles): angle = step * k + step / 2 points.append((center[0] + radius * math.cos(angle), center[1] + radius * math.sin(angle))) return points \end{lstlisting} \noindent\small \texttt{safe\_disk\_radius}：保证接收盘半径 $1000-\rho$；小于零表示没有保证接收的区域。 \begin{lstlisting}[style=corecode] def safe_disk_radius(region: CandidateRegion, config: ProblemConfig = ProblemConfig()) -> float: """保证接收盘 B(mec, safe) 的半径；<0 表示无保证（ρ 超过最小接收半径）。""" return config.field.r_eff_min_m - region.mec()[1] - 1e-3 \end{lstlisting} \noindent\small \texttt{select\_measurement\_point}：在保证接收盘内选下一检测点：主指标为最坏后验包围半径，同分取移动代价较小者。 \begin{lstlisting}[style=corecode] def select_measurement_point(region: CandidateRegion, from_position: Vec2, config: ProblemConfig = ProblemConfig()) -> SecondPointPlan: """在保证接收盘内选择下一检测点（P3 核心；对任意候选集通用）。""" field, search = config.field, config.search center, radius = region.mec() half_angle = field.bearing_error_half_deg + search.bearing_rounding_margin_deg def move_cost(point: Vec2) -> float: return move_time(from_position, point, config.cost) # 1mm 数值余量在 safe_disk_radius 内记账。 safe_radius = safe_disk_radius(region, config) if safe_radius < 0: return SecondPointPlan(tuple(center), radius, move_cost(center), "无保证接收盘（ρ>1000）：探索 MEC 圆心，可能 no_signal") coarse = region.witnesses(search.witness_spacing_m) samples = candidate_points(center, safe_radius, search.candidate_rings, search.candidate_angles) near_radius, sides = field.near_radius_m, search.polygon_sides best_key, best_point, best_worst = None, None, None for point in samples: worst = _worst_posterior_radius(region, point, coarse, half_angle, near_radius, sides) key = (round(worst, 6), round(move_cost(point), 6)) if best_key is None or key < best_key: best_key, best_point, best_worst = key, point, worst dense = _worst_posterior_radius(region, best_point, region.witnesses(search.witness_spacing_m / 4), half_angle, near_radius, sides) reason = (f"安全盘 B(centre,{safe_radius:.0f}m) 内 {len(samples)} 候选；" f"粗采样最坏后验 {best_worst:.1f}m，1/4 间距加密复核 {dense:.1f}m") return SecondPointPlan(tuple(best_point), dense, move_cost(best_point), reason) \end{lstlisting} \noindent\small \texttt{select\_second\_point}：问题二入口：由首次示向度构造候选集并给出第二检测点。 \begin{lstlisting}[style=corecode] def select_second_point(client: SimulatorClient, first_position: Vec2, first_bearing_deg: float, config: ProblemConfig = ProblemConfig()) -> SecondPointPlan: """P3 入口：由首次 direction 构造 K1，返回保证接收的第二检测点方案。""" region = CandidateRegion.initial(first_position, config) half_angle = config.field.bearing_error_half_deg \ + config.search.bearing_rounding_margin_deg region.clip_wedge(first_position, first_bearing_deg, half_angle) return select_measurement_point(region, client.state.position, config) \end{lstlisting} \subsection*{问题三：覆盖计划、检测记录与清除证书} 覆盖计划与运行期抽样校验、逐频道检测记录，以及七点解析覆盖与全向／定向两种覆盖语义的证书判据。 \noindent\small \texttt{CoverageCertificate}：覆盖证书：采样分辨率、位置与朝向目标数、缺口数与最坏命中距离。 \begin{lstlisting}[style=corecode] class CoverageCertificate: """覆盖证书：采样语义 + 缺口 + 最坏命中距离。""" kind: str # "full_disk" | "forward_half_disk" target_step_m: float headings: int r_eff_min_m: float targets: int uncovered: int worst_hit_distance_m: float nn_tour_m: float generator_sha256: str = "" verify_heading_margin_deg: float = 0.0 @property def valid(self) -> bool: return self.uncovered == 0 and self.targets > 0 def as_dict(self) -> dict: return {"kind": self.kind, "target_step_m": self.target_step_m, "headings": self.headings, "r_eff_min_m": self.r_eff_min_m, "targets": self.targets, "uncovered": self.uncovered, "worst_hit_distance_m": round(self.worst_hit_distance_m, 3), "nn_tour_m": round(self.nn_tour_m, 1), "nn_tour_virtual_s": round(self.nn_tour_m / SPEED_MPS, 1), "generator_sha256": self.generator_sha256, "verify_heading_margin_deg": self.verify_heading_margin_deg} \end{lstlisting} \noindent\small \texttt{CoveragePlan}：覆盖计划：点集、覆盖语义与证书。 \begin{lstlisting}[style=corecode] class CoveragePlan: """一份覆盖计划：点集 + 语义 + 证书（可为解析证明）。""" name: str problem: int points: tuple[Vec2, ...] kind: str # "analytic" | "certified" coverage: str # "full_disk" | "forward_half_disk" note: str = "" certificate: CoverageCertificate | None = None generator: dict = field(default_factory=dict) def as_dict(self) -> dict: return {"name": self.name, "problem": self.problem, "kind": self.kind, "coverage": self.coverage, "note": self.note, "generator": self.generator, "certificate": (self.certificate.as_dict() if self.certificate is not None else None), "points": [[p[0], p[1]] for p in self.points]} \end{lstlisting} \noindent\small \texttt{coverage\_points}：49 点覆盖网格，单元对角 $700\sqrt2\approx989.95\le1000$ m。 \begin{lstlisting}[style=corecode] def coverage_points(config: ProblemConfig = ProblemConfig()) -> list[Vec2]: """49 点覆盖网格生成（纯几何，可独立测试）。""" step = config.search.grid_step_m half = config.search.grid_half_cells return [(step * i, step * j) for j in range(half, -half - 1, -1) for i in range(-half, half + 1)] \end{lstlisting} \noindent\small \texttt{CoverageLedger}：逐频道检测记录，按（频道, 网格点）二维记账。 \begin{lstlisting}[style=corecode] class CoverageLedger: """逐频道覆盖账本：(频道, 网格点) -> 是否已测。""" config: ProblemConfig = field(default_factory=ProblemConfig) points: list[Vec2] = field(default_factory=coverage_points) measured: set = field(default_factory=set) # {(channel, point_index)} skipped: dict = field(default_factory=dict) # channel -> 原因 def mark_measured(self, channel: int, point_index: int) -> None: self.measured.add((channel, point_index)) def measured_at(self, channel: int, point_index: int) -> bool: return (channel, point_index) in self.measured def skip_channel(self, channel: int, reason: str) -> None: """记录整频道跳过原因（如 already_cleared / located_pending_clear）。""" self.skipped[channel] = reason def pending_points(self, channel: int) -> list[int]: """该频道尚未测量的网格点下标。扫描调度策略未实现（P4）。""" return [index for index in range(len(self.points)) if not self.measured_at(channel, index)] def coverage_complete(self, channel: int) -> bool: """该频道覆盖是否完成（含显式跳过）。""" return channel in self.skipped or not self.pending_points(channel) \end{lstlisting} \noindent\small \texttt{seven\_points}：原点与环上六个等分点的七站覆盖集。 \begin{lstlisting}[style=corecode] def seven_points(ring_radius_m: float = RING_RADIUS_M) -> list[tuple[float, float]]: """原点 + 环上 6 等分点的七站覆盖集。""" return [(0.0, 0.0)] + [(ring_radius_m * math.cos(k * math.pi / 3), ring_radius_m * math.sin(k * math.pi / 3)) for k in range(6)] \end{lstlisting} \noindent\small \texttt{seven\_plan}：第三问解析覆盖计划；环半径 $1159.7$ m 时最坏覆盖距离不超过 $984.539$ m。 \begin{lstlisting}[style=corecode] def seven_plan(ring_radius_m: float) -> CoveragePlan: """第三问解析计划（覆盖证明由 validate_seven_coverage 给出）。""" return CoveragePlan( name=f"seven-r{ring_radius_m:.1f}", problem=3, kind="analytic", coverage="full_disk", note=("解析证明 max(r/√3, sqrt(R²+r²−√3Rr)) < 1000；" "r∈[1122.9, 1732] 时成立，巡检半径越小巡回越短"), points=tuple(seven_points(ring_radius_m))) \end{lstlisting} \noindent\small \texttt{verify\_full\_disk}：第三问覆盖语义：圆域内每点都落在某个被访点的最小接收圆盘内。 \begin{lstlisting}[style=corecode] def verify_full_disk(points: list[Vec2], step_m: float = VERIFY_STEP_M, r_eff: float = R_EFF_MIN_M) -> CoverageCertificate: """第三问语义：Ω 内每点都在某个被访点的 r_eff 圆盘内。""" points = list(points) positions = omega_positions(step_m) table = _bucket(points, r_eff) missing, worst = 0, 0.0 for (px, py) in positions: ci, cj = int(px // r_eff), int(py // r_eff) best = None for di in (-1, 0, 1): for dj in (-1, 0, 1): for c_index in table.get((ci + di, cj + dj), ()): cx, cy = points[c_index] distance = math.hypot(cx - px, cy - py) if distance <= r_eff and (best is None or distance < best): best = distance if best is None: missing += 1 else: worst = max(worst, best) return CoverageCertificate("full_disk", step_m, 1, r_eff, len(positions), missing, worst, nn_tour_length(points)) \end{lstlisting} \noindent\small \texttt{verify\_forward\_half\_disk}：第四问覆盖语义：任意位置与朝向都存在前向半盘内的被访点。 \begin{lstlisting}[style=corecode] def verify_forward_half_disk(points: list[Vec2], step_m: float = VERIFY_STEP_M, headings: int = VERIFY_HEADINGS, r_eff: float = R_EFF_MIN_M, heading_margin_deg: float = 0.0) -> CoverageCertificate: """第四问语义：∀(p,h) 存在前向半盘内的被访点。""" points = list(points) targets = forward_targets(step_m, headings) table = _bucket(points, r_eff) missing, worst = 0, 0.0 sin_margin = math.sin(math.radians(max(heading_margin_deg, 0.0))) for (px, py, heading) in targets: ux, uy = math.cos(heading), math.sin(heading) ci, cj = int(px // r_eff), int(py // r_eff) best = None for di in (-1, 0, 1): for dj in (-1, 0, 1): for c_index in table.get((ci + di, cj + dj), ()): cx, cy = points[c_index] dx, dy = cx - px, cy - py distance_sq = dx * dx + dy * dy if distance_sq > r_eff * r_eff: continue if dx * ux + dy * uy < sin_margin * math.sqrt(distance_sq) - 1e-9: continue distance = math.sqrt(distance_sq) if best is None or distance < best: best = distance if best is None: missing += 1 else: worst = max(worst, best) return CoverageCertificate("forward_half_disk", step_m, headings, r_eff, len(targets), missing, worst, nn_tour_length(points)) \end{lstlisting} \noindent\small \texttt{plan\_instance}：按名称取覆盖计划；证书计划在运行期按同一分辨率复核。 \begin{lstlisting}[style=corecode] def plan_instance(name: str, verify: bool = True) -> CoveragePlan: """统一入口：按名称取覆盖计划（内置解析计划 / 内嵌证书计划）。""" if name == "square49": return square_49_plan() if name == "seven": return seven_plan(1200.0) if name not in EMBEDDED_PLANS: raise FileNotFoundError(f"未找到覆盖计划 {name!r}；现有：{plan_names()}") plan = embedded_plan(name) if verify and plan.kind == "certified": certificate = verify_plan(plan) if not certificate.valid: raise ValueError( f"覆盖计划 {name!r} 运行期校验失败：" f"{certificate.uncovered}/{certificate.targets} 个 (p,h) 未被覆盖") plan = CoveragePlan( name=plan.name, problem=plan.problem, kind=plan.kind, coverage=plan.coverage, note=plan.note, generator=plan.generator, certificate=CoverageCertificate( kind=certificate.kind, target_step_m=certificate.target_step_m, headings=certificate.headings, r_eff_min_m=certificate.r_eff_min_m, targets=certificate.targets, uncovered=certificate.uncovered, worst_hit_distance_m=certificate.worst_hit_distance_m, nn_tour_m=certificate.nn_tour_m, generator_sha256=(plan.certificate.generator_sha256 if plan.certificate else "")), points=plan.points) return plan \end{lstlisting} \subsection*{问题三：网格候选集} 逐频道的网格幸存集，按格心裕量 $\ell/\sqrt2$ 内缩后再删格，据此支持负观测更新、跳测规则与缺席证明。 \noindent\small \texttt{ChannelParticles}：单个频道的网格候选集：格心裕量 $\ell/\sqrt2$，所有删除判定内缩该裕量，故真源所在格永不被删。 \begin{lstlisting}[style=corecode] class ChannelParticles: """单个频道的网格幸存集。""" def cell_center(self, index: int) -> Vec2: i, j = divmod(index, self.width) return (self.x0 + (i + 0.5) * self.spacing_m, self.x0 + (j + 0.5) * self.spacing_m) def discard(self, indices: list[int]) -> None: """公开删除入口：联合信念把"朝向不可能"的格一并删除时使用。""" self._discard(indices) def _discard(self, indices: list[int]) -> None: if not indices: return coarse_w = self._coarse_width() for index in indices: if index not in self.alive: continue self.alive.discard(index) key = (index // self.width // _COARSE_FACTOR) * coarse_w \ + (index % self.width) // _COARSE_FACTOR count = self._coarse.get(key, 0) - 1 if count <= 0: self._coarse.pop(key, None) else: self._coarse[key] = count self.version += 1 self._mec = None def _kill_disk(self, center: Vec2, radius_m: float) -> None: """删除格心距 center <= radius_m 的全部幸存格（整格判据见模块注释）。""" if not self.alive: return reach = radius_m + _DIST_EPS_M ilo = max(0, int(math.floor((center[0] - reach - self.x0) / self.spacing_m))) ihi = min(self.width - 1, int(math.ceil((center[0] + reach - self.x0) / self.spacing_m))) jlo = max(0, int(math.floor((center[1] - reach - self.x0) / self.spacing_m))) jhi = min(self.width - 1, int(math.ceil((center[1] + reach - self.x0) / self.spacing_m))) limit_sq = reach * reach doomed = [] for i in range(ilo, ihi + 1): dx = self.x0 + (i + 0.5) * self.spacing_m - center[0] base = i * self.width for j in range(jlo, jhi + 1): index = base + j if index not in self.alive: continue dy = self.x0 + (j + 0.5) * self.spacing_m - center[1] if dx * dx + dy * dy <= limit_sq: doomed.append(index) self._discard(doomed) def _kill_outside_disk(self, center: Vec2, radius_m: float) -> None: """删除格心距 center > radius_m 的全部幸存格（保留盘内）。""" if not self.alive: return reach = radius_m + _DIST_EPS_M doomed = [index for index in self.alive if math.dist(self.cell_center(index), center) > reach] self._discard(doomed) def _kill_outside_wedge(self, station: Vec2, reading_deg: float) -> None: """删除角向偏离读数超出 δ_eff + 角向格裕量的幸存格。 角向格裕量 = asin(m / max(dist−m, 1))：格内任一点相对格心的方位 偏差不超过它（点距测点还可再近 m，故分母再内缩 m）。 """ if not self.alive: return sx, sy = float(station[0]), float(station[1]) doomed = [] for index in self.alive: gx, gy = self.cell_center(index) dist = math.hypot(gx - sx, gy - sy) if dist <= self.margin_m + 1.0: continue # 格贴测点：格内方位可指向任意方向，不可删 angular_margin = math.degrees( math.asin(min(0.9999, self.margin_m / (dist - self.margin_m)))) gap = angular_gap(math.degrees(math.atan2(gy - sy, gx - sx)), reading_deg) if gap > self.half_angle_deg + angular_margin + _ANG_EPS_DEG: doomed.append(index) self._discard(doomed) def apply_measure(self, station: Vec2, result: str, svd_deg: float | None = None, *, directional: bool = False) -> None: """按读数更新幸存集。真源保留性由 Voronoi 裕量保证（模块注释）。 ``directional=True``（第四问）：no_signal **不删空间候选**——定向源 的静默还可能来自朝向遮挡，单独删空间候选会误删真源；朝向维度由 :class:`~mcm26b.belief.joint.JointBelief` 维护。 """ if result == NO_SIGNAL: if directional: return # no_signal ⟹ dist > r_eff >= 1000：整格落入 B(s, 1000−m) 的删。 self._kill_disk(station, self.r_eff_min_m - self.margin_m) elif result == NEAR: self._kill_outside_disk(station, self.near_radius_m + self.margin_m) elif result == DIRECTION: if svd_deg is None: raise ValueError("direction 读数必须携带 svd_deg") # direction 隐含 dist <= r_max [题面] self._kill_outside_disk(station, self.r_eff_max_m + self.margin_m) self._kill_outside_wedge(station, svd_deg) else: raise ValueError(f"未知测向结果: {result!r}") def alive_points(self, limit: int | None = None) -> list[Vec2]: sorted_alive = sorted(self.alive) if limit is not None and len(sorted_alive) > limit: step = len(sorted_alive) / limit indices = {sorted_alive[int(k * step) % len(sorted_alive)] for k in range(limit)} indices = sorted(indices) else: indices = sorted_alive return [self.cell_center(index) for index in indices] def nearest_alive_distance(self, point: Vec2, limit_m: float) -> float | None: """最近幸存格心距离；> limit_m 时返回 None（粗格环扫 + 下界剪枝）。 返回值语义：恰为最近幸存格心距离（若 <= limit_m）；None 表示所有 幸存格心都超出 limit_m。粗格聚合半径已计入扫描范围，不漏判。 """ if not self.alive: return None coarse_cell_m = self.spacing_m * _COARSE_FACTOR coarse_reach = coarse_cell_m * math.sqrt(2) / 2 cw = self._coarse_width() ci_base = int(math.floor((point[0] - self.x0) / coarse_cell_m)) cj_base = int(math.floor((point[1] - self.x0) / coarse_cell_m)) ring_max = int(math.ceil((limit_m + coarse_reach) / coarse_cell_m)) + 1 best: float | None = None for ring in range(ring_max + 1): lower_bound = max(0.0, ring - 1) * coarse_cell_m - coarse_reach if best is not None and lower_bound >= best: break # 更外环的所有格心都 >= 下界 >= 当前最优 for ci in range(max(0, ci_base - ring), min(cw - 1, ci_base + ring) + 1): for cj in range(max(0, cj_base - ring), min(cw - 1, cj_base + ring) + 1): if max(abs(ci - ci_base), abs(cj - cj_base)) != ring: continue # 只扫当前环（Chebyshev 环） if (ci * cw + cj) not in self._coarse: continue ccx = self.x0 + (ci + 0.5) * coarse_cell_m ccy = self.x0 + (cj + 0.5) * coarse_cell_m if math.hypot(ccx - point[0], ccy - point[1]) \ > limit_m + coarse_reach: continue i0, j0 = ci * _COARSE_FACTOR, cj * _COARSE_FACTOR for i in range(i0, min(i0 + _COARSE_FACTOR, self.width)): base = i * self.width for j in range(j0, min(j0 + _COARSE_FACTOR, self.width)): index = base + j if index not in self.alive: continue dist = math.dist(self.cell_center(index), point) if dist <= limit_m and (best is None or dist < best): best = dist return best def guaranteed_no_signal(self, station: Vec2) -> bool: """真源与该测点距离必 > r_max（>= r_eff）⟹ 必 no_signal，测了也无信息。""" return self.nearest_alive_distance( station, self.r_eff_max_m + self.margin_m) is None def mec(self) -> tuple[Vec2, float]: """幸存集最小包围圆（圆心, 半径）；按 version 缓存。""" if self._mec is None: if not self.alive: raise ValueError("幸存集为空（该频道已可证明缺席）") points = [self.cell_center(index) for index in self.alive] self._mec = smallest_enclosing_circle(points) return self._mec def subsample(self, size: int) -> list[Vec2]: """按索引序均匀抽取的子采样（索引序≈空间序，代表点空间铺开）。""" return self.alive_points(limit=size) def expected_survival(self, point: Vec2, subsample, representatives) -> float: """E[log 幸存占比]：越小该测点的期望信息增益越大。 r_eff 不确定性解析积分：dist∈(1000,1500] 时按 r_eff~U(1000,1500) 取 p_ns=(dist−1000)/500 为 no_signal 分支概率。 - no_signal 分支：删 B(point, 1000−m)，幸存 = 子采样中盘外比例； - direction 分支：以代表点方位为期望读数，角域 ∩ B(point,1500+m) 内幸存； - near 分支（dist <= 5+m）：可行集收缩到 5m 盘，幸存占比趋近 0。 """ if not self.alive or len(subsample) == 0 or len(representatives) == 0: return 0.0 if _np is not None and isinstance(subsample, _np.ndarray) \ and isinstance(representatives, _np.ndarray): return self._expected_survival_vector(point, subsample, representatives) eps = 1.0 / len(subsample) keep_ns = max(1.0 - self.fraction_within(subsample, point, self.r_eff_min_m - self.margin_m), eps) log_keep_ns = math.log(keep_ns) total = 0.0 for rep in representatives: dist = math.hypot(rep[0] - point[0], rep[1] - point[1]) if dist <= self.near_radius_m + self.margin_m: total += math.log(eps) continue p_ns = min(1.0, max(0.0, (dist - self.r_eff_min_m) / (self.r_eff_max_m - self.r_eff_min_m))) reading = math.atan2(rep[1] - point[1], rep[0] - point[0]) survivors = 0 for p in subsample: if math.dist(p, point) > self.r_eff_max_m + self.margin_m: continue gap = angular_gap(math.degrees( math.atan2(p[1] - point[1], p[0] - point[0])), math.degrees(reading)) if gap <= self.half_angle_deg + _ANG_EPS_DEG: survivors += 1 keep_dir = max(survivors / len(subsample), eps) total += p_ns * log_keep_ns + (1.0 - p_ns) * math.log(keep_dir) return total / max(len(representatives), 1) def fraction_within(self, samples, point: Vec2, radius_m: float) -> float: if len(samples) == 0: return 0.0 reach = radius_m + _DIST_EPS_M if _np is not None and isinstance(samples, _np.ndarray): distance = _np.hypot(samples[:, 0] - point[0], samples[:, 1] - point[1]) return float(_np.count_nonzero(distance <= reach)) / len(samples) count = sum(1 for p in samples if math.dist(p, point) <= reach) return count / len(samples) def expected_posterior_radius(self, point: Vec2, subsample, representatives) -> float: """一次测量后期望后验半径（米）：选点主指标，直指清除证书阈值。 对每个代表假设（真值候选）按 r_eff~U(1000,1500) 解析积分 no_signal/direction 分支，near 分支（dist <= 5+m）后验为 5m 盘； 各分支半径用幸存子采样的质心最大距离近似（偏差 << 格距，排序用）。 与几何版 worst-posterior 同族，但作用于含负观测的非凸幸存集， 且以期望代替最坏（代表点覆盖幸存集主轴）。 """ if not self.alive or len(subsample) == 0 or len(representatives) == 0: return 0.0 if _np is not None and isinstance(subsample, _np.ndarray) \ and isinstance(representatives, _np.ndarray): return self._expected_posterior_radius_vector( point, subsample, representatives) eps_radius = self.spacing_m / 2.0 total = 0.0 for rep in representatives: dist = math.hypot(rep[0] - point[0], rep[1] - point[1]) if dist <= self.near_radius_m + self.margin_m: total += self.near_radius_m + self.margin_m continue p_ns = min(1.0, max(0.0, (dist - self.r_eff_min_m) / (self.r_eff_max_m - self.r_eff_min_m))) members_ns = [p for p in subsample if math.dist(p, point) > self.r_eff_min_m - self.margin_m] reading = math.atan2(rep[1] - point[1], rep[0] - point[0]) members_dir = [] for p in subsample: if math.dist(p, point) > self.r_eff_max_m + self.margin_m: continue gap = angular_gap(math.degrees( math.atan2(p[1] - point[1], p[0] - point[0])), math.degrees(reading)) if gap <= self.half_angle_deg + _ANG_EPS_DEG: members_dir.append(p) r_ns = self._branch_radius(members_ns) if members_ns else eps_radius r_dir = self._branch_radius(members_dir) if members_dir else eps_radius total += p_ns * r_ns + (1.0 - p_ns) * r_dir return total / max(len(representatives), 1) \end{lstlisting} \subsection*{问题三：信息增益选点与滚动任务排序} 统一任务表由站点扫频、定位测量与清除三类任务组成，每执行一项即重排；定位点按期望信息增益加移动代价评分。 \noindent\small \texttt{BayesLoop}：问题三闭环：网格候选集的正负观测更新、信息增益选点与滚动任务排序。 \begin{lstlisting}[style=corecode] class BayesLoop(OmnidirectionalLoop): """问题3贝叶斯闭环：负观测感知 + 信息增益选点 + 滚动任务排序。""" def _build_tasks(self) -> list[Task]: tasks: list[Task] = [] for index, point in enumerate(self.ledger.points): pending = self._fused_channels_at(index) if not pending: continue switches = len(pending) - (1 if self.client.state.channel in pending else 0) tasks.append(Task("station", tuple(point), self.config.cost.measure_s * len(pending) + self.config.cost.switch_s * max(switches, 0), station_index=index, station_channels=tuple(pending))) for channel in self.belief.channels_by_status(ChannelStatus.DETECTED): cert = self._best_certificate(channel) if cert is not None and cert.certified \ and channel not in self._clear_failed: tasks.append(Task("clear", cert.center, self.config.cost.clear_success_s, channel=channel)) continue if self._defer_to_station(channel): continue # 下一条方位线可由"本就要去的站点"顺路提供 point = self._select_point(channel) if point is None: continue tasks.append(Task("locate", point, self.config.cost.measure_s + self.config.cost.switch_s, channel=channel)) return tasks def _execute_station(self, index: int) -> None: """单站批量扫频：只测待测频道（跳测规则已过滤），扫完即返回重排。 P0-1：本站同时测量"已发现未定证"频道，使第二/第三条方位线在巡回 途中顺路取得；证书一旦成立即就地清除，不额外绕行。 """ point = self.ledger.points[index] pending = self._fused_channels_at(index) unknown = [c for c, b in self.belief.channels.items() if b.status is ChannelStatus.UNKNOWN] self.skipped_measures += len(unknown) - len( [c for c in pending if self.belief.channels[c].status is ChannelStatus.UNKNOWN]) # 当前频道优先，减少换频；同站 near 直接清除不打断扫频 pending.sort(key=lambda c: (c != self.client.state.channel, c)) for channel in pending: if self.stop_reason is not None or self._count_limit_reached() \ or self._virtual_cap_reached(): return detected_before = (self.belief.channels[channel].status is ChannelStatus.DETECTED) outcome = self.client.measure(point, channel) self.ledger.mark_measured(channel, index) self._record_measure(channel, point, outcome) if outcome.result in (DIRECTION, NEAR): self._may_clear.add(channel) if outcome.result == NEAR: region = self.belief.candidate_set(channel) cert = clear_certificate(region, channel, self.config) if cert.certified: if self._attempt_clear(channel, point, cert): self._clear_failed.discard(channel) else: self._clear_failed.add(channel) self._may_clear.discard(channel) else: self._anomaly("near_without_certificate", channel, position=point) continue if not detected_before: continue if outcome.result == NO_SIGNAL: # 该频道已被发现：保证接收盘内不可能 no_signal，属守恒性异常 region = self.belief.candidate_set(channel) safe = (safe_disk_radius(region, self.config) if region is not None else -1.0) if safe >= 0 and math.dist(point, region.mec()[0]) <= safe: self._anomaly("no_signal_in_safe_disk", channel, position=point) continue cert = self._best_certificate(channel) if cert is not None and cert.certified: if self._attempt_clear(channel, cert.center, cert): self._clear_failed.discard(channel) else: self._clear_failed.add(channel) self._may_clear.discard(channel) def _record_measure(self, channel: int, position: Vec2, outcome) -> None: """belief（凸证书）与 PF（幸存集）同步更新；PF 空即证缺席。""" self.belief.update_measure(channel, position, outcome.result, outcome.svd_deg) # 同一 0.5m 测点 × 同一频道的重复检测：历史上一次官方演练有 98.84% # 的虚拟时间耗在空复扫上，本计数是该问题的回归指标。 key = _point_key(position) seen = self._measured_points.setdefault(channel, set()) if key in seen: self.redundant_rescan_measures += 1 seen.add(key) self._clear_failed.discard(channel) # 新证据落地：允许再次尝试清除 belief = self.belief.channels[channel] if belief.status in (ChannelStatus.CLEARED, ChannelStatus.ABSENT): return # 已证毕：清除/缺席后的 no_signal 与距离无关，不进 PF pf = self.particles[channel] pf.apply_measure(position, outcome.result, outcome.svd_deg) if not pf: if belief.status is ChannelStatus.UNKNOWN: self.belief.mark_absent(channel) self.absence_evidence[channel] = "pf_empty" if channel not in self.pf_absent_channels: self.pf_absent_channels.append(channel) else: self._anomaly("pf_empty_while_detected", channel, position=position) self._apply_count_bound() def _best_certificate(self, channel: int) -> ClearCertificate | None: """凸证书优先；负观测挖洞使 PF 幸存集更紧时用网格证书补位。""" region = self.belief.candidate_set(channel) cert = clear_certificate(region, channel, self.config) \ if region is not None else None if cert is not None and cert.certified: return self._near_clear_target(cert) if self.config.bayes.grid_cert_enabled: pf = self.particles[channel] if len(pf): center, radius = pf.mec() limit = (self.config.field.clear_radius_m - pf.margin_m - self.config.search.clear_margin_m) if radius <= limit: return self._near_clear_target( ClearCertificate(channel, center, radius, limit, True)) return cert def _select_point(self, channel: int) -> Vec2 | None: """定位测量点：期望 log 幸存占比 + λ·移动时间 最小的候选。""" context = self._probe_context(channel) if context is None: return None candidates, info_score, center, position = context bayes = self.config.bayes def score(q: Vec2) -> float: # 一次定位后仍需接近源才能清除；仅计到检测点的移动会偏好绕行。 move = bayes.lambda_move_s * (move_time(position, q, self.config.cost) + move_time(q, center, self.config.cost)) return info_score(q) + move return min(candidates, key=lambda q: (score(q), _point_key(q))) def _localize_and_clear(self, channel: int) -> bool: """阻塞式定位子循环（resolve 阶段安全网）：PF 选点，失败退凸选点。""" search = self.config.search self._may_clear.add(channel) for _ in range(search.max_localize_steps): if self._virtual_cap_reached(): return False cert = self._best_certificate(channel) if cert is not None and cert.certified and channel in self._may_clear: if self._attempt_clear(channel, cert.center, cert): self._clear_failed.discard(channel) return True self._clear_failed.add(channel) self._may_clear.discard(channel) continue point = self._select_point(channel) if point is None: region = self.belief.candidate_set(channel) if region is None: return False point = select_measurement_point( region, self.client.state.position, self.config).position outcome = self.client.measure(point, channel) self._record_measure(channel, point, outcome) self._may_clear.add(channel) if outcome.result == NO_SIGNAL: region = self.belief.candidate_set(channel) if region is not None and safe_disk_radius(region, self.config) >= 0 \ and math.dist(point, region.mec()[0]) \ <= safe_disk_radius(region, self.config): self._anomaly("no_signal_in_safe_disk", channel, position=point) cert = self._best_certificate(channel) if cert is not None and cert.certified: if self._attempt_clear(channel, cert.center, cert): self._clear_failed.discard(channel) return True self._clear_failed.add(channel) self._may_clear.discard(channel) return self.belief.channels[channel].status is ChannelStatus.CLEARED \end{lstlisting} \subsection*{问题四：朝向候选、可见性选点与网格遍历清除} 每个幸存格附朝向档掩码构成联合候选集，可见性权重用于探点排序；定位停滞时以 $20$ m 方格遍历清除兜底。 \noindent\small \texttt{JointBelief}：单频道的位置与朝向联合候选集；朝向分 $24$ 档，负观测只在整格必落在最小接收半径内时才约束朝向。 \begin{lstlisting}[style=corecode] class JointBelief: """单频道 (位置 × 朝向) 联合可行集。""" def is_absent(self) -> bool: """联合集空 ⟹ 该频道不存在源（真源保留定理的二维推广）。""" return not self.position.alive def sample(self, size: int = 256) -> list[int]: """均匀抽样幸存格（评分用；证书始终用完整集合）。""" cached = self._sample_cache if cached is not None and cached[0] == self.position.version: return cached[1] alive = sorted(self.position.alive) if len(alive) > size: step = len(alive) / size alive = [alive[int(k * step) % len(alive)] for k in range(size)] self._sample_cache = (self.position.version, alive) return alive def expected_posterior_radius(self, reader: Vec2, size: int | None = None) -> float: """一次测量后的期望后验 MEC 半径（米）；D1 值函数。 与 :meth:`ChannelParticles.expected_posterior_radius` 同族，差别只在 no_signal 分支概率：定向源不仅会因距离 > r_eff 静默，还会因朝向遮挡 静默，因此 p_vis(p) = P(r_eff ≥ |p−q|) × cover_ratio(p, q) 其中 ``cover_ratio`` 是该格朝向掩码中"可能覆盖 q 方位"的档位占比。 分支半径语义： - direction（可见，非 near）：读数为 q 处的方位 ±δ，后验收缩到该方位 楔形与当前可行集之交，半径取幸存子采样的质心最大距离； - no_signal：定向源的静默**不删空间候选**，位置后验半径保持当前值； - near：后验为 5m 盘。 返回值只用于**选点排序**，不参与任何安全判定；证书仍由凸/网格口径给出。 """ pf = self.position if not pf.alive: return 0.0 if size is None: size = self.config.bayes.score_subsample indices = self.sample(max(2, size)) if not indices: return 0.0 stride = max(len(indices) // max(self.config.bayes.mc_representatives, 1), 1) representatives = indices[::stride] if _np is not None: return self._epr_vector(reader, indices, representatives) return self._epr_python(reader, indices, representatives) # pragma: no cover def prob_detect(self, reader: Vec2, size: int = 256) -> float: """在 ``reader`` 处能看见任一可行 (位置, 朝向) 的概率质量（启发式）。 r_eff ~ U[1000, 1500] 解析积分；朝向按可行档位占比计权。 只用于选点排序，不参与任何安全判定。 """ cells = self.sample(size) if not cells: return 0.0 total = 0.0 span = self.field.r_eff_max_m - self.field.r_eff_min_m for index in cells: gx, gy = self.position.cell_center(index) distance = math.hypot(reader[0] - gx, reader[1] - gy) if distance > self.field.r_eff_max_m: continue if distance <= self.field.r_eff_min_m: p_range = 1.0 else: p_range = (self.field.r_eff_max_m - distance) / span if p_range <= 0.0: continue bearing = math.degrees(math.atan2(reader[1] - gy, reader[0] - gx)) cover = _COVER_TABLE[_table_index(bearing)] bits = bin(self._mask(index) & cover).count("1") if bits: total += p_range * bits / HEADING_BINS return total / len(cells) def apply(self, station: Vec2, result: str, svd_deg: float | None = None) -> None: """按读数更新联合可行集（保守，真源 (位置, 朝向) 永不被删）。""" if result == DIRECTION: self.position.apply_measure(station, DIRECTION, svd_deg, directional=True) self._apply_heading(station, keep=True) elif result == NEAR: # near 同样要求可见（评测机先判朝向再判距离） self.position.apply_measure(station, NEAR, directional=True) self._apply_heading(station, keep=True) elif result == NO_SIGNAL: # 距离阈值必须取**保守侧**：只有"整格都必然落在 r_eff_min 圆盘内" # 时，no_signal 才一定由朝向遮挡解释。格心到格内任一点 ≤ margin， # 故阈值取 r_eff_min − margin；取 + margin 会把 1000~1010m 的模糊带 # 也当成"必须不覆盖"，从而可能删掉真源的朝向档（破坏真源保留性）。 self._apply_heading(station, keep=False, within_m=max( 0.0, self.field.r_eff_min_m - self.position.margin_m)) else: raise ValueError(f"未知测向结果: {result!r}") def summary(self) -> dict: """联合集规模快照（写入策略 summary，供审计与论文复算）。""" cells = self.sample(256) bins = [bin(self._mask(index)).count("1") for index in cells] return {"cells": len(self.position), "sampled_mean_heading_bins": (round(sum(bins) / len(bins), 2) if bins else 0.0), "heading_bins_total": HEADING_BINS} \end{lstlisting} \noindent\small \texttt{DirectionalExtension}：问题四闭环：定向覆盖计划、遮挡解释与网格遍历清除兜底。 \begin{lstlisting}[style=corecode] class DirectionalExtension(OmnidirectionalLoop): """问题4混合源闭环（P5）：全向基线 + 定向观测处理 + 后备清除。""" def _localize_and_clear(self, channel: int) -> bool: """定位子循环：选点测量收敛候选集，拿到清除证书后清除。 与 P4 的差异：保证接收盘内的 no_signal 说明探点被朝向遮挡 [观测推理]，走"对跖复核 -> 环探"升级梯子（不计入定位步数， 由 max_blocked_probes 限流），之后仍未收敛再转后备清除。 """ search = self.config.search self._may_clear.add(channel) steps = 0 while steps < search.max_localize_steps: if self._virtual_cap_reached(): return False region = self.belief.candidate_set(channel) if region is None: return False cert = clear_certificate(region, channel, self.config) if cert.certified and channel in self._may_clear: if self._attempt_clear(channel, cert.center, cert): return True self._may_clear.discard(channel) # 证书失败：补测后才可重试 continue if len(self.blocked_probes.get(channel, [])) >= search.max_blocked_probes: break # 梯子用尽仍被遮挡：定位停滞，转后备清除 point = self._next_probe(region, channel) within_range = (math.dist(point, region.mec()[0]) <= safe_disk_radius(region, self.config)) outcome = self.client.measure(point, channel) self._observe(channel, point, outcome.result, outcome.svd_deg) self._may_clear.add(channel) if outcome.result == NO_SIGNAL: if within_range and channel not in self.directional_observed: self.directional_observed.append(channel) self.blocked_probes.setdefault(channel, []).append(point) continue # 盘外无信号也可能是距离不足，不据此断言定向 self.blocked_probes.pop(channel, None) steps += 1 if self.belief.channels[channel].status is ChannelStatus.CLEARED: return True region = self.belief.candidate_set(channel) if region is None: return False cert = clear_certificate(region, channel, self.config) if cert.certified and channel in self._may_clear: if self._attempt_clear(channel, cert.center, cert): return True self._may_clear.discard(channel) return self._fallback_clear(channel, region) def _observe(self, channel: int, point: Vec2, result: str, svd_deg: float | None) -> None: """同步更新凸候选集与 (位置 × 朝向) 联合可行集。 联合集只在频道首次被看见（direction/near）时建立：从未被看见的 频道由覆盖证明负责缺席，不需要朝向维度。联合集一旦为空，直接给出 该频道的缺席证明（严格，与覆盖证明同一条件），可提前停止对其检测。 """ self.belief.update_measure(channel, point, result, svd_deg) if not self.config.search.directional_joint_belief: return joint = self.joint.get(channel) if joint is None: if result not in (DIRECTION, NEAR): return joint = JointBelief(self.config) self.joint[channel] = joint joint.apply(point, result, svd_deg) if joint.is_absent() and self.belief.channels[channel].status \ is not ChannelStatus.CLEARED: self.belief.mark_absent(channel) self.absence_evidence[channel] = "joint_empty" if channel not in self.joint_absent_channels: self.joint_absent_channels.append(channel) def _next_probe(self, region: CandidateRegion, channel: int) -> Vec2: """定向遮挡下的探点升级梯子 [工程]。 - 未被遮挡过：worst-posterior 选点器（保证接收盘内）； - 第 1 次被遮挡：对 MEC 圆心取对跖方向，偏移裁进保证接收盘； - 之后：绕 MEC 圆心、半径 = ρ 的环探，方位自首个被遮挡点起 每次 +45°。几何边界必须如实陈述：环上可见弧没有下界—— 源可位于环上（|g-c| = ρ）且朝向沿外切线，可见集退化为单个 切点，固定步进环探可以全部漏过。梯子只提高提前命中概率， 不构成接收或命中保证；完备性完全由后备清除兜底。 """ visibility = self._visibility_probe(region, channel) if visibility is not None: return visibility blocked = self.blocked_probes.get(channel, []) center, radius = region.mec() if not blocked: return select_measurement_point(region, self.client.state.position, self.config).position safe = safe_disk_radius(region, self.config) first = blocked[0] if len(blocked) == 1 and safe >= 0: dx, dy = center[0] - first[0], center[1] - first[1] norm = math.hypot(dx, dy) if norm <= 1e-9: # 首个被遮挡探点恰在圆心：方向信息为零，沿 +x 走满安全盘。 return (center[0] + safe, center[1]) scale = min(1.0, safe / norm) # 对跖方向，偏移模长 = min(norm, safe)：圆心对称点裁进保证接收盘。 return (center[0] + dx * scale, center[1] + dy * scale) base = math.atan2(first[1] - center[1], first[0] - center[0]) azimuth = base + math.radians(45.0 * (len(blocked) - 1)) return (center[0] + radius * math.cos(azimuth), center[1] + radius * math.sin(azimuth)) def _fallback_clear(self, channel: int, region: CandidateRegion) -> bool: """定位停滞的后备：候选集有限圆盘覆盖逐格清除（格数即请求上界）。""" search, session = self.config.search, self.config.session width = max(p[0] for p in region.polygon) - min(p[0] for p in region.polygon) height = max(p[1] for p in region.polygon) - min(p[1] for p in region.polygon) cells = fallback_clear_cells(region.polygon, self.config) entry = {"channel": channel, "cells": len(cells), "bound": fallback_cell_bound(width, height, self.config), "posterior_radius_m": round(region.mec()[1], 3)} self.fallback_log.append(entry) if len(cells) > search.fallback_max_cells: # 宽候选集盲扫违背预算纪律：留待定位收敛，不得硬扫 [计划§4]。 entry["outcome"] = "skipped_cell_cap" self._anomaly("fallback_cell_cap_exceeded", channel, cells=len(cells), cap=search.fallback_max_cells) return False # 现实预算核对：N 次请求 × 单请求高分位估计 + 退出余量，不足不启动。 remaining = self.client.real_remaining_s if remaining is None or len(cells) * session.real_rtt_p95_s \ + session.reserve_exit_s > remaining: entry["outcome"] = "skipped_budget" return False attempts = 0 order = self._fallback_order(cells) for cell in order: if self._virtual_cap_reached(): entry["outcome"] = "skipped_virtual_budget" entry["attempts"] = attempts return False attempts += 1 outcome = self.client.clear(cell, channel) if outcome.result == CLEAR_SUCCESS: self.belief.mark_cleared(channel) self.cleared.append(channel) entry["outcome"] = "cleared" entry["attempts"] = attempts return True # 候选集保守包含真源，格心距格内任一点 <= 10*sqrt(2) < 20m：全部 # 失败说明守恒性被打破，按异常留痕，不得据此删除候选。 entry["outcome"] = "exhausted" entry["attempts"] = attempts self._anomaly("fallback_exhausted_without_success", channel, cells=len(cells)) return False \end{lstlisting} \endgroup