算額（その1307）

算額（その1307）

百三十四 群馬県富岡市一ノ宮 貫前神社 明治20年(1887)
群馬県和算研究会:群馬の算額，上武印刷株式会社，高崎市，1987年3月31日.
キーワード：円6個，正三角形

正三角形の内外に甲円，乙円を容れ，両円を包括する円（無名）の上に正三角形の頂点を通る天円を置く。また，無名円と正三角形の間に丙円を置く。天円の直径が 6 寸，乙円の直径が 1 寸のとき，丙円の直径はいかほどか。

正三角形の一辺の長さを 2a
無名円の半径と中心座標を r0, (0, √3a - 2r1 - r0)
天円の半径と中心座標を r1, (0, √3a - r1)
甲円の半径と中心座標を r2, (0, r2)
乙円の半径と中心座標を r3, (0, -r3)
丙円の半径と中心座標を r4, (x4, r4)
とおき，以下の連立方程式を解く。
SymPy の能力上，一度に解くことができないので，逐次解いていく。

include("julia-source.txt");

using SymPy

@syms a::positive, r0::positive,
     r1::positive, r2::positive,
     r3::positive, r4::positive, x4::positive;
r0 = r2 + r3
eq1 = 2r1 + 2r2 - √Sym(3)a
eq2 = x4^2 + (√Sym(3)a - 2r1 - r0 - r4)^2 - (r0 + r4)^2
eq3 = dist2(a, 0, 0, √Sym(3)a, 0, √Sym(3)a - 2r1 - r0, r0)
eq4 = dist2(a, 0, 0, √Sym(3)a, x4, r4, r4);

ans_r2 = solve(eq3, r2)[1]
ans_r2 |> println

   2*r1 - r3

eq11 = eq1(r2 => ans_r2)
ans_a = solve(eq11, a)[1]
ans_a |> println

   2*sqrt(3)*(3*r1 - r3)/3

eq12 = eq2(r2 => ans_r2, a => ans_a) |> simplify
ans_x4 = solve(eq12, x4)[1]
ans_x4 |> println

   2*sqrt(2*r1 - r3)*sqrt(r3 + r4)

eq14 = eq4(r2 => ans_r2, a => ans_a, x4 => ans_x4) |> simplify
ans_r4 = solve(eq14, r4)[2]  # 2 of 3
ans_r4 |> println

   -4*sqrt(2)*sqrt(r1)*sqrt(2*r1 - r3)/3 + 10*r1/3 - 4*r3/3

eq1, eq2, eq3 を連立させて a, r2, x4 を同時に解くことはできる。前述の最終段階と同じく，得られた解を eq4 に代入して解けば r4 が得られる。

(ans_a, ans_r2, ans_x4) = solve([eq1, eq2, eq3], (a, r2, x4))[1]

   (2*sqrt(3)*(3*r1 - r3)/3, 2*r1 - r3, sqrt(8*r1*r3 + 8*r1*r4 - 4*r3^2 - 4*r3*r4))

ans_x4 は以下のように簡約化できる。

ans_x4 |> factor |> println

   2*sqrt(2*r1 - r3)*sqrt(r3 + r4)

最後に r4 を求める。

eq14 = eq4(r2 => ans_r2, a => ans_a, x4 => ans_x4) |> simplify
ans_r4 = solve(eq14, r4)[2]  # 2 of 3
ans_r4 |> println

   -4*sqrt(2)*sqrt(r1)*sqrt(2*r1 - r3)/3 + 10*r1/3 - 4*r3/3

いずれにせよ，天円の直径が 6 寸，乙円の直径が 1 寸のとき，丙円の直径は 3.34783294256526 寸である。

ans_r4(r1 => 6/2, r3 => 1/2).evalf() * 2 |> println

   3.34783294256526

すべてのパラメータは以下のとおりである。

   r1 = 3;  r3 = 0.5;  r2 = 5.5;  a = 9.81495;  r4 = 1.67392;  x4 = 6.91565;  r0 = 6

「答」では「丙円径三寸」となっている？
天円の直径が 6，乙円の直径が 1.55 のとき，丙円の直径は 3.00238 になる。

function draw(r1, r3, more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   r2 = 2r1 - r3
   a = 2√3(3r1 - r3)/3
   r4 = 10r1/3 - 4r3/3 - 4√2*sqrt(r1)*sqrt(2r1 - r3)/3
   x4 = 2sqrt(2r1 - r3)*sqrt(r3 + r4)
   r0 = r2 + r3
   @printf("天円の直径が %g，乙円の直径が %g のとき，丙円の直径は %g である。\n", 2r1, 2r3, 2r4)
   @printf("r1 = %g;  r3 = %g;  r2 = %g;  a = %g;  r4 = %g;  x4 = %g;  r0 = %g\n", r1, r3, r2, a, r4, x4, r0)
   plot([a, 0, -a, a], [0, √3a, 0, 0], color=:green, lw=0.5)
   circle(0, √3a - r1, r1)
   circle(0, r2, r2, :blue)
   circle(0, √3a - 2r1 - r0, r0, :magenta)
   circle(0, - r3, r3, :orange)
   circle2(x4, r4, r4, :brown)
   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, 0, "a", :green, :left, :bottom, delta=delta/2, deltax=delta/2) 
       point(0, √3a, "√3a", :green, :left, :bottom, delta=delta/2, deltax=delta/2)
       point(0, √3a - r1, "天円:r1\n(0,√3a-r1)", :red, :center, delta=-delta)
       point(0, r2, "甲円:r2,(0,r2)", :blue, :center,:bottom, delta=delta)
       point(0, -r3, "乙円:r3,(0,-r3)", :black, :left, :vcenter, deltax=2delta)
       point(x4, r4, "丙円:r4\n(x4,r4)", :brown, :center, deltax=-delta)
       point(0, √3a - 2r1 - r0, "無名円:r0,(0,√3a-2r1-r0)", :magenta, :center, delta=-delta)
   end
end;

draw(6/2, 1/2, true)

#SymPy #算額 #和算