裏 RjpWiki

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

算額(その1433)

2024年11月29日 | Julia

算額(その1433)

(福島県福島市岡部 春日神社 明治28年(1895))
街角の数学 Street Wasan 街角の数学 ~落書き帳「○△□」~ 39.正方形と四分円(4)
http://streetwasan.web.fc2.com/math15.6.20.html
キーワード:円2個,四分円2個,正方形2個
#Julia, #SymPy, #算額, #和算

外正方形(外方)と 2 つの四分円とで囲まれた上下の隙間に,小円と大円が内接している。
外正方形の一辺の長さを 1 寸として大円と小円の半径をそれぞれ求めよ。

外正方形の一辺の長さ(=四分円の半径)を R
内正方形の一辺の長さを a
大円の半径と中心座標を r1, (0, r1)
小円の半径と中心座標を r2, (R - r2)
とおき,以下の連立方程式を解く。

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

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

# r1
res[1] |> println
res[1](R => 1) |> println

    3*R/8
    3/8

# r2
res[2] |> println
res[2](R => 1) |> println

    R/16
    1/16

# a
res[3] |> println
res[3](R => 1) |> println

    3*R/5
    3/5

外方の一辺の長さが R のとき,大円の半径 r1 はその 3/8 倍,小円の半径 r2 は 1/16 倍,内方の一辺の長さ a は 3/5 倍である。
外方の一辺の長さが 1 のとき,大円の直径は 0.75,小円の直径は 0.125,, 内方の一辺の長さは 0.6 である。
外方の一辺の長さが 80 のとき,大円の直径は 60,小円の直径は 10,内方の一辺の長さは 48 である。

四分円と大円の接点は,内方の頂点と一致する。

function draw(R, more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    (a, r1, r2) = (3*R/5, 3*R/8, R/16)
    @printf("外方の一辺の長さが %g のとき,大円の直径は %g,小円の直径は %g,, 内方の一辺の長さは %g である。\n", R, 2r1, 2r2, a)
    plot([R/2, R/2, -R/2, -R/2, R/2], [0, R, R, 0, 0], color=:green, lw=0.5)
    circle(-R/2, 0, R, beginangle=0, endangle=90)
    circle(R/2, 0, R, beginangle=90, endangle=180)
    circle(0, r1, r1, :orange)
    circle(0, R - r2, r2, :blue)
    plot!([a, a, -a, -a, a] ./2, [0, a, a, 0, 0], color=:magenta, 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(a/2, 0, " a/2", :magenta, :left, :bottom, delta=delta/2)
        point(0, R - r2, " 小円:r2,(0,R-r2)", :blue, :center, delta=7delta)
        point(0, r1, "大円:r1,(0,r1)", :orange, :center, delta=-delta/2)
        point(a/2, a, " (a/2,a)", :magenta, :left, :vcenter)
        point(R/2, R, "(R/2,R) ", :red, :right, :bottom, delta=delta/2)
    end
end;

draw(1, true)

 

コメント
  • X
  • Facebookでシェアする
  • はてなブックマークに追加する
  • LINEでシェアする

算額(その1432)

2024年11月29日 | Julia

算額(その1432)

(福島県田村市船引町石森戸屋 稲荷神社 明治18年(1885))
街角の数学 Street Wasan ~落書き帳「○△□」~ 37.正方形と四分円(2)
http://streetwasan.web.fc2.com/math15.6.18.html
キーワード:四分円2個,正方形3個
#Julia, #SymPy, #算額, #和算

算額(その1431)の拡張で,正方形をもう 1 個加える。
正方形の中に四分円を 2 個を容れ,交差した部分に正方形を 2 個容れる。外側の正方形の一辺の長さが 5 寸のとき,内側の一番小さな正方形の一辺の長さはいかほどか。

外側の正方形の左下の頂点を原点とし,一辺の長さを a
内側の正方形の右下の頂点の座標を (b, 0)
内側の一番小さな正方形の右下の頂点の座標を (c, 2(b - a/2))
とおき,以下の連立方程式を解く。
内側の 2 つの正方形の右上の頂点いずれもが,原点を中心とする四分円の円周上にあるということを利用する。
内側の正方形の一辺の長さは a - 2(a - b) = 2b - a である。
内側の一番小さい正方形の一辺の長さは a - 2(a - c) = 2c - a である。

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
eq1 = b^2 + (2b - a)^2 - a^2
eq2 = c^2 + ((2b - a) + (2c - a))^2 - a^2
res = solve([eq1, eq2], (b, c))[1]
res |> println

    (4*a/5, 3*a/5)

# b
res[1](a => 5) |> println

    4

# c
res[2](a => 5) |> println

    3

a = 5 のとき b = 4, c = 3 なので,内側の一番小さな正方形の一辺の長さは 2c - a = 1 寸である。

function draw(a, more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    b = 4a/5
    c = 3a/5
    plot([0, a, a, 0, 0], [0, 0, a, a, 0], color=:green, lw=0.5)
    plot!([a - b, b, b, a - b], [0, 0, 2(b - a/2), 2(b - a/2), 0], color=:blue, lw=0.5)
    plot!([a - c, c, c, a - c], 2(b - a/2) .+ [0, 0, 2(c - a/2), 2(c - a/2), 0], color=:magenta, lw=0.5)
    circle(0, 0, a, beginangle=0, endangle=90)
    circle(a, 0, a, beginangle=90, endangle=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(a, 0, " a", :green, :left, :bottom, delta=delta/2)
        point(b, 0, " b", :blue, :left, :bottom, delta=delta/2)
        point(b, 2(b - a/2), "(b,2(b-a/2))", :blue, :left, :bottom, delta=delta/2)
        point(c, 2(b - a/2) + 2(c - a/2), "(c,2(b-a/2)+2(c-a/2))", :magenta, :left, :bottom, delta=delta/2)
    end
end;

draw(5, true)

コメント
  • X
  • Facebookでシェアする
  • はてなブックマークに追加する
  • LINEでシェアする

算額(その1431)

2024年11月29日 | Julia

算額(その1431)

福島県田村市船引町石森戸屋 稲荷神社 明治18年(1885)
街角の数学 Street Wasan ~落書き帳「○△□」~ 36.正方形と四分円

http://streetwasan.web.fc2.com/math15.6.17.html
キーワード:四分円2個,正方形2個
#Julia, #SymPy, #算額, #和算

正方形の中に四分円を 2 個を容れ,交差した部分に正方形を 1 個容れる。外側の正方形の一辺の長さが 5 寸のとき,内側の正方形の一辺の長さはいかほどか。

外側の正方形の左下の頂点を原点とし,一辺の長さを a
内側の正方形の右下の頂点の座標を (b, 0)
とおき,以下の連立方程式を解く。
内側の正方形の右上の頂点が原点を中心とする四分円の円周上にあるということを利用する。
内側の正方形の一辺の長さは a - 2(a - b) = 2b - a である。

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

using SymPy
@syms a::positive, b::positive
eq = b^2 + (2b - a)^2 - a^2
res = solve(eq, b)[1]
res |> println
res(a => 5) |> println

    4*a/5
    4

a = 5 のとき b = 4 なので,内側の正方形の一辺の長さは 2b - a = 3 寸である。

function draw(a, more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    b = 4a/5
    plot([0, a, a, 0, 0], [0, 0, a, a, 0], color=:green, lw=0.5)
    plot!([a - b, b, b, a - b], [0, 0, 2(b - a/2), 2(b - a/2), 0], color=:blue, lw=0.5)
    circle(0, 0, a, beginangle=0, endangle=90)
    circle(a, 0, a, beginangle=90, endangle=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(a, 0, " a", :green, :left, :bottom, delta=delta/2)
        point(b, 0, " b", :blue, :left, :bottom, delta=delta/2)
        point(b, 2(b - a/2), "(b,2(b-a/2))", :blue, :left, :bottom, delta=delta/2)
    end
end;

draw(5, true)

コメント
  • X
  • Facebookでシェアする
  • はてなブックマークに追加する
  • LINEでシェアする

算額(その1430)

2024年11月29日 | Julia

算額(その1430)

八十九 陸前高田市小友町 常膳寺観音堂 天保13年(1842)
山村善夫:現存 岩手の算額,昭和52年1月30日,熊谷印刷,盛岡市.

http://www.wasan.jp/yamamura/yamamura.html
キーワード:円12個,外円,正三角形3個
#Julia, #SymPy, #算額, #和算

正三角形の中に斜線を 2 本,全円 1 個,甲円 1 個,乙円 3 個,全円 1 個を容れる。甲円の直径が与えられたとき,丙円の直径を得る術を述べよ。

正三角形の一辺の長さを 2a
甲円の半径と中心座標を r1, (0, 2r0 - r1)
乙円の半径と中心座標を r2, (x2, y2), (0, 2r0 + r2)
丙円の半径と中心座標を r3, (x3, y3)
全円の半径と中心座標を r0, (0, r0)
とおき,以下の連立方程式を解く。

なお,r0 = √3a/3, r2 = √3a/9, y3 = y2 である。

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

using SymPy
@syms a, r0, r1, r2, x2, y2, r3, x3, y3, b
# eq1 = r0 - a/√Sym(3)
# eq2 = 2r2 - (√Sym(3)a - 2r0 - r2)
# res = solve([eq1, eq2], (r0, r2))
# r0 = √Sym(3)a/3
# r2 = √Sym(3)a/9

r0 = √Sym(3)a/3
r2 = √Sym(3)a/9
y3 = y2
# dist2() は使えない
eq3 = dist(-a, 0, b, √Sym(3)*(a - b), x2, y2) - r2^2
eq4 = dist(-a, 0, b, √Sym(3)*(a - b), 0, 2r0 - r1) - r1^2
eq5 = dist(a, 0, -b, √Sym(3)*(a - b), x2, y2) - r2^2
eq6 = dist(a, 0, 0, √Sym(3)a, x2, y2) - r2^2
eq7 = dist(-a, 0, b, √Sym(3)*(a - b), x3, y3) - r3^2
eq8 = (x2 - x3) - (r2 + r3);

using NLsolve

function nls(func, params...; ini = [0.0])
    if typeof(ini) <: Number
        r = nlsolve((vout, vin) -> vout[1] = func(vin[1], params..., [ini]), ftol=big"1e-40")
        v = r.zero[1]
    else
        r = nlsolve((vout, vin)->vout .= func(vin, params...), ini, ftol=big"1e-40")
        v = r.zero
    end
    return Float64.(v), r.f_converged
end;

function H(u)
    (x2, y2, r1, b, r3, x3) = u
    return [
-a^2/27 + (y2 - sqrt(3)*(a - b)*(sqrt(3)*y2*(a - b) + (a + b)*(a + x2))/(3*(a - b)^2 + (a + b)^2))^2 + (a + x2 - (a + b)*(sqrt(3)*y2*(a - b) + (a + b)*(a + x2))/(3*(a - b)^2 + (a + b)^2))^2,  # eq3
-r1^2 + (a - (a + b)*(a*(a + b) + sqrt(3)*(a - b)*(2*sqrt(3)*a/3 - r1))/(3*(a - b)^2 + (a + b)^2))^2 + (2*sqrt(3)*a/3 - r1 - sqrt(3)*(a - b)*(a*(a + b) + sqrt(3)*(a - b)*(2*sqrt(3)*a/3 - r1))/(3*(a - b)^2 + (a + b)^2))^2,  # eq4
-a^2/27 + (y2 - sqrt(3)*(a - b)*(sqrt(3)*y2*(a - b) + (-a - b)*(-a + x2))/((-a - b)^2 + 3*(a - b)^2))^2 + (-a + x2 - (-a - b)*(sqrt(3)*y2*(a - b) + (-a - b)*(-a + x2))/((-a - b)^2 + 3*(a - b)^2))^2,  # eq5
-a^2/27 + (y2 - sqrt(3)*a*(sqrt(3)*a*y2 - a*(-a + x2))/(4*a^2))^2 + (-a + x2 + a*(sqrt(3)*a*y2 - a*(-a + x2))/(4*a^2))^2,  # eq6
-r3^2 + (y2 - sqrt(3)*(a - b)*(sqrt(3)*y2*(a - b) + (a + b)*(a + x3))/(3*(a - b)^2 + (a + b)^2))^2 + (a + x3 - (a + b)*(sqrt(3)*y2*(a - b) + (a + b)*(a + x3))/(3*(a - b)^2 + (a + b)^2))^2,  # eq7
-sqrt(3)*a/9 - r3 + x2 - x3,  # eq8
    ]
end;

a = 10
r0 = √3a/3  # a*0.5773502691896257
r2 = √3a/9  # a*0.19245008972987523
iniv = BigFloat[0.575411, 0.34912, 0.387466, 0.662547, 0.095137, 0.284131] .* a
res = nls(H, ini=iniv)

    ([5.711594854986474, 3.5787338000056996, 3.8641388667743293, 6.575250125060548, 0.9544492311458473, 2.832644726541875], true)

正三角形の一辺の長さが 2a = 20 のとき,
全円の半径 r0 = 5.773502691896257
甲円の半径 r1 = 3.8641388667743293
乙円の半径 r2 = 1.9245008972987523
丙円の半径 r3 = 0.9544492311458473

丙円の直径/甲円の直径 = 0.9544492311458473/3.8641388667743293 = 0.2470017936862072
甲円の直径/丙円の直径 = 3.8641388667743293/0.9544492311458473 = 4.048553595810754

術は,「15572 の平方根から 121 を引き,甲円の直径を掛ければ全円の直径がえられる」とあるが,√15572 - 121 = 3.787819918451973 であり,題意を満たさない。

function draw(a, more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    (x2, y2, r1, b, r3, x3) = [5.711594854986474, 3.5787338000056996, 3.8641388667743293, 6.575250125060548, 0.9544492311458473, 2.832644726541875]
    y3 = y2
    r0 = √3a/3
    r2 = √3a/9
    plot([a, 0, -a, a], [0, sqrt(3)a, 0, 0], color=:green, lw=0.5)
    circle(0, a/sqrt(3), r0)
    circle(0, 2r0 + r2, r2, :blue)
    circle2(x2, y2, r2, :blue)
    circle2(x3, y3, r3, :brown)
    circle(0, 2r0 - r1, r1, :orange)
    segment(-a, 0, b, (a - b)*√3)
    segment(a, 0, -b, (a - b)*√3)
    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, √3a, " √3a", :green, :left, :vcenter)
        point(a, 0, "a", :green, :left, :bottom, delta=delta/2)
        point(0, r0, "全円:r0,(0,r0)", :red, :center, delta=-delta)
        point(0, 2r0 - r1, "甲円:r1,(0,2r0-r1)", :orange, :center, delta=-delta)
        point(x2, y2, "乙円:r2\n(x2,y2)", :blue, :left, delta=-delta, deltax=-3delta)
        point(0, 2r0 + r2, " 乙円:r2\n(x2,r2)", :blue, :center, delta=-delta)
        point(x3, y3, "丙円:r3,(x3,y3)", :brown, :right, delta=-4delta, deltax=-delta)
    end
end;

draw(10, true)

コメント
  • X
  • Facebookでシェアする
  • はてなブックマークに追加する
  • LINEでシェアする

橋本製麺所

2024年11月29日 | さぬきうどん

高松市仏生山町 橋本製麺所
旧二級国道沿いにあるが,店名表示も暖簾も表札すらもない。何百回となく店の前を通っているのに,製麺所とは知らなかった。
いわゆるうどん店ではない。うどん玉1個90円で売っており,丼と箸を持っていけば店内で食べることもできる(店内に数脚のパイプ椅子があるがテーブルはない)。
家に帰り,スーパーで買ってきた天ぷらを乗せて掛けうどんでいただく。

田舎そば(これも玉売)


コメント
  • X
  • Facebookでシェアする
  • はてなブックマークに追加する
  • LINEでシェアする

PVアクセスランキング にほんブログ村

PVアクセスランキング にほんブログ村