30 岩手県一関市赤荻字宿 正慶山観音寺 天保2年(1831)
安富有恒：和算—岩手の現存算額のすべて，青磁社，東京都，1987.
http://www.wasan.jp/iwatenosangaku_yasutomi.pdf
キーワード：円5個，直線上
#Julia #SymPy #算額 #和算 #数学
直線の上に互いに外接する大円 2 個が載っており，その大円 2 個に下側から外接する大円が 1 個ある。
左右の大円と直線および下の大円が作る隙間に小円が左右に 2 個ある。
小円の直径が 1 寸のとき，大円の直径はいかほどか。
大円の半径と中心座標を \(r_1,\ (r_1,\ 0),\ (0,\ \sqrt{3}r_1)\)
小円の半径と中心座標を \(r_2,\ (x_2,\ r_2 - r_1)\)
とおき，以下の連立方程式を解く。
include("julia-source.txt"); # julia-source.txt ソース
using SymPy
@syms r1::positive, r2::positive, x2::positive
eq1 = x2^2 + ( (r2 - r1) + √Sym(3)*r1)^2 - (r1 + r2)^2
# eq1 = x2^2 + ( (r2 - r1) + √3*r1)^2 - (r1 + r2)^2 # 記号計算を行わない場合
eq2 = (x2 - r1)^2 + (r2 - r1)^2 - (r1 + r2)^2
res = solve([eq1, eq2], (r1, x2))[1];
大円の半径は以下のようになる。桁数の大きい整数係数を持つ，二重根号を含む長い式になる。
res[1] |> println
(-25413135*sqrt(3)*r2 + 44016841*r2 + (-28103451 + 16225535*sqrt(3))*(4*sqrt(3)*r2 + 7*r2 - 4*sqrt(2)*r2*sqrt(1939428460839052729369227845122062114242 - 1119729543939448663603914987253430016841*sqrt(3))/(164355571104683741015 - 94890733220103835857*sqrt(3))))*(4*sqrt(3)*r2 + 7*r2 - 4*sqrt(2)*r2*sqrt(1939428460839052729369227845122062114242 - 1119729543939448663603914987253430016841*sqrt(3))/(164355571104683741015 - 94890733220103835857*sqrt(3)))^2/(2*r2^2*(36060146 - 20819335*sqrt(3)))
SymPy では自動的に簡約化できないので，SymPy の手を借りながら手作業で簡約化する。
まず二重根号を外すと，4 つの項に分解される。
ans_r1 = res[1] |> sympy.sqrtdenest;
それぞれの項を個別に簡約化すると以下のようになる。
t1 = ans_r1.args[1] |> simplify
t1 |> println
3008680048456320*sqrt(3)*r2/121 + 5211186707645184*r2/121
t2 = ans_r1.args[2] |> simplify
t2 |> println
-3008680048441074*sqrt(3)*r2/121 - 5211186707618685*r2/121
t3 = ans_r1.args[3] |> simplify
t3 |> println
52055161703399520*r2/121 + 30054061622167230*sqrt(3)*r2/121
t4 = ans_r1.args[4] |> simplify
t4 |> println
-52055161703423478*r2/121 - 30054061622181024*sqrt(3)*r2/121
4 項を足し合わせ，更に簡約化すると，\(r_1 = 3r_2(4\sqrt{3} + 7) = 3r_2(4\sqrt{3} + 7) = r_2(12\sqrt{3} + 21) = r_2(\sqrt{432} + 21)\) となる。
最後の式は，術に述べられているものである。すなわち，「432 の平方根に 21 を加え，小円の直径を掛けると大円の直径が得られる」。
t1 + t2 + t3 + t4 |> simplify |> println
3*r2*(4*sqrt(3) + 7)
大円の半径 \(r_1\) は，小円の半径 \(r_2\) の \(\sqrt{432} + 21\) 倍である。
小円の直径が 1 寸のとき，大円の直径は \(\sqrt{432} + 21 = 41.78460969082653\) 寸である。
もちろん方程式を解いて得られた解に直接 \(r_2 = 1/2\) を代入して 2 倍すれば，上のように迂遠な手順を踏まなくても大円の直径が得られる。
2res[1](r2 => 1/2).evalf() |> println
41.7846096908265
なお，当初の方程式 eq1 を，eq1 = x2^2 + ( (r2 - r1) + √3*r1)^2 - (r1 + r2)^2 とすれば，記号計算は行われず 41.7846096908267*r2 という解が得られる。
描画関数プログラムのソースを見る
function draw(r2, more=false)
pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
r1 = 41.7846096908267*r2 # = r2*(√432 + 21)
x2 = 28.8564064605511*r2
@printf("小円の直径が %g のとき，大円の直径は %g である。\n", 2r2, 2r1)
plot()
circle2(r1, 0, r1, :blue)
circle(0, -√3r1, r1, :blue)
circle2(x2, r2 - r1, r2, :red)
segment(-2r1, -r1, 2r1, -r1, :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(r1, 0, "大円:r1,(r1,0)", :blue, :center, :bottom, delta=delta)
point(0, -√3r1, "大円:r1,(0,-√3r1)", :blue, :center, delta=-delta)
point(x2, r2 - r1, "小円:r2,(x2,r2-r1)", :red, :left, delta=-delta, deltax=3delta, mark=false)
end
end;
draw(1/2, true)
以下のアイコンをクリックして応援してください
