算額（その1485）

算額（その1485）

四十四 岩手県一関市真滝 熊野白山滝神社 文久元年(1861)
山村善夫：現存 岩手の算額，昭和52年1月30日，熊谷印刷，盛岡市.
http://www.wasan.jp/yamamura/yamamura.html
キーワード：円5個，二等辺三角形，弦2本
#Julia, #SymPy, #算額, #和算

問題文，答，術ともに不明であるが，以下のようなものでもあろう。

外円の中に二等辺三角形と 2 本の弦を容れ，隙間に大円 1 個，中円 1 個，小円 2 個を容れる。外円の直径が与えられたとき，大円，中円，小円の直径を求める術を述べよ。

注：これだけの条件では大円の大きさは不定である（中円・小円とは無関係に，どのような大きさにもなれる）。これを解消するための制約条件はいくつか考えられるが，ここでは，外円の円周上にある ５ 点が，外円に内接する正五角形の頂点であるということにする。

外円の半径と中心座標を R, (0, 0)
大円の半径と中心座標を r1, (0, r1 - R)
中円の半径と中心座標を r2, (0, R*sind(18) + r2)
小円の半径と中心座標を r3, (R*sind(18) + r3)
とおき，以下の連立方程式を解く。
なお，∠EAD = 54°，∠BAC = 18°，∠ADE = 36° である。

include("julia-source.txt");

using SymPy

@syms R, r1, r2, r3, x3;
# eq1 = dist2(R*sind(Sym(36)), -R*cosd(Sym(36)), 0, R, 0, r1 - R, r1)
# eq2 = dist2(R*cosd(Sym(18)), R*sind(Sym(18)), 0, R, 0, R*sind(Sym(18)) + r2, r2)
# eq3 = dist2(R*cosd(Sym(18)), R*sind(Sym(18)), 0, R, x3, R*sind(Sym(18)) + r3, r3)
# eq4 = x3^2 + (r2 - r3)^2 - (r2 + r3)^2

eq11 = r1/(2R - r1) - sind(Sym(18));
eq12 = r2/(R - (r2 + R*sind(Sym(18)))) - sind(Sym(54)) |> simplify;
eq13 = r2*(R*cosd(Sym(18)) - x3) - r3*(R*cosd(Sym(18)));
eq14 = x3 - 2sqrt(r2*r3);
res = solve([eq11, eq12, eq13, eq14], (r1, r2, r3, x3))[2]

    (-4*R + 2*sqrt(5)*R, R*(-1 + sqrt(5))/4, R*(-35 - 6*sqrt(30 - 10*sqrt(5)) + 10*sqrt(6 - 2*sqrt(5)) + 19*sqrt(5))/20, R*((-5*sqrt(2) + 2*sqrt(10))*sqrt(sqrt(5) + 5) + 4*sqrt(25 - 10*sqrt(5)))/10)

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

    2*R*(-2 + sqrt(5))

# r2
res[2] |> println

    R*(-1 + sqrt(5))/4

# r3
res[3] |> sympy.sqrtdenest |> factor |> println

    R*(-15 + 7*sqrt(5))/4

# x3
res[4] |>sympy.sqrtdenest |> simplify |> println

    R*(-5*sqrt(2*sqrt(5) + 10) + 4*sqrt(25 - 10*sqrt(5)) + 2*sqrt(10*sqrt(5) + 50))/10

大円，中円，小円の半径は，外円の半径の 2(√5 - 2) 倍，(√5 - 1)/4 倍，(7√5 - 15)/4 倍である。

外円の直径が 1 のとき，大円，中円，小円の直径は 0.4721359549995796, 0.30901699437494745, 0.16311896062463216 である。

function draw(R, more=false)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    (r1, r2, r3, x3) = (-4*R + 2*sqrt(5)*R, R*(-1 + sqrt(5))/4, R*(-35 - 6*sqrt(30 - 10*sqrt(5)) + 10*sqrt(6 - 2*sqrt(5)) + 19*sqrt(5))/20, R*((-5*sqrt(2) + 2*sqrt(10))*sqrt(sqrt(5) + 5) + 4*sqrt(25 - 10*sqrt(5)))/10)
    @printf("外円の直径が %g のとき，大円，中円，小円の直径は %g, %g, %g である。\n", 2R, 2r1, 2r2, 2r3)
    r1 = 2*R*(-2 + sqrt(5))
    r2 = R*(-1 + sqrt(5))/4
    r3 = R*(-15 + 7*sqrt(5))/4
    x3 = R*(-5*sqrt(2*sqrt(5) + 10) + 4*sqrt(25 - 10*sqrt(5)) + 2*sqrt(10*sqrt(5) + 50))/10
    y = R*sind(18)
    plot([0, -R*cosd(18), R*cosd(18), 0], [R, y, y, R], color=:green, lw=0.5)
    circle(0, 0, R, :blue)
    circle(0, r1 - R, r1)
    circle(0, R*sind(18) + r2, r2, :magenta)
    circle2(x3, R*sind(18) + r3, r3, :orange)
    segment(0, R, R*sind(36), -R*cosd(36))
    segment(0, R, -R*sind(36), -R*cosd(36))
    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 - R, "大円:r1,(0,r1-R)", :red, :center, delta=-delta)
        point(0, R - r2, "中円:r2,(0,R-r2)", :magenta, :center, delta=-delta)
        point(x3, R*sind(18) + r3, "小円:r3,(x3,R*sind(18)+r3)", :black, :center, delta=-delta)
        point(0, R, "A", :blue, :center, :bottom, delta=delta/2)
        point(0, -R, "B", :blue, :center, delta=-delta)
        point(R*sind(36), -R*cosd(36), " C", :blue, :left)
        point(R*cosd(18), R*sind(18), "  D", :blue, :left, :vcenter)
        point(0, R*sind(18), "E ", :blue, :right, :vcenter)
    end
end;

draw(1/2, true)

算額（その1484）

算額（その1484）

四十四 岩手県一関市真滝 熊野白山滝神社 文久元年(1861)
山村善夫：現存 岩手の算額，昭和52年1月30日，熊谷印刷，盛岡市.
http://www.wasan.jp/yamamura/yamamura.html
キーワード：円7個，半円
#Julia, #SymPy, #算額, #和算

問題文，答，術ともに不明であるが，以下のようなものでもあろう。

半円の中に，全円 1 個，大円 4 個，小円 2 個を容れる。小円の直径が与えられたとき，大円の直径を求める術を述べよ。

注：算額の図（山村が描いたもの）はいつもの通り不適切なものである。

半円ではなく円弧（弓形）が入れ物になっているが，左右に余白がある。このような状態では円弧（弓形）が特定されない。左右端にある大円・小円は円弧に内接していなければいけない。そして，そのような場合には，半円になるのが適切である。

半円の半径と中心座標を 4r1, (0, 0)
全円の半径と中心座標を 2r1, (0, 2r1)
大円の半径と中心座標を r1, (r1, 2r1), ((x1, r1)
小円の半径と中心座標を r2, (x2, y2)
とおき，以下の連立方程式を解く。

include("julia-source.txt");

using SymPy

@syms r1::positive, x1::positive, r2::positive, x2::positive, y2::positive
eq1 = x2^2 + y2^2 - (4r1 - r2)^2
# eq2 = x1^2 + r1^2 - (2r1 + r1)^2
eq2 = x2^2 + (y2 - 2r1)^2 - (2r1 + r2)^2
eq3 = (x1 - x2)^2 + (y2 - r1)^2 - (r1 + r2)^2
eq4 = x1^2 + r1^2 - (4r1 - r1)^2
res = solve([eq1, eq2, eq3, eq4], (r1, x1, x2, y2))[1]

    (sqrt(2)*r2/2 + 5*r2/4, 2*r2 + 5*sqrt(2)*r2/2, 2*r2*sqrt(2*sqrt(2) + 3), 2*r2*(1 + sqrt(2)))

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

    r2*(2*sqrt(2) + 5)/4

大円の半径 r1 は，小円の半径 r2 の (2√2 + 5)/4 倍である。
小円の直径が 373 寸のとき，大円の直径は 730.000有奇 寸である（730.0008293825822 寸）。
きれいな数になる。

function draw(r2, more=false)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    (r1, x1, x2, y2) = (sqrt(2)*r2/2 + 5*r2/4, 2*r2 + 5*sqrt(2)*r2/2, 2*r2*sqrt(2*sqrt(2) + 3), 2*r2*(1 + sqrt(2)))
    plot()
    circle(0, 0, 4r1, :magenta, beginangle=0, endangle=180)
    circle(0, 2r1, 2r1, :green)
    circle2(r1, 2r1, r1)
    circle2(x1, r1, r1)
    circle2(x2, y2, 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, 4r1, "4r1", :magenta, :center, :bottom, delta=delta)
        point(0, 0, "半円:4r1,(0,0)", :magenta, :center, :bottom, delta=delta)
        point(0, 2r1, "全円:2r1\n(0,2r1)", :green, :right, :vcenter, deltax=-2delta)
        point(0, 2r1, "全円:2r1\n(0,2r1)", :green, :right, :vcenter, deltax=-2delta)
        point(r1, 2r1, "大円:r1\n(r1,2r1)", :red, :center, delta=-delta)
        point(x1, r1, "大円:r1\n(x1,r1)", :red, :center, delta=-delta)
        point(x2, y2, "小円:r2\n(x2,r2)", :blue, :center, delta=-delta)
    end
end;

draw(1/2, true)

 

算額（その1483）

算額（その1483）

四十四 岩手県一関市真滝 熊野白山滝神社 文久元年(1861)
山村善夫：現存 岩手の算額，昭和52年1月30日，熊谷印刷，盛岡市.
http://www.wasan.jp/yamamura/yamamura.html
キーワード：円6個，正方形
#Julia, #SymPy, #算額, #和算

問題文，答，術ともに不明であるが，以下のようなものでもあろう。

正方形と大円 2 個が交差しており，正方形に内接し大円に外接する小円 4 個を容れる。正方形の一辺の長さが与えられたとき，小円の直径を求める術を述べよ。

正方形の一辺の長さを 2a
大円の半径と中心座標を r1, (a/2, a/2); r1 = a/√2
小円の半径と中心座標を r2, (a - r2, y2), (-y2, r2 - a)
とおき，以下の連立方程式を解く。

r2 に具体的な数値を与えれば連立方程式は解けるが，数値解になる。
解析解を求めるにはまず eq2 を解いて得られる a を eq1 に代入し，その式を解いて y2 を求める。

include("julia-source.txt");

using SymPy

@syms a::positive, r1::positive, r2::positive, y2::negative
eq1 = (a - r2 - a/2)^2 + (a/2 - y2)^2 - (a/√Sym(2) + r2)^2 |> simplify
eq2 = dist2(0, 0, a, -a, a - r2, y2, r2);
# res = solve([eq1, eq2], (a, y2))

ans_a = solve(eq2, a)[1]
ans_a |> println

    r2 + sqrt(2)*r2 - y2

# a
eq11 = eq1(a => ans_a) |> simplify
eq11 |> println

    -3*r2^2 - 2*sqrt(2)*r2^2 + 2*y2^2

# y2
ans_y2 = solve(eq11, y2)[1] |> sympy.sqrtdenest
ans_y2 |> println

    -r2*(sqrt(2) + 2)/2

ans_a の y2 に ans_y2 を代入する。
正方形の一辺の長さ 2a は，小円の半径 r2 の (4 + 3√2) 倍である。
小円の直径が 1 寸のとき，正方形の一辺の長さは 4.121320343559643 寸である。

# a
ans_a(y2 => ans_y2) |> simplify |> println

    r2*(4 + 3*sqrt(2))/2

function draw(r2, more=false)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    y2 = -r2*(sqrt(2) + 2)/2
    a = r2 + sqrt(2)*r2 - y2
    r1 = a/√2
    println("a = $a; r1 = $r1; r2 = $r2;  y2 = $y2")
    plot([a, a, -a, -a, a], [-a, a, a, -a, -a], color=:green, lw=0.5)
    circle(a/2, a/2, r1, :blue)
    circle(-a/2, -a/2, r1, :blue)
    circle(a - r2, y2, r2)
    circle(-y2, r2 - a, r2)
    circle(r2 - a, -y2, r2)
    circle(y2, a - r2, 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(0, a, "a", :green, :right, :bottom, delta=delta/2)
        point(a, 0, "a ", :green, :right, :bottom, delta=delta/2)
        point(a - r2, y2, "小円:r2\n(a-r2,y2)", :red, :center, :bottom, delta=delta)
        point(-y2, r2 - a, "小円:r2\n(-y2,r2-a)", :red, :center, :bottom, delta=delta)
        point(a/2, a/2, "大円:r1,(a/2,a/2)", :blue, :center, delta=-delta)
    end
end;

draw(1/2, true)

算額（その1482）

算額（その1482）

四十四 岩手県一関市真滝 熊野白山滝神社 文久元年(1861)
山村善夫：現存 岩手の算額，昭和52年1月30日，熊谷印刷，盛岡市.
http://www.wasan.jp/yamamura/yamamura.html
キーワード：円11個，外円
#Julia, #SymPy, #算額, #和算

問題文，答，術ともに不明であるが，以下のようなものでもあろう。

算額（その446），算額（その24）と同じ図形である。

外円の中に白円 2 個，赤円 4 個，青1円 2 個，青2円 2 個を容れる。それぞれの円は外円に内接するほか，互いに外接し合っている。各円の直径を求めよ。

注：算額の図（山村が描いたもの）はいつもの通り不適切なものである。

上の青と赤の位置関係と下のそれとの位置関係が違うのは，実は上のように白円と 青赤の 4 個の円が互いに外接し合うような図は描けない。無理に描こうとすれば「赤円が青・白円に外接するが外円には内接しない」か「赤円は外円に内接し外側の青円と白円には外接するが内側の青円には外接しない」図にしかならない。

山村の図で，上の円のグループは白，赤，青，外円が互いに接するという条件を満たしているが，よく見れば青の円の大きさが違うのがわかる。つまり，山村は同じ青で塗っているがこの 2 つの円は直径が違う別物である。そこで，上の問題では青1円，青2円と書き分けたのである。

外円の半径と中心座標を R, (0, 0)
白円の半径と中心座標を r1, (r1, 0)
赤円の半径と中心座標を r2, (x2, y2)
青1円の半径と中心座標を r3, (0, R - r3)
青2円の半径と中心座標を r4, (0, R - 2r3 - r4)
とおき，以下の連立方程式を解く。

include("julia-source.txt");

using SymPy
@syms R::positive, r1::positive,
      r2::positive, x2::positive, y2::positive,
      r3::positive, r4::positive
r1 = R/2
eq1 = x2^2 + y2^2 - (R - r2)^2
eq2 = r1^2 + (R - 2r3 - r4)^2 - (r1 + r4)^2
eq3 = (r1 - x2)^2 + y2^2 - (r1 + r2)^2
eq4 = x2^2 + (R - r3 - y2)^2 - (r2 + r3)^2
eq5 =  x2^2 + (y2 - R + 2r3 + r4)^2 - (r2 + r4)^2
# res = solve([eq1, eq2, eq3, eq4, eq5], (r1, r2, x2, y2, r3))

using NLsolve

function nls(func, params...; ini = [0.0])
    if typeof(ini) <: Number
        r = nlsolve((vout, vin) -> vout[1] = func(vin[1], params..., [ini]), ftol=big"1e-40")
        v = r.zero[1]
    else
        r = nlsolve((vout, vin)->vout .= func(vin, params...), ini, ftol=big"1e-40")
        v = r.zero
    end
    return Float64.(v), r.f_converged
end;

function H(u)
    (r2, x2, y2, r3, r4) = u
    return [
       x2^2 + y2^2 - (R - r2)^2,  # eq1
       R^2/4 - (R/2 + r4)^2 + (R - 2*r3 - r4)^2,  # eq2
       y2^2 - (R/2 + r2)^2 + (R/2 - x2)^2,  # eq3
       x2^2 - (r2 + r3)^2 + (R - r3 - y2)^2,  # eq4
       x2^2 - (r2 + r4)^2 + (-R + 2*r3 + r4 + y2)^2,  # eq5
    ]
end;
R = 10/2
iniv = BigFloat[2, 3.5, 7, 2, 2]
res = nls(H, ini=iniv)

    ([1.0827118232955024, 1.751864530113493, 3.503729060226986, 0.8017872829546411, 0.9781948282296248], true)

外円の直径が 10 のとき，白円，赤円，青1円，青2円の直径は 5, 2.16542, 1.60357, 1.95639 である。

function draw(R, more=false)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    #(r1, r2, x2, y2) = (r3*(sqrt(17) + 7)/4, r3*(-7*sqrt(2*sqrt(17) + 46) - sqrt(34*sqrt(17) + 782) + 22*sqrt(17) + 122)/112, r3*(-10*sqrt(17) + 26 + 3*sqrt(34*sqrt(17) + 782) + 21*sqrt(2*sqrt(17) + 46))/112, r3*(sqrt(2*sqrt(17)/49 + 46/49) + 19/14 + 5*sqrt(17)/14))
    #(r1, r2, x2, y2) = (r3*(sqrt(17) + 7)/4, r3*(-42 - sqrt(82*sqrt(17) + 350) - 6*sqrt(17) + sqrt(1394*sqrt(17) + 5950))/16, -sqrt(2)*r3*sqrt(41*sqrt(17) + 175)/8 + r3*(19 + 5*sqrt(17))/8, r3*(3/2 + sqrt(17)/2))
    (r2, x2, y2, r3, r4) = res[1]#[6.23606797749979, 2.700745812045397, 4.369898518863388, 8.739797037726776, 2.440035777631468]
    r1 = R/2
    @printf("外円の直径が %g のとき，白円，赤円，青1円，青2円の直径は %g, %g, %g, %g である。\n", 2R, 2r1, 2r2, 2r3, 2r4)
    plot()
    circle(0, 0, R, :green)
    circle2(r1, 0, r1, :gray50)
    circle4(x2, y2, r2, :red)
    circle22(0, R - r3, r3, :blue)
    circle22(0, R - 2r3 - r4, r4, :blue4)
    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, 0, "白円:r1,(r1,0)", :gray50, :center, delta=-delta)
        point(x2, y2, "赤円:r2\n(x2,y2)", :red, :center, delta=-delta)
        point(0, R - r3, "青1円:r3,(0,R-r3)", :blue, :center, delta=-delta)
        point(0, R - 2r3 - r4, "青2円:r4,(0,R-2r3-r4)", :blue4, :center, delta=-delta)
    end
end;

draw(10/2, true)