百三十 現存するが奉納場所不明 明治14年(1881)
群馬県和算研究会:群馬の算額，上武印刷株式会社，高崎市，1987年3月31日.
キーワード：円3個，正三角形，正方形
#Julia #SymPy #算額 #和算 #数学
正三角形と正方形が交わっており，その区分領域に大円，中円，小円を容れる。正方形の一辺の長さが与えられたとき小円の直径を求めよ。
正三角形の一辺の長さを \(a\)
正方形の一辺の長さを \(b\)
大円の半径と中心座標を \(r_1,\ (r_1,\ r_1)\)
中円の半径と中心座標を \(r_2,\ (x_2,\ r_2)\)
小円の半径と中心座標を \(r_3,\ (x_3,\ r_3)\)
とおき，以下の連立方程式を解く。
1. a, r2, x2 を求める
include("julia-source.txt"); # julia-source.txt ソース
using SymPy
@syms a::positive, b::positive,
r1::positive, r2::positive, x2::positive;
# r1 = b/2
b = 2r1
eq1 = b/(a - b) - √Sym(3)
eq2 = (a - x2)/r2 - √Sym(3)
eq3 = (x2 - r1)^2 + (r1 - r2)^2 - (r1 + r2)^2 |> expand;
ans_a = solve(eq1, a)[1]
ans_a |> println
2*r1*(sqrt(3) + 3)/3
eq12 = eq2(a => ans_a) |> simplify;
eq13 = eq3(a => ans_a);
res1 = solve([eq12, eq13], (r2, x2))[1]
(r1*(-2*sqrt(sqrt(3) + 3) + sqrt(3) + 4)/3, r1*(-2*sqrt(3) + 3 + 2*sqrt(3*sqrt(3) + 9))/3)
2. r3, x3 を求める
前節で求めた \(a,\ r_2,\ x_2\) を既知として \(r_3,\ x_3\) を求める。
@syms a::positive, b::positive,
r1::positive, r2::positive, x2::positive,
r3::positive, x3::positive;
b = 2r1
# r1 = b/2
a = 2r1*(√Sym(3) + 3)/3
r2 = r1*(-2sqrt(√Sym(3) + 3) + √Sym(3) + 4)/3
x2 = r1*(-2√Sym(3) + 3 + 2sqrt(3√Sym(3) + 9))/3
eq4 = (x3 - r1)^2 + (r1 - r3)^2 - (r1 + r3)^2 |> expand
eq5 = (x2 - x3)^2 + (r2 - r3)^2 - (r2 + r3)^2 |> expand;
res2 = solve([eq4, eq5], r3, x3)[2]; # 2 of 2
@syms d
ans_r3 = apart(res2[1]) |> simplify
ans_r3 |> println
r1*(-204*sqrt(sqrt(3) + 3) - 92*sqrt(3*sqrt(3) + 9) + 229*sqrt(3) + 449)/338
\(r_1 = b/2\) として，\(r_3\) は \(r_1\) の \( (-204\sqrt{\sqrt{3} + 3} - 92\sqrt{3\sqrt{3} + 9} + 229\sqrt{3} + 449)/338\) 倍である。
たとえば，\(b = 10\) のとき，\(r_3 = 0.8171180745896389\)，直径は \(1.6342361491792778\) である。
ans_x3 = apart(res2[2]) |> simplify
ans_x3 |> println
r1*(-9*sqrt(sqrt(3) + 3) - sqrt(3*sqrt(3) + 9) + 8*sqrt(3) + 33)/13
\(b = 10\) のとき，各パラメータは以下のとおりである。
\( b = 10;\ a = 15.7735;\ r_1 = 5;\ r_2 = 2.30233;\ x_2 = 11.7858\)
\(r_3 = 0.817118;\ x_3 = 9.04257\)
描画関数プログラムのソースを見る
function draw(b, more)
pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
r1 = b/2
a = 2r1*(√3 + 3)/3
r2 = r1*(√3 + 4 - 2sqrt(√3 + 3))/3
x2 = r1*(2sqrt(3√3 + 9) -2√3 + 3)/3
r3 = r1*(229√3 + 449 - 204sqrt(√3 + 3) - 92sqrt(3√3 + 9))/338
x3 = r1*(8√3 + 33 - 9sqrt(√3 + 3) - sqrt(3√3 + 9))/13
@printf("b = %g; a = %g; r1 = %g; r2 = %g; x2 = %g; r3 = %g; x3 = %g\n", b, a, r1, r2, x2, r3, x3)
plot([0, a, a/2, 0], [0, 0, √3a/2, 0], color=:blue, lw=0.5)
plot!([0, b, b, 0, 0], [0, 0, b, b, 0], color=:green, lw=0.5)
circle(r1, r1, r1)
circle(x2, r2, r2, :orange)
circle(x3, r3, r3, :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(a, 0, "a", :blue, :center, delta=-delta)
point(0, b, "b ", :green, :right, :vcenter)
point(b, 0, "b", :green, :left, delta=-delta)
point(a/2, √3a/2, "(a/2,√3a/2)", :blue, :left, :bottom, delta=delta/2)
point(r1, r1, "大円:r1,(r1,r1)", :red, :center, delta=-delta)
point(x2, r2, "中円:r2\n(x2,r2)", :orange, :center, delta=-delta)
point(x3, r3, "小円:r3,(x3,r3)", :magenta, :right, :bottom, delta=2delta, deltax=-4delta)
plot!(xlims=(-5delta, a + 4delta), ylims=(-5delta, √3a/2 + 5delta))
end
end;
draw(10, true)
以下のアイコンをクリックして応援してください
