算額（その1486）

算額（その1486）

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

大円と斜線 6 本で区分された領域に，小円 6 個，楕円 3 個を容れる。大円，小円の直径が 3 寸，1 寸のとき，楕円の長径はいかほどか。

算法助術の公式 97 を適用する。
⊿ABC に楕円が内接し，⊿ABC の底辺 BC と楕円の長軸が平行である。楕円の長軸，短軸をそれぞれ p, q, BC = a, AC = b, AB = c, AH = h とすると
-(b^2 - c^2)^2*h*q^2 + (b^2 - c^2)^2*q^3 + a^4*h*(2h - q)^2 - a^4*q*(2h - q)^2 - a^2*h*p^2*(2h -  q)^2 = 0 
が成り立つ。
大円の半径と中心座標を r1, (0, 0)
小円の半径と中心座標を r2, (0, r1 + r2), (0, r2 - r1)
q = 2r2, h = r1 + 2r2, a = b = c = 2h/√3 である。

include("julia-source.txt");

using SymPy

@syms p, q, a, b, c, h
eq97 = -(b^2 - c^2)^2*h*q^2 + (b^2 - c^2)^2*q^3 + a^4*h*(2h - q)^2 - a^4*q*(2h - q)^2 - a^2*h*p^2*(2h -  q)^2

q = 1, a = b = c = 5/√3, h = 2.5 が既知で，p が未知なので，以下の方程式を解いて p = √5 = 2.23606797749979 を得る。

eq97_2 = eq97(q => 1, a =>5/√Sym(3), b => 5/√Sym(3), c => 5/√Sym(3), h => 25//10)
eq97_2 |> println

    5000/3 - 1000*p^2/3

ans_p = solve(eq97_2,  p)[2] # 2 of 2
ans_p |> println
ans_p.evalf() |> println

    sqrt(5)
    2.23606797749979

function draw(more=false)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    r1 = 3/2
    r2 = 1/2
    q = 1
    a = b = c = 5/√3
    h = 2.5
    θ = 0:60:360
    x = 5/√3 .*cosd.(θ)
    y = 5/√3 .*sind.(θ)
    plot()
    # circle(0, 0, 5/√3)
    ellipse(0, r1 + r2, √5/2, 1/2, color=:orange)
    ellipse((r1 + r2)*cosd(30), -(r1 + r2)*sind(30), √5/2, 1/2, φ=60, color=:orange)
    ellipse(-(r1 + r2)*cosd(30), -(r1 + r2)*sind(30), √5/2, 1/2, φ=120, color=:orange)
    circle(0, 0, r1, :blue)
    rotate(0, r2 - r1, r2)
    rotate(0, r1 + r2, r2)
    for i = 1:3
        j = i + 3
        segment(x[i], y[i], x[i + 1], y[i + 1], :green)
        segment(x[j], y[j], x[j + 1], y[j + 1], :green)
        segment(x[i], y[i], x[i + 3], y[i + 3], :green)
    end
    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, "r1", :blue, :center, delta=-delta)
        point(0, r1 + r2, "r1+r2", :red, :center, delta=-delta)
    end
end;

draw(true)