続神壁算法 付録〇一 東都麹町平川天満宮 寛政八年(1796)
丸山因平良玄門人 参州苅屋 林政右衛門盛保
藤田嘉言(1807)：続神壁算法
http://www.wasan.jp/jinpeki/jinpeki.html
キーワード：円5個
#Julia #SymPy #算額 #和算 #数学
十字線を隔てて乾円，坤円，巽円，艮円の4円と，中央に容円を置く。艮円，坤円，巽円の直径がそれぞれ 15 寸，10 寸，6 寸のとき，乾円の直径を求めよ。
容円の半径と中心座標を \(r_0,\ (x_0,\ y_0)\)
坤円の半径と中心座標を \(r_1,\ (r_1,\ r_1)\)
乾円の半径と中心座標を \(r_2,\ (r_2,\ -r_2)\)
巽円の半径と中心座標を \(r_3,\ (-r_3,\ r_3)\)
艮円の半径と中心座標を \(r_4,\ (-r_4,\ -r_4)\)
とおき，以下の連立方程式を解く。
include("julia-source.txt"); # julia-source.txt ソース
using SymPy
@syms r0::positive, x0::negative, y0::positive,
r1::positive, r2::positive, r3::positive, r4::positive
eq1 = (r1 - x0)^2 + (r1 - y0)^2 - (r0 + r1)^2
eq2 = (r2 - x0)^2 + (-r2 - y0)^2 - (r0 + r2)^2
eq3 = (-r3 - x0)^2 + (r3 - y0)^2 - (r0 + r3)^2
eq4 = (-r4 - x0)^2 + (-r4 - y0)^2 - (r0 + r4)^2
res = solve([eq1, eq2, eq3, eq4], (r2, r0, x0, y0))[4]; # 4 of 4
4 番目の組が適解であるが，SymPy では簡約化できない長い式になる。
res[1](r1 => 5, r3 => 3, r4 => 7.5).evalf() |> println
res[2](r1 => 5, r3 => 3, r4 => 7.5).evalf() |> println
res[3](r1 => 5, r3 => 3, r4 => 7.5).evalf() |> println
res[4](r1 => 5, r3 => 3, r4 => 7.5).evalf() |> println
10.0000000000000
2.70833333333333
0.666666666666667
-1.37500000000000
まえもって \(r_1,\ r_3,\ r_4\) に定数を代入しておいてから解くと数値解が求まる。
@syms r0::positive, x0::negative, y0::positive,
r1::positive, r2::positive, r3::positive, r4::positive
@syms r0, x0, y0, r1, r2, r3, r4
(r1, r3, r4) = (10, 6, 15) .// 2
eq1 = (r1 - x0)^2 + (r1 - y0)^2 - (r0 + r1)^2
eq2 = (r2 - x0)^2 + (-r2 - y0)^2 - (r0 + r2)^2
eq3 = (-r3 - x0)^2 + (r3 - y0)^2 - (r0 + r3)^2
eq4 = (-r4 - x0)^2 + (-r4 - y0)^2 - (r0 + r4)^2
res = solve([eq1, eq2, eq3, eq4], (r2, r0, x0, y0))[4] # 4 of 4
(10, 65/24, 2/3, -11/8)
艮円，坤円，巽円の直径がそれぞれ 15 寸，10 寸，6 寸のとき，乾円の直径は 20 寸である。
描画関数プログラムのソースを見る
function draw(more)
pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
(r1, r3, r4) = (10, 6, 15) .// 2
(r2, r0, x0, y0) = [10, 65/24, 2/3, -11/8]
plot()
circle(r1, r1, r1)
circle(r2, -r2, r2, :blue)
circle(-r3, r3, r3, :green)
circle(-r4, -r4, r4, :magenta)
circle(x0, y0, r0, :orange)
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=:black, lw=0.5)
point(x0, y0, "容円:r0,(x0,y0)", :black, :left, delta=-delta/2)
point(r1, r1, "坤円:r1,(r1,r1)", :red, :center, delta=-delta/2)
point(r2, -r2, "乾円:r2,(r2,-r2)", :blue, :center, delta=-delta/2)
point(-r3, r3, "巽円:r3\n(-r3,r3)", :green, :center, :bottom, delta=delta/2)
point(-r4, -r4, "艮円:r4\n(-r4,-r4)", :magenta, :center, :bottom, delta=delta/2)
end
end;
以下のアイコンをクリックして応援してください
