Tenka1 Programmer Contest E - CARtesian Coodinate

問題ページ
E - CARtesian Coodinate

解法

解説を見た。
https://img.atcoder.jp/tenka1-2017/editorial.pdf
www.youtube.com

まず、N点からのマンハッタン距離の和が最小になる点について考える。マンハッタン距離はx座標、y座標それぞれ独立に足すことができるのx方向とy方向でそれぞれ独立に考えることができる。また、距離の和が最小になる点は中央値になる。中央値から右にΔt動かしたとすると(右の頂点の数)Δtの分距離の総和が減るが、(左の頂点の数)Δtの分距離の総和が増えるため総和が減ることはありえない。左に動かすのも同様。よって中央値から動かして得になるケースは絶対に存在しない。

(Δt*2減るがΔt*3増えるので得しない)

よって、交点の座標の中央値を求めることができれば解ける。交点はO(N^2)個存在するので愚直に全列挙するのは間に合わない(はず)。二分探索することを考えると座標sより左側にある交点を高速に列挙できれば二分探索が可能。
座標sの地点で直線が並ぶ順番と左側に交点が存在しない十分左の座標で直線が並ぶ順番はそれぞれO(N)で求めることができる。このとき、交点は順番が入れ替わっている直線の数となり、これは転倒数となる。

転倒数はBITを用いてO(NlogN)で求めることができるので二分探索の判定部分はO(NlogN)で可能。したがって求めることができる。

学び

マンハッタン距離の総和は次元ごとに独立に考えられる。総和が最小になる地点は中央値。
直線の交点の数は転倒数に帰着できる。

#include <bits/stdc++.h>

using namespace std;
typedef long long ll;
// #define int ll
typedef vector<int> VI;
typedef vector<VI> VVI;
typedef vector<ll> VL;
typedef vector<VL> VVL;
typedef pair<int, int> PII;

#define FOR(i, a, n) for (ll i = (ll)a; i < (ll)n; ++i)
#define REP(i, n) FOR(i, 0, n)
#define ALL(x) x.begin(), x.end()
#define IN(a, b, x) (a<=x&&x<b)
#define MP make_pair
#define PB push_back
#ifdef int
const int INF = (1LL<<30);
#else
const ll INF = (1LL<<60);
#endif
const double PI = 3.14159265359;
const double EPS = 1e-12;
const int MOD = 1000000007;

template <typename T> T &chmin(T &a, const T &b) { return a = min(a, b); }
template <typename T> T &chmax(T &a, const T &b) { return a = max(a, b); }

int dx[] = {0, 1, 0, -1}, dy[] = {1, 0, -1, 0};

int n;
int a[40010], b[40010], c[40010];
pair<double, int> p[40010], p2[40010], p3[40010];

ll con;

int bit[40010], m = 40010;

//0からi-1番目までの和
int sum(int i) {
  i++;
  int s = 0;
  while(i > 0) {
    s += bit[i];
    i -= i&-i;
  }
  return s;
}

//i番目(0-index)にxを加える
void add(int i, int x) {
  i++;
  while(i <= m) {
    bit[i] += x;
    i += i&-i;
  }
}

bool check(double mid) {
  // 数列を求める
  // 画像と数列の並べ方が反対になってる
  REP(i, n) {
    int idx = p[i].second;
    p2[i] = {-(double)a[idx]/b[idx]*mid+(double)c[idx]/b[idx], i};
  }
  sort(p2, p2+n);
  // p2.secondの転倒数を求める
  memset(bit, 0, sizeof(bit));
  ll ans = 0;
  REP(i, n) {
    ans += i - sum(p2[i].second);
    add(p2[i].second, 1);
  }
  return ans >= con;
}

bool check2(double mid) {
  // 数列を求める
  REP(i, n) {
    int idx = p3[i].second;
    p2[i] = {-(double)b[idx]/a[idx]*mid+(double)c[idx]/a[idx], i};
  }
  sort(p2, p2+n);
  // p2.secondの転倒数を求める
  memset(bit, 0, sizeof(bit));
  ll ans = 0;
  REP(i, n) {
    ans += i - sum(p2[i].second);
    add(p2[i].second, 1);
  }
  return ans >= con;
}

signed main(void)
{
  cin >> n;
  REP(i, n) {
    cin >> a[i] >> b[i] >> c[i];
    p[i] = {-(double)a[i]/b[i], i};
    p3[i] = {(double)b[i]/a[i]*INF+(double)c[i]/a[i], i};
  }
  sort(p, p+n); sort(p3, p3+n);
  reverse(p, p+n);

  // 中央値はcon番目
  con = n*(n-1)/4 + !!(n*(n-1)%4);

  // x座標
  // 小数の二分探索は定数回で固定して回す
  double low = -1e10, high = 1e10;
  REP(_, 300) {
    double mid = (low+high)/2;
    if(check(mid)) {
      high = mid;
    } else {
      low = mid;
    }
  }
  cout << fixed << setprecision(15) << low << " ";

  // y座標
  low = -1e10, high = 1e10;
  REP(_, 300) {
    double mid = (low+high)/2;
    if(check2(mid)) {
      high = mid;
    } else {
      low = mid;
    }
  }
  cout << low << endl;

  return 0;
}