POJ 1031: Fence

poj.org

問題概要

座標平面上があり、原点にはランプが置かれています。

また、 N 頂点の多角形があり、その各辺に沿って高さ  h のフェンスが個置かれています。多角形の  i \ (1 \leq i \leq N) 番目の頂点の座標は点  (x _ i, y _ i) です。フェンスはランプの光を反射させたり通過させたりすることはありません。

原点からの距離が  r の点における光の強さを  \displaystyle{\frac{k}{r}} とします。ここで、 k は点の位置によらない定数です。

幅が微小な幅  dl、高さ  h の長方形に当たる光の量  dI I _ 0 | \cos \alpha | dl h (ただし  dl はその点における光の強さ、 \alpha はその平面に垂直な直線と光とその点を結ぶ直線のなす角とします。

全てのフェンスに対する当たる光の量の合計を求め、小数第 3 位以下を四捨五入して小数第 2 位までで答えてください。

制約

  •  N は整数
  •  3 \leq N \leq 100
  • 多角形の辺は原点を通らない

解法

問題文が少々曖昧ですが、頑張ってエスパーしましょう。厳密な議論は行いません。

以下、原点を  O とします。

 A(a, d),  B(b, d) ( 0 \lt d, a \lt b) の間にフェンスが張られているとき、これに当たる光の量を求めることを考えます。これができれば、適切に移動や回転を行うことで全てのフェンスをこの場合に帰着できるので、答えが得られます。

 P(x, d) ( a \leq x \leq b) に当たる光の強さは  OP = \sqrt{x ^ 2 + d ^ 2} で、 |\cos\alpha| y 軸と直線  OP のなす角に等しいので、 \displaystyle \frac {d} {\sqrt{x ^ 2 + d ^ 2}} となります。よって、当たる光の量の合計は  \displaystyle {kh \int _ {a} ^ {b} \frac {d}{x ^ 2 + d ^ 2} dx} です。

この積分の値を  I とします。 x = d \tan \theta と置換積分を行い、 \displaystyle {\arctan \frac{a}{d}, \arctan \frac {b}{d}} をそれぞれ  \theta _ 1, \theta _ 2 とすると

 \displaystyle {I = kh \int _ {\theta _ 1} ^ {\theta _ 2} \frac {d}{d ^ 2 (\tan ^ 2 \theta + 1)} \cdot \frac {d}{\cos ^ 2 \theta} d \theta = kh \int _ {\theta _ 1} ^ {\theta _ 2} d \theta = kh(\theta _ 2 - \theta _ 1)} となります。

また、 \theta _ 1, \theta _ 2 はそれぞれ直線  OA, OB y 軸となす角なので、 \theta _ 2 - \theta _ 1 は角  \angle AOB の大きさになります。

よって、各フェンスについてその両端と原点の間の角 (符号も考慮) を求め、その累積和の最大値と最小値の差 (ただし  2 \pi より大きいときは  2 \pi とする) の  kh 倍が答えになります。

実装

#include <iostream>
#include <iomanip>
#include <vector>
#include <cmath>
using namespace std;
const double PI = acos((double) -1);
int main(){
  double k, h;
  int N;
  cin >> k >> h >> N;
  vector<double> x(N + 1), y(N + 1);
  for (int i = 0; i < N; i++){
    cin >> x[i] >> y[i];
  }
  x[N] = x[0];
  y[N] = y[0];
  vector<double> t(N + 1);
  for (int i = 0; i <= N; i++){
    t[i] = atan2(y[i], x[i]);
  }
  vector<double> d(N);
  for (int i = 0; i < N; i++){
    d[i] = t[i + 1] - t[i];
    if (d[i] > PI){
      d[i] -= PI * 2;
    }
    if (d[i] < -PI){
      d[i] += PI * 2;
    }
  }
  double mx = 0, mn = 0, sum = 0;
  for (int i = 0; i < N; i++){
    sum += d[i];
    mx = max(mx, sum);
    mn = min(mn, sum);
  }
  double ans = min(mx - mn, PI * 2) * k * h;
  cout << fixed << setprecision(2) << ans << endl;
}