AGC012 B - Holes - くれなゐの雑記

問題

agc021.contest.atcoder.jp

問題概要

めちゃめちゃでかい $R$ を持つ円が与えられる。
その真ん中 $\pm 10^6$ の範囲内に穴を置く。
めちゃめちゃでかい円 $R$ 内に点を起き、その点は最も近くにある穴に落ちる。
円の内部すべてに点を起き、各穴に落ちる面積の割合をそれぞれ求める。
同じ座標に2つ置くようなケースはない
ちなみに与えられる座標は整数値(今回はこの条件は使わないゴリ押し解法だと必要かも？)

考察

入力例2に注目すると0のケースが有る。点の凸包を作り、その内部にある面積は円に比べて微小なので0だと考えることが出来る。
結果、凸包の外側のみ考えればよいということになる。
凸包の形が例えば正六角形のような形の時、それぞれの穴に落ちる部分は図のようになる気がする。(実際は円はもっとでかいので、六角形の外側のほうがだいぶ面積は大きい)
この図から察するに、(中心が一致している時は)各面積は点と点の間の垂直二等分線で円を区切った面積になりそう
形が正n角系からずれると、中心点は当然ずれて面積の比率は変わってくるが円が非常に大きいのでそれは微々たる差
結局、ある点に注目した時の面積はその両隣の点で作ることができる2本の垂直二等分線の角度のみが面積の比率に影響する
角度は全部足すと $2 \pi$ になるので、両隣の線の垂直二等分線の角度を $2 \pi$ で割ればよさそう
2本の垂直二等分線の角度って結局それは2本の線分を90度回しただけなので、2本の線分の角度を求めれば良い。

f:id:kurenaif:20180227020709p:plain

実装方針

まずはconvex_hullを使っていらない点を削除
convex_hullの重心が0になるように移動させる(この処理は不要だった)
convex_hullで点をいい感じにソートするので、両隣の点を比較し、線分を作って角度を測る
$2 \pi$ で割り、各点がどの確率になるかを求める。
各点情報からindexを復元し、出力する

ソースコード

agc021.contest.atcoder.jp

#include <iostream>
#include <queue>
#include <map>
#include <list>
#include <vector>
#include <string>
#include <stack>
#include <limits>
#include <climits>
#include <cassert>
#include <fstream>
#include <cstring>
#include <cmath>
#include <bitset>
#include <iomanip>
#include <algorithm>
#include <functional>
#include <cstdio>
#include <ciso646>
#include <set>
#include <array>
#include <unordered_map>
#include <complex>

using namespace std;

#define FOR(i,a,b) for (int i=(a);i<(b);i++)
#define RFOR(i,a,b) for (int i=(b)-1;i>=(a);i--)
#define REP(i,n) for (int i=0;i<(n);i++)
#define RREP(i,n) for (int i=(n)-1;i>=0;i--)

#define inf 0x3f3f3f3f
#define PB push_back
#define MP make_pair
#define ALL(a) (a).begin(),(a).end()
#define SET(a,c) memset(a,c,sizeof a)
#define CLR(a) memset(a,0,sizeof a)
#define VS vector<string>
#define VI vector<ll>
#define DEBUG(x) cout<<#x<<": "<<x<<endl
#define MIN(a,b) (a>b?b:a)
#define MAX(a,b) (a>b?a:b)
#define pi 2*acos(0.0)
#define INFILE() freopen("in0.txt","r",stdin)
#define OUTFILE()freopen("out0.txt","w",stdout)
#define ll long long
#define ull unsigned long long
#define pii pair<ll,ll>
#define pcc pair<char,char>
#define pic pair<ll,char>
#define pci pair<char,ll>
#define eps 1e-14
#define FST first
#define SEC second
#define SETUP cin.tie(0), ios::sync_with_stdio(false), cout << setprecision(15)

namespace {
	struct input_returnner {
		ll N; input_returnner(ll N_ = 0) :N(N_) {}
		template<typename T> operator vector<T>() const { vector<T> res(N); for (auto &a : res) cin >> a; return std::move(res); }
		template<typename T> operator T() const { T res; cin >> res; return res; }
		template<typename T> T operator - (T right) { return T(input_returnner()) - right; }
		template<typename T> T operator + (T right) { return T(input_returnner()) + right; }
		template<typename T> T operator * (T right) { return T(input_returnner()) * right; }
		template<typename T> T operator / (T right) { return T(input_returnner()) / right; }
		template<typename T> T operator << (T right) { return T(input_returnner()) << right; }
		template<typename T> T operator >> (T right) { return T(input_returnner()) >> right; }
	};
	template<typename T> input_returnner in() { return in<T>(); }
	input_returnner in() { return input_returnner(); }
	input_returnner in(ll N) { return std::move(input_returnner(N)); }
}

const ll MOD = 1e9 + 7;

void solve();

signed main() {
	SETUP;
	solve();
#ifdef _DEBUG
	system("pause");
#endif
	return 0;
}

typedef complex<double> P;
namespace std {
	bool operator < (const P& a, const P& b) {
		return real(a) != real(b) ? real(a) < real(b) : imag(a) < imag(b);
	}
}
double cross(const P& a, const P& b) {
	return imag(conj(a)*b);
}
double dot(const P& a, const P& b) {
	return real(conj(a)*b);
}

int ccw(P a, P b, P c) {
	b -= a; c -= a;
	if (cross(b, c) > 0)   return +1;       // counter clockwise
	if (cross(b, c) < 0)   return -1;       // clockwise
	if (dot(b, c) < 0)     return +2;       // c--a--b on line
	if (norm(b) < norm(c)) return -2;       // a--b--c on line
	return 0;
}

vector<P> convex_hull(vector<P> ps) {
	int n = ps.size(), k = 0;
	sort(ps.begin(), ps.end());
	vector<P> ch(2 * n);
	for (int i = 0; i < n; ch[k++] = ps[i++]) // lower-hull
		while (k >= 2 && ccw(ch[k - 2], ch[k - 1], ps[i]) <= 0) --k;
	for (int i = n - 2, t = k + 1; i >= 0; ch[k++] = ps[i--]) // upper-hull
		while (k >= t && ccw(ch[k - 2], ch[k - 1], ps[i]) <= 0) --k;
	ch.resize(k - 1);
	return ch;
}

void solve() {
	int N; cin >> N;
	vector<P> ps;
	REP(i, N) {
		double x, y; cin >> x >> y;
		ps.push_back(P(x,y));
	}
	if (N == 2) {
		cout << 0.5 << endl << 0.5 << endl;
		return;
	}
	vector<P> res = convex_hull(ps);
	P sum;
	REP(i, res.size()) {
		sum += res[i];
	}
	P center = P(sum.real() / res.size(), sum.imag() / res.size());
	double dsum = 0.0;
	vector<double> test(res.size());
	double summm = 0;
	REP(i, res.size()) {
		int prv = i;
		int tar = (i + 1) % res.size();
		int nxt = (i + 2) % res.size();
		P p1 = res[tar] - res[prv];
		P p2 = res[nxt] - res[tar];
		double a = abs(arg(p1) - arg(p2));
		if (a > pi) a = 2 * pi - a;
		test[tar]  = (a / (2 * pi));
	}
	vector<pair<P, double> > pd;
	REP(i, res.size()) {
		pd.push_back(make_pair(P(res[i]), test[i]));
	}
	vector<double> output(ps.size());
	REP(i, pd.size()) {
		auto result = std::find(ALL(ps), pd[i].first);
		size_t index = std::distance(ps.begin(), result);
		output[index] = pd[i].second;
	}
	REP(i, N) {
		cout << output[i] << endl;
	}
}

感想

中心にずらしてから復元して、indexを検索したのだが、doubleの計算のズレが発生してエラーが発生した反省
垂直二等分線の中心がconvex_hullの中心からずれていても大丈夫ということに気づくのにすごい時間がかかった
角度が $\pi$ の時にコーナーケースだった

References

http://www.prefield.com/algorithm/geometry/convex_hull.html