算額（その1357）

算額（その1357）

十八 群馬県高崎市八幡町 八幡宮 文化7年(1810)
群馬県和算研究会:群馬の算額，上武印刷株式会社，高崎市，1987年3月31日.
キーワード：円5個，正方形
#Julia, #SymPy, #算額, #和算

正方形の中に，大円 2 個，小円 3 個を容れる。正方形の一辺の長さが 5 寸のとき，小円の直径はいかほどか。

正方形の一辺の長さを 2a
大円の半径と中心座標を r1, (a - r1, a - r1)
小円の半径と中心座標を r2, (a - r2, r2 - a)
とおき，以下の連立方程式を解く。

include("julia-source.txt");

using SymPy
@syms a::positive, r1::positive, r2::positive
eq1 = 2(r1 - a)^2 - (r1 - r2)^2
eq2 = (a - r1 - (a - r2))^2 + (r2 - a - (a - r1))^2 - (r1 + r2)^2
res = solve([eq1, eq2], (r2, r1))[4];  # 4 of 4

r2 = res[1] |> simplify
r2 |> println

   a*(-sqrt(2) + 2*sqrt(sqrt(2) + 2))*(-2 - sqrt(2) + 2*sqrt(sqrt(2) + 2))/2

小円の半径は，a*(2sqrt(√2 + 2) - √2)*(2sqrt(√2 + 2) - 2 - √2)/2
a = 5 のとき，小円の直径は 1.60435348784316

r1 = res[2] |> simplify
r1 |> println

   a*(-2*sqrt(sqrt(2) + 2) + sqrt(2) + 3)

大円の半径は，a*(√2 + 3 - 2sqrt(√2 + 2))
a = 5 のとき，大円の直径は 3.59347716163974

正方形の一辺の長さが 5 のとき，小円の直径は 1.60435348784316，大円の直径は 3.59347716163974 である。

2*res[1](a => 5/2).evalf() |> println
2*res[2](a => 5/2).evalf() |> println

   1.60435348784316
   3.59347716163974

function draw(a, more)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   (r1, r2) = (-(-10*a^3 - 20*a^3*(-sqrt(28*sqrt(2) + 40) + 3*sqrt(2) + 5)^2 + a^3*(-sqrt(28*sqrt(2) + 40) + 3*sqrt(2) + 5)^3 + 33*a^3*(-sqrt(28*sqrt(2) + 40) + 3*sqrt(2) + 5))/(2*a^2), a*(-sqrt(28*sqrt(2) + 40) + 3*sqrt(2) + 5))
   @printf("正方形の一辺の長さが %g のとき，小円の直径は %g，大円の直径は %g である。\n", 2a, 2r2, 2r1)
   plot([a, a, -a, -a, a], [-a, a, a, -a, -a], color=:green, lw=0.5)
   circle(a - r1, a - r1, r1)
   circle(r1 - a, r1 - a, r1)
   circle(0, 0, r2, :blue)
   circle(a - r2, r2 - a, r2, :blue)
   circle(r2 - a, a - r2, r2, :blue)
   if more
       delta = (fontheight = (ylims()[2]- ylims()[1]) / 500 * 10 * 2) /3  # size[2] * fontsize * 2
       hline!([0], color=:gray80, lw=0.5)
       vline!([0], color=:gray80, lw=0.5)
       point(a - r1, a - r1, "大円:r1,(a-r1,a-r1)", :red, :left, :bottom, delta=delta, deltax=-4delta)
       point(a - r2, r2 - a, "小円:r2,(a-r2,r2-a)", :blue, :center, delta=-delta/2)
       point(a, a, "(a,a)", :green, :right, :bottom, delta=delta)
   end  
end;

draw(5/2, true)