裏 RjpWiki

Julia ときどき R, Python によるコンピュータプログラム,コンピュータ・サイエンス,統計学

算額(その2053)

2024年08月30日 | Julia

算額(その2053)

七十八 群馬県甘楽郡下仁田町上小坂 中之嶽神社 安政3年(1856)
群馬県和算研究会:群馬の算額,上武印刷株式会社,高崎市,1987年3月31日.
キーワード:円8個

外円 2 個が交わり,区画された領域に大円 2 個,小円 4 個を容れる。外円の直径が 1 寸のとき,小円の直径が最大になるときの大円の直径はいかほどか。

外円の半径と中心座標を R, (x1, 0); x1 = r1
大円の半径と中心座標を r1, (x1 + R - r1, 0); x1 + R - r1 = R
小円の半径と中心座標を r2, (x2, y2)
とおき,以下の連立方程式を解く。

include("julia-source.txt");

using SymPy
@syms R::positive, r1::positive, x1::positive,
     r2::positive, x2::positive, y2::positive;
x1 = r1
eq1 = (x2 - x1)^2 + y2^2 - (R - r2)^2
eq2 = (x1 + R - r1 - x2)^2 + y2^2 - (r1 + r2)^2
eq3 = (x1 + x2)^2 + y2^2 - (R + r2)^2;
res = solve([eq1, eq2, eq3], (r2, x2, y2))[1]

   (-r1*(-R + r1)*(R + r1)/(R^2 + r1^2), -R*(-R + r1)*(R + r1)/(R^2 + r1^2), 2*R*r1*sqrt(R - r1)*sqrt(R + r1)/(R^2 + r1^2))

ans_r2 = res[1]
ans_r2 |> println

   -r1*(-R + r1)*(R + r1)/(R^2 + r1^2)

小円の半径 r2 は,大円の半径 r1 と外円の半径 R の関数である。

たとえば,R = 5/2 のときには,下図のように,r1 が 0.2 〜 0.3 の間で r2 が最大になることがわかる。

pyplot(size=(300, 150), grid=false, aspectratio=:none, label="")
plot(ans_r2(R => 1/2), xlims=(0, 0.5), xlabel="r1", ylabel="r2")

r1 がどのような値をとると r2 が最大になるかを知るためには,ans_r2 の導関数をとり,それが 0 になるときの r1 を求めればよい。

ans_r1 = solve(diff(ans_r2, r1), r1)[1]
ans_r1 |> println

   R*sqrt(-2 + sqrt(5))

ans_r2 = ans_r2(r1 => ans_r1) |> simplify |> factor
ans_r2 |> println

   R*sqrt(-2 + sqrt(5))*(-1 + sqrt(5))/2

r1 が R*sqrt(-2 + sqrt(5)) のとき,r2 が最大値 R*sqrt(-2 + sqrt(5))*(-1 + sqrt(5))/2 になる。

外円の直径が 1 寸のとき,大円の直径が 0.485868271756646 寸のときに,小円の直径は最大値 0.300283106000778 寸になる。

2ans_r1(R => 1/2).evalf() |> println
2ans_r2(R => 1/2).evalf() |> println

   0.485868271756646
   0.300283106000778

function draw(R, more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   r1 = R*sqrt(-2 + sqrt(5))
   (r2, x2, y2) = (-r1*(-R + r1)*(R + r1)/(R^2 + r1^2), -R*(-R + r1)*(R + r1)/(R^2 + r1^2), 2*R*r1*sqrt(R - r1)*sqrt(R + r1)/(R^2 + r1^2))
   @printf("外円の直径が %g のとき,大円の直径が %g のときに小円の直径は最大値 %g になる。\n", 2R, 2r1, 2r2)
   plot()
   circle2(r1, 0, R)
   circle2(R, 0, r1, :blue)
   circle4(x2, y2, 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(r1, 0, "外円:R,(r1,0)", :red, :center, delta=-delta/2)
       point(R, 0, "大円:r1,(R,0)", :blue, :center, :bottom, delta=delta/2)
       point(x2, y2, "小円:r2\n(x2,y2)", :green, :center, delta=-delta/2)
   end
end;

draw(1/2, true)


コメント    この記事についてブログを書く
  • X
  • Facebookでシェアする
  • はてなブックマークに追加する
  • LINEでシェアする
« 算額(その2052) | トップ | 算額(その2054) »
最新の画像もっと見る

コメントを投稿

Julia」カテゴリの最新記事