裏 RjpWiki

Julia ときどき R, Python によるコンピュータプログラム,コンピュータ・サイエンス,統計学

算額(その1572)

2025年01月31日 | Julia

算額(その1572)

四十 岩手県一関市牧沢 牧沢八幡神社 明治5年(1872)
山村善夫:現存 岩手の算額,昭和52年1月30日,熊谷印刷,盛岡市.

http://www.wasan.jp/yamamura/yamamura.html

一関市博物館 第23回 令和6年度 和算に挑戦!
https://www.city.ichinoseki.iwate.jp/index.cfm/7,177295,c,html/177295/20241125-090413.pdf

キーワード:円,長方形,折り紙
#Julia, #SymPy, #算額, #和算

長辺,短辺が 4,3 の長方形の紙を折ってできる直角三角形に内接する円の直径はいかほどか。

山村の図は相変わらず不明瞭・不正確で(そもそも何を描いているのかさえわからない),算額の穴埋めも失敗している。「一関博物館」は流石である。

長辺,短辺の長さを a, b
折り線と長辺の交点座標を (c, 0)
元の頂点が対角線上に重なる座標を (x0, y0)
円の半径と中心座標を r, (x, y)
とおき,以下の連立方程式を解く。

円の半径だけを求めるならば eq1, eq2 だけの連立方程式を解けばよい。

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

using SymPy
@syms a, b, c, r, diag, x0, y0, x, y
eq1 = ((sqrt(a^2 + b^2) - b)^2 + c^2) - (a - c)^2
eq2 = c + sqrt(a^2 + b^2) - b - (a - c) - 2r
eq3 = y0/(a - x0) - b/a
eq4 = (x0 - c)^2 + y0^2 - c^2;

res = solve([eq1, eq2, eq3, eq4], (c, r, x0, y0))[1]

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

対角線の長さを d = sqrt(a^2 + b^2) とすれば,
折り線と長方形の頂点の交点の x 座標は c = b*(d - b)/a
円の半径は (d - a - b)/2 + c である。

以下は図を描くために円の中心座標を求める連立方程式であるが,a, b を記号のままで解くと長い式になるので,eq5, eq6 の a, b に特定の値を与えて解くほうが現実的である。

eq5 = dist2(a, 0, 0, b, x, y, r)(r => res[2], x0 => res[3], y0 => res[4])(a => 4, b => 3)
eq6 = dist2(c, 0, x0, y0, x, y, r)(c => res[1], r => res[2], x0 => res[3], y0 => res[4])(a => 4, b => 3);

solve([eq5, eq6], (x, y))[3]  # 3 of 4

    (5/2, 1/2)

function draw(a, b, more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    d = sqrt(a^2 + b^2)
    c = b*(d - b)/a
    r = (d - a - b)/2 + c
    (x0, y0) = (a*b/d, -b^2/d + b)
    (x, y) = (5/2, 1/2)
    @printf("a = %g;  b = %g;  d = %g;  c = %g;  r = %g\n", a, b, d, c, r)
    plot([0, a, a, 0, 0], [0, 0, b, b, 0], color=:green, lw=0.5, seriestype=:shape, fillcolor=:lemonchiffon)
    plot!([c, a, x0, c], [0, 0, y0, 0], color=:gray, seriestype=:shape, fillcolor=:orange)
    plot!(a .- [c, a, x0, c], b .- [0, 0, y0, 0], color=:gray, seriestype=:shape, fillcolor=:orange)
    plot!([0, c, x0, 0], [b, 0, y0, b], color=:red, seriestype=:shape, fillcolor=:white)
    plot!(a .- [0, c, x0, 0], b .- [b, 0, y0, b], color=:red, seriestype=:shape, fillcolor=:white)
    circlef(x, y, r)
    circlef(a - x, b - y, 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, b, " b", :blue, :left, :bottom, delta=delta/2) 
        point(a, b, "(a,b)", :blue, :right, :bottom, delta=delta/2) 
        point(a - c, b, "(a-c,b)", :blue, :center, :bottom, delta=delta/2) 
        point(c, 0, "c ", :blue, :right, :bottom, delta=delta/2) 
        point(x0, y0, "(x0,y0)", :blue, :left, :bottom, delta=delta/2)
        point(x, y, "(x,y)", :blue, :left, :bottom, delta=delta/2)
    end  
end;

draw(4, 3, true)

コメント    この記事についてブログを書く
  • X
  • Facebookでシェアする
  • はてなブックマークに追加する
  • LINEでシェアする
« 算額(その1571) | トップ | Julia の歪度と尖度 »
最新の画像もっと見る

コメントを投稿

ブログ作成者から承認されるまでコメントは反映されません。

Julia」カテゴリの最新記事