一〇〇 桶川市小針領家 氷川諏訪神社 明治30年(1897)
埼玉県立図書館：埼玉県史料集 第二集『埼玉の算額』，昭和44年，誠美堂印刷所，埼玉県与野市.
キーワード：円3個，正方形
#Julia #SymPy #算額 #和算 #数学
正方形の中に甲円 1 個，乙円 2 個を容れる。正方形の一辺の長さが 2 寸 5 分，乙円の直径が 7 分のとき，甲円の直径はいかほどか。
正方形の一辺の長さと頂点座標を \(a，(a/\sqrt{2}，0),\ (0,\ a/\sqrt{2}),\ (-a/\sqrt{2},\ 0),\ (0,\ -a/\sqrt{2})\)
甲円の半径と中心座標を \(r_1,\ (0,\ y_1)\)
乙円の半径と中心座標を \(r_2,\ (r_2,\ y_2);\ y_2 < 0\)
とおき，以下の連立方程式を立てる。
include("julia-source.txt"); # julia-source.txt ソース
using SymPy
@syms a::positive, a2::positive, r1::positive,
y1::positive, r2::positive, y2::negative
a2 = a/√Sym(2)
eq1 = r2^2 + (y1 - y2)^2 - (r1 + r2)^2
eq2 = dist2(a2, 0, 0, a2, 0, y1, r1)
eq3 = dist2(a2, 0, 0, -a2, r2, y2, r2);
eq1 |> println
eq2 |> println
eq3 |> println
r2^2 - (r1 + r2)^2 + (y1 - y2)^2
a^2/4 - sqrt(2)*a*y1/2 - r1^2 + y1^2/2
a^2/4 - sqrt(2)*a*r2/2 + sqrt(2)*a*y2/2 - r2^2/2 - r2*y2 + y2^2/2
SymPy の性能上，一度に解くことができないので，逐次解いてゆく。
まず eq2 を解き，\(y_1\) を求める。
ans_y1 = solve(eq2, y1)[1]
ans_y1 |> println
sqrt(2)*(a/2 - r1)
eq1, eq3 に代入する。
eq11 = eq1(y1 => ans_y1) |> expand
eq11 |> println
eq13 = eq3(y1 => ans_y1) |> expand
eq13 |> println
a^2/2 - 2*a*r1 - sqrt(2)*a*y2 + r1^2 - 2*r1*r2 + 2*sqrt(2)*r1*y2 + y2^2
a^2/4 - sqrt(2)*a*r2/2 + sqrt(2)*a*y2/2 - r2^2/2 - r2*y2 + y2^2/2
eq13 を解き，\(y_2\) を求める。この式は \(a,\ r_2\) だけを含むので \(y_2\) は決定される。
ans_y2 = solve(eq13, y2)[1]
ans_y2 |> println
-sqrt(2)*a/2 + r2 + sqrt(2)*r2
eq11 に \(y_1,\ y_2\) を代入する。
eq21 = eq11(y1 => ans_y1, y2 => ans_y2) |> simplify
eq21 |> println
2*a^2 - 4*a*r1 - 4*a*r2 - 2*sqrt(2)*a*r2 + r1^2 + 2*r1*r2 + 2*sqrt(2)*r1*r2 + 2*sqrt(2)*r2^2 + 3*r2^2
eq21 は \(a,\ r_2\) だけを含むので \(r_1\) は決定される。
ans_r1 = solve(eq21, r1)[1]
ans_r1 |> println
-sqrt(2)*sqrt(a)*sqrt(a - sqrt(2)*r2) + 2*a - sqrt(2)*r2 - r2
\(r_1\) が決まれば，逆順にたどり，\(y_1\) が決定される。
r2 = 0.7/2 # 半径を与える
a = 2.5
r1 = -sqrt(2)*sqrt(a)*sqrt(a - sqrt(2)*r2) + 2*a - sqrt(2)*r2 - r2
y2 = -sqrt(2)*a/2 + r2 + sqrt(2)*r2
y1 = sqrt(2)*(a/2 - r1)
(r1, y1, y2)
# (0.9887772739600998, 0.3694247219656983, -0.9227922061357855)
(0.9887772739600998, 0.3694247219656983, -0.9227922061357855)
甲円の直径は \(2\left(2a - (1 + \sqrt{2})r_2 - \sqrt{2}\sqrt{a^2 - \sqrt{2}r_2 a}\right)\)
正方形の一辺の長さ \(a\) が 2 寸 5 分，乙円の直径 \(r_2\) が 7 分のとき，甲円の直径は 1.9775545479202004 寸である。
r2 = 0.7/2 # 半径を与える
a = 2.5
2(2a - (1 + √2)r2 -√2sqrt(a^2 - √2r2*a))
1.9775545479202004
その他のパラメータは以下のとおりである。
\(a = 2.5;\ r_2 = 0.35;\ r_1 = 0.988777;\ y_1 = 0.369425;\ y_2 = -0.922792\)
描画関数プログラムのソースを見る
function draw(more=false)
pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
r2 = 0.7/2
a = 2.5
r1 = -sqrt(2)*sqrt(a)*sqrt(a - sqrt(2)*r2) + 2*a - sqrt(2)*r2 - r2
y2 = -sqrt(2)*a/2 + r2 + sqrt(2)*r2
y1 = sqrt(2)*(a/2 - r1)
a2 = a/√2
@printf("正方形の一辺の長さ = %g，乙円の直径 = %g のとき，甲円の直径は %g である。\n", a, 2r2, 2r1)
@printf("a = %g; r2 = %g; r1 = %g; y1 = %g; y2 = %g\n", a, r2, r1, y1, y2)
plot([a2, 0, -a2, 0, a2], [0, a2, 0, -a2, 0], color=:blue, lw=0.5, xlims=(-2, 2.2))
circle(0, y1, r1, :green)
circle2(r2, y2, r2)
if more
delta = (fontheight = (ylims()[2]- ylims()[1]) / 500 * 10 * 2) /3 # size[2] * fontsize * 2
hline!([0], color=:black, lw=0.5)
vline!([0], color=:gray80, lw=0.5)
point(0, a/√2, " a/√2", :blue, :left, :bottom, delta=delta/2)
point(a/√2, 0, " a/√2", :blue, :left, :bottom, delta=delta/2)
point(0, y1, "甲円:r1,(0,y1)", :green, :center, delta=-delta/2)
point(r2, y2, "乙円:r2\n(r2,y2)", :red, :center, delta=-delta/2)
end
end;
以下のアイコンをクリックして応援してください
