93 岩手県大船渡市盛町館下 根城八幡宮 天保12年(1841)
安富有恒：和算—岩手の現存算額のすべて，青磁社，東京都，1987.
http://www.wasan.jp/iwatenosangaku_yasutomi.pdf
キーワード：楕円，正六角形
#Julia #SymPy #算額 #和算 #数学
正六角形の中に等楕円を 6 個容れる。楕円の長径・短径が与えられたとき，正六角形の一辺の長さを求めよ。
計算式を簡単にするため，図形を時計方向に 30° 回転させたもので考える。\(x\) 軸上に長径を持つ楕円について以下のように記号を定める。
正六角形の一辺の長さを \(R\)
楕円の長半径，短半径，中心座標を \(a,\ b,\ (x_0,\ 0)\)
隣の楕円との接点，正三角形の一辺との接点の座標を \( (x_1,\ y_1),\ (x_2,\ y_2);\ y_1 = \sqrt{3}(R - x_1),\ y2 = \sqrt{3}x_2/3\)
とおき，以下の連立方程式を解く。
include("julia-source.txt"); # julia-source.txt ソース
using SymPy
@syms a::positive, b::positive, R::positive, x0::positive,
x1::positive, y1::positive, x2::positive, y2::positive
y1 = √Sym(3)*(R - x1)
y2 = √Sym(3)*x2/3
eq1 = (x1 - x0)^2/a^2 + y1^2/b^2 - 1
eq2 = -b^2*(x1 - x0)/(a^2*y1) + √Sym(3)
eq3 = (x2 - x0)^2/a^2 + y2^2/b^2 - 1
eq4 = -b^2*(x2 - x0)/(a^2*y2) - 1/√Sym(3);
res = solve([eq1, eq2, eq3, eq4], (R, x0, x1, x2))[2]
( (a^2 + 3*b^2 + sqrt(9*a^4 + 30*a^2*b^2 + 9*b^4)/3)/sqrt(a^2 + 3*b^2), sqrt(a^2 + 3*b^2), (a^2*sqrt(9*a^4 + 30*a^2*b^2 + 9*b^4) + (a^2 + 3*b^2)*(3*a^2 + b^2))/(sqrt(a^2 + 3*b^2)*(3*a^2 + b^2)), 3*b^2/sqrt(a^2 + 3*b^2))
各パラメータは以下のように簡約化して求めることができる。
t = a^2 + 3b^2
u = 3a^2 + b^2
x0 = sqrt(t)
R = x0 + sqrt(u/3)
x1 = (a^2*sqrt(3t*u) + t*u)/(x0*u)
x2 = 3b^2/x0;
楕円の長径・短径が 10，7 のとき，正六角形の一辺の長さは \(\sqrt{ (10/2)^2 + 3(7/2)^2} + \sqrt{ (10/2)^2 + (7/2)^2/3} = 13.251013385205335\) である。
sqrt( (10/2)^2 + 3(7/2)^2) + sqrt( (10/2)^2 + (7/2)^2/3) |> println
13.251013385205335
その他のパラメータは以下のとおりである。
\(a = 5;\ b = 3.5;\ x_0 = 7.85812;\ R = 13.251;\ x_1 = 12.4938;\ x_2 = 4.67669\)
描画関数プログラムのソースを見る
function draw(more=false)
pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
(a, b) = (10, 7)./2
t = a^2 + 3b^2
u = 3a^2 + b^2
x0 = sqrt(t)
R = x0 + sqrt(u/3)
R = sqrt(a^2 + 3b^2) + sqrt(a^2 + b^2/3)
x1 = (a^2*sqrt(3t*u) + t*u)/(x0*u)
x2 = 3b^2/x0
y1 = √3(R - x1)
y2 = √3x2/3
@printf("楕円の長径・短径が %g, %g のとき，正六角形の一辺の長さは %g である。\n", 2a, 2b, R)
@printf("a = %g; b = %g; x0 = %g; R = %g; x1 = %g; x2 = %g\n", a, b, x0, R, x1, x2)
θ = 0:60:420
x = @. R*cosd(θ)
y = @. R*sind(θ)
plot()
circle(0, 0, R)
for i = 1:6
ellipse(x0*cosd(60(i)), x0*sind(60(i)), a, b, φ=60(i), color=:green)
segment(x[i], y[i], x[i + 1], y[i + 1], :blue)
end
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(x1, y1, "(x1,y1) ", :green, :right, :vcenter)
point(x2, y2, "(x2,y2)", :green, :left, delta=-delta/2)
point(x0, 0, "楕円:a, b,(x0, 0)", :green, :center, delta=-delta/2)
end
end;
以下のアイコンをクリックして応援してください
