算額（その1423）

算額（その1423）

九十七 大船渡市猪川町 雨宝堂 現雨宝山竜宝院 文政7年(1824)
山村善夫：現存 岩手の算額，昭和52年1月30日，熊谷印刷，盛岡市.
http://www.wasan.jp/yamamura/yamamura.html
キーワード：円12個，外円，1/3円3個
#Julia, #SymPy, #算額, #和算

全円の中に全円と同じ直径の 1/3 円が 6 個と，甲円および乙円を 6 個ずつ容れる。乙円の直径が 1 寸のとき，全円の直径はいかほどか。

全円の半径と中心座標を R, (0, 0)
条件式に使う 1/3 円の中心座標を (0, R)
甲円の半径と中心座標を r1, (R - r1, 0)
乙円の半径と中心座標を r2, (R - 2r1 - r2, 0)
として以下の連立方程式の解を求める。

include("julia-source.txt");

using SymPy
@syms R, r1, r2
eq1 = ( R - r1)^2 + R^2 - (R + r1)^2
eq2 = (R - 2r1 - r2)^2 + R^2 - (R + r2)^2
res = solve([eq1, eq2], (R, r1))[3]

    (12*r2, 3*r2)

全円の半径 R は，乙円の半径 r2 の 12 倍である。
乙円の直径が 1 寸のとき，全円の直径は 12 寸である。

function draw(r2, more=false)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    (R, r1) = (12*r2, 3*r2)
    plot()
    circle(0, 0, R, :magenta)
    for i = 90:60:390
        circle(R*cosd(i), R*sind(i), R, :blue, beginangle=i + 120, endangle=i + 240)
    end
    rotate(R - r1, 0, r1, angle=60)
    rotate(R - 2r1 - r2, 0, r2, :green, angle=60)
    if more
        delta = (fontheight = (ylims()[2]- ylims()[1]) / 500 * 10 * 2) /  3  # size[2] * fontsize * 2
        vline!([0], color=:gray80, lw=0.5)
        hline!([0], color=:gray80, lw=0.5)
        point(R - r1, 0, "甲円:r1\n(R-r1,0)", :red, :center, delta=-delta)
        point(R - 2r1 - r2, 0, "乙円:r2\n(R-2r1-r2,0)", :green, :left, delta=-6delta/2, deltax=-9delta)
        point(0, R, "R", :magenta, :center, :bottom, delta=delta/2)
    end
end;

draw(1/2, true)

算額（その1422）

算額（その1422）

百一 大船渡市猪川町 田茂山神社 奉納年不明 現存せず
山村善夫：現存 岩手の算額，昭和52年1月30日，熊谷印刷，盛岡市.
http://www.wasan.jp/yamamura/yamamura.html
キーワード：円16個，外円
#Julia, #SymPy, #算額, #和算

全円の中に中円 1 個，甲円，乙円，丙円を 5 個ずつ入れる。中円の直径が 1 寸のとき，全円の直径はいかほどか。

全円の半径と中心座標を R, (0, 0)
中円の半径と中心座標を r4, (0, 0)
甲円の半径と中心座標を r1, (0, R - r1)
乙円の半径と中心座標を r2, (x2, y2)
丙円の半径と中心座標を r3, (x3, y3)
ただし，
x2 = (R - r2)*cosd(54)
y2 = (R - r2)*sind(54)
x3 = (R - 2r2 - r3)*cosd(54)
y3 = (R - 2r2 - r3)*sind(54)
として以下の連立方程式の解を求める。

include("julia-source.txt");

using SymPy
@syms R, r1, r2, x2, y2, r3, x3, y3, r4
x2 = (R - r2)*cosd(Sym(54))
y2 = (R - r2)*sind(Sym(54))
x3 = (R - 2r2 - r3)*cosd(Sym(54))
y3 = (R - 2r2 - r3)*sind(Sym(54))
eq1 = 2r1 + r4 - R
eq2 = r4 + 2r2 + 2r3 - R
eq3 = x2^2 + (R - r1 - y2)^2 - (r1 + r2)^2
eq4 = x3^2 + (R - r1 - y3)^2 - (r1 + r3)^2;

res = solve([eq1, eq2, eq3, eq4], (R, r1, r2, r3))[3]  # 3 of 3

    (-56*sqrt(5)*r4*sqrt(34 - 10*sqrt(5))/5 - 24*sqrt(5)*r4*sqrt(170 - 50*sqrt(5))/5 - 16*sqrt(5)*r4/5 + 43*r4/5 + 54*r4*sqrt(170 - 50*sqrt(5))/5 + 126*r4*sqrt(34 - 10*sqrt(5))/5, -28*sqrt(5)*r4*sqrt(34 - 10*sqrt(5))/5 - 12*sqrt(5)*r4*sqrt(170 - 50*sqrt(5))/5 - 8*sqrt(5)*r4/5 + 19*r4/5 + 27*r4*sqrt(170 - 50*sqrt(5))/5 + 63*r4*sqrt(34 - 10*sqrt(5))/5, -31*sqrt(5)*r4*sqrt(34 - 10*sqrt(5))/10 - 13*sqrt(5)*r4*sqrt(170 - 50*sqrt(5))/10 - 7*sqrt(5)*r4/5 + 33*r4/10 + 117*r4*sqrt(170 - 50*sqrt(5))/40 + 279*r4*sqrt(34 - 10*sqrt(5))/40, (-sqrt(5)*r4 + 3*r4 + sqrt(2)*r4*sqrt(17 - 5*sqrt(5)))/(2*sqrt(5) + 10))

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

    -r4*(-43 - 6*sqrt(2)*sqrt(17 - 5*sqrt(5)) + 2*sqrt(10)*sqrt(17 - 5*sqrt(5)) + 16*sqrt(5))/5

全円の半径 R は，中円の半径 r4 の (43 + 6√2sqrt(17 - 5√5) - 2√10*sqrt(17 - 5√5) - 16√5)/5 = 2.487088356114319 倍である。
中円の直径が 1 寸のとき，全円の直径は 2.487088356114319 寸である。

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

    -r4*(-19 - 3*sqrt(2)*sqrt(17 - 5*sqrt(5)) + sqrt(10)*sqrt(17 - 5*sqrt(5)) + 8*sqrt(5))/5

# r2
res[3] |> factor |> println

    -r4*(-132 - 19*sqrt(2)*sqrt(17 - 5*sqrt(5)) + 7*sqrt(10)*sqrt(17 - 5*sqrt(5)) + 56*sqrt(5))/40

# r3
@syms d
res[4]/r4 |> factor |> x -> apart(x, d)*r4 |> factor |> println

    -r4*(-20 - 5*sqrt(2)*sqrt(17 - 5*sqrt(5)) + sqrt(10)*sqrt(17 - 5*sqrt(5)) + 8*sqrt(5))/40

R, r1, r2, r3 は，共通部分式 sqrt(17 - 5√5) の括りだしと，r3 については分母の有理化などで，以下のように規則性のある簡潔な計算式になる。

r4 = 1/2
s = sqrt(17 - 5√5)
t = √2s
u = √10s
R  = r4*( 43 +  6t - 2u - 16√5)/5
r1 = r4*( 19 +  3t -  u -  8√5)/5
r2 = r4*(132 + 19t - 7u - 56√5)/40
r3 = r4*( 20 +  5t -  u -  8√5)/40
(R, r1, r2, r3)

    (1.2435441780571594, 0.3717720890285797, 0.22750945795996708, 0.14426263106861265)

function draw(r4, more=false)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    (R, r1, r2, r3) = (-56*sqrt(5)*r4*sqrt(34 - 10*sqrt(5))/5 - 24*sqrt(5)*r4*sqrt(170 - 50*sqrt(5))/5 - 16*sqrt(5)*r4/5 + 43*r4/5 + 54*r4*sqrt(170 - 50*sqrt(5))/5 + 126*r4*sqrt(34 - 10*sqrt(5))/5, -28*sqrt(5)*r4*sqrt(34 - 10*sqrt(5))/5 - 12*sqrt(5)*r4*sqrt(170 - 50*sqrt(5))/5 - 8*sqrt(5)*r4/5 + 19*r4/5 + 27*r4*sqrt(170 - 50*sqrt(5))/5 + 63*r4*sqrt(34 - 10*sqrt(5))/5, -31*sqrt(5)*r4*sqrt(34 - 10*sqrt(5))/10 - 13*sqrt(5)*r4*sqrt(170 - 50*sqrt(5))/10 - 7*sqrt(5)*r4/5 + 33*r4/10 + 117*r4*sqrt(170 - 50*sqrt(5))/40 + 279*r4*sqrt(34 - 10*sqrt(5))/40, (-sqrt(5)*r4 + 3*r4 + sqrt(2)*r4*sqrt(17 - 5*sqrt(5)))/(2*sqrt(5) + 10))
    s = sqrt(17 - 5√5)
    t = √2s
    u = √10s
    R  = r4*( 43 +  6t - 2u - 16√5)/5
    r1 = r4*( 19 +  3t -  u -  8√5)/5
    r2 = r4*(132 + 19t - 7u - 56√5)/40
    r3 = r4*( 20 +  5t -  u -  8√5)/40
    x2 = (R - r2)*cosd(54)
    y2 = (R - r2)*sind(54)
    x3 = (R - 2r2 - r3)*cosd(54)
    y3 = (R - 2r2 - r3)*sind(54)
    @printf("外円の直径 = %g\n", 2R)
    @printf("r4 = %g;  R = %g;  r1 = %g;  r2 = %g;  r3 = %g;  x2 = %g;  y2 = %g;  x3 = %g;  y3 = %g\n", r4, R, r1, r2, r3, x2, y2, x3, y3)
    plot()
    circle(0, 0, R, :magenta)
    rotate(0, R - r1, r1, angle=72)
    rotate(x2, y2, r2, :green, angle=72)
    rotate(x3, y3, r3, :blue, angle=72)
    circle(0, 0, r4, :orange)
    if more
        delta = (fontheight = (ylims()[2]- ylims()[1]) / 500 * 10 * 2) /  3  # size[2] * fontsize * 2
        vline!([0], color=:gray80, lw=0.5)
        hline!([0], color=:gray80, lw=0.5)
        point(0, R - r1, "甲円:r1,(0,R-r1)", :red, :center, delta=-delta/2)
        point(x2, y2, "乙円:r2\n(x2,y2)", :green, :center, delta=-delta/2)
        point(x3, y3, "丙円:r3,(x3,y3)", :green, :left, delta=-delta/2, deltax=-2delta)
        point(0, 0, "中円:r4,(0,0)", :orange, :center, delta=-delta/2)
        point(0, R, "R", :magenta, :center, :bottom, delta=delta/2)
    end
end;

draw(1/2, true)