算額（その1514）

算額（その1514）

七十三 岩手県一関市川崎町門崎宮畑 伊吹神社 明治17年(1884)
山村善夫：現存 岩手の算額，昭和52年1月30日，熊谷印刷，盛岡市.
http://www.wasan.jp/yamamura/yamamura.html

今有如図 03080
https://w.atwiki.jp/sangaku/pages/320.html

キーワード：球8個，3次元
#Julia, #SymPy, #算額, #和算

大球の中に小球 5 個，中球 2 個を容れる。小球は互いに外接しあい，すべてが中球とも外接し，大球に内接している。中球は大球に内接している（中球同士は外接していない）。小球の直径が与えられたとき，大球の直径はいかほどか。

大球の半径と中心座標を r1, (0, 0, 0)
中球の半径と中心座標を r2, (0, 0, z)
小球の半径と中心座標を r3, (l, 0, 0); l = r3/sin(π/10)
とおき，以下の方程式を解く。

include("julia-source.txt")
# julia-source.txt ソース　https://blog.goo.ne.jp/r-de-r/e/ad3a427b84bb416c4f5b73089ae813cf

using SymPy
@syms r1, r2, r3, z
l = r3/sind(Sym(36))
eq1 = z^2 + l^2 - (r2 + r3)^2
eq2 = z + r2 - r1
eq3 = l + r3 - r1;
res = solve([eq1, eq2, eq3], (r1, r2, z))[1]

    (r3 + 2*sqrt(2)*r3/sqrt(5 - sqrt(5)), r3*(-15*sqrt(2) - 10*sqrt(5 - sqrt(5)) + 2*sqrt(25 - 5*sqrt(5)) + 5*sqrt(10))/(5*(-3*sqrt(5 - sqrt(5)) - 3*sqrt(2) + sqrt(10) + sqrt(5)*sqrt(5 - sqrt(5)))), 2*r3*(-5*sqrt(5) - sqrt(50 - 10*sqrt(5)) + 5*sqrt(10 - 2*sqrt(5)) + 15)/(-10*sqrt(5) - sqrt(10)*sqrt(5 - sqrt(5)) + 5*sqrt(2)*sqrt(5 - sqrt(5)) + 30))

# r1
res[1] |> factor |> println

    r3*(sqrt(5 - sqrt(5)) + 2*sqrt(2))/sqrt(5 - sqrt(5))  

大球の半径 r1 は，小球の半径 r3 の (sqrt(5 - √5) + 2√2)/sqrt(5 - √5) = 2.70130161670408 倍である。

術では単に「3 倍」とあるが，理解不能である。

# r2
@syms d
apart(res[2]/r3, d) |> factor |> println

    -(-30 - 10*sqrt(5) + 5*sqrt(2)*sqrt(5 - sqrt(5)) + 3*sqrt(10)*sqrt(5 - sqrt(5)))/20

# z
res[3] |> simplify |> println

    2*r3*(-5*sqrt(5) - sqrt(50 - 10*sqrt(5)) + 5*sqrt(10 - 2*sqrt(5)) + 15)/(-10*sqrt(5) - sqrt(10)*sqrt(5 - sqrt(5)) + 5*sqrt(2)*sqrt(5 - sqrt(5)) + 30)

function draw(r3, more=false)
    pyplot(size=(500, 250), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    (r1, r2, z) = (r3 + 2*sqrt(2)*r3/sqrt(5 - sqrt(5)), r3*(-15*sqrt(2) - 10*sqrt(5 - sqrt(5)) + 2*sqrt(25 - 5*sqrt(5)) + 5*sqrt(10))/(5*(-3*sqrt(5 - sqrt(5)) - 3*sqrt(2) + sqrt(10) + sqrt(5)*sqrt(5 - sqrt(5)))), 2*r3*(-5*sqrt(5) - sqrt(50 - 10*sqrt(5)) + 5*sqrt(10 - 2*sqrt(5)) + 15)/(-10*sqrt(5) - sqrt(10)*sqrt(5 - sqrt(5)) + 5*sqrt(2)*sqrt(5 - sqrt(5)) + 30))
    l = r3/sind(Sym(36))
    p1 = plot()
    circle(0, 0, r1, :green)
    circle(0, 0, r2, :blue)
    rotate(l, 0, r3, angle=72)
    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)
    end
    p2 = plot()
    circle(0, 0, r1, :green)
    circle22(0, z, r2, :blue)
    circle(l, 0, r3)
    circle(l*cosd(72), 0, r3)
    circle(l*cosd(144), 0, r3)
    plot(p1, p2)
end;

draw(1, true)

算額（その1513）

算額（その1513）

六十六 岩手県一関市花泉町金沢大柳 金沢八幡宮 明治29年(1896)
山村善夫：現存 岩手の算額，昭和52年1月30日，熊谷印刷，盛岡市.
http://www.wasan.jp/yamamura/yamamura.html

今有如図 03080
https://w.atwiki.jp/sangaku/pages/318.html

キーワード：球6個，3次元
#Julia, #SymPy, #算額, #和算

大球 3 個が互いに接し合って盤上に載っている。大球 3 個の上部の凹みに，互いに接し合って 3 個の 小球が載っている。3 個の大球の頂点と 3 個の小球の頂点は同じ水平面に接している（小球の頂点の z 座標は大球の直径と等しい）。大球の直径が与えられたとき，小球の直径はいかほどか。

山村の図はひどいが，「今有如図」も似たりよったりだ。3次元の物体の図を描くのは難しいが，球の間の関係式を解くのは円の間の関係式を解くのと同じだ。

互いに接する大球と小球の位置関係を考える。図の A と B
大球の半径と中心座標を r1, (r1, -r1/√3, r1)
小球の半径と中心座標を r2, (r2, r2/√3, 2r1 - r2)
とおき，以下の方程式を解く。

include("julia-source.txt")
# julia-source.txt ソース　https://blog.goo.ne.jp/r-de-r/e/ad3a427b84bb416c4f5b73089ae813cf

using SymPy
@syms r1, r2
eq = (r2 - r1)^2 + (r2 + r1)^2/3 + (r1 - r2)^2 - (r1 + r2)^2
ans_r2 = solve(eq, r2)[1]
ans_r2 |> println

    r1*(2 - sqrt(3))

小球の半径 r2 は，大球の半径 r1 の (2 - √3) 倍である。
大球の直径が 1 寸のとき，小球の直径は (2 - √3) = 0.2679491924311228

function draw(r1, more=false)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    r2 = r1*(2 - √3)
    plot()
    rotate(0, 2r1/√3, r1, :blue)
    rotate(0, -2r2/√3, r2)
    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/√3, "A:大球,(r1,-r1/√3,r1)", :blue, :center, delta=-delta)
        point(r2, r2/√3, "B:小球\n(r2/r2/√3,2r1-r2)", :red, :left, delta=delta, deltax=7delta)
    end
end;

draw(1, true)

算額（その1512）

算額（その1512）

六十六 岩手県一関市花泉町金沢大柳 金沢八幡宮 明治29年(1896)
山村善夫：現存 岩手の算額，昭和52年1月30日，熊谷印刷，盛岡市.
http://www.wasan.jp/yamamura/yamamura.html

今有如図 03080
https://w.atwiki.jp/sangaku/pages/318.html

キーワード：円2個，四分円2個，正三角形
#Julia, #SymPy, #算額, #和算

正三角形と四分円2個を組み合わせ，斜線を 2 本引き，隙間に円を容れる。正三角形の一辺の長さが与えられたとき，円の直径を求める術を述べよ。

正三角形の一辺の長さを 2a
円の半径と中心座標を r,(x, y)
とおき，以下の連立方程式を解く。

include("julia-source.txt")
# julia-source.txt ソース　https://blog.goo.ne.jp/r-de-r/e/ad3a427b84bb416c4f5b73089ae813cf

using SymPy
@syms a, r, x, y
eq1 = dist(-a, 0, 2a*cosd(Sym(30)), 2a*sind(Sym(30)), x, y) - r^2
eq2 = dist(0, √Sym(3)*a, a, 0, x, y) - r^2
eq3 = x^2 + (√Sym(3)*a - y)^2 - (2a - r)^2;
# res = solve([eq1, eq2, eq3], (r, x, y))

function H(u)
    (r, x, y) = u
    return [
        dist(-a, 0, 2a*cosd(30), √3a + 2a*sind(30), x, y) - r^2,
        dist(0, √3a, a, 0, x, y) - r^2,
        x^2 + (√3a - y)^2 - (2a - r)^2
    ]
end;
a = 1/2
iniv = BigFloat[0.3, 0.6, 0.6]
res = nls(H, ini=iniv)

    ([0.35430062553187885, 0.576740154417105, 0.5756834046326507], true)

正三角形の一辺の長さが 1 のとき，円の直径は res[1][1]*2 = 0.7086012510637577 である。

function draw(more=false)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    (r, x, y) = res[1]
    plot([a, 0, -a, a], [0, √3a, 0, 0], color=:magenta, lw=0.5)
    println("正三角形の一辺の長さが $(2a) のとき，円の直径は $(2r) である。")
    circle(0, √3a, 2a, beginangle=-60, endangle=30)
    circle(x, y, r, :blue)
    segment(-a, 0, 2a*cosd(30), √3a + 2a*sind(30), :green)
    segment(0, √3a, 2a*cosd(30), √3a + 2a*sind(30), :green)
    circle(-a, 0, 2a, beginangle=60, endangle=150)
    segment(-a, 0, -a - 2a*cosd(30), 2a*sind(30), :green)
    segment(a, 0, -a - 2a*cosd(30), 2a*sind(30), :green)
    ox = 0
    oy = √3a/3
    len = sqrt(x^2 + (y - oy)^2)
    θ = atand(y - √3a/3, x) + 120
    x2 = len*cosd(θ)
    y2 = len*sind(θ) + oy
    circle(x2, y2, r, :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(x, y, "r,(x,y)", :blue, :center, delta=-delta)
        point(0, √3a, "√3a", :magenta, :center, :bottom, delta=2delta)
        point(a, 0, "a", :magenta, :center, delta=-1.5delta)
        ylims!(-7delta, √3a + 2a*sind(30) + 5delta)
    end
end;

draw(true)

術は以下のようである。SymPy では，無理だ。

三角面 = 1
天 = 2 - √3
地 = √天/2 + 1
人 = 1 - √天/2
A = (1 - √0.5)/人
円径 = 4(2天 - A)*地*三角面

    0.7086012510637592