算額（その76）

算額（その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