宮城県栗原市瀬峰泉谷 瀬峰泉谷熊野神社 天保6年(1835) 昭和61年(1986) 復元奉納
徳竹亜紀子，谷垣美保，萬伸介：瀬峰泉谷熊野神社奉納算額をめぐる諸問題，仙台高等専門学校名取キャンパス 研究紀要 第60号(2024)
https://www.jstage.jst.go.jp/article/nitsendainatori/60/0/60_1/_pdf/-char/ja
キーワード：円2個，楕円，斜線
#Julia #SymPy #算額 #和算 #数学
楕円内に等円 2 個が入っている。楕円の短径が 1 のとき，原点を点対称とする円と楕円の 2 接点間の距離を求めよ。
条件が「楕円の短径」のみなので，解けるかどうか不安であったが，「楕円の長径」は無関係であることがわかった。
楕円の長半径と短半径，中心座標を \(a,\ b,\ (0,\ 0)\)
円の半径と中心座標を \(r,\ (r,\ 0)\)
円と，楕円の接点を \( (x_0,\ y_0)\)
とおき，以下の連立方程式を解く。
include("julia-source.txt"); # julia-source.txt ソース
using SymPy
@syms a::positive, b::positive, x0::positive, y0::positive, r::positive
eq1 = x0^2/a^2 + y0^2/b^2 - 1
eq2 = -b^2*x0/(a^2*y0) + (x0 - r)/y0
eq3 = (x0 - r)^2 + y0^2 - r^2
res = solve([eq1, eq2, eq3], (x0, y0, r))
2-element Vector{Tuple{Sym{PyCall.PyObject}, Sym{PyCall.PyObject}, Sym{PyCall.PyObject}}}:
(a*b/sqrt(a^2 - b^2), -b*sqrt( (a^2 - 2*b^2)/(a - b))/sqrt(a + b), b*sqrt(a^2 - b^2)/a)
(a*b/sqrt(a^2 - b^2), b*sqrt( (a^2 - 2*b^2)/(a - b))/sqrt(a + b), b*sqrt(a^2 - b^2)/a)
2 組の解が得られるが，符号の違いのみなのでどちらでもよい。
接点間の距離 \(diag\) は以下のようになり，点対称な接点間の距離は，楕円の短径のみで決まる。常に，楕円の短径の \(\sqrt{2}\) 倍である。
すなわち，問のように短径が 1 寸ならば，接点間の距離は 1.41421356... である。
@syms d
diag = 2sqrt(res[2][1]^2 + res[2][2]^2)
diag |> println
2*sqrt(a^2*b^2/(a^2 - b^2) + b^2*( (a^2 - 2*b^2)/(a - b))/(a + b))
長い式は，以下のように簡約化される。
apart(diag, d) |> simplify |> println
2*sqrt(2)*b
\(b\) に特定の値を代入して数値解を求める。
apart(diag, d)(b => 1/2) |> simplify |> println
1.41421356237310
その他のパラメータは以下のとおり。
\(b = 1.234;\ x_0 = 1.35383;\ y_0 = 1.1012;\ r = 1.12477\)
描画関数プログラムのソースを見る
function draw(more)
pyplot(size=(500, 500), grid=false, showaxis=true, aspectratio=1, label="", fontfamily="IPAMincho")
(a, b) = (3, 1.234)
(x0, y0, r) = (a*b/sqrt(a^2 - b^2), b*sqrt( (a^2 - 2*b^2)/(a - b))/sqrt(a + b), b*sqrt(a^2 - b^2)/a)
@printf("点対称な接点間の距離 = %.15g = %.15g", 2sqrt(x0^2 + y0^2), 2√2b)
@printf("b = %g; x0 = %g; y0 = %g; r = %g\n", b, x0, y0, r)
plot()
ellipse(0, 0, a, b, color=:red)
circle2(r, 0, r, :blue)
segment(x0, y0, -x0, -y0, :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(0, b, " b", :red, :left, :bottom, delta=delta/2)
point(x0, y0, "(x0,y0)", :red, :left, :bottom, delta=delta)
point(-x0, -y0, "(-x0,-y0)", :red, :left, :bottom, delta=delta)
#plot!(xlims=(-0.1,0.1), ylims=(-0.1,0.1))
end
end;
以下のアイコンをクリックして応援してください
