算額（その258） - 算額あれこれ

長野県長野市久保寺観音堂享和3年(1803)

中村信弥「改訂増補長野県の算額」県内の算額(P.70)
http://www.wasan.jp/zoho/zoho.html
キーワード：円3個，菱形
#Julia #SymPy #算額 #和算 #数学

菱形の中に大円 1 個と小円 2 個が入っている。菱形の辺の長さの和と菱形の横の対角線の長さの積が 160 平方寸，菱形の辺の長さの和と大円の直径の積が 96 平方寸のとき，小円の直径を求めよ。

菱形の中心を原点とし，長い方と短い方の対角線の長さを \(2x,\ 2y\) とおく。
大円の半径，中心座標を \(r_1,\ (0,\ 0)\)
小円の半径，中心座標を \(r_2,\ (r_2,\ 0)\)
とおき，以下の連立方程式を解く。

include("julia-source.txt"); # julia-source.txt ソース

using SymPy

@syms r1::positive, r2::positive, x::positive, y::positive, sinθ::positive;

xy = sqrt(x^2 + y^2)
sinθ = y / xy
eq1 = 2x * 4xy - 160
eq2 = 2r1 * 4xy - 96
eq3 = x*sinθ - r1
eq4 = (x - r2)sinθ - r2;

# res = solve([eq1, eq2, eq3, eq4], (r1, r2, x, y))

なぜか solve() では解けないので nlsolve() を用いる。

function H(u)
(r1, r2, x, y) = u
return [
8*x*sqrt(x^2 + y^2) - 160, # eq1
8*r1*sqrt(x^2 + y^2) - 96, # eq2
-r1 + x*y/sqrt(x^2 + y^2), # eq3
-r2 + y*(-r2 + x)/sqrt(x^2 + y^2), # eq4
]
end;

iniv = [1.0, 1, 1, 1]
res = nls(H, ini=iniv);
println(res);

([2.4, 1.5, 4.0, 2.9999999999999996], true)

小円の直径は 3 寸である。

描画関数プログラムのソースを見る

function draw(more)
pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
(r1, r2, x, y) = res[1]
@printf("r1 = %.5f; r2 = %.5f; x = %.5f; y = %.5f\n", r1, r2, x, y)
println("小円の直径 = ", 2r2)
plot([x, 0, -x, 0, x], [0, y, 0, -y, 0], linecolor=:black, linewidth=0.5)
circlef(0, 0, r1, :red)
circlef(r2, 0, r2, :yellow)
circlef(-r2, 0, r2, :yellow)
circle(0, 0, r1, :black)
circle(r2, 0, r2, :black)
circle(-r2, 0, r2, :black)
if more
point(0, y, " y", :black, :left, :bottom)
point(x, 0, " y", :black, :left, :bottom)
point(0.2, 1.8, "大円:r1", :black, mark=false)
point(r2, 0, "小円:r2\n", :black, :center, :bottom)
hline!([0], color=:black, lw=0.5)
vline!([0], color=:black, lw=0.5)
end
end;