算額（その1129）改訂版

算額（その1129）改訂版

算額（その1129）は依拠した図が間違っているので，改訂版を書いた。

四十五 岩手県一関市真滝 熊野白山滝神社 明治13年(1880)
山村善夫：現存 岩手の算額，昭和52年1月30日，熊谷印刷，盛岡市. http://www.wasan.jp/yamamura/yamamura.html

今有如図 03061
https://w.atwiki.jp/sangaku/pages/327.html

キーワード：円4個，長方形，斜線2本
#Julia, #SymPy, #算額, #和算

山村の図では，条件不足で解が一意に定まらない。「今有如図」の図に基づいて得られた解は「答」，「術」と一致する。

長方形の中に 2 本の斜線を引き，区画された領域に大円 1 個，中円 1 個，小円 2 個を容れる。小円の直径が 1 寸のとき，大円の直径はいかほどか。

長方形の長辺，短辺を a, b
斜線と長方形の辺の交点座標を (c, 0), (a, d)
大円の半径と中心座標を r1, (r1, r1)
中円の半径と中心座標を r2, (a - r2, r2)
小円の半径と中心座標を r3, (r3, y31), (a - r3, y32)
とおき，以下の連立方程式を解く。

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

using SymPy
@syms a::positive, b::positive, c::positive, d::positive,
      r1::positive, r2::positive,
      r3::positive, y31::positive, y32::positive
eq1 = dist(0, b, a, d, r1, r1) - r1^2
eq2 = dist(0, b, a, d, a - r2, r2) - r2^2
eq3 = dist(0, b, a, d, r3, y31) - r3^2
eq4 = dist(0, b, a, d, a - r3, y32) - r3^2
eq5 = dist(c, 0, a, b, r1, r1) - r1^2
eq6 = dist(c, 0, a, b, a - r2, r2) - r2^2
eq7 = dist(c, 0, a, b, a - r3, y32) - r3^2
eq8 = (r1 - r3)^2 + (y31 - r1)^2 - (r1 + r3)^2;
# res = solve([eq1, eq2, eq3, eq4, eq5, eq6, eq7, eq8], (a, b, c, d, r1, r2, y31, y32))

function H(u)
    (a, b, c, d, r1, r2, y31, y32) = u
    return [
        dist(0, b, a, d, r1, r1) - r1^2,
        dist(0, b, a, d, a - r2, r2) - r2^2,
        dist(0, b, a, d, r3, y31) - r3^2,
        dist(0, b, a, d, a - r3, y32) - r3^2,
        dist(c, 0, a, b, r1, r1) - r1^2,
        dist(c, 0, a, b, a - r2, r2) - r2^2,
        dist(c, 0, a, b, a - r3, y32) - r3^2,
        (r1 - r3)^2 + (y31 - r1)^2 - (r1 + r3)^2
    ]
end;

r3 = 1/2
iniv = BigFloat[4.6, 4.1, 2.4, 1.4, 1.5, 0.8, 3.2, 2.2]
res = nls(H, ini=iniv)

    ([4.732050807568878, 4.098076211353316, 2.366025403784439, 1.3660254037844386, 1.5, 0.8660254037844386, 3.232050807568877, 2.232050807568877], true)

小円の直径が 1 のとき，大円の直径は 3, 中円の直径は 1.73205 である。

全てのパラメータは，以下のとおり。

r3 = 0.5;  a = 4.73205;  b = 4.09808;  c = 2.36603;  d = 1.36603;  r1 = 1.5;  r2 = 0.866025;  y31 = 3.23205;  y32 = 2.23205

function draw(r3, more=false)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    (a, b, c, d, r1, r2, y31, y32) = res[1]
    @printf("小円の直径が %g のとき，大円の直径は %g, 中円の直径は %g である。\n", 2r3, 2r1, 2r2)
    @printf("r3 = %g;  a = %g;  b = %g;  c = %g;  d = %g;  r1 = %g;  r2 = %g;  y31 = %g;  y32 = %g\n", r3, a, b, c, d, r1, r2, y31, y32)
    plot([0, a, a, 0, 0], [0, 0, b, b, 0], color=:green, lw=0.5)
    circle(r1, r1, r1)
    circle(a - r2, r2, r2, :blue)
    circle(r3, y31, r3, :magenta)
    circle(a - r3, y32, r3, :magenta)
    segment(0, b, a, d, :orange)
    segment(c, 0, a, b, :orange)
    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(a, b, "(a,b)", :green, :right, :bottom, delta=delta/2)
        point(r1, r1, "大円:r1,(r1,r1)", :red, :center, delta=-delta)
        point(a - r2, r2, "中円:r2,(a-r2,r2)", :blue, :center, delta=-delta)
        point(r3, y31, "小円:r3\n(r3,y31)", :magenta, :center, delta=-delta)
        point(a - r3, y32, "小円:r3\n(a-r3,y32)", :magenta, :center, delta=-delta)
    end
end;

draw(1/2, true)