算額（その1317）

算額（その1317）

百四十七 群馬県甘楽郡妙義町下高田 高太神社 大正12年(1923)
群馬県和算研究会:群馬の算額，上武印刷株式会社，高崎市，1987年3月31日.
キーワード：円4個，等脚台形
#Julia, #SymPy, #算額, #和算

等脚台形の中に大円 2 個，小円 2 個を容れる。台形の高さが 20 寸，大円の直径が 14.4 寸のとき，小円の直径はいかほどか。

等脚台形の上底，下底の長さをそれぞれ 2b, 2a 
等脚台形の斜辺を延長してできる二等辺三角形の高さを h0
大円の半径と中心座標を r1, (r1, r1)
小円の半径と中心座標を r2, (r2, h - r2)
とおき，以下の連立方程式を解く。

include("julia-source.txt");

using SymPy
@syms a::positive, b::positive, h::positive, h0::positive,
     r1::positive, r2::positive, y2::positive;
eq1 = r1/(h0 - r1) - r2/(h0 - (h - r2))
eq2 = b/(h0 - h) - a/h0
eq3 = dist2(0, h0, a, 0, r1, r1, r1)
eq4 = (r1 - r2)^2 + (h - r2 - r1)^2 - (r1 + r2)^2
res = solve([eq1, eq2, eq3, eq4], (r2, a, b, h0))[2]

   ((-2*h^(5/2)*sqrt(r1) + 12*h^(3/2)*r1^(3/2) - 16*sqrt(h)*r1^(5/2) + h^3 - 5*h^2*r1 + 2*h*r1^2 + 8*r1^3)/(h^2 - 6*h*r1 + 8*r1^2), 4*(h^(17/2)*r1^(3/2) - 26*h^(15/2)*r1^(5/2) + 249*h^(13/2)*r1^(7/2) - 1192*h^(11/2)*r1^(9/2) + 3160*h^(9/2)*r1^(11/2) - 4768*h^(7/2)*r1^(13/2) + 3984*h^(5/2)*r1^(15/2) - 1664*h^(3/2)*r1^(17/2) + 256*sqrt(h)*r1^(19/2) - h^8*r1^2 + 26*h^7*r1^3 - 249*h^6*r1^4 + 1192*h^5*r1^5 - 3160*h^4*r1^6 + 4768*h^3*r1^7 - 3984*h^2*r1^8 + 1664*h*r1^9 - 256*r1^10)/(h^9 - 28*h^8*r1 + 301*h^7*r1^2 - 1690*h^6*r1^3 + 5544*h^5*r1^4 - 11088*h^4*r1^5 + 13520*h^3*r1^6 - 9632*h^2*r1^7 + 3584*h*r1^8 - 512*r1^9), 4*(-h^(9/2)*sqrt(r1) + 13*h^(7/2)*r1^(3/2) - 52*h^(5/2)*r1^(5/2) + 68*h^(3/2)*r1^(7/2) - 16*sqrt(h)*r1^(9/2) - h^4*r1 + h^3*r1^2 + 24*h^2*r1^3 - 52*h*r1^4 + 16*r1^5)/(h^4 - 18*h^3*r1 + 84*h^2*r1^2 - 120*h*r1^3 + 32*r1^4), (sqrt(h)*r1*(h - 6*r1) + 2*r1^(3/2)*(-h + 2*r1))/(sqrt(h)*(h - 4*r1)))

res[1] |> println

   (-2*h^(5/2)*sqrt(r1) + 12*h^(3/2)*r1^(3/2) - 16*sqrt(h)*r1^(5/2) + h^3 - 5*h^2*r1 + 2*h*r1^2 + 8*r1^3)/(h^2 - 6*h*r1 + 8*r1^2)

小円の半径 r2 は，h, r1 を用いてい以下のように計算できる。

f(h, r1) = (-2*h^(5/2)*sqrt(r1) + 12*h^(3/2)*r1^(3/2) - 16*sqrt(h)*r1^(5/2) + h^3 - 5*h^2*r1 + 2*h*r1^2 + 8*r1^3)/(h^2 - 6*h*r1 + 8*r1^2)
f(20, 14.4/2) * 2

   6.399999999999852

台形の高さが 20 寸，大円の直径が 14.4 寸のとき，小円の直径は 6.4 寸である。

すべてのパラメータは以下のとおりである。

   h = 20;  r1 = 7.2;  r2 = 3.2;  a = 24.6857;  b = 4.51765;  h0 = 24.48

function draw(h, r1, more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   (r2, a, b, h0) = ((-2*h^(5/2)*sqrt(r1) + 12*h^(3/2)*r1^(3/2) - 16*sqrt(h)*r1^(5/2) + h^3 - 5*h^2*r1 + 2*h*r1^2 + 8*r1^3)/(h^2 - 6*h*r1 + 8*r1^2), 4*(h^(17/2)*r1^(3/2) - 26*h^(15/2)*r1^(5/2) + 249*h^(13/2)*r1^(7/2) - 1192*h^(11/2)*r1^(9/2) + 3160*h^(9/2)*r1^(11/2) - 4768*h^(7/2)*r1^(13/2) + 3984*h^(5/2)*r1^(15/2) - 1664*h^(3/2)*r1^(17/2) + 256*sqrt(h)*r1^(19/2) - h^8*r1^2 + 26*h^7*r1^3 - 249*h^6*r1^4 + 1192*h^5*r1^5 - 3160*h^4*r1^6 + 4768*h^3*r1^7 - 3984*h^2*r1^8 + 1664*h*r1^9 - 256*r1^10)/(h^9 - 28*h^8*r1 + 301*h^7*r1^2 - 1690*h^6*r1^3 + 5544*h^5*r1^4 - 11088*h^4*r1^5 + 13520*h^3*r1^6 - 9632*h^2*r1^7 + 3584*h*r1^8 - 512*r1^9), 4*(-h^(9/2)*sqrt(r1) + 13*h^(7/2)*r1^(3/2) - 52*h^(5/2)*r1^(5/2) + 68*h^(3/2)*r1^(7/2) - 16*sqrt(h)*r1^(9/2) - h^4*r1 + h^3*r1^2 + 24*h^2*r1^3 - 52*h*r1^4 + 16*r1^5)/(h^4 - 18*h^3*r1 + 84*h^2*r1^2 - 120*h*r1^3 + 32*r1^4), (sqrt(h)*r1*(h - 6*r1) + 2*r1^(3/2)*(-h + 2*r1))/(sqrt(h)*(h - 4*r1)))
   @printf("台形の高さが %g，大円の直径が %g のとき，小円の直径は %g である。\n", h, 2r1, 2r2)
   @printf("h = %g;  r1 = %g;  r2 = %g;  a = %g;  b = %g;  h0 = %g\n", h, r1, r2, a, b, h0)
   plot([a, 0, -a, a], [0, h0, 0, 0], color=:green, lw=0.5)
   segment(-b, h, b, h, :green)
   circle2(r1, r1, r1)
   circle2(r2, h - r2, r2, :blue)
   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, r1, "大円:r1,(r1,r1)", :red, :center, delta=-delta/2)
       point(r2, h - r2, " 小円:r2,(r2,h-r2)", :blue, :left, :vcenter)
       point(a, 0, " a", :green, :left, :bottom, delta=delta/2)
       point(0, h, " h", :green, :left, :bottom, delta=delta/2)
       point(0, h0, " h0", :green, :left, :bottom, delta=delta/2)
       point(b, h, " (b,h)", :green, :left, :vcenter)
   end
end;

draw(20, 14.4/2, true)

 

算額（その1316）

算額（その1316）

百四十七 群馬県甘楽郡妙義町下高田 高太神社 大正12年(1923)
群馬県和算研究会:群馬の算額，上武印刷株式会社，高崎市，1987年3月31日.
キーワード：円4個，外円，正三角形
#Julia, #SymPy, #算額, #和算

外円の中に大円 2 個，小円 1 個，直角三角形 1 個を容れる。正三角形の一辺の長さが 7.5 寸，小円の直径が 4.755 寸のとき，大円の直径はいかほどか。

正三角形の一辺の長さを 2a，頂点と円の接点座標を (a, x0, y0);x0 = a, y0 = -sqrt(R^2 - a^2)
外円の半径と中心座標を R, (0, 0)
大円の半径と中心座標を r1, (x1, y1); x1 = r1, y1 = √Sym(3)a + y0
小円の半径と中心座標を r2, (0, R - r2)
とおき，以下の連立方程式を解く。

include("julia-source.txt");

using SymPy
@syms a::positive, R::positive, r1::positive,
     x1::positive, y1::positive, r2::positive;
x0 = a
y0 = -sqrt(R^2 - a^2)
x1 = r1
y1 = √Sym(3)a + y0
eq1 = x1^2 + y1^2 - (R - r1)^2
eq2 = R - 2r2 - √Sym(3)a - y0
res = solve([eq1, eq2], (r1, R))[1]

   (a*(-3*a^3 - 4*sqrt(3)*a^2*r2 - 4*a*r2^2 + 3*a*sqrt(a^4 + 4*sqrt(3)*a^3*r2 + 16*a^2*r2^2 + 8*sqrt(3)*a*r2^3 + 4*r2^4) + 2*r2*sqrt(3*a^4 + 12*sqrt(3)*a^3*r2 + 48*a^2*r2^2 + 24*sqrt(3)*a*r2^3 + 12*r2^4))/(2*(sqrt(3)*a^3 + 5*a^2*r2 + 3*sqrt(3)*a*r2^2 + 2*r2^3)), 2*(a^2 + sqrt(3)*a*r2 + r2^2)/(sqrt(3)*a + 2*r2))

正三角形の一辺の長さが 7.5 寸，小円の直径が 4.755 寸のとき，大円の直径は 5.89244501321687 寸である。

「答」には「大円径六寸二分五厘」と書いているが，それは多分，外円の半径である。

res[1](a => 7.5/2, r2 => 4.755/2).evalf() * 2 |> println

   5.89244501321687

f(a, r2) = a*(-3*a^3 - 4*sqrt(3)*a^2*r2 - 4*a*r2^2 + 3*a*sqrt(a^4 + 4*sqrt(3)*a^3*r2 + 16*a^2*r2^2 + 8*sqrt(3)*a*r2^3 + 4*r2^4) + 2*r2*sqrt(3*a^4 + 12*sqrt(3)*a^3*r2 + 48*a^2*r2^2 + 24*sqrt(3)*a*r2^3 + 12*r2^4))/(2*(sqrt(3)*a^3 + 5*a^2*r2 + 3*sqrt(3)*a*r2^2 + 2*r2^3))
f(7.5/2, 4.755/2) * 2

   5.892445013216866

function draw(a, r2, more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   (r1, R) = (a*(-3*a^3 - 4*sqrt(3)*a^2*r2 - 4*a*r2^2 + 3*a*sqrt(a^4 + 4*sqrt(3)*a^3*r2 + 16*a^2*r2^2 + 8*sqrt(3)*a*r2^3 + 4*r2^4) + 2*r2*sqrt(3*a^4 + 12*sqrt(3)*a^3*r2 + 48*a^2*r2^2 + 24*sqrt(3)*a*r2^3 + 12*r2^4))/(2*(sqrt(3)*a^3 + 5*a^2*r2 + 3*sqrt(3)*a*r2^2 + 2*r2^3)), 2*(a^2 + sqrt(3)*a*r2 + r2^2)/(sqrt(3)*a + 2*r2))
   x1 = r1
   y1 = √3a - sqrt(R^2 - a^2)
   x0 = a
   y0 = -sqrt(R^2 - x0^2)
   @printf("正三角形の一辺の長さが %g，小円の直径が %g のとき，大円の直径は %g である。\n", 2a, 2r2, 2r1)
   @printf("a = %g;  r2 = %g;  R = %g;  r1 = %g;  x1 = %g;  y1 = %g;  x0 = %g;  y0 = %g\n", a, r2, R, r1, x1, y1, x0, y0)
   plot([a, 0, -a, a], [y0, √3a + y0, y0, y0], color=:green, lw=0.5)
   circle(0, 0, R, :magenta)
   circle2(x1, y1, r1)
   circle(0, R - r2, r2, :blue)
   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(0, R, "R", :red, :center, :bottom, delta=delta/2)
       point(0, √3a + y0, "√3a+y0  ", :green, :right, delta=-delta/2)
       point(x1, y1, "大円:r1,(x1,y1)", :red, :center, delta=-delta)
       point(0, R - r2, "小円:r2,(0,R-r2)", :blue, :center, :bottom, delta=2delta)
       point(x0, y0, " (x0,y0)", :green, :left, :vcenter)
   end
end;

draw(7.5/2, 4.755/2, true)