算額（その1503）

算額（その1503）

六十三 花泉町金沢 大門神社・大門観世音菩薩堂 明治13年(1880)
山村善夫：現存 岩手の算額，昭和52年1月30日，熊谷印刷，盛岡市.
http://www.wasan.jp/yamamura/yamamura.html

今有如図 03064
https://w.atwiki.jp/sangaku/pages/292.html

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

大円の一部である円弧（弓形）の中に小円を 2 個と界斜を容れる。大円の径が 8 寸，矢が 2 寸，小円の直径が 1 寸のとき，界斜はいかほどか。

山村の図は右側の小円と界斜の端点の位置が違う。

「今有如図」の仮定で解を導き，答と一致することを確認した。

大円の半径と中心座標を R, (0, 0)
矢をそのまま「矢」
小円の半径と中心座標を r, (x, R - 矢 + r)
界斜の端点を (x1, R - 矢), (x2, sqrt(R^2 - x2^2))
とおき，以下の連立方程式を解く。

include("julia-source.txt");

using SymPy

@syms R::positive, 矢::positive, x::positive, r::positive,
      x1::positive, x2::negative
y = R - 矢
eq1 = x^2 + (y + r)^2 - (R - r)^2
eq2 = dist(x2, sqrt(R^2 - x2^2), x1, y, x, y + r) - r^2
eq3 = (sqrt(R^2 - x2^2) - y)/(x1 - x2) - r/x1;
# res = solve([eq1, eq2, eq3], (x, x1, x2))

function H(u)
    (x, x1, x2) = u
    return [
        x^2 + (y + r)^2 - (R - r)^2,
        dist(x2, sqrt(R^2 - x2^2), x1, y, x, y + r) - r^2,
        (sqrt(R^2 - x2^2) - y)/(x1 - x2) - r/x1
    ]
end;
R = 8//2
矢 = 2
y = R - 矢
r = 1//2
iniv = BigFloat[2.5, 2.5, -2.6]
res = nls(H, ini=iniv)

    ([2.449489742783178, 2.3979157616563596, -2.597742075127723], true)

界斜は sqrt((x1 - x2)^2 + (sqrt(R^2 - x2^2) - (R - 矢))^2) により求めることができる。
大円の直径が 8 寸，矢が 2 寸，円の直径が 1 寸のとき，5.103103630798287 寸である。

(x, x1, x2) = res[1]
sqrt((x1 - x2)^2 + (sqrt(R^2 - x2^2) - (R - 矢))^2)

    5.103103630798287

function draw(r1, r2, more=false)
    pyplot(size=(500, 500), showaxis=true, grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    R = 8//2
    矢 = 2
    r = 1//2
    y = R - 矢
    (x, x1, x2) = res[1]
    界斜 = sqrt((x1 - x2)^2 + (sqrt(R^2 - x2^2) - (R - 矢))^2)
    println("界斜 = $界斜")
    plot()
    circle(0, 0, R, :purple)
    segment(x2, sqrt(R^2 - x2^2), x1, R - 矢, :orange)
    circle2(x, y + r, r, :blue)
    segment(-sqrt(R^2 - y^2), y, sqrt(R^2 - y^2), y)
    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(x2, sqrt(R^2 - x2^2), "(x2,sqrt(R^2-x2^2))", :green, :left, :bottom, deltax=3delta)
        point(x1, y, "(x1,y)", :green, :center, delta=-delta)
        point(R, 0, " R", :purple, :left, :bottom, delta=delta)
        point(x, y + r, "(x,y+r)", :blue, :center, delta=-delta)
        ylims!(0, R)
    end
end;

draw(12/2, 8/2, true)