算額（その172）
岐阜県大垣市外野釜笛 釜笛八幡神社 慶応元年(1865)
http://ryugen3.sakura.ne.jp/toukou3/wasankibousya.PDF
第4問: 囲円を 8 個の等円で囲む。等円 8 個は互いに外接し，囲円に接している。囲円の径を知って等円の径を求めよ。
この算額は「算額（その60）」の簡単版である。https://blog.goo.ne.jp/r-de-r/e/06710f85cfba141184ddd8f21529fae0
囲円の径 r1，等円の径 r2，右上の等円の中心座標を (x, r2) とする。図のように記号を定め方程式を解く。方程式は 2 本なので，x, r2 を r1 で表す解を求める。
using SymPy
@syms r1::positive, r2::positive, x::positive;
eq1 = x^2 + r2^2 - (r1 + r2)^2 # 等円と囲円が外接する
eq2 = 2(x- r2)^2 - 4r2^2; # 等円同士が外接する
res = solve([eq1, eq2], (x, r2)) # x, r2 は r1 を含む式で得られる
3-element Vector{Tuple{Sym, Sym}}:
(-(-12*r1^3*(-2*sqrt(2) - sqrt(20 - 14*sqrt(2)) + 3)^2 + r1^3*(-2*sqrt(2) - sqrt(20 - 14*sqrt(2)) + 3)^3 - r1^3*(-2*sqrt(2) - sqrt(20 - 14*sqrt(2)) + 3) + 2*r1^3)/(2*r1^2), r1*(-2*sqrt(2) - sqrt(20 - 14*sqrt(2)) + 3))
(-(-12*r1^3*(-2*sqrt(2) + sqrt(20 - 14*sqrt(2)) + 3)^2 - r1^3*(-2*sqrt(2) + sqrt(20 - 14*sqrt(2)) + 3) + r1^3*(-2*sqrt(2) + sqrt(20 - 14*sqrt(2)) + 3)^3 + 2*r1^3)/(2*r1^2), r1*(-2*sqrt(2) + sqrt(20 - 14*sqrt(2)) + 3))
(-(-12*r1^3*(2*sqrt(2) + 3 + sqrt(14*sqrt(2) + 20))^2 - r1^3*(2*sqrt(2) + 3 + sqrt(14*sqrt(2) + 20)) + 2*r1^3 + r1^3*(2*sqrt(2) + 3 + sqrt(14*sqrt(2) + 20))^3)/(2*r1^2), r1*(2*sqrt(2) + 3 + sqrt(14*sqrt(2) + 20)))
r2 は正で，r1 より小さくなければならないので，2組目のものが適解である。
res[1][2](r1 => 1).evalf() |> println
-0.276768653914155
res[2][2](r1 => 1).evalf() |> println
0.619914404421775
res[3][2](r1 => 1).evalf() |> println
12.1370711845441
3 -2√2 + sqrt(20 - 14√2)
0.6199144044217748
等円の径は囲円の径の 3 -2√2 + sqrt(20 - 14√2) = 0.619914404421775 倍である。
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")
r1 = 1
(x, r2) = (-(-12*r1^3*(-2*sqrt(2) + sqrt(20 - 14*sqrt(2)) + 3)^2 - r1^3*(-2*sqrt(2) + sqrt(20 - 14*sqrt(2)) + 3) + r1^3*(-2*sqrt(2) + sqrt(20 - 14*sqrt(2)) + 3)^3 + 2*r1^3)/(2*r1^2), r1*(-2*sqrt(2) + sqrt(20 - 14*sqrt(2)) + 3))
@printf("囲円（赤）の径は等円（青）の径の 3 -2√2 + sqrt(20 - 14√2) = %.6f 倍である。", r2/r1)
plot()
circle(0, 0, r1, :indianred1, fill=true)
circle(x, r2, r2, :steelblue1, fill=true)
circle(-x, r2, r2, :steelblue1, fill=true)
circle(r2, x, r2, :steelblue1, fill=true)
circle(-r2, x, r2, :steelblue1, fill=true)
circle(x, -r2, r2, :steelblue1, fill=true)
circle(-x, -r2, r2, :steelblue1, fill=true)
circle(r2, -x, r2, :steelblue1, fill=true)
circle(-r2, -x, r2, :steelblue1, fill=true)
if more == true
point(x, r2, "(x,r2)", :black)
point(r2, x, "(r2,x)", :black)
point(r1, 0, " r1", :black, :left, :bottom)
point(0, 0, " 0", :black, :left, :bottom)
hline!([0], color=:black, lw=0.5)
vline!([0], color=:black, lw=0.5)
else
plot!(showaxis=false)
end
end;
囲円（赤）の径は等円（青）の径の 3 -2√2 + sqrt(20 - 14√2) = 0.619914 倍である。
