算額（その1546）
神壁算法 越後州長岡 蒼柴大明神 寛政8年(1796)
關流藤田貞資門人 長岡 石垣作右衛門光隆
藤田貞資(1789)：神壁算法巻上
http://www.wasan.jp/jinpeki/jinpekisanpo1.pdf
キーワード：円2個，等脚台形，対角線
#Julia, #SymPy, #算額, #和算
外円の中に等脚台形と上円，下円の 2 個の円を容れる。内斜（対角線）は旁斜（斜辺）より 4 寸長い。上円の直径が 1 寸のとき，下円の直径はいかほどか。
外円の半径と中心座標を R, (0, 0)
上円の半径と中心座標を r1, (0, R - r1)
下円の半径と中心座標を r2,（0, r2 - R)
上頭（上底），下頭（下底） を 2b，2a
内斜（対角線）と 旁斜（斜辺）の差を K; 内斜 > 傍斜
とおき，以下の方程式を解く。
include("julia-source.txt");
# # julia-source.txt ソース https://blog.goo.ne.jp/r-de-r/e/ad3a427b84bb416c4f5b73089ae813cf
@syms R, a, b, r1, r2, K
a = sqrt(R^2 - (2r2 - R)^2)
b = sqrt(R^2 - (R - 2r1)^2)
eq1 = sqrt((a + b)^2 + (2R - 2r1 - 2r2)^2) - sqrt((a - b)^2 + (2R - 2r1 - 2r2)^2) - K
res = solve(eq1, r2)[1] # 1 of 2
res |> println
K^2/(16*r1)
不思議なことに（？），解の中に外円の半径 R が出てこない。つまり，内斜（対角線）と 旁斜（斜辺）の差 K と，上円の半径 r1 が決まれば r2 は自動的に決まり，R はなんでもよい（R により a, b が決まるから）ということだ。
それはさておき，下円の半径は『「内斜 - 旁斜」の二乗を上円の半径の16倍で割る』。
上円の半径が 1/2 のとき，下円の半径は 4^2/(16*(1/2)) = 2 である。術と一致する。
function draw(K, r1, R, more=false)
pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
r1 = 1/2
r2 = K^2/(16r1)
@printf("「旁斜-内斜」が %g，上円の直径が %g のとき，下円の直径は %g である。外円の直径は %g である。\n", K, 2r1, 2r2, 2R)
b = sqrt(R^2 - (R - 2r1)^2)
a = sqrt(R^2 - (2r2 - R)^2)
println("a = $a, b = $b")
l1 = sqrt((a + b)^2 + (2R - 2r1 - 2r2)^2)
l2 = sqrt((a - b)^2 + (2R - 2r1 - 2r2)^2)
println((l1, l2, l1 - l2))
plot([a , b, -b, -a, a], [2r2 - R, R - 2r1, R - 2r1, 2r2 - R, 2r2 - R], color=:magenta, lw=0.5)
circle(0, 0, R, :blue)
circle(0, R - r1, r1)
circle(0, r2 - R, 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, r2 - R, "下円:r1,(a,r2-R)", :green, :center, delta=-delta/2)
point(0, R - r1, "上円:r2,(b,R-r1)", :red, :center, delta=-delta/2)
point(a, 2r2 - R, "(a,2r2-R) ", :magenta, :right, delta=-delta/2)
point(b, R - 2r1, "(b,R-2r1) ", :magenta, :right, delta=-delta/2)
dimension_line(-b, R - 2r1, a, 2r2 - R, "内斜", :orange, :left, :vcenter, deltax=2delta)
dimension_line(b, R - 2r1, a, 2r2 - R, "旁斜", :tomato, :left, :vcenter, deltax=delta)
end
end;
draw(4, 1/2, 6, true)
