算額（その76）
大阪府茨木市 井於神社 弘化3年(1846)
http://www.wasan.jp/osaka/iyo.html
一般の三角形において，辺の長い順に大斜，中斜，小斜という。それぞれが 4 寸，2 寸 7 分，1 寸 8 分であるとき，内接する円の径を求めよ。
図のように記号を定め方程式を解けば，半径は求まる。単位は「分」で表しておく。
using SymPy
@syms r::positive, s::positive, S::positive;
(大斜, 中斜, 小斜) = (40, 27, 18)
s = (大斜 + 中斜 + 小斜)//2
eq1 = (大斜 + 中斜 + 小斜)*r//2 - S
eq2 = sqrt(Sym(s) * (s - 大斜) * (s - 中斜) * (s - 小斜)) - S # ヘロンの公式
solve([eq1, eq2])
Dict{Any, Any} with 2 entries:
S => 35*sqrt(527)/4
r => 7*sqrt(527)/34
7*sqrt(527)/34 * 2 # 直径
9.452668468557999
もとの単位では 9分4厘5毛である。算額では，径は 1寸7厘5毛となっている。
どちらが正しいか，図を描いて確かめると，得られた解で正しいようである。算額における三角形とは随分異なる。
BF = BE = a, CD = CF = b, AE = AD = c, ∠CAB = α, ∠ABC = β として，プログラム中で必要な座標を独自に計算する。
using Plots
function circle(ox, oy, r, color=:red; beginangle=0, endangle=360)
θ = beginangle:0.1:endangle
x = r.*cosd.(θ)
y = r.*sind.(θ)
plot!(ox .+ x, oy .+ y, color=color, linewidth=0.5)
end;
function point(x, y, string="", color=:green, position=:left, vertical=:top; fontsize=10, mark=true)
mark && scatter!([x], [y], color=color, markerstrokewidth=0)
annotate!(x, y, text(string, fontsize, vertical, position, color))
end;
function segment(x1, y1, x2, y2, color=:black; linestyle=:solid, linewidth=0.5)
plot!([x1, x2], [y1, y2], color=color, linestyle=linestyle, linewidth=linewidth)
end;
function draw(more=false)
pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
(大斜, 中斜, 小斜) = (40, 27, 18)
α = acos((小斜^2 + 大斜^2 - 中斜^2) / (2*小斜*大斜))
β = acos((中斜^2 + 大斜^2 - 小斜^2) / (2*中斜*大斜))
a = 中斜/2 + 大斜/2 - 小斜/2
b = 中斜/2 - 大斜/2 + 小斜/2
c = -中斜/2 + 大斜/2 + 小斜/2
s = (大斜 + 中斜 + 小斜)/2
S = sqrt(s * (s - 大斜) * (s - 中斜) * (s - 小斜)) # ヘロンの公式
r = 2S/(大斜 + 中斜 + 小斜)
plot([0, 大斜, 小斜*cos(α), 0], [0, 0, 小斜*sin(α), 0],
lwd=0.25, linestyle=more ? :dot : :solid,
xlims=(-3, 大斜+3), ylims=(-2, 小斜*sin(α)+1))
circle(c, r, r)
if more
println("直径 = $(2r)")
point(0, 0, " α", :red, :bottom, :left, mark=false)
point(大斜, 0, "β ", :red, :bottom, :right, mark=false)
segment(c, 0, c, r)
point(c, r, " O", :black)
point(c*cos(α), c*sin(α), "D ", :red, :center, :right)
segment(c*cos(α), c*sin(α), c, r)
point(c, 0, " E", :red, :top, :center)
point(大斜 - a*cos(β), a*sin(β), " F", :red, :center)
segment(大斜 - a*cos(β), a*sin(β), c, r)
point(0, 0, "A ", :black, :right)
point(大斜, 0, " B", :black)
point(小斜*cos(α), 小斜*sin(α), " C", :black, :bottom)
segment(0, 0, c, 0, :green)
point(c/2, 0, "c", :green, mark=false)
segment(大斜, 0, 大斜 - a*cos(β), a*sin(β), :magenta)
point(大斜 - a*cos(β)/2, a*sin(β)/2, "a ", :magenta, :right, mark=false)
segment(小斜*cos(α), 小斜*sin(α), c*cos(α), c*sin(α), :blue)
point((小斜 + c)*cos(α)/2, (小斜+c)*sin(α)/2, "b", :blue, :bottom, :right, mark=false)
end
end;
直径 = 9.452668468557997
