算額（その1145） - 算額あれこれ

一〇一埼玉県大宮市高鼻町氷川神社明治31年(1898)

埼玉県立図書館：埼玉県史料集第二集『埼玉の算額』，昭和44年，誠美堂印刷所，埼玉県与野市.
キーワード：円2個，楕円，長方形，斜線
#Julia #SymPy #算額 #和算 #数学

長方形の中に2本の斜線を引き，大円 2 個と楕円 2 個を容れる。楕円の長径が 4 寸，大円の直径が 8 寸のとき，楕円の短径はいかほどか。

長方形の長辺の長さを \(2c\), （短辺の長さは \(2r\)）
斜線と長方形の長辺との交点の座標を \( (d,\ 0)\)
大円の半径と中心座標を \(r,\ (c - r,\ r)\)
楕円の長半径，短半径，中心座標を \(a,\ b,\ (0,\ b),\ (0,\ 2b + a)\)
楕円と斜線の接点座標を \( (x_1,\ y_1),\ (x_2,\ y_2)\)
とおき，以下の連立方程式を解く。

include("julia-source.txt"); # julia-source.txt ソース

using SymPy

@syms a::positive, b::positive, r::positivr,
c::positive, d::positive,
x1::positive, y1::positive,
x2::positive, y2::positive
@syms a, b, r, c, d, x1, y1, x2, y2
slope = -2r/d
eq1 = x1^2/a^2 + (y1 - b)^2/b^2 - 1
eq2 = -b^2*x1/(a^2*(y1 - b)) - slope
eq3 = y1/(d - x1) + slope

eq4 = x2^2/b^2 + (y2 - a - 2b)^2/a^2 - 1
eq5 = -a^2*x2/(b^2*(y2 - a - 2b)) - slope
eq6 = y2/(d - x2) + slope

eq7 = dist2(0, 2r, d, 0, c - r, r, r);

まず，eq1, eq2, eq3 を解き，\(b,\ x_1,\ x_2\) を求める。

res1 = solve([eq1, eq2, eq3], (b, x1, y1))[1]

(-a^2*r/d^2 + r, 2*a^2*d/(a^2 + d^2), 2*r*(-a^2 + d^2)/(a^2 + d^2))

eq4, eq5 に b を代入し，eq6, eq7 とともに連立方程式を解き，\(x_2,\ y_2,\ c,\ d\) を求める。

eq4 = eq4(b => res1[1]);
eq5 = eq5(b => res1[1]);
res2 = solve([eq4, eq5, eq6, eq7], (x2, y2, c, d))[8] # 8 of 8

(2*sqrt(a^2*r^2/(a^2 - a*r + r^2))*(a^4 - 3*a^3*r + 4*a^2*r^2 - 3*a*r^3 + r^4)/(r^2*(2*a^2 - 3*a*r + 2*r^2)), 2*a*(-2*a^3 + 6*a^2*r - 6*a*r^2 + 3*r^3)/(r*(2*a^2 - 3*a*r + 2*r^2)), r*sqrt( (5*a^2 - 4*a*r + 4*r^2)/(a^2 - a*r + r^2))/2 + r + sqrt(a^2*r^2/(a^2 - a*r + r^2))/2, sqrt(a^2*r^2/(a^2 - a*r + r^2)))

短半径は最初の連立方程式の解で \(-a^2 r/d^2 + r\) となっており，まだ未知数 \(d\) を含んでいる。
\(d\) は次の連立方程式の解で \(\sqrt{a^2 r^2/(a^2 - a r + r^2)}\) であることがわかる。
そこで，\(d\) を短半径の式に代入すれば \(a,\ r\) のみを含む式 \(a(r - a)/r = (a r - a^2)/r \) になる。これは「術」と同じである。
長半径が \(a = 4/2\)，大円の半径が \(r = 8/2\) のとき，短半径は \(a(-a + r)/r = 1\) である（短径は 2）。

res1[1](d => res2[4]) |> simplify |> println

a*(-a + r)/r

術は以下のようになっている。

長径 = 4
大径 = 8
(長径*大径 - 長径^2)/大径

2.0

その他のパラメータは以下の通りである。

\(a = 2;\ r = 4;\ b = 1;\ x_1 = 1.97949;\ y_1 = 1.14286\)
\(x_2 = 0.866025;\ y_2 = 5;\ c = 9.31803;\ d = 2.3094\)

描画関数プログラムのソースを見る

function draw(a, r, more=false)
pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
(x2, y2, c, d) = (2*sqrt(a^2*r^2/(a^2 - a*r + r^2))*(a^4 - 3*a^3*r + 4*a^2*r^2 - 3*a*r^3 + r^4)/(r^2*(2*a^2 - 3*a*r + 2*r^2)), 2*a*(-2*a^3 + 6*a^2*r - 6*a*r^2 + 3*r^3)/(r*(2*a^2 - 3*a*r + 2*r^2)), r*sqrt( (5*a^2 - 4*a*r + 4*r^2)/(a^2 - a*r + r^2))/2 + r + sqrt(a^2*r^2/(a^2 - a*r + r^2))/2, sqrt(a^2*r^2/(a^2 - a*r + r^2)))
(b, x1, y1) = (-a^2*r/d^2 + r, 2*a^2*d/(a^2 + d^2), 2*r*(-a^2 + d^2)/(a^2 + d^2))
@printf("楕円の長径が %g，大円の直径が %g のとき，楕円の短径は %g である。\n", 2a, 2r, 2b)
@printf("a = %g; r = %g; b = %g; x1 = %g; y1 = %g; x2 = %g; y2 = %g; c = %g; d = %g\n", a, r, b, x1, y1, x2, y2, c, d)
plot([c, c, -c, -c, c], [0, 2r, 2r, 0, 0], color=:blue, lw=0.5)
plot!([-d, 0, d], [0, 2r, 0], color=:magenta, lw=0.5)
circle2(c - r, r, r)
ellipse(0, b, a, b, color=:green)
ellipse(0, a + 2b, b, a, color=: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(c, 0, " c", :blue, :left, :bottom, delta=delta/2)
point(d, 0, " d", :magenta, :left, :bottom, delta=delta/2)
point(c - r, r, "大円:r,(c-r,r)", :red, :center, delta=-delta/2)
point(0, 2b, "2b", :green, :center, delta=-2delta)
point(0, 2b + 2a, "2b+2a", :black, :center, :bottom, delta=delta)
end
end;