以下の内容はhttps://sangaku0418.hatenablog.com/entry/2024/01/30/171815より取得しました。


算額(その0664)

北海道函館市谷地頭町 函館八幡宮 文化2年(1805)

中村信弥「幻の算額」
http://www.wasan.jp/maborosi/maborosi.html
キーワード:円1個,半円,正方形,面積
#Julia #SymPy #算額 #和算 #数学


正方形の中に半円 2 個と容円 1 個が入っている容円の直径を「黒積」で表わせ。

黒積とは,名前とは相反して「図形の背景の空白領域の面積」を表すことが常識的であるようだ。この問題の場合は,正方形の内部の半円,容円以外の部分,図で灰色で示した部分の面積を指す。

半円の半径と中心座標を \(r_1,\ (2r_1,\ r_1)\)
容円の半径と中心座標を \(r_2,\ (r_2,\ r_2)\)
黒積を \(A\) とする。
正方形の一辺の長さは \(2r_1\) である。

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

using SymPy
@syms r1::positive, r2::positive, A::positive
s3 = sqrt(Sym(3));

黒積 A は
\( (2r_1)^2 -(\pi r_1^2 - 2(\pi r_1^2/4 - r_1^2/2)) - \pi r_2^2\)
\(\ \ \ \ = -\pi r_1^2/2 + 3r_1^2 - \pi r_2^2\)
である。

(2r1)^2 -(PI*r1^2 - 2(PI*r1^2/4 - r1^2/2)) - PI*r2^2 |> simplify |> println

   -pi*r1^2/2 + 3*r1^2 - pi*r2^2

eq1 = A - (-PI*r1^2/2 + 3*r1^2 - PI*r2^2)
eq1 |> println

   A - 3*r1^2 + pi*r1^2/2 + pi*r2^2

eq1 から solve() により \(r_1\) を求める。

res = solve(eq1, r1)[1]
res |> println  # r1

   sqrt(2*A + 2*pi*r2^2)/sqrt(6 - pi)

\(r_1\) を表すこの式には \(r_2\) が含まれるので,\(r_2\) を消去する。

そこで,半円と容円が外接することに基づき,半円,容円の半径の関係を求める。
最初の解が適解である。すなわち,\(r_2 = 2(2 - \sqrt{3})r_1\) である。

eq2 = (2r1 - r2)^2 + (r1 - r2)^2 - (r1 + r2)^2
res = solve(eq2, r2)[1]
res |> println

   2*r1*(2 - sqrt(3)), 2*r1*(sqrt(3) + 2)

この \(r_1\) に上で求めた式を代入し,方程式 eq2 を解いて \(r_2\) を求めれば,\(A\) を含む式で解が得られる。

eq2 = 2*(sqrt(2*A + 2*PI*r2^2)/sqrt(6 - PI))*(2 - s3) - r2 |> simplify
eq2 |> println

   (-r2*sqrt(6 - pi) + (4 - 2*sqrt(3))*sqrt(2*A + 2*pi*r2^2))/sqrt(6 - pi)

ans_r2 = solve(eq2, r2)[1] |> sympy.sqrtdenest |> simplify;
ans_r2 |> println

   2*sqrt(2)*sqrt(A)*(2 - sqrt(3))/sqrt(-57*pi + 6 + 32*sqrt(3)*pi)

容円の半径 \(r_2\) は
\(\displaystyle  \frac{2\sqrt{2}\sqrt{A}(2 - \sqrt{3})}{\sqrt{-57\pi + 6 + 32\sqrt{3}\pi}}\)
である。

この式が正しいことを逆算して確かめる。
黒積 \(A\) が 2 のとき,容円の半径 \(r_2\) は 1.0440008465400472 である。

r2 = (2*sqrt(2)*sqrt(2)*(2 - sqrt(3))/sqrt(-57*pi + 6 + 32*sqrt(3)*pi))

   1.0440008465400472

\(r_2 = 2r_1(2 - \sqrt{3})\) なので,半円の半径 \(r_1\) は 1.9481321012161867 である。

r1 = r2/(2(2 - sqrt(3)))

   1.9481321012161867

黒積は,以下の通り 1.9999999999999942 ≒ 2 になる。

-pi*r1^2/2 + 3*r1^2 - pi*r2^2

   1.9999999999999942

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

function draw(more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   r1 = 1
   r2 = 2*r1*(2 - sqrt(3))
   @printf("r1 = %g;  r2 = %g\n", r1, r2)
   plot([0, 2r1, 2r1, 0, 0], [0, 0, 2r1, 2r1, 0], color=:black, lw=0.5,seriestype=:shape, fillcolor=:gray)
   circlef(2r1, r1, r1, beginangle=90, endangle=270)
   circlef(r1, 2r1, r1, beginangle=180, endangle=360)
   circle(2r1, r1, r1, :black, beginangle=90, endangle=270)
   circle(r1, 2r1, r1, :black, beginangle=180, endangle=360)
   circlef(r2, r2, r2, :blue)
   plot!([0, 2r1, 2r1, 0, 0], [0, 0, 2r1, 2r1, 0], color=:black)
   if more
       delta = (fontheight = (ylims()[2]- ylims()[1]) / 500 * 10 * 2) /3  # size[2] * fontsize * 2
       vline!([0], color=:black, lw=0.5)
       hline!([0], color=:black, lw=0.5)
       point(2r1, 0, "2r1 ", :black, :right, :bottom, delta=delta/2)
       point(2r1, r1, "(2r1,r1) ", :black, :right, :vcenter)
       point(r1, 2r1, "(r1,2r1) ", :black, :center, delta=-delta/2)
       point(r2, r2, "容円:r2,(r2, r2)", :black, :center, delta=-delta/2)
   end
end;


以下のアイコンをクリックして応援してください




以上の内容はhttps://sangaku0418.hatenablog.com/entry/2024/01/30/171815より取得しました。
このページはhttp://font.textar.tv/のウェブフォントを使用してます

不具合報告/要望等はこちらへお願いします。
モバイルやる夫Viewer Ver0.14