算額(その1242)
(2) 兵庫県三田市藍本大丸 酒垂神社 文化8年(1811)
近畿数学史学会:近畿の算額「数学の絵馬を訪ねて」,平成4年5月16日 初版第一刷,大阪教育図書株式会社,大阪市.
キーワード:円2個,菱形,扇
地紙長が 6.4 寸,菱形の横の対角線の長さが 14.664 寸のとき,等円の直径はいかほどか。
扇長(要から先端までの長さ)を R
地紙長(紙の貼られている部分の長さ)を k
菱形の横の対角線の長さを 2a
等円の半径と中心座標を r, (x, y)
とおき,以下の連立方程式を解く。
include("julia-source.txt");
using SymPy
@syms R::positive, r::positive, x::positive, y::positive,
a::positive, k::positive
eq1 = x^2 + y^2 - ((R - k) + r)^2
eq2 = dist2(0, R - k, a, sqrt(R^2 -a^2), x, y, r)
eq3 = dist2(0, 0, a, sqrt(R^2 -a^2), x, y, r)
eq4 = sqrt(R^2 -a^2) - (R - k/2)
using NLsolve
function nls(func, params...; ini = [0.0])
if typeof(ini) <: Number
r = nlsolve((vout, vin) -> vout[1] = func(vin[1], params..., [ini]), ftol=big"1e-40")
v = r.zero[1]
else
r = nlsolve((vout, vin)->vout .= func(vin, params...), ini, ftol=big"1e-40")
v = r.zero
end
return Float64.(v), r.f_converged
end;
function H(u)
(R, r, x, y) = u
return [
x^2 + y^2 - (R - k + r)^2, # eq1
4*R^2*a^2 - 8*R*a^2*k - 8*R*a^2*y + 4*R*a*k*x + 4*a^2*k^2 + 8*a^2*k*y - 4*a^2*r^2 + 4*a^2*y^2 - 4*a*k^2*x - 4*a*k*x*y - k^2*r^2 + k^2*x^2, # eq2
-4*R^2*r^2 + 4*R^2*x^2 - 8*R*a*x*y + 4*R*k*r^2 - 4*R*k*x^2 - 4*a^2*r^2 + 4*a^2*y^2 + 4*a*k*x*y - k^2*r^2 + k^2*x^2, # eq3
-R + k/2 + sqrt(R^2 - a^2), # eq4
]
end;
a = 14.664/2
k = 6.4
iniv = BigFloat[10, 1, 2.6, 3.7]
res = nls(H, ini=iniv)
([9.999722499999999, 0.9429406670508707, 2.61703497874062, 3.7130737360479014], true)
地紙長が 6.4,菱横が 14.664 のとき,等円の直径は 1.88588 である。
その他のパラメータは以下のとおりである
a = 7.332; k = 6.4; R = 9.99972; r = 0.942941; x = 2.61703; y = 3.71307
function draw(a, k, more=false)
pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
(R, r, x, y) = res[1]
@printf("地紙長が %g,菱横が %g のとき,等円の直径は %g である。\n", k, 2a, 2r)
@printf("a = %g; k = %g; R = %g; r = %g; x = %g; y = %g\n", a, k, R, r, x, y)
θ = atand(R - k/2, a)
plot([-a, 0, a, 0, -a], [R - k/2, R, R - k/2, R - k, R - k/2], color=:green, lw=0.5)
plot!([-a, 0, a], [R - k/2, 0, R - k/2], color=:red, lw=0.5)
circle(0, 0, R, beginangle=θ, endangle=180 - θ)
circle(0, 0, R - k, :blue, beginangle=θ, endangle=180 - θ)
circle2(x, y, r, :magenta)
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, R, "R", :red, :center, :bottom, delta=delta)
point(0, R - k/2, "R-k/2", :green, :center, :bottom, delta=delta)
point(a, R - k/2, "(a,R-k/2) ", :green, :right, :vcenter)
point(0, R - k, "R-k", :green, :center, :bottom, delta=delta)
point(x, y, "r,(x,y)", :magenta, :center, delta=-delta/2)
dimension_line(2delta, -2delta, a + 2delta, R - k/2 - 2delta, " 扇長", :red, :left)
dimension_line(a*(R - k)/R + 4delta, (R - k/2)*(R - k)/R - 4delta, a + 4delta, R - k/2 - 4delta, " 地紙長", :green, :left)
end
end;
draw(14.664/2, 6.4, true)