Aizu Online Judge 1132 - Circle and Points

問題

http://judge.u-aizu.ac.jp/onlinejudge/description.jsp?id=1132

昔といた問題です(また後で新しく解き直したコードを置きます)。

解法

答えを導き出せるような円の位置は無数にあるので、円を固定して考えることにします。
ある2点を通ることがわかれば円の位置は2つに定まります。また含む点を円周上に置くことによって、他の点を含むための面積を最大限使えると考えることができます。後は通る2点を全通り試して、他の点を含むかどうかを調べればokです。

ソースコード

#include <iostream>
#include <cstdio>
#include <math.h>
using namespace std;
 
const int MAX_N = 300;
 
int N;
double px[MAX_N], py[MAX_N];
 
int main() {
    while (cin >> N, N) {
        for (int i = 0; i < N; i++) {
            cin >> px[i] >> py[i];
        }
        int ans = 1;
        for (int i = 0; i < N; i++) {
            for (int j = i + 1; j < N; j++) {
                // 傾きのベクトル => [dx, dy]
                double dx = px[i] - px[j],
                       dy = py[i] - py[j];
                if (dx * dx + dy * dy > 2.0001 * 2.0001) continue;
                // 法線ベクトル => [dy, -dx]
                // 弦を二等分する点 P から中心に向かうベクトルをqとする。
                double norm = 1.0 / sqrt(dy * dy + dx * dx) *
                            sqrt(1.0 * 1.0 - (dx * dx + dy * dy) / 4.0);
                double qx = dy * norm,
                       qy = -dx * norm;
                for (int k = 0; k < 2; k++) {
                    qx *= -1;
                    qy *= -1;
                    double rx = (px[i] + px[j]) / 2.0 + qx,
                           ry = (py[i] + py[j]) / 2.0 + qy;
                    int tmp_ans = 0;
                    for (int l = 0; l < N; l++) {
                        if ((rx - px[l]) * (rx - px[l]) + (ry - py[l]) * (ry - py[l]) <= 1.0001) {
                            tmp_ans++;
                        }
                    }
                    ans = max(ans, tmp_ans);
                }
            }
        }
        cout << ans << endl;
    }
    return 0;
}

感想

Nが1のケースなどでバグった記憶があります...