算額（その1004）

算額（その1004）

二十 武州不動岡村 不動堂 文政元年(1818)

埼玉県立図書館：埼玉県史料集 第二集『埼玉の算額』，昭和44年，誠美堂印刷所，埼玉県与野市.

楕円の中に甲円，乙円が 2 個ずつ入っている。楕円の長径，短径がそれぞれ 25 寸，15 寸のとき，甲円の直径はいかほどか。

楕円の長半径，短半径，中心座標を a, b, (0, 0)
甲円の半径と中心座標を r1, (x1, 0)
乙円の半径と中心座標を r2, (0, R - r2)
とおき，以下の連立方程式を解く。

include("julia-source.txt");

using SymPy

@syms a::poitive, b::poitive,
     r1::poitive, x1::poitive, r2::poitive
r2 = b/2
eq1 = x1^2 + r2^2 - (r1 + r2)^2
eq2 = (a^2 - b^2)*(b^2 - r1^2)/b^2 - x1^2  # 「算法助術」公式84
res = solve([eq1, eq2], (r1, x1))[3]

   ((-a^2*b^2 + b^4 + b^2*(a - b)*(a + b)*(2*a^2 - b^2)/a^2)/(b*(a - b)*(a + b)), b*sqrt((a - b)*(a + b)*(2*a^2 - b^2))/a^2)

r1 の式は簡略化できる。

楕円の長径，短径がそれぞれ a, b のとき甲円の直径は (b - b^3/a^2) である。

res[1] |> simplify |> println

   b - b^3/a^2

楕円の長径，短径がそれぞれ 25 寸，15 寸のとき甲円の直径は 9.6 寸である。

(a, b) = (25, 15) ./ 2
(b - b^3/a^2, b*sqrt((a - b)*(a + b)*(2*a^2 - b^2))/a^2)

   (4.8, 7.683749084919418)

なお，例としては「楕円の長径，短径がそれぞれ 8, 4 のとき，甲円の直径は 3 である（乙円の直径は 2)」などがきれいである（上の図はこのときのもの）。

function draw(more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   (a, b) = (25, 15) ./ 2
   r2 = b/2
   (r1, x1) = (b - b^3/a^2, b*sqrt((a - b)*(a + b)*(2*a^2 - b^2))/a^2)
   @printf("楕円の長径，短径がそれぞれ %g, %g のとき，甲円の直径は %g である（乙円の直径は %g)。\n", 2a, 2b, 2r1, b)
   @printf("a = %g;  b = %g;  r2 = %g;  r1 = %g;  x1 = %g\n", a, b, r2, r1, x1)
   plot()
   ellipse(0, 0, a, b, color=:red)
   circle2(x1, 0, r1, :blue)
   circle22(0, r2, r2, :green)
   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", :red, :left, :bottom, delta=delta)
       point(0, b, " b", :red, :left, :bottom, delta=delta)
       point(x1, 0, "甲円:r1,(x1,0)", :blue, :center, delta=-delta)
       point(0, r2, "乙円:r2,(0,r2)", :green, :center, delta=-delta)
   end
end;