算額あれこれ

算額問題をコンピュータで解きます

算額(その414)

愛媛県松山市桜谷町 伊佐爾波神社 嘉永3年(1850)

早坂平四郎: 算額の一考察,苫小牧工業専門学校紀要
https://www.tomakomai-ct.ac.jp/uploads/files/wp/kiyou/old/kiyou5-8.pdf
キーワード:円7個,斜線
#Julia #SymPy #算額 #和算 #数学


直線の上に等円 3 個と大円が互いに接している。さらに,大円の中に4本の弦を引き分割された領域に,等円と同じ半径を持ち弦と接する 2 個の円と,中央の菱形に内接する中円を容れる(菱形の長い方の対角線は,上述の直線と並行である)。

大円の半径と中心座標を \(a,\ (0,\ a)\)
等円の半径と中心座標を \(r_1,\ (x_1,\ r_1),\ ( x_1 - r_1,\ y_1)\) ただし,\(y_1 = r_1(1 + \sqrt{3})\)
菱形の中の円の半径と中心座標を \(r_2,\ (0,\ a)\)
弦の 1 つの端点の座標を \( (b_x,\ b_y)\)
とおき,以下の連立方程式を解く。

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

using SymPy

@syms a::positive, r1::positive, x1::negative, y1::positive,
     r2::positive, bx::positive;
a = 10
eq1 = x1^2 + (a - r1)^2 - (a + r1)^2
eq2 = (x1 - r1)^2 + (a - r1*(1+sqrt(Sym(3))))^2 - (a + r1)^2
eq3 = distance(-bx, a - sqrt(a^2 - bx^2), a, a, 0, a) - r2^2
eq4 = distance(-bx, a - sqrt(a^2 - bx^2), a, a, 0, r1) - r1^2;
res = solve([eq1, eq2, eq3, eq4], (r1, x1, bx, r2))[1]

    (210 - 120*sqrt(3), 60 - 40*sqrt(3), -5891176044722613657600*sqrt(-154227 + 89043*sqrt(3))/52283326179841 - 3401272075258036224000*sqrt(3)*sqrt(-154227 + 89043*sqrt(3))/52283326179841 - 18210720*sqrt(3)/7230721 + 47757270/7230721 + 5891187146667514675200*sqrt(3)*sqrt(-51409 + 29681*sqrt(3))/52283326179841 + 10203835454942977476480*sqrt(-51409 + 29681*sqrt(3))/52283326179841, 10*sqrt(-185511737623 - 132968544*sqrt(-51409 + 29681*sqrt(3)) + 76769280*sqrt(3)*sqrt(-51409 + 29681*sqrt(3)) + 107105251656*sqrt(3))/sqrt(3836935845121 - 2215255943040*sqrt(3)))

   \(r_1 = 2.1539;\  x_1 = -9.28203\)
   \(b_x = 5.48676;\  b_y = 1.63965;\  r_2 = 4.75039\)

菱形の中の円の半径は 4.75038687805787 である。

res[4] |> simplify

res[4].evalf()

   4.75038687805787

術では以下の式が示されているが,計算される値は不適切なもののようだ。

47/(16*sqrt(57√3 - 73) - 141(7 - 4√3)) * 10

   6.616796128691446

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

function draw(more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   a = 10
   (r1, x1, bx, r2) = res[1]
   by = a - sqrt(a^2 - bx^2)
   @printf("r1 = %g;  x1 = %g\n", r1, x1)
   @printf("bx = %g;  by = %g;  r2 = %g\n", bx, by, r2)
   plot()
   circle(0, a, a, :blue)
   circle(x1, r1, r1)
   circle(x1 - 2r1, r1, r1)
   circle(x1 - r1, (1 + √3)r1, r1)
   circle(0, r1, r1)
   circle(0, 2a - r1, r1)
   circle(0, a, r2, :green)
   segment(a, a, -bx, by)
   segment(-a, a, bx, by)
   segment(a, a, -bx, 2a-by)
   segment(-a, a, bx, 2a-by)
   if more
       delta = (fontheight = (ylims()[2]- ylims()[1]) / 500 * 10 * 2) /  3  # size[2] * fontsize * 2
       point(x1, r1, "(x1,r1)", :red, :center, delta=-delta)
       point(x1 - r1, r1*(1 + √3), "(x1-r1,(1+√3)r1)", :red, :center, delta=-delta)
       point(0, a, " a", :green)
       point(float(bx), float(by), " (bx,by)", :black)
       vline!([0], color=:black, lw=0.5)
       hline!([0], color=:black, lw=0.5)
   end
end;


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