算額（その1269）
百三十三 群馬県高崎市山名町 八幡宮 明治18年(1885)
群馬県和算研究会:群馬の算額，上武印刷株式会社，高崎市，1987年3月31日.
キーワード：円3個，円弧，扇
扇面の，地紙の部分に大円 2 個，骨の部分に小円 1 個を容れる。扇長（扇の要から先端まで）と大円の直径が与えられたとき，小円の直径はいかほどか。
扇長を R，大円の中心と要と扇の端がなす角をθ（扇の中心角が 4θ）
大円の半径と中心座標を r1, (r1, y1)
小円の半径と中心座標を r2, (0, R - 2r1 - r2)
とおき，以下の連立方程式を解く。
include("julia-source.txt");
using SymPy
@syms R::positive, θ::positive,
r1::positive, y1::positive, r2::positive;
eq1 = r1^2 + y1^2 - (R - r1)^2
eq2 = r1^2 + y1^2 - (r1 + (R - 2r1))^2
eq3 = (R - r1)*sind(θ) - r1
eq4 = (R - 2r1 - r2)*sind(2θ) - r2
res = solve([eq1, eq3, eq4], (r2, y1, θ))[2] # 2 of 2
(2*sqrt(R)*r1*sqrt((R - 2*r1)/(R^2 - 2*R*r1 + r1^2))*(R - 2*r1)/(2*sqrt(R)*r1*sqrt((R - 2*r1)/(R^2 - 2*R*r1 + r1^2)) + R - r1), sqrt(R)*sqrt(R - 2*r1), 180*asin(r1/(R - r1))/pi)
小円の半径を表す式は複雑で SymPy では簡約化できない。
res[1]
たとえば，扇長 R = 1，大円の直径 0.6 のとき，小円の直径は 0.349152353114678 である。
2res[1](R => 1, r1 => 0.6/2) |> println
0.349152353114678
すべてのパラメータは以下のとおりである。
R = 1; r1 = 0.3; y1 = 0.632456; r2 = 0.174576; θ = 25.3769
function draw(R, r1, more)
pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
(r2, y1, θ) = (2*sqrt(R)*r1*sqrt((R - 2*r1)/(R^2 - 2*R*r1 + r1^2))*(R - 2*r1)/(2*sqrt(R)*r1*sqrt((R - 2*r1)/(R^2 - 2*R*r1 + r1^2)) + R - r1), sqrt(R)*sqrt(R - 2*r1), 180*asin(r1/(R - r1))/pi)
@printf("扇長が %g，大円の直径が %g のとき，小円の直径は %g である。\n", R, 2r1, 2r2)
@printf("R = %g; r1 = %g; y1 = %g; r2 = %g; θ = %g\n", R, r1, y1, r2, θ)
plot()
circle(0, 0, R, beginangle=90 - 2θ, endangle=90 + 2θ, n=500)
circle(0, 0, R - 2r1, :magenta, beginangle=90 - 2θ, endangle=90 + 2θ, n=500)
x = R*sind(2θ)
y = R*cosd(2θ)
segment(0, 0, x, y)
segment(0, 0, -x, y)
circle2(r1, y1, r1, :blue)
circle(0, R - 2r1 - r2, r2, :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, R, "R", :red, :center, :bottom, delta=delta)
point(0, R - 2r1, "R-2r1", :magenta, :center, :bottom, delta=delta)
point(r1, y1, "大円:r1,(r1,y1)", :blue, :center, :bottom, delta=delta)
point(0, R - 2r1 - r2, "小円:r2\n(0,R-2r1-r2)", :green, :center, :bottom, delta=delta)
end
end;
draw(1, 0.6/2, true)
