算額（その132）

算額（その132）

岐阜県大垣市西外側町 大垣八幡神社 天保年間
http://www.wasan.jp/gifu/ogakihatiman.html

外円の中に 2 個の甲円，1 個の乙円と 4 個の丙円が入っている。両脇の半円の直径は 2 つの甲円の中心間の距離に等しい。甲円，乙円，丙円の径を求めよ。

外円，甲円，乙円，丙円の半径をそれぞれ 1，r1, r2, r4 として，方程式を解く。

using SymPy

@syms r1::positive, r2::positive, r4::positive, x4::positive, y4::positive;

eq1 = (r1 + 2r2)^2 + (r1 + r2)^2 - (2r1 + r2)^2  # 甲円と半円が接する
eq2 = 2r1 + r2 - 1  # 外円の半径との関係
eq3 = x4^2 + (1 - r1 - y4)^2 - (r1 + r4)^2;  # 甲円と丙円が接する
eq4 = (r1 + 2r2 - x4)^2 + y4^2 - (r4 + r1 + r2)^2  # 半円と丙円が接する
eq5 = x4^2 + y4^2 - (r2 + r4)^2;  # 乙円と丙円が接する

res = solve([eq1, eq2, eq3, eq4, eq5], (r1, r2, r4, x4, y4))

   1-element Vector{NTuple{5, Sym}}:
    (2/5, 1/5, 6/115, 4/23, 21/115)

算額では「乙円の径を知って丙円の径を求めよ」とある。乙円の径を a とすれば，丙円の径は 6/115 × 5a = 6a / 23 である。

using Plots
using Printf

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; mark=true)
   mark && scatter!([x], [y], color=color, markerstrokewidth=0)
   annotate!(x, y, text(string, 10, position, color, vertical))
end;

function draw(more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   (r1, r2, r4, x4, y4) = float.(res[1])
   @printf("甲円の半径 = %.5f,  乙円の半径 = %.5f,  丙円の半径 = %.5f\n", r1, r2, r4)
   @printf("丙円の中心座標 = (%.5f, %.5f)\n", x4, y4)
   plot()
   θ = asin((r1 + r2) / (2r1 + r2))*180.0/pi
   (ba, ea) = round.(Int, [θ, 180 - θ])
   circle(0, 0, 1, :black, beginangle=ba, endangle=ea)
   circle(0, 0, 1, :black, beginangle=180+ba, endangle=180+ea)
   circle(0, 1 - r1, r1)
   circle(0, r1 - 1, r1)
   circle(0, 0, r2, :blue)
   circle(x4, y4, r4, :orange)
   circle(x4, -y4, r4, :orange)
   circle(-x4, y4, r4, :orange)
   circle(-x4, -y4, r4, :orange)
   circle(r1 + 2r2, 0, r1 + r2, :green, beginangle=90, endangle=270)
   circle(-r1 - 2r2, 0, r1 + r2, :green, beginangle=270, endangle=450)
   if more == true
       point(0, 1-r1, "甲:(0,1-r1) ", :red, :right)
       point(0, 0, "乙:(0,0,r2) ", :blue, :right)
       point(r1 + r2, 0, "(r1+2r2,0,r1+r2)", :green, :center, :bottom)
       point(x4, y4, "  (x4,y4,r4)", :orange, :left)
       hline!([0], color=:black, lw=0.5)
       vline!([0], color=:black, lw=0.5)
   else
       plot!(showaxis=false)
   end
end;

   甲円の半径 = 0.40000,  乙円の半径 = 0.20000,  丙円の半径 = 0.05217
   右上の丙円の中心座標 = (0.17391, 0.18261)