41 岩手県一関市滝沢字寺田下 熊野白山滝神社 文久元年(1861)
安富有恒：和算—岩手の現存算額のすべて，青磁社，東京都，1987.
http://www.wasan.jp/iwatenosangaku_yasutomi.pdf
キーワード：円11個，外円
#Julia #SymPy #算額 #和算 #数学
問題文，答，術ともに不明であるが，以下のようなものでもあろう。
山村の図より複雑なものである。改訂版を書くべきである。
算額（その446），算額（その24）と同じ図形である。
外円の中に白円 2 個，赤円 4 個，青1円 2 個，青2円 2 個を容れる。それぞれの円は外円に内接するほか，互いに外接し合っている。各円の直径を求めよ。
上の青と赤の位置関係と下のそれとの位置関係が違うのは，実は上のように白円と 青赤の 4 個の円が互いに外接し合うような図は描けない。無理に描こうとすれば「赤円が青・白円に外接するが外円には内接しない」か「赤円は外円に内接し外側の青円と白円には外接するが内側の青円には外接しない」図にしかならない。
外円の半径と中心座標を \(R,\ (0,\ 0)\)
白円の半径と中心座標を \(r_1,\ (r_1,\ 0)\)
赤円の半径と中心座標を \(r_2,\ (x_2,\ y_2)\)
青1円の半径と中心座標を \(r_3,\ (0,\ R - r_3)\)
青2円の半径と中心座標を \(r_4,\ (0,\ R - 2r_3 - r_4)\)
とおき，以下の連立方程式を解く。
include("julia-source.txt"); # julia-source.txt ソース
using SymPy
@syms R::positive, r1::positive,
r2::positive, x2::positive, y2::positive,
r3::positive, r4::positive
r1 = R/2
eq1 = x2^2 + y2^2 - (R - r2)^2
eq2 = r1^2 + (R - 2r3 - r4)^2 - (r1 + r4)^2
eq3 = (r1 - x2)^2 + y2^2 - (r1 + r2)^2
eq4 = x2^2 + (R - r3 - y2)^2 - (r2 + r3)^2
eq5 = x2^2 + (y2 - R + 2r3 + r4)^2 - (r2 + r4)^2
# res = solve([eq1, eq2, eq3, eq4, eq5], (r1, r2, x2, y2, r3))
function H(u)
(r2, x2, y2, r3, r4) = u
return [
x2^2 + y2^2 - (R - r2)^2, # eq1
R^2/4 - (R/2 + r4)^2 + (R - 2*r3 - r4)^2, # eq2
y2^2 - (R/2 + r2)^2 + (R/2 - x2)^2, # eq3
x2^2 - (r2 + r3)^2 + (R - r3 - y2)^2, # eq4
x2^2 - (r2 + r4)^2 + (-R + 2*r3 + r4 + y2)^2, # eq5
]
end;
R = 10/2
iniv = BigFloat[2, 3.5, 7, 2, 2]
res = nls(H, ini=iniv)
([1.0827118232955024, 1.751864530113493, 3.503729060226986, 0.8017872829546411, 0.9781948282296248], true)
外円の直径が 10 のとき，白円，赤円，青1円，青2円の直径は 5, 2.16542, 1.60357, 1.95639 である。
描画関数プログラムのソースを見る
function draw(R, more=false)
pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
#(r1, r2, x2, y2) = (r3*(sqrt(17) + 7)/4, r3*(-7*sqrt(2*sqrt(17) + 46) - sqrt(34*sqrt(17) + 782) + 22*sqrt(17) + 122)/112, r3*(-10*sqrt(17) + 26 + 3*sqrt(34*sqrt(17) + 782) + 21*sqrt(2*sqrt(17) + 46))/112, r3*(sqrt(2*sqrt(17)/49 + 46/49) + 19/14 + 5*sqrt(17)/14))
#(r1, r2, x2, y2) = (r3*(sqrt(17) + 7)/4, r3*(-42 - sqrt(82*sqrt(17) + 350) - 6*sqrt(17) + sqrt(1394*sqrt(17) + 5950))/16, -sqrt(2)*r3*sqrt(41*sqrt(17) + 175)/8 + r3*(19 + 5*sqrt(17))/8, r3*(3/2 + sqrt(17)/2))
(r2, x2, y2, r3, r4) = res[1]#[6.23606797749979, 2.700745812045397, 4.369898518863388, 8.739797037726776, 2.440035777631468]
r1 = R/2
@printf("外円の直径が %g のとき，白円，赤円，青1円，青2円の直径は %g, %g, %g, %g である。\n", 2R, 2r1, 2r2, 2r3, 2r4)
plot()
circle(0, 0, R, :green)
circle2(r1, 0, r1, :gray50)
circle4(x2, y2, r2, :red)
circle22(0, R - r3, r3, :blue)
circle22(0, R - 2r3 - r4, r4, :blue4)
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)", :gray50, :center, delta=-delta)
point(x2, y2, "赤円:r2\n(x2,y2)", :red, :center, delta=-delta)
point(0, R - r3, "青1円:r3,(0,R-r3)", :blue, :center, delta=-delta)
point(0, R - 2r3 - r4, "青2円:r4,(0,R-2r3-r4)", :blue4, :center, delta=-delta)
end
end;
draw(10/2, true)
以下のアイコンをクリックして応援してください
