AOJ1190 Anchored Balloon

問題
http://judge.u-aizu.ac.jp/onlinejudge/description.jsp?id=1190&lang=jp

Z=0の平面上にN個の整数座標が与えられる。各々に杭を打ち、それぞれ異なる長さのロープで一つのバルーンに繋ぐ。
上空に飛ばしたとき、バルーンの上昇する最大の高さを求めよ。
1 <= N <= 10
−100 <= X[i] <= 100, −100 <= Y[i] <= 100, 1 <= L[i] <= 300

解答
それぞれのロープがL[i]の長さを持つ。1つのロープの自由度を想像してみると、半径L[i]の半球をなすことが分かる。
よってN個の半球が作れる。バルーンが飛ぶことのできる位置は、N個全ての半球の交わる領域に限られる。
その領域のうち、高さが最大となる位置はどこかを調べたい。

最大となる高さを二分探索でZ = Hと決め打ちしよう。すると、全ての半球はZ = Hの平面で切断されるか、平面に接触するかとなる。すなわち、Z = Hの平面上には半径0以上のN個の半球の切断面となる円が描かれる。各々の円の半径は、ロープの長さと切断する平面の高さZ = Hにより三平方の定理を用いれば求まる。中心座標は高さが変化しただけなので、杭の位置と変わらない。

Z = H平面上で切断面の円が交差するならば、3次元平面上でも領域が交差するので、問題を3次元から2次元平面上の問題に帰着できた。あとは、全ての円が交わる領域が存在するか調べればよい。その境界となる高さが答えとなる。

全ての円が交わるためには、まず任意の2円が、交わっている or 一方が他方を包含している必要がある。言い換えると、2円が外側に離れているケースがあれば不適となる。

しかし、これは交差領域の存在判定の必要条件であるが、十分条件は満たさない。実際にサンプルの一部が合うが、サンプル出力より少し高い位置を示す出力が得られる。

3円が互いに交わっていても、間に隙間ができる場合がある。この場合はすべての円が交わらない領域が生じてしまう。

よって、任意の2円の交点のうち少なくともどちらか一方を、異なる任意の円が包含している必要がある。
このとき2円が（互いに外側に離れていることはないが）包含関係にあることによって交点を生成しない場合があるので、そのような場合を省いて交わる場合のみ考えることに注意する。

任意の3円が交わる領域が存在すれば、4円以上を考慮する必要はないので、これで十分となる。

int N;
double x[11], y[11], l[11];
 
bool compare(double h) {
  vector<circle> cs;
  rep(i, N) {
    if(l[i] < h - EPS) return false;
    cs.push_back({{x[i], y[i]}, sqrt(sq(l[i]) - sq(h))});
  }
  rep(i, N) REP(j, i+1, N) {
    if(out(cs[i], cs[j])) return false;
  }
 
  rep(i, N) REP(j, i+1, N) REP(k, j+1, N) {
    auto p = crosspoint(cs[i], cs[j]);
    auto q = crosspoint(cs[j], cs[k]);
    auto r = crosspoint(cs[k], cs[i]);
    if(!contains(cs[k], p.first) && !contains(cs[k], p.second) &&
       !contains(cs[i], q.first) && !contains(cs[i], q.second) &&
       !contains(cs[j], r.first) && !contains(cs[j], r.second) &&
       !contains(cs[i], cs[j]) && !contains(cs[j], cs[k]) && !contains(cs[k], cs[i])) return false;
  }
   
  return true;
}
 
double solve() {
  double l = 0.0, r = 1e10;
  rep(_, 150) {
    double m = (l + r) / 2.0;
    (compare(m) ? l : r) = m;
  }
  return l;
}
 
int main() {
 
  for(; cin >> N && N;) {
    rep(i, N) {
      cin >> x[i] >> y[i] >> l[i];
    }
    printf("%.10f\n", solve());
  }
   
  return 0;
}