Ex - Disk and Segments Editorial
by 
MMNMM
この問題は、関数 \(f(x,y)\coloneqq\) 点 \((x,y)\) を中心とした閉円盤がすべての線分と共有点を持つような最小の半径 について、次の値を求めるものだと考えることができます。 \[\min _ {(x,y)\in\mathbb R ^ 2}f(x,y)\]
重要な事実として、\(f(x,y)\) は凸関数です。
証明
補題 1
\(p\in\mathbb R^2\) について、点 \(p\) と点 \((x,y)\) の距離を \(d _ p(x,y)\) とする。\(d _ p\colon\mathbb R ^ 2\to\mathbb R _ {{}\geq0}\) は凸関数である。
補題 1 の証明
エピグラフ \(\left\lbrace(x,y,z)\mid z\geq d _ p(x,y)\right\rbrace\) は円錐になるため、これは凸です。よって、\(d _ p(x,y)\) は凸です。
補題 2
\(a,b\in\mathbb R^2\) を結ぶ線分 \(s\) について、線分 \(s\) と点 \((x,y)\) の距離を \(d _ s(x,y)\) とする。\(d _ s\colon\mathbb R^2\to\mathbb R _ {{}\geq0}\) は凸関数である。
補題 2 の証明
\(c,d\in\mathbb R^2\) および \(0\leq t\leq1\) について、\(d _ s((1-t)c+td)\leq(1-t)d _ s(c )+td _ s(d)\) を示します。
\(d _ s(c )=d _ {p _ c}(c ),d _ s(d)=d _ {p _ d}(d)\) を満たすような \(s\) 上の点 \(p _ c,p _ d\) が存在します。 \(p=(1-t)p _ c+tp _ d\) とすると、\(d _ s\) の定義から \(d _ s((1-t)c+td)\leq d _ p((1-t)c+td)\) です。
ここで、原点 \(O=(0,0)\) と \(q,r\in\mathbb R^2\) について \(d _ q(r )=d _ O(r-q)\) であることに注意すると、 \[d _ p((1-t)c+td)=d _ O((1-t)(c-p _ c)+t(d-p _ d))\] および \[d _ {p _ c}(c )=d _ O(c-p _ c),d _ {p _ d}(d)=d _ O(d-p _ d)\] が成り立ちます。 補題 1 より、\(d _ O((1-t)(c-p _ c)+t(d-p _ d))\leq(1-t)d _ O(c-p _ c)+td _ O(d-p _ d)\) なので、 \[d _ s((1-t)c+td)\leq d _ p((1-t)c+td)\leq(1-t)d _ s(c )+td _ s(d)\] となり、示されました。
補題 3
\(g _ i\colon\mathbb R ^ 2\to\mathbb R\ (i=1,2,\ldots,n)\) を凸関数とする。\(\displaystyle g(x,y)\coloneqq\max _ ig _ i(x,y)\) は凸関数である。
補題 3 の証明
\(a,b\in\mathbb R^2\) および \(0\leq t\leq1\) について、\(g((1-t)a+tb)\leq(1-t)g(a)+tg(b)\) を示します。
\(g\) の定義より、\(g((1-t)a+tb)=g _ i((1-t)a+tb)\) なる \(i\) が存在します。 \(g _ i\) は凸なので、\(g _ i((1-t)a+tb)\leq(1-t)g _ i(a)+tg _ i(b)\) です。 \(g\) の定義より \({} ^ \forall i.\;g\geq g _ i\) なので、\((1-t)g _ i(a)+tg _ i(b)\leq(1-t)g(a)+tg(b)\) が成り立ち、示されました。
与えられた線分 \(s _ 1,s _ 2,\ldots,s _ N\) について、\[f(x,y)=\max _ {i=1,2,\ldots,N}d _ {s _ i}(x,y)\]です。
補題 2 より各 \(d _ {s _ i}\) は凸関数なので、補題 3 より \(f(x,y)\) は凸関数です。
一変数凸（あるいは単峰）関数の最小値がある範囲に存在するとわかっているとき、三分探索でその最小値を近似できることが知られています。
ここで、次の関数を考えます。 \[F(x):=\min _ {y\in\mathbb R}f(x,y)\] すると、\(F(x)\) は一変数凸関数であることが示せます。
証明
\(F\) のエピグラフ \(\left\lbrace(x,z)\mid z\geq F(x)\right\rbrace\) が凸であることを示します。\(F\) のエピグラフ内の任意の \(2\) 点 \(a,b\) を取ると、\(x\) 座標および \(z\) 座標がそれぞれ等しいような \(f\) のエピグラフ \(\left\lbrace(x,y,z)\mid z\geq f(x,y)\right\rbrace\) の点 \(a ^ \prime,b ^ \prime\) が存在します。
\(f\) の凸性より、\(a ^ \prime,b ^ \prime\) を結んだ線分は \(f\) のエピグラフに含まれます。つまり、線分上の任意の点 \((x,y,z)\) は \(z\geq f(x,y)\) を満たします。 \(F\) の定義より \(F(x)\leq f(x,y)\) なので \(z\geq F(x)\) となり、\(a,b\) を結んだ線分が \(F\) のエピグラフに含まれることが示されました。
\(x\) を固定したとき \(f(x,y)\) が一変数凸関数であることはすぐわかるため、三分探索の中で三分探索を行うことでこの問題を解くことができました。
実装例は以下のようになります。
#include <iostream>
#include <tuple>
#include <vector>
#include <cmath>
#include <ranges>
int main() {
using namespace std;
unsigned N;
cin >> N;
using real = double;
using segment_type = tuple<real, real, real, real>;
using point_type = tuple<real, real>;
vector<segment_type> segments(N);
for (auto&&[a, b, c, d] : segments)
cin >> a >> b >> c >> d;
// 2 点の距離
const auto distance_point_and_point{[](const point_type &point0, const point_type &point1) {
const auto&[x, y]{point0};
const auto&[z, w]{point1};
return hypot(x - z, y - w);
}};
// 点と線分の距離
const auto distance_segment_and_point{[&](const segment_type &segment, const point_type &point) {
const auto&[a, b, c, d]{segment};
const auto&[x, y]{point};
if ((a - c) * (a - x) + (b - d) * (b - y) < 0)
return distance_point_and_point({a, b}, {x, y});
if ((c - a) * (c - x) + (d - b) * (d - y) < 0)
return distance_point_and_point({c, d}, {x, y});
return abs((a - c) * (b - y) - (b - d) * (a - x)) / distance_point_and_point({a, b}, {c, d});
}};
// 与えられた点を中心とする閉円盤であってすべての線分と共有点をもつもののうち最小の半径
const auto minimum_crossing_circle{[&](const point_type &point) {
real ret{};
for(const auto& segment : segments)
ret = max(ret, distance_segment_and_point(segment, point));
return ret;
}};
// 三分探索で単峰関数の最小値を求める
const auto minimize_unimodal_function{[](real L, real R, auto &&function) {
for (const auto _ : ranges::views::iota(0, 100)) {
auto M1{(L * 2 + R) / 3};
auto M2{(L + 2 * R) / 3};
if (function(M1) < function(M2))
R = M2;
else
L = M1;
}
return function(L);
}};
cout << minimize_unimodal_function(0, 1000, [&](auto x) {
// x を固定したときの f(x, y) の最小値を求める
return minimize_unimodal_function(0, 1000, [&](auto y) {
// x, y を決めると半径の最小値がわかる
return minimum_crossing_circle({x, y});
});
}) << endl;
return 0;
}
また、\(f(x,y)\) の最小値を scipy.optimize.minimize 関数で求めることもできます。
from scipy.optimize import minimize
from math import sqrt
N = int(input())
segments = [tuple(map(int, input().split())) for i in range(N)]
def f(param):
x, y = param
ret = 0
for a, b, c, d in segments:
if (a - c) * (a - x) + (b - d) * (b - y) < 0:
ret = max(ret, sqrt((a - x) ** 2 + (b - y) ** 2))
elif (c - a) * (c - x) + (d - b) * (d - y) < 0:
ret = max(ret, sqrt((c - x) ** 2 + (d - y) ** 2))
else:
ret = max(ret, abs((a - c) * (b - y) - (b - d) * (a - x)) / sqrt((a - c) ** 2 + (b - d) ** 2))
return ret
print(minimize(f, (500, 500), args=(), method='Nelder-Mead').fun)
posted:
last update: