算額（その1148）

算額（その1148）

番外三 川口市峯後 峯が岡八幡宮 （昭和10年代よりは前）
埼玉県立図書館：埼玉県史料集 第二集『埼玉の算額』，昭和44年，誠美堂印刷所，埼玉県与野市.
キーワード：円2個，直角三角形，正三角形

直角三角形内に大円，小円，正三角形を容れる。鈎（直角を挟む二辺の短い方）が 5 寸のとき，図形の各寸法を得よ。

大円の半径と中心座標を r1, (r1, y1)
小円の半径と中心座標を r2, (r2, r2)
正三角形の頂点の座標を (a, 0), (b, 0), ((a + b)/2, √3(b - a)/2)
とおき，以下の連立方程式を解く。

include("julia-source.txt")

using SymPy

@syms 鈎::positive, a::positive, b::positive,
     r1::positive, y1::positive, r2::positive
s3 = √Sym(3)
eq1 = dist2(a, 0, (a + b)/2, s3*(b - a)/2, r1, y1, r1)
eq2 = dist2(a, 0, (a + b)/2, s3*(b - a)/2, r2, r2, r2)
eq3 = dist2(0, 鈎, 2b - a, 0, r1, y1, r1)
eq4 = (r1 - r2)^2 + (y1 - r2)^2 - (r1 + r2)^2
eq5 =  s3*鈎 - (2b - a);

まず，eq1, eq2, eq4, eq5 を解いて a, b, y1, y2 を求める。

(ans_a, ans_b, ans_y1, ans_r2) = solve([eq1, eq2, eq4, eq5], (a, b, y1, r2))[1]

   (r1*(-9*sqrt(3) - 12*sqrt(26 - 15*sqrt(3)) - 4*sqrt(78 - 45*sqrt(3)) + 21)/3, -3*sqrt(3)*r1/2 - 2*r1*sqrt(26 - 15*sqrt(3)) - 2*r1*sqrt(78 - 45*sqrt(3))/3 + 7*r1/2 + sqrt(3)*鈎/2, r1*(-6*sqrt(3) + 4*sqrt(26 - 15*sqrt(3)) + 4*sqrt(78 - 45*sqrt(3)) + 11), r1*(-8*sqrt(3) - 4*sqrt(26 - 15*sqrt(3)) + 15))

ans_a = ans_a |> sympy.sqrtdenest |> factor
ans_a |> println

   r1*(-8*sqrt(6) - 9*sqrt(3) + 12*sqrt(2) + 21)/3

ans_b = ans_b |> sympy.sqrtdenest |> factor
ans_b |> println

   (-8*sqrt(6)*r1 - 9*sqrt(3)*r1 + 12*sqrt(2)*r1 + 21*r1 + 3*sqrt(3)*鈎)/6

ans_y1 = ans_y1 |> sympy.sqrtdenest |> factor
ans_y1 |> println

   r1*(-6*sqrt(3) - 4*sqrt(6) + 11 + 8*sqrt(2))

ans_r2 = ans_r2 |> sympy.sqrtdenest |> simplify
ans_r2 |> println

   r1*(-6*sqrt(6) - 8*sqrt(3) + 10*sqrt(2) + 15)

これらを eq13 に代入する。

eq13 = eq3(a => ans_a, b => ans_b, y1 => ans_y1, r2 => ans_r2) |> simplify
eq13 |> println

   鈎^2*(-536*sqrt(6)*r1^2 - 758*sqrt(3)*r1^2 + 1320*r1^2 + 936*sqrt(2)*r1^2 - 48*sqrt(2)*r1*鈎 - 66*r1*鈎 + 24*sqrt(6)*r1*鈎 + 34*sqrt(3)*r1*鈎 + 3*鈎^2)

eq13 を解いて r1 を求める。

ans_r1 = solve(eq13, r1)[1];
ans_r1 |> println

   鈎*(-19*sqrt(3) - 12*sqrt(6) + 33 + 24*sqrt(2))/(2*(-268*sqrt(6) - 379*sqrt(3) + 660 + 468*sqrt(2)))

den = denominator(ans_r1) |> sympy.sqrtdenest |> simplify
ans_r1 = numerator(ans_r1)/den
ans_r1 = apart(ans_r1) |> simplify
ans_r1 |> println

   鈎*(-16*sqrt(6) - 4*sqrt(2) + 19 + 29*sqrt(3))/94

鈎 = 5 を代入すると数値解が得られる。

ans_r1(鈎 => 5).evalf() |> println

   1.29685017475927

全てのパラメータは以下のように計算される。

鈎 = 5
r1 = 鈎*(-16*sqrt(6) - 4*sqrt(2) + 19 + 29*sqrt(3))/94
r2 = r1*(-6*sqrt(6) - 8*sqrt(3) + 10*sqrt(2) + 15)
y1 = r1*(-6*sqrt(3) - 4*sqrt(6) + 11 + 8*sqrt(2))
b = (-8*sqrt(6)*r1 - 9*sqrt(3)*r1 + 12*sqrt(2)*r1 + 21*r1 + 3*sqrt(3)*鈎)/6
a = r1*(-8*sqrt(6) - 9*sqrt(3) + 12*sqrt(2) + 21)/3
println("正三角形の一辺の長さ = $(b - a)")
println("大円の直径 = $(2r1)")
println("小円の直径 = $(2r2)")

   正三角形の一辺の長さ = 3.727915719641115
   大円の直径 = 2.593700349518533
   小円の直径 = 1.5271466611926818

「答」では以下のようになっている。
正三角形の一辺の長さ = 3.66
大円の直径 = 2.68
小円の直径 = 1.55
小数点以下2位までしか書かれていないので比較しても意味がないが，少なくとも，このパラメータで図を描くと赤のようになり，上で求めた正確な図（黒）とは違う。大円と小円の共通接線が水辺線となす角は 60° にならない。



2位までしか書かれていないから不正確なのではなく，「術」に疑問がある。たとえば正三角形の一辺の長さは

鈎 = 5
甲 = sqrt(3) + 1
三角面 = (鈎 / 甲)*2

で求めるとしているが，これで計算すると三角面は 3.6602540378443864 となり，そもそも正しい値 3.727915719641115 とは異なった術式である。

function draw(鈎, more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   r1 = 鈎*(-16*sqrt(6) - 4*sqrt(2) + 19 + 29*sqrt(3))/94
   r2 = r1*(-6*sqrt(6) - 8*sqrt(3) + 10*sqrt(2) + 15)
   y1 = r1*(-6*sqrt(3) - 4*sqrt(6) + 11 + 8*sqrt(2))
   b = (-8*sqrt(6)*r1 - 9*sqrt(3)*r1 + 12*sqrt(2)*r1 + 21*r1 + 3*sqrt(3)*鈎)/6
   a = r1*(-8*sqrt(6) - 9*sqrt(3) + 12*sqrt(2) + 21)/3
   @printf("大正三角形の一辺の長さが %g のとき，正方形の一辺の長さは %g である。\n", 2a, 2b)
   plot([0, 2b - a, 0, 0], [0, 0, 鈎, 0], color=:blue, lw=0.5)
   plot!([a, b, (a + b)/2, a], [0, 0, √3(b - a)/2, 0], color=:magenta, lw=0.5)
   circle(r1, y1, r1)
   circle(r2, r2, 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(a, 0, " a", :magenta, :left, :bottom, delta=delta/2)
       point(b, 0, " b", :magenta, :left, :bottom, delta=delta/2)
       point(2b - a, 0, " 2b-a", :blue, :left, :bottom, delta=delta/2)
       point(0, 鈎, " 鈎", :blue, :left, :bottom, delta=delta/2)
       point((a + b)/2, √3*(b - a)/2, "((a+b)/2,√3(b-a)/2)", :magenta, :left, :bottom, delta=delta/2)
       point(r1, y1, "大円:r1,(r1,y1)", :red, :center, delta=-delta/2)
       point(r2, r2, "小円:r2\n(r2,r2)", :green, :center, :bottom, delta=delta)
   end
end;