Help us understand the problem. What is going on with this article?

三角形の外接円の半径と中心座標(外心)を求める

More than 1 year has passed since last update.

出典

Aizu Online Judge 0010
http://judge.u-aizu.ac.jp/onlinejudge/description.jsp?id=0010&lang=jp

平面上の点$(x_1, y_1),(x_2, y_2), (x_3, y_3)$ を頂点とした三角形の外接円の中心座標$(p_x, p_y)$と半径rを出力するプログラムを作成して下さい。

※問題文では三角形の頂点の座標を$(x_1, y_1),(x_2, y_2), (x_3, y_3)$と定義していますが、式変形する時に訳がわからなくなるので$(a, b), (c, d), (e, f)$と読み替えます。

外心座標の計算

Wikipediaに「外心の位置」という項目がありますのでこちらが理解できる方はそのほうがたぶん早いです。私は理解できなかったので別の方法を取りました。

まず、3つの頂点からまでの外心までの距離は外接円の半径であり、全て同じであることに注目しました(下図の赤線は全て同じ長さ)。

すると、各頂点を中心とした3つの円の交点が外心座標であることに気付きます。

つまり、3つの円の方程式が並んだ3連立方程式を解けば外心座標を求められそうです。

\begin{eqnarray}
  \left\{
    \begin{array}{l}
      (p_x - a)^2 + (p_y - b)^2 = r^2 (1)\\
      (p_x - c)^2 + (p_y - d)^2 = r^2 (2)\\ 
      (p_x - e)^2 + (p_y - f)^2 = r^2 (3)\\
    \end{array}
  \right.
\end{eqnarray}

あとはゴリ押しで解くだけです。とりあえず左辺を展開します。

\begin{eqnarray}
  \left\{
    \begin{array}{l}
      p_x^2 - 2ap_x + a^2 + p_y^2 - 2bp_y + b^2 = r^2 (1')\\
      p_x^2 - 2cp_x + c^2 + p_y^2 - 2dp_y + d^2 = r^2 (2')\\
      p_x^2 - 2ep_x + e^2 + p_y^2 - 2fp_y + f^2 = r^2 (3')\\
    \end{array}
  \right.
\end{eqnarray}

$(1') - (2')$を計算すると$p_x$が求まります。

\begin{align}
  -2ap_x + 2cp_x + a^2 - c^2 - 2bp_y + 2dp_y + b^2 - d^2 &= 0 \\
  2(c - a)p_x + a^2 - c^2 + 2(d - b)p_y + b^2 -d^2   &= 0 \\
  2(c - a)p_x &= -2(d - b)p_y - a^2 - b^2 + c^2 + d^2     \\
  p_x &= \frac{2(b - d)p_y - a^2 - b^2 + c^2 + d^2}{2(c - a)} (4)
\end{align}

同様に$(1') - (3')$を計算すると$p_x$がもう1パターン求まります。

\begin{align}
  p_x &= \frac{2(b - f)p_y - a^2 - b^2 + e^2 + f^2}{2(e - a)} (5)
\end{align}

$(4)=(5)$なので$p_y$が求まります。

\begin{align}
  \frac{2(b - d)p_y - a^2 - b^2 + c^2 + d^2}{2(c - a)} &= \frac{2(b - f)p_y - a^2 - b^2 + e^2 + f^2}{2(e - a)} \\
  2(e - a)(b - d)p_y - (e - a)(a^2 + b^2 - c^2 - d^2) &= 2(c - a)(b - f)p_y - (c - a)(a^2 + b^2 - e^2 - f^2) \\
  2(e - a)(b - d)p_y - 2(c - a)(b - f)p_y &= (e - a)(a^2 + b^2 - c^2 - d^2) - (c - a)(a^2 + b^2 - e^2 - f^2) \\
  (2(e - a)(b - d) - 2(c - a)(b - f))p_y &= (e - a)(a^2 + b^2 - c^2 - d^2) - (c - a)(a^2 + b^2 - e^2 - f^2) \\
  p_y &= \frac{(e - a)(a^2 + b^2 - c^2 - d^2) - (c - a)(a^2 + b^2 - e^2 - f^2)}{2(e - a)(b - d) - 2(c - a)(b - f)} (6)
\end{align}

$p_y$が求まりましたので、$(6)$を$(4)$あるいは$(5)$に代入すれば$p_x$も求まります。

まとめると以下になります。

\begin{align}
  p_y &= \frac{(e - a)(a^2 + b^2 - c^2 - d^2) - (c - a)(a^2 + b^2 - e^2 - f^2)}{2(e - a)(b - d) - 2(c - a)(b - f)} \\
  p_x &= \left\{
    \begin{array}{ll}
      \frac{2(b - d)p_y - a^2 - b^2 + c^2 + d^2}{2(c - a)} & ((c-a) \neq 0) \\
      \frac{2(b - f)p_y - a^2 - b^2 + e^2 + f^2}{2(e - a)} & ((c-a) = 0)
    \end{array}
  \right.
\end{align}

外接円の半径

外心座標$p_x$と$p_y$がわかりましたので各頂点からのユーグリッド距離を求めれば良いです。

\begin{align}
  r = \sqrt{(p_x - a)^2 + (p_y - b)^2}
\end{align}

プログラム例(Python3)

import math

n = int(input())

for i in range(n):
    points = list(map(float, input().split()))
    a = points[0]
    b = points[1]
    c = points[2]
    d = points[3]
    e = points[4]
    f = points[5]

    aa = a * a
    bb = b * b
    cc = c * c
    dd = d * d
    ee = e * e
    ff = f * f

    py = ((e - a) * (aa + bb - cc - dd) - (c - a) * (aa + bb - ee- ff)) / (2 * (e - a)*(b - d) - 2 * (c - a) * (b - f))

    px = (2 * (b - f) * py - aa - bb + ee + ff) / (2 * (e - a)) \
        if (c == a) else (2 * (b - d) * py - aa - bb + cc + dd) / (2 * (c - a))

    r = math.hypot(px - a, py - b)

    print("%.3f %.3f %.3f" % (px, py, r))

蛇足

なんだかゴリ押し過ぎる。$p_x$を求めるにも条件分岐が必要だし美しくないですね・・・他にもっときれいな解法がありそう。
まあ解けたのでよしとするか。

Why not register and get more from Qiita?
  1. We will deliver articles that match you
    By following users and tags, you can catch up information on technical fields that you are interested in as a whole
  2. you can read useful information later efficiently
    By "stocking" the articles you like, you can search right away
Comments
Sign up for free and join this conversation.
If you already have a Qiita account
Why do not you register as a user and use Qiita more conveniently?
You need to log in to use this function. Qiita can be used more conveniently after logging in.
You seem to be reading articles frequently this month. Qiita can be used more conveniently after logging in.
  1. We will deliver articles that match you
    By following users and tags, you can catch up information on technical fields that you are interested in as a whole
  2. you can read useful information later efficiently
    By "stocking" the articles you like, you can search right away