算額（その222）
中村信弥「改訂増補 長野県の算額」
http://www.wasan.jp/zoho/zoho.html
県内の算額3（200）
長野県北佐久郡軽井沢町峠 熊野神社 安政4年(1857)
外円の中に大円，中円，小円とそれぞれに外接する正三角形が入っている。大円は中円と小円にも外接している。大円と正三角形の1辺の長さが等しく 625 寸であるとき，中円と小円の径はいかほどか。
正三角形の底辺と，大円の中心は外円の直径上にある。
外円の半径を r0 とする。
大円の中心座標と半径が (-r0, 0), r1
中円の中心座標と半径が (x2, -r2), r2
小円の中心座標と半径が (0, r0 - r3), r3
とする。
以下の連立方程式を解く。
using SymPy
function distance(x1, y1, x2, y2, x0, y0)
p1, p2 = sympy.Point(x1, y1), sympy.Point(x2, y2)
l = sympy.Line(p1, p2)
l.distance(sympy.Point(x0, y0))^2 # 二乗距離を返す！！
end;
@syms r0::positive, r1::positive, r2::positive, r3::positive, x2::positive;
r0 = 625
r1 = r0 // 2
eq1 = (x2 + r1)^2 + r2^2 - (r1 + r2)^2
eq2 = r1^2 + (r0 -r3)^2 - (r1 + r3)^2
eq3 = distance(0, 0, r0/2, sqrt(Sym(3))/2*r0, 0, r0 - r3) - r3^2
eq4 = x2^2 + r2^2 - (r0 - r2)^2;
res = solve([eq1, eq4, eq3], (x2, r2, r3))
1-element Vector{Tuple{Sym, Sym, Sym}}:
(625/3, 2500/9, 625/3)
中円の直径 = 555.556, 小円の直径 = 416.667; 中円の中心座標 = (208.333, -277.778)
中円の径は 2\*2500/9 = 555 + 5/9 である。術では (505 + 2/9)寸となっているが，転記ミスであろう。
小円の径は 2\*625/3 = 416 + 2/3 である。術では正しく (416 + 2/3)寸となっている。
using Plots
using Printf
function circle(ox, oy, r, color=:red; beginangle=0, endangle=360, fill=false)
θ = beginangle:0.1:endangle
x = r.*cosd.(θ)
y = r.*sind.(θ)
if fill
plot!(ox .+ x, oy .+ y, linecolor=color, linewidth=0.5, seriestype=:shape, fillcolor=color)
else
plot!(ox .+ x, oy .+ y, color=color, linewidth=0.5)
end
end;
function point(x, y, string="", color=:green, position=:left, vertical=:top; mark=true)
mark && scatter!([x], [y], color=color, markerstrokewidth=0)
annotate!(x, y, text(string, 10, position, color, vertical))
end;
function draw(more)
pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
r0 = 625
r1 = r0 // 2
(x2, r2, r3) = res[1]
@printf("中円の直径 = %.3f, 小円の直径 = %.3f; 中円の中心座標 = (%.3f, %.3f)\n", 2r2, 2r3, x2, -r2)
plot([r0, r0/2, 0, r0], [0, √3r0/2, 0, 0], color=:orange, lw=0.5)
circle(0, 0, r0, :blue)
circle(-r1, 0, r1, :magenta)
circle(x2, -r2, r2)
circle(0, r0 - r3, r3, :green)
if more == true
point(r0, 0, " r0", :blue, :left, :bottom)
point(0, r0 - r3, " ro-r3")
point(0, r0, " r0", :blue, :left, :bottom)
point(-r1, 0, " -r1", :magenta, :left, :bottom)
point(r0/2, √3*r0/2, "(r0/2,√3/2*r0)", :orange, :left, :bottom)
point(x2, -r2, "(x2,-r2)", :red)
hline!([0], color=:black, lw=0.5)
vline!([0], color=:black, lw=0.5)
else
plot!(showaxis=false)
end
end;
