算額（その1155）

算額（その1155）

九八 鴻巣市三ツ木山王 三木神社 明治28年(1895)
埼玉県立図書館：埼玉県史料集 第二集『埼玉の算額』，昭和44年，誠美堂印刷所，埼玉県与野市.
キーワード：円4個，楕円，正方形，菱形

正方形の中に楕円，菱形を 1 個ずつ，等円を 4 個容れる。楕円の長径，短径が 4 寸，3 寸のとき，等円の直径はいかほどか。

正方形の一辺の長さを 2s
楕円の長半径，短半径を a, b
楕円と菱形の一辺との接点座標を (x0, y0)
等円の半径と中心座標を r, (s - r, s - r)
とおき，以下の連立方程式を解く。

include("julia-source.txt")

using SymPy

@syms s::positive, r::positive,
     x0::positive, y0::positive, a::positive, b::positive
eq1 = dist2(s, 0, 0, s, s - r, s - r, r)
eq2 = x0^2/a^2 + y0^2/b^2 - 1
eq3 = -b^2*x0/(a^2*y0^2) + 1
eq4 = y0/(s - x0) - 1
res = solve([eq1, eq2, eq3, eq4], (r, s, x0, y0))[1]

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

@syms t
res[1](sqrt(4a^2 +  1) => t) |> println
res[2](sqrt(4a^2 +  1) => t) |> println
res[3](sqrt(4a^2 +  1) => t) |> println
res[4](sqrt(4a^2 +  1) => t) |> println

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

楕円の長半径 a，短半径 b が与えられたとき，等円の半径は
(a*t - a + sqrt(2)*b*sqrt(t - 1) - sqrt(2*a^4 - a^2*t + a^2 + sqrt(2)*a*b*t*sqrt(t - 1) - sqrt(2)*a*b*sqrt(t - 1) + b^2*t - b^2))/(2*a)
ただし，t = sqrt(4a^2 +  1) と，長々しい式になる。

楕円の長半径，短半径が 4/2, 3/2 のとき，等円の直径は 1.463744764599768 である。

「答」は 「1.461 余り」としている。

a = 4/2
b = 3/2
t = sqrt(4a^2 +  1)
r = (a*t - a + sqrt(2)*b*sqrt(t - 1) - sqrt(2*a^4 - a^2*t + a^2 + sqrt(2)*a*b*t*sqrt(t - 1) - sqrt(2)*a*b*sqrt(t - 1) + b^2*t - b^2))/(2*a)
2r

   1.463744764599768

function draw(a, b, more=false)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   t = sqrt(4a^2 +  1)
   r = (a*t - a + sqrt(2)*b*sqrt(t - 1) - sqrt(2*a^4 - a^2*t + a^2 + sqrt(2)*a*b*t*sqrt(t - 1) - sqrt(2)*a*b*sqrt(t - 1) + b^2*t - b^2))/(2*a)
   s = (a*t - a + sqrt(2)*b*sqrt(t - 1))/(2*a)
   x0 = t/2 - 1/2
   y0 = b*sqrt(2*t - 2)/(2*a)
   @printf("等円の直径は %g である。\n", 2r)
   plot([s, s, -s, -s, s], [-s, s, s, -s, -s], color=:magenta, lw=0.5)
   plot!([s, 0, -s, 0, s], [0, s, 0, -s, 0], color=:green, lw=0.5)
   circle4(s - r, s - r, r)
   ellipse(0, 0, a, b, color=: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(s - r, s - r, "等円:r,(s-r,s-r)", :red, :center, delta=-delta/2)
       point(x0, y0)
   end
end;