算額（その543）

算額（その543）

群馬の算額 43-1 清水寺 文政7年
http://takasakiwasan.web.fc2.com/gunnsann/g043-1.html

外円内に大円，中円，小円が入っている。外円の直径が 4 寸のとき，小円の大きさを最大にしたときの直径を求めよ。

本問は，「算額（その527）」と基本的に同じである（r3 が最も大きくなるときの r1 を求める）。

外円の半径と中心座標を R, (0, 0)
大円の半径と中心座標を r1, (0, R - r1)
中円の半径と中心座標を r2, (r2, y2)
小円の半径と中心座標を r3, (x3, y3)
とおき，以下の連立方程式を r1 について解く。

include("julia-source.txt");

using SymPy

@syms R, r1, r2, y2::negative, r3, x3, y3::negative
R = 4//2
eq1 = r2^2 + y2^2 - (R - r2)^2
eq2 = x3^2 + y3^2 - (R - r3)^2
eq3 = r2^2 + (R - r1 - y2)^2 - (r1 + r2)^2
eq4 = x3^2 + (R - r1 - y3)^2 - (r1 + r3)^2
eq5 = (x3 - r2)^2 + (y3 - y2)^2 - (r3 + r2)^2;
res = solve([eq1, eq2, eq3, eq4, eq5], (r2, y2, r3, x3, y3))

   2-element Vector{NTuple{5, Sym}}:
    (-8*r1*(r1 - 2)/(r1 + 2)^2, -2*(3*r1 - 2)/(r1 + 2), -8*r1*(r1 - 2)/(r1 + 2)^2, 8*r1*(r1 - 2)/(r1 + 2)^2, (4 - 6*r1)/(r1 + 2))
    (-8*r1*(r1 - 2)/(r1 + 2)^2, -2*(3*r1 - 2)/(r1 + 2), -8*r1*(r1 - 2)/(9*r1^2 - 28*r1 + 36), -24*r1*(r1 - 2)/(9*r1^2 - 28*r1 + 36), (10*r1^2 - 72*r1 + 72)/(9*r1^2 - 28*r1 + 36))

2 組の解が得られるが，2 番目の解が適解である。その解で r3 は r1 の式で表される。

r3 = res[2][3]
r3 |> println

   -8*r1*(r1 - 2)/(9*r1^2 - 28*r1 + 36)

図で表すと以下のようになり，r1 が 1.0 〜 1.5 のあたりで最大値になる。

pyplot(size=(300, 200), grid=false, aspectratio=:none, label="")
plot(r3, xlims=(0, 2), xlabel="r1", ylabel="r3")

数値的に解くには，r3 式を微分しその式の値が 0 になるときの r1 を求めればよい。

solve(diff(r3)) |> println

   Sym[6/5, 6]

解は 2 通り得られるが r1 = 6/5 が適解である。

r1 = 6/5 のときの r3 は，式に r1 を代入すれば，r3 の最大値が求まる。

r3(r1 => 6//5) |> println

   1/2

大円の半径が 6/5 のとき，小円の半径は最大値 1/2 となる。
大円の直径が 12/5 のとき，小円の直径は最大値 1 となる。

そのときのその他のパラメータは以下の通り。

  r2 = 0.75;  y2 = -1;  r3 = 0.5;  x3 = 1.5;  y3 = 0

function draw(more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   R = 4//2
   r1 = (12/5)/2
   (r2, y2, r3, x3, y3) = (-8*r1*(r1 - 2)/(r1 + 2)^2, -2*(3*r1 - 2)/(r1 + 2), -8*r1*(r1 - 2)/(r1 + 2)^2, 8*r1*(r1 - 2)/(r1 + 2)^2, (4 - 6*r1)/(r1 + 2))
   (r2, y2, r3, x3, y3) = (-8*r1*(r1 - 2)/(r1 + 2)^2, -2*(3*r1 - 2)/(r1 + 2), -8*r1*(r1 - 2)/(9*r1^2 - 28*r1 + 36), -24*r1*(r1 - 2)/(9*r1^2 - 28*r1 + 36), (10*r1^2 - 72*r1 + 72)/(9*r1^2 - 28*r1 + 36))
   @printf("r2 = %g;  y2 = %g;  r3 = %g;  x3 = %g;  y3 = %g\n", r2, y2, r3, x3, y3)
   plot()
   circle(0, 0, R)
   circle(0, R - r1, r1, :blue)
   circle(r2, y2, r2, :green)
   circle(-r2, y2, r2, :green)
   circle(x3, y3, r3, :orange)
   circle(-x3, y3, r3, :orange)
   if more
       delta = (fontheight = (ylims()[2]- ylims()[1]) / 500 * 10 * 2) /3  # size[2] * fontsize * 2
       hline!([0], color=:black, lw=0.5)
       vline!([0], color=:black, lw=0.5)
       point(R, 0, "R ", :red, :right, :bottom, delta=delta/2)
       point(0, R - r1, " 大円:r1,(0,R-r1)", :blue, :left, :vcenter)
       point(r2, y2, "中円:r2,(r2,y2)", :green, :center, delta=-delta/2)
       point(x3, y3, "小円:r3,(x3,y3)", :orange, :center, delta=-delta/2)
   end
end;

 

算額（その542）

算額（その542）

『福島の算額』 矢祭神社の算額 明治21年(1888)
http://takasakiwasan.web.fc2.com/gunnsann/yamaturi.html

正三角形内に全円，大円，小円を入れる。大円，小円の直径が 8 寸，6 寸のとき，全円の直径はいかほどか。

正三角形の一辺の長さを 2a とする。
全円の半径と中心座標を r1, (0, r1)
大円の半径と中心座標を r2, (x2, 0)
小円の半径と中心座標を r3, (x3, r3)
として，以下の連立方程式を解く。

include("julia-source.txt");

using SymPy

@syms a, r1, r2, x2::negative, r3, x3
eq1 = (x3 - x2)^2 + (r3 - r2)^2 - (r3 + r2)^2
eq2 = (a - x3)/r3 - sqrt(Sym(3))
eq3 = (a + x2)/r2 - sqrt(Sym(3))
eq4 = a/r1 - sqrt(Sym(3))
res = solve([eq1, eq2, eq3, eq4], (a, r1, x2, x3))

   2-element Vector{NTuple{4, Sym}}:
    (sqrt(3)*r2/2 + sqrt(3)*r3/2 - sqrt(r2*r3), r2/2 + r3/2 - sqrt(3)*sqrt(r2*r3)/3, sqrt(3)*r2/2 - sqrt(3)*r3/2 + sqrt(r2*r3), -sqrt(r2*r3) + sqrt(3)*(r2 - r3)/2)
    (sqrt(3)*r2/2 + sqrt(3)*r3/2 + sqrt(r2*r3), r2/2 + r3/2 + sqrt(3)*sqrt(r2*r3)/3, sqrt(3)*r2/2 - sqrt(3)*r3/2 - sqrt(r2*r3), sqrt(r2*r3) + sqrt(3)*(r2 - r3)/2)

2 組の解が得られるが，2 番目のものが適解である。

全円の半径は 11/2 である。

res[2][2] |> println
res[2][2](r2 => 8//2, r3 => 6//2) |> println

   r2/2 + r3/2 + sqrt(3)*sqrt(r2*r3)/3
   11/2

その他のパラメータは以下の通りである。

a = 9.52628;  r1 = 5.5;  x2 = -2.59808;  x3 = 4.33013

function draw(more=false)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    (r2, r3) = (8, 6) .// 2
    (a, r1, x2, x3) = (
        (r2 + r3)√3/2 + sqrt(r2*r3),
        (r2 + r3)/2 + sqrt(r2*r3)√3/3,
        (r2 - r3)√3/2 - sqrt(r2*r3),
        sqrt(r2*r3) + (r2 - r3)√3/2)
    @printf("a = %g;  r1 = %g;  x2 = %g;  x3 = %g\n", a, r1, x2, x3)
    plot([a, 0, -a, a], [0, √3a, 0, 0], color=:brown, lw=0.5)
    circle(0, r1, r1)
    circle(x2, r2, r2, :green)
    circle(x3, r3, r3, :orange)
    if more
        delta = (fontheight = (ylims()[2]- ylims()[1]) / 500 * 10 * 2) /3  # size[2] * fontsize * 2
        hline!([0], color=:black, lw=0.5)
        vline!([0], color=:black, lw=0.5)
        point(0, r1, " 全円:r1,(0,r1)", :red, :left, :vcenter)
        point(x2, r2, "大円:r2,(x2,r2)", :green, :center, :top, delta=-delta)
        point(x3, r3, "小円:r3,(x3,r3)", :black, :center, :top, delta=-delta)
        point(a, 0, "a ", :brown, :right, :bottom, delta=delta/2)
        point(0, √3a, " √3a", :brown, :left, :vcenter)
    end
end;

算額（その541）

算額（その541）

群馬県吉岡町漆原 鎮守社 弘化3年(1846)
http://takasakiwasan.web.fc2.com/gunnsann/pdf/070_2016_08_03.pdf

弓形内に長さの等しい弦を 2 本，黃円 1 個，青円 2 個を入れる。黃円，青円の直径がそれぞれ 12 寸，5 寸のとき，弦（赤斜と呼ぶ）の長さはいかほどか。

なお，算額の図はほぼ上の図に近いものであるが，この図は「黃円，青円の直径がそれぞれ 12 寸，3.6 寸のとき」のものである。
算額の問の通り「黃円，青円の直径がそれぞれ 12 寸，5 寸のとき」のものは下図に示すもののようになる。

弓形を作る円の半径と中心座標を R, (0, 0)
黃円の半径と中心座標を r1, (0, y + r1)
青円の半径と中心座標を r2, (x2, y2)
黃円と接する x 軸に平行な弦と y 軸の交点の y 座標を y とする。弦と円の交点座標(sqrt(R^2 - y^2), y) である。
赤斜は sqrt((R - y)^2 + x^2) である。
以下の連立方程式を解く。

include("julia-source.txt");

using SymPy

@syms R, x, y, r1, r2, x2::positive, y2::positive, a
a = sqrt((R - y)^2 + x^2)/2
eq1 = R^2 - (R - 2r2)^2 - a^2
eq2 = R*(R - y) - 2a^2
eq3 = R*r1 - (R - y - r1)*(R - 2r2)
eq4 = x2^2 + y2^2 - (R - r2)^2
eq5 = y2/x2 - (R - 2r2)/a;

res = solve([eq1, eq2, eq3, eq4, eq5], (R, x, y, x2, y2));

8 組の解が得られるが，4 番目のものが適解である。

for i in 1:5
   println(res[4][i])
end

   -8*r2^2/(r1 - 4*r2)
   r1*sqrt((-r1 - 4*r2)/(r1 - 4*r2))
   (-r1^2 + 8*r2^2)/(r1 - 4*r2)
   -sqrt(-(r1^6 - 48*r1^4*r2^2 + 768*r1^2*r2^4 - 4096*r2^6)/(r1^2 - 8*r1*r2 + 16*r2^2))/(4*r1 - 16*r2)
   -r1*(r1 + 4*r2)/(4*r1 - 16*r2)

r1 = 6;  r2 = 2.5
R = 12.5;  x = 12;  y = -3.5;  x2 = 8;  y2 = 6

このとき，赤斜は 20 寸である。

R = 12.5;  x = 12;  y = -3.5
sqrt((R - y)^2 + x^2)

   20.0

function draw(r1, r2, more=false)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   (R, x, y, x2, y2) = (
       -8*r2^2/(r1 - 4*r2),
       r1*sqrt((-r1 - 4*r2)/(r1 - 4*r2)),
       (-r1^2 + 8*r2^2)/(r1 - 4*r2),
       -sqrt(-(r1^6 - 48*r1^4*r2^2 + 768*r1^2*r2^4 - 4096*r2^6)/(r1^2 - 8*r1*r2 + 16*r2^2))/(4*r1 - 16*r2),
       -r1*(r1 + 4*r2)/(4*r1 - 16*r2))
   a = sqrt((R - y)^2 + x^2)
   θ = atand(y, x)
   @printf("r1 = %g;  r2 = %g;  R = %g;  x = %g;  y = %g;  x2 = %g;  y2 = %g\n", r1, r2, R, x, y, x2, y2)
   @printf("赤斜 = %g\n", a)
   plot()
   circle(0, 0, R, :black, beginangle=round(θ), endangle=round(180-θ))
   circle(0, y + r1, r1, :orange)
   circle(x2, y2, r2, :blue)
   circle(-x2, y2, r2, :blue)
   plot!([x, -x, 0], [y, y, R], color=:green, lw=0.5)
   segment(0, R, x, y, :red)
    if more
       delta = (fontheight = (ylims()[2]- ylims()[1]) / 500 * 10 * 2) /3  # size[2] * fontsize * 2
       hline!([0], color=:black, lw=0.5)
       vline!([0], color=:black, lw=0.5)
       point(x, y, "(x,y)", :blue, :right, :top, delta=-delta/2)
       point(0, y + r1, " 黃円:r1\n (0,y+r1)", :orange, :left, :vcenter)
       point(x2, y2, " 青円:r1\n (0,y+r1)", :black, :left, :vcenter)
       point(0, y, " y", :green, :left, :top, delta=-delta/2)
       point(0, R, " R", :black, :left, :bottom, delta=delta/2)
    end
end;