UVa : 361 - Cops and Robbers

問題概要

点集合C, R, Oが与えられる。
それぞれの集合には、点が最大200個含まれている。
Oの中の各点について、以下の3つのパターンの内どれになるか答えよ

  • この点が、Cの任意の3点により作られる三角形に包含される : safeと呼ぶ
  • この点が、Rの任意の3点により作られる三角形に包含される、かつsafeではない : robbedと呼ぶ
  • safeでもrobbedでもない場合 : neitherと呼ぶ

ただし、「点が三角形に包含される」とは、点が三角形の中に完全に入っている、もしくは点が三角形の辺上に乗っていることを意味する。

解法

さて、この問題で難しいところは「点集合の任意の3点から成る三角形に、ある点が含まれるかどうか」を計算することです。
一番単純な方法は、3点を三重ループで選んで点がその三角形に含まれるか判定する方法です。しかし、そうするとその三重ループに、さらにOの点を選ぶループがかかってくるので、O(200^4)になってしまいます。

もっと効率のいい方法で計算するには、凸包を使いましょう。
とりあえず、CとRの集合を凸包して凸多角形を作ります。
もし、この凸多角形の中に点が包含されていれば、「点集合の任意の3点から成る三角形に、点が含まれている」ということになります。

点集合の任意の3点から成る三角形というのは、必ず凸包の多角形に含まれているため、この計算法で判定することができます。

コーナーケース

以下のような感じのセットがあるっぽいので注意です。

0 2 1
    0 0  0 0
    0 0
0 3 1
    0 0  0 0  0 0
    0 0
2 3 1
    -5 0  5 0
    -5 0  5 0  5 0
    0 0
3 0 1
    0 0  10 0  -10 0
    20 0
0 0 0
Data set 1:
     Citizen at (0,0) is neither.

Data set 2:
     Citizen at (0,0) is robbed.

Data set 3:
     Citizen at (0,0) is robbed.

Data set 4:
     Citizen at (20,0) is neither.

プログラム

幾何のコードは、スパゲッティソースのコードを使用しています。

  • 凸包:O(n log n)
  • 凸多角形包含判定 : O(log n)
    • ただし、今回は点が辺上に乗っているかの判定を入れているので O(n) になります
#include <iostream>
#include <complex>
#include <vector>
#include <algorithm>
#include <cstring>
using namespace std;

#define EPS (1e-8)

typedef complex<double> P;

bool cmp(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;
  if(cross(b, c) < 0)   return -1;
  if(dot(b, c) < 0)     return +2;
  if(norm(b) < norm(c)) return -2;
  return 0;
}

bool intersectSP(const P &a, const P &b, const P &p) {
  return abs(a - p) + abs(b - p) - abs(b - a) < EPS;
}

vector<P> convexHull(vector<P> ps){
  if(ps.size() < 3) return vector<P>();

  int sz = ps.size();
  int k = 0;
  vector<P> ch(2 * sz);

  sort(ps.begin(), ps.end(), cmp);

  for(int i = 0; i < sz; ch[k++] = ps[i++]){
    while(k >= 2 && ccw(ch[k - 2], ch[k - 1], ps[i]) <= 0){
      --k;
    }
  }
  for(int i = sz - 2, j = k + 1; i >= 0; ch[k++] = ps[i--]){
    while(k >= j && ccw(ch[k - 2], ch[k - 1], ps[i]) <= 0){
      --k;
    }
  }
  ch.resize(k - 1);

  vector<P> res;

  for(int i = 0; i < ch.size(); i++){
    bool flg = true;

    for(int j = 0; j < i; j++){
      if(ch[i] == ch[j]){
        flg = false;
        break;
      }
    }
    if(flg) res.push_back(ch[i]);
  }

  return res;
}

enum{OUT, ON, IN};

int convexContains(const vector<P> &pol, const P &p){
  const int n = pol.size();
  P g = (pol[0] + pol[n / 3] + pol[2 * n / 3]) / 3.0;
  int a = 0, b = n;

  while(a + 1 < b){
    int c = (a + b) / 2;
    if(cross(pol[a] - g, pol[c] - g) > 0){
      if(cross(pol[a] - g, p - g) > 0 && cross(pol[c] - g, p - g) < 0) b = c;
      else                                                             a = c;
    }
    else{
      if(cross(pol[a] - g, p - g) < 0 && cross(pol[c] - g, p - g) > 0) a = c;
      else                                                             b = c;
    }
  }
  b %= n;
  if(cross(pol[a] - p, pol[b] - p) < 0) return OUT;
  if(cross(pol[a] - p, pol[b] - p) > 0) return IN;

  for(int i = 0; i < pol.size(); i++){
    if(intersectSP(pol[i], pol[(i + 1) % pol.size()], p)){
      return ON;
    }
  }

  return OUT;
}

enum{NEITHER, SAFE, ROBBED};
string output[] = {"neither", "safe", "robbed"};

int c, r, o;
vector<P> C, R, O;
int kind[202];

void input(int n, vector<P> &v){
  P p;
  v.clear();

  for(int i = 0; i < n; i++){
    cin >> p.real() >> p.imag();
    v.push_back(p);
  }
}

int main(){
  for(int SET = 1; cin >> c >> r >> o, c || r || o; SET++){
    cout << "Data set " << SET << ":" << endl;

    input(c, C);
    input(r, R);
    input(o, O);

    C = convexHull(C);
    R = convexHull(R);

    memset(kind, NEITHER, sizeof(kind));

    for(int i = 0; i < o; i++){
      if(C.size() >= 1 && convexContains(C, O[i])){
        kind[i] = SAFE;
      }
    }

    for(int i = 0; i < o; i++){
      if(kind[i] != 0) continue;
      if(R.size() >= 1 && convexContains(R, O[i])){
        kind[i] = ROBBED;
      }
    }

    for(int i = 0; i < o; i++){
      cout << "     Citizen at ("
           << O[i].real() << "," << O[i].imag() << ") is "
           << output[kind[i]] << "." << endl;
    }

    cout << endl;
  }
}