算額（その19）

算額（その19）

円の中に大小 2 種類の円が入っている。外円の径を 17寸としたとき，小円の半径を求めよ。

福島県楢葉町 波倉稲荷神社
http://www.wasan.jp/fukusima/hakura.html

必要とする座標は図に示すとおりである。
大円の半径: r1, 小円の半径: r2, 大円・小円の中心座標 (x, y1), (x, y2)

using SymPy

@syms r1::positive, r2::positive, x::positive, y1::negative, y2::positive;

以下の 5 式を立て，解く。簡単な式もあるがそのまま書く。

eq1 = x^2 + (1 - r1 - y2)^2 - (r1 + r2)^2;
eq2 = x^2 +(1-3r1)^2 - (1 - r1)^2;
eq3 = x^2 + (1-3r1-(r1-1))^2 - 4r1^2;
eq4 = y1 - 1 + 3r1;
eq5 = y1 + r1 + r2 - y2;

res = solve([eq1, eq2, eq3, eq4, eq5], (r1, r2, x, y1, y2))

   1-element Vector{NTuple{5, Sym}}:
    (3/2 - sqrt(5)/2, -2 + sqrt(5), sqrt(-22 + 10*sqrt(5)), -7/2 + 3*sqrt(5)/2, -4 + 2*sqrt(5))

for i = 1:5
   simplify(res[1][i]) |> println
end

   3/2 - sqrt(5)/2
   -2 + sqrt(5)
   sqrt(-22 + 10*sqrt(5))
   -7/2 + 3*sqrt(5)/2
   -4 + 2*sqrt(5)

r1 = (3 - √5)/2
r2 = -2 + √5
x  = √(-22 + 10√5)
y1 = (-7 + 3√5)/2
y2 = -4 + 2√5

小円の径は外円の径を 17 寸とすると，17*(√5 - 2) = 4.013155617496427 寸である。

算額では，「答 4 寸」と書いている。

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

function draw(more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   (r1, r2, x, y1, y2) = (3/2 - sqrt(5)/2, -2 + sqrt(5), sqrt(-22 + 10*sqrt(5)), -7/2 + 3*sqrt(5)/2, -4 + 2*sqrt(5))
   println("r1 = $r1,  r2 = $r2,  y1 = $y1,  y2 = $y2")
   plot()
   circle(0, 0, 1)
   circle(0, 1 - r1, r1, :blue)
   circle( x, y1, r1, :brown)
   circle(-x, y1, r1, :brown)
   circle( x, y2, r2, :green)
   circle(-x, y2, r2, :green)
   circle(0, r1 - 1, r1, :brown)
   hline!([1 - 2r1], color=:black, linewidth=0.5)
   vline!([0], color=:black, linewidth=0.5)
   if more
       point(0, 1 - r1, "(0,1-r1)", :blue, :right)
       point(x, y2, "(x,y2)", :green)
       point(x, y1, "(x,y1)", :brown)
       point(0, r1 - 1, "(0,-1+r1)", :brown)
       plot!([x, x], [y2, y2 - r2], color=:green)
       point(x, y2 - r2, "(x,y2-r2)", :green, :center)
       hline!([0, 1 - 2r1], linestyle=:dot, linewidth=0.25)
   end
end;

算額（その18）

算額（その18）

大きな 2 つの円が同じ直線に接し，なおかつ互いに接している。更に 3 種類の半径を持つ円が図のように配置されている。
それぞれの円の半径を求めよ。

福島県伊達郡梁川町 八幡神社 文化14年
http://www.wasan.jp/fukusima/fyahata.html

一番大きい 2 つの円の半径を1 とする。3 種類の円の半径の小さい方から r1, r2, r3 とし，半径 r2 の円の中心座標を (x2, r2) とする。

using SymPy

@syms r1::positive, r2::positive, r3::positive, x2::positive;

以下の 4 式を立て，解く。

eq1 = x2^2 + (r2 - r1)^2 - (r1 + r2)^2;
eq2 = (1 - x2)^2 + (1 - r2)^2 - (1 + r2)^2;
eq3 = 1^2 + (1 - (r3 + 2r1))^2 - (1 + r3)^2;
eq4 = x2^2 + ((r3 + 2r1) - r2)^2 - (r2 + r3)^2;

res = solve([eq1, eq2, eq3, eq4], (r1, r2, r3, x2))

   1-element Vector{NTuple{4, Sym}}:
    (3/8 - sqrt(5)/8, 7/2 - 3*sqrt(5)/2, sqrt(5)/40 + 1/8, -2 + sqrt(5))

r1 = (3 - √5)/8
r2 = (7 - 3√5)/2
r3 = (5 + √5)/40
x2 = -2 + √5

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

function draw(more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   (r1, r2, r3, x2) = (3/8 - sqrt(5)/8, 7/2 - 3*sqrt(5)/2, sqrt(5)/40 + 1/8, -2 + sqrt(5))
   println("r1 = $r1,  r2 = $r2,  r3 = $r3,  x2 = $x2")
   plot()
   circle(1, 1, 1)#, beginangle=170, endangle=280)
   circle(-1, 1, 1)#, beginangle=260, endangle=370)
   circle(0, r1, r1, :blue)
   circle(x2, r2, r2, :brown)
   circle(-x2, r2, r2, :brown)
   circle(0, r3 + 2r1, r3, :green)
   hline!([0])
   if more
       point(1, 1, "(1,1)")
       point(0, r1, "r1", :blue, :right)
       point(0, r3 + 2r1, "r3+2r1", :green, :right)
       point(x2, r2, "(r2,x2)", :brown, :top)
       vline!([0], xlims=(-0.5, 1.05), ylims=(-0.05, 1.05))
   else
       plot!(xlims=(-2.05, 2.05), ylims=(-0.05, 2.05))
   end
end;