5-5. SymPyで３点を通る円を求める

import matplotlib.pyplot as plt
import matplotlib.patches as pat
import sympy

#赤点（動かす点）
x = 120
y =  150

#黒点（固定する２点）
x_fix = [-100, 100]
y_fix = [20, -20]

#表示関数
def show(center, r):
    plt.figure()
    ax = plt.axes()
    #動かす点の描画
    ax.plot(x, y, 'or')
    #固定点の描画
    ax.plot(x_fix, y_fix, 'ok')
    #円の描画
    e = pat.Circle(xy=center, radius=r, color='k', alpha=0.3)
    ax.add_patch(e)
    #軸の設定
    ax.set_aspect('equal')
    ax.set_xlim(-200, 200)
    ax.set_ylim(-100, 300)
    ax.spines['bottom'].set_position(('data', 0))
    ax.spines['left'].set_position(('data', 0))
    ax.spines['right'].set_visible(False)
    ax.spines['top'].set_visible(False)
    plt.show()

#確認
show((40, 50), 20)

#３点を通る円の方程式を求める
"""
中心(s, t), 半径rの円の方程式は
(x - s)**2 + (y - t)**2 = r**2

３点を通るという情報から未知数(s, t, r)を求める
"""
s, t, r = sympy.symbols('s t r') #変数を定義
eq1 = (x_fix[0] - s)**2 + (y_fix[0] - t)**2 - r**2 #１つ目の固定点を代入した式（＝０の形）
eq2 = (x_fix[1] - s)**2 + (y_fix[1] - t)**2 - r**2 #２つ目の固定点を代入した式
eq3 = (x - s)**2 + (y - t)**2 - r**2 #動く点Pを代入した式

#方程式を解く
ans = sympy.solve([eq1, eq2, eq3], [s, t, r])
print('ans:{}'.format(ans))
s, t, r = ans[0]
print('s, t, r = {}, {}, {}'.format(s, t, r))

#表示
show((s, t), r)