算額（その2119）

算額（その2119）

八 武周足立郡桶川宿 不動堂 文化3年(1806)
埼玉県立図書館：埼玉県史料集 第二集『埼玉の算額』，昭和44年，誠美堂印刷所，埼玉県与野市.
キーワード：円8個
#Julia, #SymPy, #算額, #和算

大円 4 個が交わっている隙間に天円 2 個，等円 2 個を容れる。天円の直径が 15 寸のとき，等円の直径はいかほどか。

大円の半径と中心座標を r1, (0, r1/2), (√3r1/2, 0)
天円の半径と中心座標を r2, (0, 3r1/2 - r2)
等円の半径と中心座標を r3, (√3r1/2 + r1 - r3, 0)
とおき，以下の方程式を解く。
まず eq1 を解いて r1 を求め，解を eq2 に代入して解いて r3 を求める。
（連立して解くと SymPy では解を求めることができない）

include("julia-source.txt");

using SymPy
@syms r1::positive, r2::positive, r3::positive;
eq1 = (√Sym(3)*r1/2)^2 + (3r1/2 - r2)^2 - (r1 + r2)^2
eq2 = (√Sym(3)*r1/2 + r1 - r3)^2 + (r1/2)^2 - (r1 + r3)^2;

大円の半径 r1 は，天円の半径 r2 の 5/2 倍である。
天円の直径が 15 寸のとき，大円の直径は 15*5/2 = 57.5 である。

ans_r1 = solve(eq1, r1)[1]
ans_r1 |> println
ans_r1(r2 => 15/2) * 2 |> println

   5*r2/2
   37.5000000000000

等円の半径 r3 は，天円の半径 r2 の 5(1 + 3√3)/26 倍である。
天円の直径が 15 寸のとき，等円の直径は 15*5(1 + 3√3)/26 = 17.873516603961438 である。

eq12 = eq2(r1 => ans_r1)
ans_r3 = solve(eq12, r3)[1] |> factor
ans_r3 |> println
ans_r3(r2 => 15/2).evalf() * 2 |> println

   5*r2*(1 + 3*sqrt(3))/26
   17.8735166039614

「術」は以下のようになっており，簡約化すると同じ関係式になる。

@syms 天径, d
子 = √Sym(3) + 1
等円径 = 子/2(子 + 3)*5天径
等円径 |> println
等円径(天径 => 15).evalf() |> println

   5*天径*(1 + sqrt(3))/(2*sqrt(3) + 8)
   17.8735166039614

apart(等円径, d) |> simplify |> println

   5*天径/26 + 15*sqrt(3)*天径/26

function draw(r2, more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   r1 = 5r2/2
   r3 = 5r2*(1 + 3√3)/26
   @printf("天円の直径が %g のとき，等円の直径は %g，大円の直径は %g である。\n", 2r2, 2r3, 2r1)
   plot()
   circle2(√3r1/2, 0, r1)
   circle22(0, r1/2, r1)
   circle22(0, 3r1/2 - r2, r2, :green)
   circle2(√3r1/2 + r1 - r3, 0, r3, :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, r1/2, " r1/2", :red, :left, :bottom, delta=delta/4)
       point(√3r1/2, 0, "√3r1/2 ", :red, :right, delta=-delta/2)
       point(0, 3r1/2 - r2, "天円:r2\n(0,3r1/2-r2)", :green, :center, delta=delta/2)
       point(√3r1/2 + r1 - r3, 0, "等円:r3\n(√3r1/2+r1-r3, 0)", :blue, :center, :bottom, delta=delta/2)
   end
end;

draw(15/2, true)

算額（その2118）

算額（その2118）

八六 加須市多聞寺 愛宕神社 明治13年(1880)
埼玉県立図書館：埼玉県史料集 第二集『埼玉の算額』，昭和44年，誠美堂印刷所，埼玉県与野市.
キーワード：円7個，外円，円弧，最大値
#Julia, #SymPy, #算額, #和算

外円の中に，外円と同じ半径を持つ円弧 1 個と円弧に接する弦と，甲円，乙円，丙円を 2 個ずつ容れる。乙円の直径が 1 寸のとき，丙円の直径が取りうる最大値はいかほどか。
注：弦の位置（弦と y 軸の交点の y 座標）により丙円の大きさが変化する。

弦と y 軸の交点の y 座標を y
外円の半径と中心座標を R, (0, 0)
甲円の半径と中心座標を r1, (r1, 0); r1 = R/2
乙円の半径と中心座標を r2, (0, R - r2)
丙円の半径と中心座標を r3, (x3, y - r3)
とおき，以下の連立方程式を解く。

include("julia-source.txt");

using SymPy
@syms y, R::positive, r1::positive,
     r2::positive, r3::positive, x3::positive;
R = 2r1
eq1 = r1^2 + (R - r2)^2 - (r1 + r2)^2
eq2 = x3^2 + (R - r3)^2 - (r3 + R)^2
eq3 = x3^2 + (y - r3)^2 - (R - r3)^2 |> expand
res = solve([eq1, eq2, eq3], (r3, x3, r1))[1]

   ((-3*r2 + y)*(3*r2 + y)/(2*(-9*r2 + y)), -sqrt(6)*sqrt(r2)*sqrt((-3*r2 + y)*(3*r2 + y)/(-9*r2 + y)), 3*r2/2)

丙円の半径は r3 = (-3*r2 + y)*(3*r2 + y)/(2*(-9*r2 + y)) となり，y と r2 の関数である。
r2 は事前に与えられるので，y の取る値で r3 が変化する。

たとえば，乙円の直径が 1 寸のとき(注2)，-1.5 ≦ y ≦ 1.5 で，y が 0.0 〜 0.5 の区間でr3 が最大値を取る。

pyplot(size=(300, 150), grid=false, aspectratio=:none, label="")
plot(res[1](r2 => 1/2), xlims=(-1.5,1.5), xlabel="y", ylabel="r3")

r3 が最大値を取るときの y，および，そのときの最大値を求める。

r3_dash = diff(res[1], y) |> factor;  # 導関数
r3_dash |> println

   (9*r2^2 - 18*r2*y + y^2)/(2*(-9*r2 + y)^2)

r3_dash = 0 を解いて y を求める。

ans_y = solve(r3_dash, y)[1]
ans_y |> println

   3*r2*(3 - 2*sqrt(2))

甲円の半径 r2 が 1/2 のときの丙円の半径 r3 を求める。

ans_r3 = res[1](y => ans_y) |> simplify
ans_r3 |> println
ans_r3(r2 => 1/2).evalf() |> println

   3*r2*(3 - 2*sqrt(2))
   0.257359312880715

丙円の中心の x 座標を求める。
res[2] は負の値を取るが，左右いずれの甲円でも問題を解く上では無関係である。

ans_x3 = res[2](y => ans_y) |> simplify
ans_x3 |> println
ans_x3(r2 => 1/2).evalf() |> println

   -3*2^(3/4)*r2*sqrt(-4 + 3*sqrt(2))
   -1.24264068711929

最後に甲円の半径を求める。

ans_r1 = res[3]
ans_r1 |> println
ans_r1(r2 => 1/2) |>  println

   3*r2/2
   0.750000000000000

図を描くために，甲円と円弧の交点座標を求める。

@syms, x0, y0
eq4 = x0^2 +  y0^2 - R^2
eq5 = x0^2 + (y0 - y + R)^2 - R^2
res2 = solve([eq4, eq5], (x0, y0))[1]

   (-sqrt(-(-6*r1 + y)*(2*r1 + y))/2, (-2*r1 + y)/2)

まとめ

乙円の半径 r2 が与えられれば，乙円の半径が最大になるときの弦の位置 y が決まる。y = 3r2*(3 - 2√2)
丙円の半径 r3 は r3 = 3r2*(3 - 2√2) である。丙円の半径だけを求めるのならば，y は求める必要はない。

乙円の直径が 1 寸のとき，丙円の直径は 3(3 - 2√2) = 0.5147186257614291 寸である。

なお，下に示す算額の図は，丙円が最大になったときの図ではない。問題を解こうとする人を欺くためにわざとこのような図を描いたのではないかと邪推する。

function draw(r2, more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   y = 3*r2*(3 - 2*sqrt(2))
   r3 = 3*r2*(3 - 2*sqrt(2))
   x3 = 3*2^(3/4)*r2*sqrt(-4 + 3*sqrt(2))
   r1 = 3*r2/2
   R = 2r1
   x = sqrt(R^2 - y^2)
   (x0, y0) = (-sqrt(-(-6*r1 + y)*(2*r1 + y))/2, (-2*r1 + y)/2)
   @printf("乙円の直径が %g，y が %g のとき，丙円の直径は最大値 %g になる。\n", 2r2, y, 2r3)
   plot()
   circle(0, 0, R, :green)
   circle2(r1, 0, r1)
   circle22(0, R - r2, r2, :blue)
   circle2(x3, y - r3, r3, :magenta)
   segment(-x, y, x, y)
   θ = atand(y0, x0)
   circle(0, y - R, R, :green, beginangle=-θ, endangle=180+θ)
   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", :green, :center, :bottom, delta=delta)
       point(0, y, "y", :black, :center, :bottom, delta=delta)
       point(r1, 0, "甲円:r1,(r1,0)", :red, :right, delta=-delta, deltax=3delta)
       point(0, R - r2, "乙円:r2\n(0,R-r2)", :blue, :center, delta=-delta)
       point(0, r2 - R, "乙円:r2\n(0,r2-R)", :blue, :center, :bottom, delta=delta)
       point(x3, y - r3, "丙円:r3,(y-r3,0)", :magenta, :right, :bottom, delta=delta, deltax=3delta)
       point(0, y - R, "円弧:R,(0,y-R)", :green, :center, delta=-delta)
   end
end;

draw(1/2, true)