AtCoder Beginner Contest 151 F - Enclose All

問題へのリンク

問題文

平面上の $N$ 個の点 $(x_i, y_i)$ が与えられます。

これら全てを内部または周上に含む円の半径の最小値を求めてください。

自分の解法

条件を満たす円の中心の座標を $(s, t)$ とする。中心 $(s, t)$ の円のうち、与えられた点を全て満たすものの半径の2乗は$\max_{i=1,\cdots,N} (x_i-s)^{2} + (y_i - t)^{2} $ となる。
この半径の2乗を最小にする $(s, t)$ を求めれば良いので、以下の最適化問題を解けば良い。

$$ \min_{s, t} \max_{i=1, \cdots, N} (x_i - s)^{2} + (y_i - t)^{2} $$

ここで、各 $i$ について $(x_i - s)^{2} + (y_i - t)^{2}$ は $(s, t)$ に関する凸関数であり、2つの凸関数の $\max$ もまた凸関数となるので、$ \max_{i=1, \cdots, N} (x_i - s)^{2} + (y_i - t)^{2}$ も$(s, t)$ に関する凸関数となる。
よって $(s, t)$ それぞれについて勾配法でパラメータを更新していけば最適解が求まる。ただし、ステップサイズを減少させる必要があることに留意する。

提出コード

#include <iostream>
#include <fstream>
#include <cstdio>
#include <cmath>
#include <vector>
#include <string>
#include <set>
#include <map>
#include <stack>
#include <queue>
#include <bitset>
#include <algorithm>
#include <complex>
#include <array>
#include <iomanip>
using namespace std;

#define REP(i,n) for(int i=0; i<n; ++i)
#define FOR(i,a,b) for(int i=a; i<=b; ++i)
#define FORR(i,a,b) for (int i=a; i>=b; --i)
#define ALL(c) (c).begin(), (c).end()

typedef long long ll;
typedef vector<int> VI;
typedef vector<ll> VL;
typedef vector<long double> VD;
typedef vector<VI> VVI;
typedef vector<VL> VVL;
typedef vector<VD> VVD;
typedef pair<int,int> P;
typedef pair<ll,ll> PL;

template<typename T> void chmin(T &a, T b) { if (a > b) a = b; }
template<typename T> void chmax(T &a, T b) { if (a < b) a = b; }

int in() { int x; scanf("%d", &x); return x; }
ll lin() { ll x; scanf("%lld", &x); return x; }
#define INF 1LL<<60

int N;
vector< pair<double, double> > points;

pair<int, double> calc_max(double s, double t) {
    double r_max = 0;
    int idx = -1;
    REP(i, N) {
        if(r_max < sqrt( pow(points[i].first - s, 2) + pow(points[i].second - t, 2) )) {
            r_max = sqrt( pow(points[i].first - s, 2) + pow(points[i].second - t, 2) );
            idx = i;
        }
    }
    return make_pair(idx, r_max);
}

int main() {
    cin >> N;
    points.resize(N);
    REP(i, N) {
        cin >> points[i].first >> points[i].second;
    }

    double s = 50.0, t = 50.0, lr;
    REP(itr_outer, 50000) {
        if(itr_outer < 1000)
            lr = 1.0 / sqrt(1.0 + itr_outer);
        else if(itr_outer < 5000)
            lr = 1.0 / pow(1 + itr_outer, 1);
        else
            lr = 1.0 / pow(1.1, sqrt(1 + itr_outer));
        pair<int, double> ret = calc_max(s, t);
        double gs = 2 * (s - points[ret.first].first);
        double gt = 2 * (t - points[ret.first].second);
        s -= lr * gs;
        t -= lr * gt;
    }
    double ans = calc_max(s, t).second;
    printf("%.12f\n", ans);

    return 0;
}

感想

解説の解き方がどうすれば思い付くのか謎。
ステップ幅の調整で魔チューニングを行ってるが、もう少しシンプルに出来るかもしれない。