算額（その1586）

算額（その1586）

岡山県瀬戸内市長船町土師宮森 片山日子神社 明治6年(1873)
深川英俊，トニー・ロスマン：聖なる数学：算額，p. 104，森北出版株式会社，2010年４月22日.

キーワード：直角三角形の区分法
#Julia, #SymPy, #算額, #和算, #数学

AC = 30m，BC = 40m の直角三角形状の畑がある。この畑に，幅が 2 メートルのあぜ道 DEFGHIJ を残りの面積が等しくなるように作りたい。BE, DE, HC, JC, AI, FG の長さを求めよ。

A の座標を (Ax, Ay) などとし，A ~ J までの合計 20 変数と，面積を S とおき，以下の連立方程式を解く。
いくつかの変数値は自明であるが，方程式の本数としては 7 本，変数は Dy, Ex, Fy, Hx, Iy, Jy, S の 7 個なので，過不足ない。

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

using SymPy
@syms Ax, Bx, Cx, Dx, Ex, Fx, Gx, Hx, Ix, Jx,
      Ay, By, Cy, Dy, Ey, Fy, Gy, Hy, Iy, Jy,
      S
Ax = Ix = Jx = Cx = 40
Ay = 30
Dx = Ex
Fx = Gx = Hx
Bx = By = Ey = Hy = Cy = 0
Gy = Iy
eq1 = Hx - Ex - 2
eq2 = Iy - Jy - 2
eq3 = (30*40//2 - 2(Ix - Gx) - ((Dy - Ey) + (Fy - Hy)))/3 - S
eq4 = (Dy - Ey)/(Ex - Bx) - 3//4
eq5 = (Ay - Fy)/(Cx - Hx) - 3//4
eq6 = ((Fy - Gy)+(Ay - Iy))*(Ix - Gx)/2 - S
eq7 = (Ex - Bx)*(Dy - Ey)/2 - S;
res = solve([eq1, eq2, eq3, eq4, eq5, eq6, eq7],
        (Dy, Ex, Fy, Hx, Iy, Jy, S))[2]  # 2 of 2

    (1/6 + sqrt(9406)/6, 2/9 + 2*sqrt(9406)/9, 5/3 + sqrt(9406)/6, 20/9 + 2*sqrt(9406)/9, sqrt(9406)/24 + 107/12, sqrt(9406)/24 + 83/12, sqrt(9406)/27 + 9407/54)

このあと，該当する線分の長さは座標点から計算される。

3 区画の面積 = 177.79572351294942
BE = Ex - Bx = 21.77434107769659
DE = Dy - By = 16.330755808272443
HC = Cx - Hx = 16.22565892230341
JC = Jy - Cy = 10.95768895206811
AI = Ay - Iy = 17.04231104793189
FG = Fy - Gy = 4.873066856204332

深川は以下のように解いている。
小数点以下4桁までで示しているが本文(P.104)と解答編(P.125)で四捨五入の誤りなどがあるが，正確な値は上に示したものと同じである。

@syms t
eq = 3//8*t^2 - (38 - t)*((3t + 110)/16)
BE = t = solve(eq, t)[2].evalf()
DE = (3/4)*t
CH = x = 38 - t
JC = y = (3t +110)/16
AI = 28 - y
FG = (3/4)t + (3/4)*2 - (y + 2);

(BE, 21.7743, 21.7743) |> println
(DE, 16.331, 16.331) |> println
(CH, 16.2255, 16.2257) |> println
(JC, 10.9577, 10.9576) |> println
(AI, 17.0423, 17.0424) |> println
(FG, 4.873, 4.873125) |> println

    (21.7743410776966, 21.7743, 21.7743)
    (16.3307558082724, 16.331, 16.331)
    (16.2256589223034, 16.2255, 16.2257)
    (10.9576889520681, 10.9577, 10.9576)
    (17.0423110479319, 17.0423, 17.0424)
    (4.87306685620433, 4.873, 4.873125)

function draw(more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    (Dy, Ex, Fy, Hx, Iy, Jy, S) = (1/6 + sqrt(9406)/6, 2/9 + 2*sqrt(9406)/9, 5/3 + sqrt(9406)/6, 20/9 + 2*sqrt(9406)/9, sqrt(9406)/24 + 107/12, sqrt(9406)/24 + 83/12, sqrt(9406)/27 + 9407/54)
    Ax = Ix = Jx = Cx = 40
    Ay = 30
    Dx = Ex
    Fx = Gx = Hx
    Bx = By = Ey = Hy = Cy = 0
    Gy = Iy
    println("3 区画の面積 = $S")
    println("BE = Ex - Bx = $(Ex - Bx)")
    println("DE = Dy - By = $(Dy - By)")
    println("HC = Cx - Hx = $(Cx - Hx)")
    println("JC = Jy - Cy = $(Jy - Cy)")
    println("AI = Ay - Iy = $(Ay - Iy)")
    println("FG = Fy - Gy = $(Fy - Gy)")
    plot([Bx, Cx, Ax, Bx], [By, Cy, Ay, By], seriestype=:shape, color=:goldenrod, fillcolor=:goldenrod, lw=0.5)
    plot!([Bx, Ex, Dx, Bx], [By, Ey, Dy, By], seriestype=:shape, color=:green, fillcolor=:green, lw=0.5)
    plot!([Hx, Cx, Jx, Hx, Hx], [Hy, Cy, Jy, Jy, Hy], seriestype=:shape, color=:green, fillcolor=:green, lw=0.5)
    plot!([Gx, Ix, Ax, Fx, Gx], [Gy, Iy, Ay, Fy, Gy], seriestype=:shape, color=:green, fillcolor=:green, lw=0.5)
    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(Ax, Ay, " A")
        point(Bx, By, "B", :green, :left, delta=-delta)
        point(Cx, Cy, " C", :green, :left, delta=-delta)
        point(Dx, Dy, " D", :green, :right, :bottom, delta=delta)
        point(Ex, Ey, "E", :green, :left, delta=-delta)
        point(Fx, Fy, "F", :green, :right, :bottom, delta=delta)
        point(Gx, Gy, "G", :green, :left, delta=-delta/2)
        point(Hx, Hy, "H", :green, :left, delta=-delta)
        point(Ix, Iy, " I")
        point(Jx, Jy, " J")
        ylims!(-5delta, Ay + 5delta)
    end
end;

draw(true)

算額（その1585）

算額（その1585）

愛知県名古屋市 熱田神宮 天保13年(1842)
深川英俊，トニー・ロスマン：聖なる数学：算額，p. 103，森北出版株式会社，2010年４月22日.

キーワード：円4個，接線4本
#Julia, #SymPy, #算額, #和算, #数学

左円と右円が外接しており，左円内には左小円が内接し，右円内には右小円が内接している。左円と右小円，右円と左小円の共通接線を 2 本ずつ引く。左小円と右小円の直径は，左円と右円の直径でどのように表されるか。

左円の半径と中心座標を r1, (r1, 0)
右円の半径と中心座標を r2, (2r1 + r2, 0)
左小円の半径と中心座標を r31, (r1 + r31, 0)
右小円の半径と中心座標を r32, (2r1 + r32, 0)
とおき，以下の連立方程式を解き，r31, r32 を求める。

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

using SymPy
@syms r1, r2, r31, r32
eq1 = r31/(2r1 - r31) - r2/(2r1 + r2)
eq2 = r1/(2r2 + r1) - r32/(2r2 - r32)
res = solve([eq1, eq2], (r31, r32));
res[r31] |> println
res[r32] |> println

    r1*r2/(r1 + r2)
    r1*r2/(r1 + r2)

左円と右円がどんな値をとっても，左小円の半径と，右小円の半径は等しい値 r1*r2/(r1 + r2) になる。

たとえば，
左円，右円の直径が 18, 14 のとき，左小円，右小円の直径は 7.875, 7.875 である。
左円，右円の直径が 18, 16 のとき，左小円，右小円の直径は 8.47059, 8.47059 である。
左円，右円の直径が 18, 6 のとき，左小円，右小円の直径は 4.5, 4.5 である。
左円，右円の直径が 8, 16 のとき，左小円，右小円の直径は 5.33333, 5.33333 である。
左円，右円の直径が 10, 18 のとき，左小円，右小円の直径は 6.42857, 6.42857 である。

function draw(r1, r2, more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    r31 = r32 = r1*r2/(r1 + r2)
    @printf("左円，右円の直径が %g, %g のとき，左小円，右小円の直径は %g, %g である。\n", 2r1, 2r2, 2r31, 2r32)
    p = plot()
    circle(r1, 0, r1)
    circle(2r1 + r2, 0, r2, :magenta)
    circle(2r1 - r31, 0, r31, :blue)
    circle(2r1 + r32, 0, r32, :green)
    θ1 = asind(r2/(2r1 + r2))
    abline(0, 0, tand(θ1), 0, 2r1 + 2r2)
    abline(0, 0, -tand(θ1), 0, 2r1 + 2r2)
    θ2 = asind(-r1/(r1 + 2r2))
    abline(2r1 + 2r2, 0, tand(θ2), 0, 2r1 + 2r2)
    abline(2r1 + 2r2, 0, -tand(θ2), 0, 2r1 + 2r2)
    delta = (fontheight = (ylims()[2]- ylims()[1]) / 500 * 10 * 2) /3  # size[2] * fontsize * 2
    str = @sprintf("r1=%d, r2=%d", r1, r2)
    #point(r1 + r2, max(r1, r2), str, :black, :center, :bottom, delta=delta, mark=false)
    point(2r1, max(r1, r2), str, :black, :center, :bottom, delta=delta, mark=false)
    if more
        hline!([0], color=:gray80, lw=0.5)
        vline!([0], color=:gray80, lw=0.5)
        point(r1, 0, "左円:r1,(r1,0)", :red, :center, delta=-delta)
        point(2r1 + r2, 0, " 右円:r2,(2r1+r2,0)", :magenta, :center, delta=-delta)
        point(2r1 - r31, 0, " 左小円", :blue, :center, :bottom, delta=delta)
        point(2r1 + r32, 0, " 右小円", :green, :center, :bottom, delta=delta)
    end
    return p
end;

draw(9, 7, true)

算額（その1584）

算額（その1584）

福島県田村市 安倍文殊菩薩堂 明治10年(1877)
深川英俊，トニー・ロスマン：聖なる数学：算額，p. 103，森北出版株式会社，2010年４月22日.

キーワード：円6個，接線2本
#Julia, #SymPy, #算額, #和算, #数学

大円 2 個が外接しており，上の大円内に中円 2 個が内接している。2 個の中円の共通接線は上の大円の中心を通る。中円と共通接線に接する小円を容れる。また，上の大円の中心を通り，下の大円が内接する外円を描く。小円の直径を 1 としたとき，大円，中円，外円の直径はいかほどか。

外円の半径と中心座標を R, (0, r3 - R); R = 3r3/2
大円の半径と中心座標を r3, (0, r3), 0, -r3)
中円の半径と中心座標を r2, (0, r2), (0, 2r3 - r2)
小円の半径と中心座標を r1, (0, 2r3, - 2r2, -r1)
とおき，以下の連立方程式を解く。

eq2 の中に，"r3/2r3" があるのでわかるように，二本の共通接線のなす角は 60° である。

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

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

    (3*r1, 9*r1)

中円の直径は小円の直径の 3 倍，大円の直径は小円の直径の 9 倍である。
外円の直径は大円の直径の 3/2 倍なので，小円の 27/2 倍である。
小円の直径が 1 寸のとき，外円，大円，中円の直径はそれぞれ 27/2 寸，9 寸，3 寸 である。

function draw(r1, more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    (r2, r3) = (3*r1, 9*r1)
    R = 3r3/2
    @printf("小円の直径が %g のとき，外円，大円，中円の直径は %g, %g, %g である。\n", 2r1, 2R, 2r3, 2r2)
    plot()
    circle22(0, r3, r3)
    circle(0, 2r3 - r2, r2, :blue)
    circle(0, r2, r2, :blue)
    circle(0, 2r3 - 2r2 - r1, r1, :green)
    circle(0, r3 - R, R, :magenta)
    segment(R*cosd(30), r3-R-R*sind(30), -r3*sind(30), r3+r3*cosd(30))
    segment(-R*cosd(30), r3-R-R*sind(30), r3*sind(30), r3+r3*cosd(30))
    if more
        delta = (fontheight = (ylims()[2]- ylims()[1]) / 500 * 10 * 2) /3  # size[2] * fontsize * 2
        hline!([0], color=:gray80, lw=0.3)
        vline!([0], color=:gray80, lw=0.3)
        point(0, -r3, "大円:r3,(0,-r3)", :red, :center, delta=-delta)
        point(0, r3, " 大円:r3,(0,r3)", :red, :left, :vcenter)
        point(0, 2r3 - r2, " 中円:r2,(0,2r3-r2)", :blue, :left, :vcenter)
        point(0, r2, " 中円:r2,(0,r2)", :blue, :left, :vcenter)
        point(0, 2r3 - 2r2 - r1, "   小円:r1,(0,2r3-2r2-r1)", :green, :left, :vcenter)
        point(0, r3 - R, "外円,R,(0,r3-R)", :magenta, :center, delta=-delta)
        point(R*cosd(30), r3-R-R*sind(30))
        point(-r3*sind(30), r3+r3*cosd(30))
        point(-R*cosd(30), r3-R-R*sind(30))
        point(r3*sind(30), r3+r3*cosd(30))
    end
end;

draw(1/2, true)

算額（その21）改訂版

算額（その21）改訂版

本来の題意に沿った解を求めるように改訂した。

福島県田村郡三春町御木沢 厳島神社 明治18年(1885)
http://www.wasan.jp/fukusima/miharuitukusima3.html

キーワード：円5個，半円2個，長方形
#Julia, #SymPy, #算額, #和算, #数学

長方形の中に，半円 2 個，甲円 1 個，乙円 2 個，丙円 2 個を容れる。丙円の直径が与えられたとき，黒積を求める術を述べよ。

半円の半径と中心座標を R, (0, R/2)
甲円の半径と中心座標を r1, (0, 0)
丙円の半径と中心座標を r3, (R - r3, 0)
乙円の半径と中心座標を r2, (r2, 0)
とおく。

黒積は，「半径 R の円の面積の 1/8 の扇形の面積 ACB から，三角形の面積 OBC と乙円の面積の 1/2 を引いたものの 4 倍である」。

以下の連立方程式を解く。

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

using SymPy
@syms R, r1, r2, r3, 黒積

eq1 = r2^2 + r1^2 - (R - r2)^2
eq2 = (R - r3)^2 + r1^2 - (R + r3)^2
eq3 = 2r1^2 - R^2
eq4 = 4(PI*R^2/8 - r1^2/2 - PI*r2^2/2) - 黒積
res = solve([eq1, eq2, eq3, eq4], (r1, r2, R, 黒積))[3] # 3 of 3

    (4*sqrt(2)*r3, 2*r3, 8*r3, r3^2*(-64 + 24*pi))

# r1
res[1] |> println

    4*sqrt(2)*r3

# r2
res[2] |> println

    2*r3

# r3
res[3] |> println

    8*r3

# 黒積
res[4] |> println

    r3^2*(-64 + 24*pi)

黒積は，円周率の 24 倍から 64 を引き，丙円の半径の二乗を掛けることで得られる。

乙円の直径が 1 のとき，黒積は 2.84955592153876 である。

res[4](r3 => 1/2).evalf() |> println

    2.84955592153876

術は以下のようで，一致する。

円周率 = PI
丙円径 = 1
((6円周率 - 16)*丙円径^2).evalf() |> println

    2.84955592153876

function draw(r3, more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    (r1, r2, R, 黒積) = (4*sqrt(2)*r3, 2*r3, 8*r3, r3^2*(-64 + 24*pi))
    @printf("半円の半径が %g のとき，甲円，乙円，丙円の半径は %g, %g, %g，黒積は %g である。\n", R, r1, r2, r3, 黒積)
    plot([R, R, -R, -R, R], [-r1, r1, r1, -r1, -r1], seriestype=:shape, color=:powderblue, lw=0.5, fillcolor=:powderblue)
    circlef(0, r1, R, :tomato, beginangle=180, endangle=360)
    circlef(0, -r1, R, :tomato, beginangle=0, endangle=180)
    circlef(0, 0, r1, :moccasin)
    x = [R*cosd(z) for z in 45:0.01:135]
    y = [R*sind(z) - r1 for z in 45:0.01:135]
    append!(x, reverse(x))
    append!(y, -reverse(y))
    plot!(x, y, seriestype=:shape, color=:gray60, fcolor=:gray60, lw=0.5, )
    circle2f(R - r3, 0, r3, :white)
    circle2f(r2, 0, r2, :moccasin)
    circle(0, r1, R, :black, beginangle=180, endangle=360)
    circle(0, -r1, R, :black, beginangle=0, endangle=180)
    circle(0, 0, r1, :black)
    circle2(R - r3, 0, r3, :black)
    circle2(r2, 0, r2, :black)    
    rect(-R, -r1, R, r1, :black)
    plot!([0, 0, r1, 0], [R - r1, -r1, 0, 0], color=:black, lw=0.5)
    if more
        delta = (fontheight = (ylims()[2]- ylims()[1]) / 500 * 10 * 2) /3  # size[2] * fontsize * 2
        hline!([0], color=:gray80, lw=0.3)
        vline!([0], color=:gray80, lw=0.3)
        point(0, r1, "半円:R,(0,R/2)", :black, :center, :bottom, delta=delta)
        point(R - r3, 0, "丙円:r2,(R-r3,0)", :black, :right, :bottom, delta=2delta, deltax=4delta)
        point(r2, 0, "乙円:r2,(r2,0)", :black, :center, delta=-delta)
        point(0, R - r1, "A", :black, :center, :bottom, delta=delta)
        point(r1, 0, "B ", :black, :right, :vcenter)
        point(0, 0, "O ", :black, :right, :vcenter)
        point(0, -r1, "C ", :black, :right, :bottom, delta=delta)
        ylims!(-r1 - 2delta, r1 + 2delta)
    end
end;

draw(1/2, true)

算額（その1583）

算額（その1583）

福島県田村郡三春町御木沢 厳島神社
三春まちなか寺子屋：厳島神社新算額
https://miharukoma.com/wp-content/uploads/2018/11/12.22厳島神社-新算額.pdf
https://miharukoma.com/wp-content/uploads/2018/12/寺子屋12.22　問題解説.pdf

キーワード：円5個，半円2個，長方形
#Julia, #SymPy, #算額, #和算, #数学

長方形の中に，半円 2 個，甲円 1 個，乙円 2 個，丙円 2 個を容れる。半円の半径が 48 寸のとき，乙円，丙円の半径はいかほどか。

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

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

using SymPy
@syms R, r1, r2, r3

eq1 = r2^2 + r1^2 - (R - r2)^2
eq2 = (R - r3)^2 + r1^2 - (R + r3)^2
eq3 = 2r1^2 - R^2
res = solve([eq1, eq2, eq3], (r1, r2, r3))[2]  # 2 of 2

    (sqrt(2)*R/2, R/4, R/8)

# r1
res[1] |> println

    sqrt(2)*R/2

# r2
res[2] |> println

    R/4

# r3
res[3] |> println

    R/8

甲円，乙円，丙円の半径は，半円の直径の 1/√2, 1/4, 1/8 倍である。

半円の半径が 48 寸のとき，甲円，乙円，丙円の半径は 33.94112549695428 寸， 12 寸，6 寸である。

using Colors
function draw(R, more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    (r1, r2, r3) = R .* (1/√2, 1/4, 1/8)
    @printf("半円の半径が %g のとき，甲円，乙円，丙円の半径は %g, %g, %g である。\n", R, r1, r2, r3)
    plot()
    circlef(0, r1, R, parse(Colorant, "#FFF565"), beginangle=180, endangle=360)
    circlef(0, -r1, R, parse(Colorant, "#FFF565"), beginangle=0, endangle=180)
    circlef(0, 0, r1, parse(Colorant, "#FF8B3E"))
    x = [R*cosd(z) for z in 45:0.01:135]
    y = [R*sind(z) - r1 for z in 45:0.01:135]
    append!(x, reverse(x))
    append!(y, -reverse(y))
    plot!(x, y, seriestype=:shape, color=parse(Colorant, "#FFFCC4"), fcolor=parse(Colorant, "#FFFCC4"), lw=0.5)
    circle2f(R - r3, 0, r3, parse(Colorant, "#F02F22"))
    circle2f(r2, 0, r2, parse(Colorant, "#757CF7"))
    circle(0, r1, R, :black, beginangle=180, endangle=360)
    circle(0, -r1, R, :black, beginangle=0, endangle=180)
    circle(0, 0, r1, :black)
    circle2(R - r3, 0, r3, :black)
    circle2(r2, 0, r2, :black)    
    rect(-R, -r1, R, r1, :black)
    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,(0,R/2)", :black, :center, :bottom, delta=delta)
        point(R - r3, 0, "丙円:r2,(R-r3,0)", :black, :right, delta=-2delta, deltax=4delta)
        point(r2, 0, "乙円:r2,(r2,0)", :black, :center, :bottom, delta=delta)
        ylims!(-r1 - 2delta, r1 + 2delta)
    end
end;

draw(48, true)

算額（その1582）

算額（その1582）

福島県田村郡三春町御木沢 厳島神社
三春まちなか寺子屋：厳島神社新算額
https://miharukoma.com/wp-content/uploads/2018/11/12.22厳島神社-新算額.pdf
https://miharukoma.com/wp-content/uploads/2018/12/寺子屋12.22　問題解説.pdf

キーワード：円4個，半円2個，長方形
#Julia, #SymPy, #算額, #和算, #数学

長方形の中に，半円 2 個，乙円 2 個，丙円 2 個を容れる。半円の半径が 48 寸のとき，乙円，丙円の半径はいかほどか。

半円の半径と中心座標を R, (0, R/2)
丙円の半径と中心座標を r3, (R - r3, 0)
乙円の半径と中心座標を r2, (r2, 0)
とおき，以下の連立方程式を解く。

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

using SymPy
@syms R, r2, r3
eq1 = r2^2 + (R/2)^2 - (R - r2)^2
eq2 = (R - r3)^2 + (R/2)^2 - (R + r3)^2;

res = solve([eq1, eq2], (r2, r3));

# r2
res[r2] |> println

    3*R/8

# r3
res[r3] |> println

    R/16

乙円，丙円の半径は，半円の直径の 3/8, 1/16 倍である。

半円の半径が 48 寸のとき，乙円，丙円の半径は 18 寸，3 寸である。

using Colors
function draw(R, more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    (r2, r3) = (3R/8, R/16)
    @printf("半円の半径が %g のとき，乙円，丙円の半径は %g, %g である。\n", R, r2, r3)
    plot()
    circlef(0, R/2, R, parse(Colorant, "#FFF565"), beginangle=180, endangle=360)
    circlef(0, -R/2, R, parse(Colorant, "#FFF565"), beginangle=0, endangle=180)
    circle2f(R - r3, 0, r3, parse(Colorant, "#F02F22"))
    circle2f(r2, 0, r2, parse(Colorant, "#757CF7"))    
    circle(0, R/2, R, :black, beginangle=180, endangle=360)
    circle(0, -R/2, R, :black, beginangle=0, endangle=180)
    circle2(R - r3, 0, r3, :black)
    circle2(r2, 0, r2, :black)    
    rect(-R, -R/2, R, R/2, :black)
    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/2, "半円:R,(0,R/2)", :black, :center, :bottom, delta=delta/2)
        point(R - r3, 0, "丙円:r2,(R-r3,0)", :black, :right, delta=-2delta, deltax=-2delta)
        point(r2, 0, "乙円:r2,(r2,0)", :black, :center, :bottom, delta=delta/2)
        ylims!(-R/2 - 2delta, R/2 + 2delta)
    end
end;

draw(48, true)

 

 

算額（その1581）

算額（その1581）

福島県田村郡三春町御木沢 厳島神社
三春まちなか寺子屋：厳島神社新算額
https://miharukoma.com/wp-content/uploads/2018/11/12.22厳島神社-新算額.pdf
https://miharukoma.com/wp-content/uploads/2018/12/寺子屋12.22　問題解説.pdf

キーワード：円4個，半円2個，長方形
#Julia, #SymPy, #算額, #和算, #数学

長方形の中に，半円 2 個，正三角形 4 個，乙円 2 個，丙円 2 個を容れる。半円の半径が 48 寸のとき，乙円，丙円の半径はいかほどか。

半円の半径と中心座標を R, (0, √3R/2)
丙円の半径と中心座標を r3, (R - r3, 0)
乙円の半径と中心座標を r2, (r2, 0)
とおき，以下の連立方程式を解く。

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

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

# r2
ans_r2 = solve(eq2, r2)[1]
ans_r2 |> println

    R/8

# r3
ans_r3 = solve(eq1, r3)[1]
ans_r3 |> println

    3*R/16

乙円，丙円の半径は，半円の直径の 1/8, 3/16 倍である。

半円の半径が 48 寸のとき，乙円，丙円の半径は 6 寸，9 寸である。

using Colors
function draw(R, more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    (r2, r3) = (R/8, 3R/16)
    @printf("半円の半径が %g のとき，乙円，丙円の半径は %g, %g である。\n", R, r2, r3)
    plot()
    circlef(0, √3R/2, R, parse(Colorant, "#FF8B3E"), beginangle=180, endangle=360)
    circlef(0, -√3R/2, R, parse(Colorant, "#FF8B3E"), beginangle=0, endangle=180)
    polygon(R/2, R/√3, R/√3, 3, φ=30, color=parse(Colorant, "#FFF565"), fcolor=parse(Colorant, "#FFF565"))
    polygon(-R/2, R/√3, R/√3, 3, φ=30, color=parse(Colorant, "#FFF565"), fcolor=parse(Colorant, "#FFF565"))
    polygon(R/2, -R/√3, R/√3, 3, φ=210, color=parse(Colorant, "#FFF565"), fcolor=parse(Colorant, "#FFF565"))
    polygon(-R/2, -R/√3, R/√3, 3, φ=210, color=parse(Colorant, "#FFF565"), fcolor=parse(Colorant, "#FFF565"))
    x = [R*cosd(z) for z in 60:0.01:120]
    y = [R*sind(z) - √3R/2 for z in 60:0.01:120]
    append!(x, reverse(x))
    append!(y, -reverse(y))
    circle2f(R - r3, 0, r3, parse(Colorant, "#F02F22"))
    plot!(x, y, seriestype=:shape, color=parse(Colorant, "#FFFCC4"), fcolor=parse(Colorant, "#FFFCC4"), lw=0.5)
    circle2f(r2, 0, r2, parse(Colorant, "#757CF7"))    
    circle(0, √3R/2, R, :black, beginangle=180, endangle=360)
    circle(0, -√3R/2, R, :black, beginangle=0, endangle=180)
    circle2(R - r3, 0, r3, :black)
    circle2(r2, 0, r2, :black)    
    rect(-R, -√3R/2, R, √3R/2, :black)
    plot!([-R, 0, R], √3R/2 .* [-1, 1, -1], color=:black, lw=0.5)
    plot!([-R, 0, R], √3R/2 .* [1, -1, 1], color=:black, lw=0.5)
    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, √3R/2, "半円:R,(0,√3R/2)", :black, :center, :bottom, delta=delta/2)
        point(R - r3, 0, "丙円:r2\n(R-r3,0)", :black, :center, delta=-delta/2)
        point(r2, 0, "乙円:r2,(r2,0)", :black, :left, delta=-5delta)
    end
end;

draw(48, true)

算額（その1580）

算額（その1580）

福島県田村郡三春町御木沢 厳島神社
三春まちなか寺子屋：厳島神社新算額
https://miharukoma.com/wp-content/uploads/2018/11/12.22厳島神社-新算額.pdf
https://miharukoma.com/wp-content/uploads/2018/12/寺子屋12.22　問題解説.pdf

キーワード：円4個，半円2個，正方形
#Julia, #SymPy, #算額, #和算, #数学

正方形の中に半円 2 個，乙円 2 個，丙円 2 個を容れる。半円の半径が 48 寸のとき，乙円，丙円の半径はいかほどか。

半円の半径と中心座標を R, (0, R)
丙円の半径と中心座標を r3, (R - r3, 0)
乙円の半径と中心座標を r2, (R - 2r3 - r2, 0)
とおき，以下の連立方程式を解く。

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

using SymPy
@syms R::positive, r2::positive, r3::positive
eq1 = (R - r3)^2 + R^2 - (R + r3)^2
eq2 = (R - 2r3 - r2)^2 + R^2 - (R + r2)^2
res = solve([eq1, eq2], (r2, r3))[1]

    (R/12, R/4)

乙円，丙円の半径は，半円の半径の 1/12, 1/4 倍である。

半円の半径が 48 寸のとき，乙円，丙円の半径は 4 寸，12 寸である。

using Colors
function draw(R, more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    (r2, r3) = (R/12, R/4)
    @printf("半円の半径が %g のとき，乙円，丙円の半径は %g, %g である。\n", R, r2, r3)
    plot()
    circlef(0, R, R, parse(Colorant, "#FFF565"), beginangle=180, endangle=360)
    circlef(0, -R, R, parse(Colorant, "#FFF565"), beginangle=0, endangle=180)
    circle2f(R - r3, 0, r3, parse(Colorant, "#F02F22"))
    circle2f(R - 2r3 - r2, 0, r2, parse(Colorant, "#757CF7"))    
    circle(0, R, R, :black, beginangle=180, endangle=360)
    circle(0, -R, R, :black, beginangle=0, endangle=180)
    circle2(R - r3, 0, r3, :black)
    circle2(R - 2r3 - r2, 0, r2, :black)    
    rect(-R, -R, R, R, :black)
    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,(0,R)", :black, :center, :bottom, delta=delta/2)
        point(R - r3, 0, "丙円:r2\n(R-r3,0)", :black, :center, delta=-delta/2)
        point(R - 2r3 - r2, 0, "乙円:r2,(R-2r3-r2,0)", :black, :right, delta=-3delta)
    end
end;

draw(48, true)

算額（その1579）

算額（その1579）

福島県田村市 日渡神社 明治21年(1888)
https://www.city.tamura.lg.jp/soshiki/30/bunkazai.html
キーワード：円9個，外円
#Julia, #SymPy, #算額, #和算, #数学

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

外円の半径と中心座標を R, (0, 0)
大円の半径と中心座標を r1, (r1, 0)
小円の半径と中心座標を r2, (0, y21), (r2, y22)
とおき，以下の連立方程式を解く。

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

using SymPy
@syms R::positive, r1::positive,
      r2::positive, y21::positive, y22::positive
r1 = R/2
eq1 = r1^2 + y21^2 - (r1 + r2)^2
eq2 = r2^2 + (y22 - y21)^2 - 4r2^2
eq3 = r2^2 + y22^2 - (R - r2)^2
res = solve([eq1, eq2, eq3], (r2, y21, y22))[1]

    (R*(-5 + sqrt(33))/4, R*(-5*sqrt(2) + sqrt(66))/(2*sqrt(7 - sqrt(33))), R*sqrt(7/2 - sqrt(33)/2))

小円の半径 r2 は，大円の半径 R の (√33 - 5)/4 倍である。
大円の直径が 3954 のとき，736.0001761028412 である。

3954*(√33 - 5)/4

    736.0001761028412

using Colors
function draw(R, more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    r1 = R/2
    (r2, y21, y22) = (R*(-5 + sqrt(33))/4, R*(-5*sqrt(2) + sqrt(66))/(2*sqrt(7 - sqrt(33))), R*sqrt(7/2 - sqrt(33)/2))
    @printf("外円の直径が %g のとき，小円の直径は %g である。\n", 2R, 2r2)
    plot()
    circlef(0, 0, R, RGB(194/255, 214/255, 236/255))
    circle2f(r1, 0, r1, RGB(234/255, 51/255, 35/255))
    circle22f(0, y21, r2, RGB(229/255, 74/255, 233/255))
    circle4f(r2, y22, r2,  RGB(229/255, 74/255, 233/255))
    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/2)
        point(r1, 0, "大円:r1,(r1,0)", :white, :center, delta=-delta/2)
        point(0, y21, "小円:r2,(0,y21)", :black, :left, :bottom, delta=delta)
        point(r2, y22, " 小円:r2,(r2,y22)", :black, :left, :bottom, delta=delta)
    end
end;

draw(1977, true)

 

算額（その1578）

算額（その1578）

『算法續淺問答 巻之七』
terasima3.pdf
http://www.wasan.jp/terasima/terasima.html
http://www.wasan.jp/terasima/terasima3.pdf

キーワード：円環，外円
#Julia, #SymPy, #算額, #和算

図のように 8 個の等円が互いに外接し，外円に内接している。外円の直径が 10 寸のとき，等円の直径を求めよ。

外円の半径と中心座標を R, (0, 0)
等円の半径と中心座標を r, (R - r, 0)
等円の個数を n
とおき，以下の方程式を解く。

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

using SymPy
@syms R::positive, r::positive, n::positive
eq = (R - r)*sind(180/n) - r
res = solve(eq, r)[1]
res |> println

    R*sin(pi/n)/(sin(pi/n) + 1)

n = 8, R = 10/2 を代入し，得られた半径を 2 倍すれば直径は 2.76768653914155 である。

2*res(n => 8, R => Sym(10)/2) |> println
2*res(n => 8, R => Sym(10)/2).evalf() |> println

    10*sqrt(1/2 - sqrt(2)/4)/(sqrt(1/2 - sqrt(2)/4) + 1)
    2.76768653914155

function draw(R, n, more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    r = R*sin(pi/n)/(sin(pi/n) + 1)
    p = plot()
    circle(0, 0, R, :green)
    circle(0, 0, R - r, :gray80)
    rotate(0, R - r, r, angle=360/n)
    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/2)
        point(0, R - r, "R-r", :gray, :center, :bottom, delta=delta/2)
    end
    return p
end;

draw(10, 8, true)

算額（その1577）

算額（その1577）

愛知県名古屋市中央区 大須観音（北野山真福寺寶生院） 天保8年(1837)
藤安順：所掲大須観音堂解義
山田潤：動的幾何ソフトウエアを利用した平面図形の作図についての一考察 - 和算所にある平面図形問題の利用 -
2023018-Yamada.pdf
https://www.kurims.kyoto-u.ac.jp/~kyodo/kokyuroku/contents/pdf/2273-19.pdf

キーワード：円10個，外円
#Julia, #SymPy, #算額, #和算

外円の中に大円 4 個，甲円 3 個，乙円を 2 個容れる。甲円の直径が 1 寸のとき，乙円および外円の直径はいかほどか。

算額（その1575）は，深川が意図的に大円が交差する部分にある円を省略し，甲円と乙円の関係を述べさせるものであったようだ。

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

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

using SymPy
@syms R::positive, r1::positive, x1::positive,
      r2::positive, r3::positive
x1 = r1 - r2
eq1 = R - (2r1 + r2)
eq2 = x1^2 + r1^2 - (r1 + r2)^2  # 大円と甲円が接する
eq3 = x1^2 + (R - r3 - r1)^2 - (r1 + r3)^2;  # 大円と乙円が接する

res = solve([eq1, eq2, eq3], (R, r1, r3))[1]

    (9*r2, 4*r2, r2)

外円の直径は甲円の直径の 9 倍，大円の直径は甲円の直径の 4 倍，乙円の直径は甲円の直径に等しい。

考察：大円の交差するところに甲円が入っているのが，算額（その1575）での y1 の制約条件になっている。

function draw(r2, more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    (R, r1, r3) = (9*r2, 4*r2, r2)
    plot()
    circle(0, 0, R, :green)
    circle4(r1 - r2, r1, r1, :blue)
    circle(0, 0, r2)
    circle22(0, r1, r3)
    circle22(0, R - r3, r3, :magenta)
    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/2)
        point(r1 - r2, r1, "大円:r1,(r1-r2,r1)", :blue, :center, delta=-delta, deltax=5delta)
        point(0, R - r3, "乙円:r3,(0,R-r3)", :magenta, :left, :bottom, delta=delta/2, deltax=4delta)
        point(0, r1, "甲円:r2,(0,r1)", :red, :right, :vcenter, deltax=-4delta)
    end
end;

draw(1/2, true)

算額（その1576）

算額（その1576）

佐久間纉：算法起源集,1877.
https://kokusho.nijl.ac.jp/biblio/100234582/69?ln=ja

山田潤：動的幾何ソフトウエアを利用した平面図形の作図についての一考察 - 和算所にある平面図形問題の利用 -
2023018-Yamada.pdf
https://www.kurims.kyoto-u.ac.jp/~kyodo/kokyuroku/contents/pdf/2273-19.pdf

キーワード：円7個，外円
#Julia, #SymPy, #算額, #和算

外円の中に甲円，乙円，丙円を 2 個ずつ容れる。甲円，乙円，丙円の直径が与えられたとき，外円の直径はいかほどか。

外円の半径と中心座標を R, (0, 0)
甲円の半径と中心座標を r1, (0, R - r1)
乙円の半径と中心座標を r2, (x2, y2)
丙円の半径と中心座標を r3, (x3, y3)
とおき，以下の連立方程式を解く。

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

using SymPy
@syms R::positive, r1::positive,
      r2::positive, x2::positive, y2::negative,
      r3::positive, x3::positive, y3::positive
eq1 = x2^2 + y2^2 - (R - r2)^2 |> expand
eq2 = x3^2 + y3^2 - (R - r3)^2 |> expand
eq3 = x2^2 + (y2 - r1 + R)^2 - (r1 + r2)^2 |> expand
eq4 = x3^2 + (R - r1 - y3)^2 - (r1 + r3)^2 |> expand
eq5 = (x3 - x2)^2 + (y3 - y2)^2 - (r2 + r3)^2 |> expand;

まず，eq1, eq3 の連立方程式から x2, y2 を求める。

res2 = solve([eq1, eq3], (x2, y2))[1]

    (-2*sqrt(R)*sqrt(r1)*sqrt(r2)*sqrt(R - r1 - r2)/(-R + r1), -(-R^2 + R*r1 + R*r2 + r1*r2)/(-R + r1))

次に，eq2, eq4 の連立方程式から x3, y3 を求める。

res3 = solve([eq2, eq4], (x3, y3))[1]

    (-2*sqrt(R)*sqrt(r1)*sqrt(r3)*sqrt(R - r1 - r3)/(-R + r1), (-R^2 + R*r1 + R*r3 + r1*r3)/(-R + r1))

eq5 の x2, y2, x3, y3 に求められた解を代入し，整理する。

eq15 = eq5(x2 => res2[1], y2 => res2[2], x3 => res3[1], y3 => res3[2])/R |> simplify |> numerator;

しかし，このままでは SymPy では解けない。

eq15 |> println

    4*R^3 - 8*R^2*r1 - 4*R^2*r2 - 4*R^2*r3 + 4*R*r1^2 + 4*R*r1*r2 + 4*R*r1*r3 - 8*r1*sqrt(r2)*sqrt(r3)*sqrt(R - r1 - r2)*sqrt(R - r1 - r3) + 8*r1*r2*r3

8*r1*sqrt(r2)*sqrt(r3)*sqrt(R - r1 - r2)*sqrt(R - r1 - r3) の項とそれ以外の項をそれぞれ二乗し，方程式を再構成する。

eq25 = (4*R^3 - 8*R^2*r1 - 4*R^2*r2 - 4*R^2*r3 + 4*R*r1^2 + 4*R*r1*r2 + 4*R*r1*r3 + 8*r1*r2*r3)^2 - ( - 8*r1*sqrt(r2)*sqrt(r3)*sqrt(R - r1 - r2)*sqrt(R - r1 - r3))^2;

ここまで来ると SymPy でも解ける。

res = solve(eq25)[2];  # 2 of 5
res[R] |> println

    r1 + r2 + r3

外円の半径は，甲円，乙円，丙円の半径の合計である。

なお，x2 = x3 であることがわかる。

# x2
res2[1](R => r1 + r2 + r3) |> println
# x3
res3[1](R => r1 + r2 + r3) |> println

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

function draw(r1, r2, r3, more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    R = r1 + r2 + r3
    (x2, y2) = (-2*sqrt(R)*sqrt(r1)*sqrt(r2)*sqrt(R - r1 - r2)/(-R + r1), -(-R^2 + R*r1 + R*r2 + r1*r2)/(-R + r1))
    (x3, y3) = (-2*sqrt(R)*sqrt(r1)*sqrt(r3)*sqrt(R - r1 - r3)/(-R + r1), (-R^2 + R*r1 + R*r3 + r1*r3)/(-R + r1))
    @printf("r1 = %g;  r2 = %g;  r3 = %g;  R = %g;  x2 = %g;  y2 = %g;  x3 = %g;  y3 = %g\n",
        r1, r2, r3, R, x2, y2, x3, y3)
    plot()
    circle(0, 0, R, :green)
    rotate(0, R - r1, r1, angle=180)
    rotate(x2, y2, r2, :blue, angle=180)
    rotate(x3, y3, r3, :magenta, angle=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/2)
        point(0, R - r1, "甲円:r1,(0,R-r1)", :red, :center, delta=-delta/2)
        point(x2, y2, "乙円:r2,(x2,y2)", :blue, :center, delta=-delta/2)
        point(x3, y3, "丙円:r3\n(x3,y3)", :magenta, :center, delta=-delta/2)
    end
end;

draw(0.5, 0.3, 0.2, true)

算額（その1575）

算額（その1575）

愛知県名古屋市中央区 大須観音（北野山真福寺寶生院） 天保8年(1837)
深川英俊，トニー・ロスマン：聖なる数学：算額，p. 102，森北出版株式会社，2010年４月22日.

キーワード：円8個，外円
#Julia, #SymPy, #算額, #和算

外円の中に，大円 4 個，小円 3 個を容れる。外円と小円の直径が与えられたとき，大円の直径はいかほどか。

深川は，「中央の小円1と左右の小円2の直径の関係を示せ」という問題としている。
図を見ただけで直径は等しいとわかるが，ちゃんと式で示せということではある。
大円の中心の y 座標はどんな非負の値も取れるが，3個の小円の直径は同じであるということで，算額問題としては，「大円の直径を外円と小円の直径であらわせ」という問題にした。

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

using SymPy
@syms R::positive, r1::positive, y1::positive, r2::positive, r3::positive
eq1 = r1^2 + y1^2 - (R - r1)^2
eq2 = r1^2 + y1^2 - (r1 + r2)^2
res = solve([eq1, eq2], (r1, y1))[1]

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

大円の直径は，外円の直径から小円の直径を差し引いたものの半分である。

外円の直径が同じでも，小円の直径により印象の異なる図が得られる。

function draw(R, r2, more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    (r1, y1) = ((R - r2)/2, sqrt(R)*sqrt(r2))
    r3 = r2
    p = plot()
    circle(0, 0, R, :green)
    circle4(r1, y1, r1, :magenta)
    circle2(R - r2, 0, r2, :blue)
    circle(0, 0, r3)
    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(R, 0, "R", :green, :left, :bottom, delta=delta/2, deltax=delta/2)
        point(R - r2, 0, "R-r2", :blue, :center, delta=-delta)
        point(r3, 0, "r3", :red, :right, :vcenter, deltax=-delta/2)
        point(r1, y1, "大円:r1,(r1,y1)", :magenta, :center, :bottom, delta=delta/2)
        pos = max(R, r2)
        str = @sprintf(" R=%.3g, r2=%.3g", R, r2)
        point(0, pos, str, :black, :center, :bottom, delta=delta, mark=false)
    end
    return p
end;

draw(1, 0.1, true)

算額（その1574）

算額（その1574）

三重県鳥羽市 観世音 寛政12年(1800)
深川英俊，トニー・ロスマン：聖なる数学：算額，p. 97，森北出版株式会社，2010年４月22日.

キーワード：円1個，二等辺三角形2個
#Julia, #SymPy, #算額, #和算

外円の中に大小の正三角形を容れる。外円の直径が与えられたとき，小正三角形の一辺の長さはいかほどか。

外円の半径と中心座標を R, (0, 0)
小正三角形の一辺の長さを q
とおき，以下の方程式を解く。

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

using SymPy
@syms R, q
eq = (q/2)^2 + (-R/2 - √Sym(3)*q/2)^2 - R^2
res = solve(eq, q)[1]
res |> println

    R*(-sqrt(3) + sqrt(15))/4

小正三角形の一辺の長さ q は，外円の半径 R の (√15 - √3)/4 倍である。
外円の直径が 100 のとき，小正三角形の一辺の長さは 26.76165673298175 である。

100/2*(√15 - √3)/4

    26.76165673298175

function draw(R, more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    q = R*(√15 - √3)/4
    @printf("外円の直径が = %g のとき，小正三角形の一辺の長さは %g である。\n", 2R, q)
    plot()
    circle(0, 0, R)
    polygon(0, 0, R, 3, color=:blue)
    polygon(0, -R/2 - q/√3, q/√3, 3, color=: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, R, "R", :red, :center, :bottom, delta=delta/2)
        point(0, -R/2, "-R/2", :blue, :center, :bottom, delta=delta/2)
        point(√3R/2, -R/2, "(√3R/2,-R/2)  ", :blue, :right, :bottom, delta=delta/2)
        point(q/2, -R/2 - √3q/2, "(q/2,-R/2-√3q/2)", :blue, :left, delta=-delta/2)
    end  
end;

draw(100/2, true)

 

算額（その1573）

算額（その1573）

岡山県瀬戸内市長船町土師宮森 片山日子神社 明治6年(1873)
http://www.wasan.jp/okayama/katayamahiko.html
深川英俊，トニー・ロスマン：聖なる数学：算額，p. 96，森北出版株式会社，2010年４月22日.

キーワード：円1個，二等辺三角形
#Julia, #SymPy, #算額, #和算

二等辺三角形に円を容れる。底辺と斜辺が与えられたとき，円の直径はいかほどか。

底辺と斜辺を a, b
円の半径と中心座標を r, (0, r)
とおき，以下の連立方程式を解く。

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

using SymPy
@syms a, b, r, h
h = sqrt(b^2 - (a/2)^2)
eq1 = ((a + 2b)*r/2) - (h*a/2)
res = solve(eq1, r)[1]
res |> println

    a*sqrt(-a^2 + 4*b^2)/(2*(a + 2*b))

円の半径は a*sqrt(4b^2 - a^2)/(2a + 4b) である。
底辺が 12，斜辺が 10 のとき，内接円の直径は 6 である。

2*res(a => 12, b => 10) |> println

    6

以下のほうが簡単。

@syms a, b, r, h
h = sqrt(b^2 - (a/2)^2)
eq2 = (a/2)/b - r/(h - r)
res2 = solve(eq2, r)[1]
res2 |> println

    a*sqrt(-a^2 + 4*b^2)/(2*(a + 2*b))

function draw(a, b, more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    r = a*sqrt(4b^2 - a^2)/(2a + 4b)
    @printf("a = %g;  b = %g;  r = %g\n", a, b, r)
    plot([a/2, 0, -a/2, a/2], [0, sqrt(b^2 - (a/2)^2), 0, 0], color=:blue, lw=0.5)
    circle(0, r, r)
    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, sqrt(b^2 - a^2/4), "sqrt(b^2-a^2/4)", :blue, :center, :bottom, delta=delta/2)
        point(a/2, 0, "a/2", :blue, :left, :bottom, delta=delta/2)
        point(0, r, " r", :red, :left, :vcenter)
    end  
end;

draw(12, 10, true)