裏 RjpWiki

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

算額(その1596)

2025年02月06日 | Julia

算額(その1596)

三重県伊賀市上野東町 菅原神社(上野天神宮) 嘉永7年(1854)
深川英俊,トニー・ロスマン:聖なる数学:算額,p. 116,森北出版株式会社,2010年4月22日.
キーワード:体積,級数
#Julia, #SymPy, #算額, #和算, #数学

底面が辺の長さ a の正三角形で高さも a である大きな角錐を逆さにして酒を満たし,その体積を V0 とする。
これより 1 石  2斗 5 升(125升)をくみ出して残りを 2 倍にする。次に同じように 125 升をくみ出して残りを 2 倍にする。こうして 10 回目に 125 升をくみ出すと残りが無くなった。さて,この容器の体積 V を求めよう。また 1 升は一辺が 4.0172 寸の立方体としたとき,a を寸で求めよ。

深川英俊,トニー・ロスマン:聖なる数学:算額,p. 114,森北出版株式会社,2010年4月22日.
キーワード:等比級数
#Julia, #SymPy, #算額, #和算, #数学

酒の量は,
V0 : 初期値
V1 = 2(V0 - 125) : 1 回目の操作の後の酒の量
V2 = 2(V1 - 125) : 2 回目の操作の後の酒の量
    :
V10 = 2(V9 - 125) : 10 回目の操作の後の酒の量 = 0
となる。最後の V10 の中の V9 に直前の V9 を代入して...を繰り返し,V10 を V0 で表した方程式を解けばよいが,面倒。

繰り返し作業はコンピュータの最も得意とする所。以下のようにして簡単に解くべき方程式を得ることができる。
等比級数だとか等比級数の和だとか,なんにも考える必要がない。

using SymPy
@syms V[1:11]::positive, V0::positive
V[1] = V0  # 初期状態
for i in 2:11  # 操作を繰り返す。結果は配列の中に代入されていく。
    V[i] = 2(V[i - 1] - 125)
end

V[1] |> println

    V0

V[2] |> println

    2*V0 - 250

V[3] |> println

    4*V0 - 750

V[11]  |> println  # 10 回目の操作が終わった段階での量。

    1024*V0 - 255750

ans_V0 = solve(V[11], V0)[1]   # 1024*V0 = 255750 を解くということ
ans_V0 |> println  # これが「0」であるということなので,方程式を解く。

    127875/512

体積は 249.755859375 升である。

ans_V0.evalf() |> println

    249.755859375000

(1) 底面が一辺 a 寸の正三角形で高さが a の三角錐の体積は (a/2 * √3a/2 * a / 3) で単位は「立法寸」
(2) V0 の単位は「升」なので,「立法寸」にするためには 4.0172 を掛ける。(4.0172^3 立法寸 = 1 升)
(1),(2) が等しいとして,以下の方程式を解いて a を求める。

@syms a
eq = (a/2 * √3a/2 * a / 3) - ans_V0*4.0172^3  # (1) - (2) = 0
solve(eq, a)[1] |> println

    48.2283298085359

a = 48.2283298085359 寸である。

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

算額(その1595)

2025年02月06日 | Julia

算額(その1595)

秋田県仙北郡角館町 熊野神社 安政5年(1858)
深川英俊,トニー・ロスマン:聖なる数学:算額,p. 114,森北出版株式会社,2010年4月22日.
キーワード:俵積み
#Julia, #SymPy, #算額, #和算, #数学

N 個の俵がある。これを 2 通りの方法で積むことができた。
1 つは上が 19 個で下が m 個の台形状になり,次は上が 6 個で下が n 個の台形状である。このとき N, m, n を求めてみよう。

以下の方程式を立てる。

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

using SymPy
@syms m::positive, n::positive, N::positive,
      h1::positive, h2::positive
eq1 = m*(m + 1)/2 - 18(18 + 1)//2 - N
eq2 = n*(n + 1)/2 - 5(5 + 1)//2 - N;

辺々引いて,N を消去する。

eq1 - eq2 |> factor |> println

    (m^2 + m - n^2 - n - 312)/2

定数項を右辺にして残りを因数分解する。

m^2 + m - n^2 - n |> factor |> println

    (m - n)*(m + n + 1)

深川は
(m - n)*(m + n + 1) = 312
で,312 の素因数分解で
(m - n)*(m + n + 1) = 8×39
(m - n)*(m + n + 1) = 1×312
の 2 通りあげているが,もう一組ある。
(m - n)*(m + n + 1) = 3×104

そもそも,因数の組は (a, b) = (1, 312), (2, 156), (3, 104), (4, 78), (6, 52), (8, 39), (12, 26), (13, 24) のすべての組み合わせについて
eq1 = (m - n) - a
eq2 = (m + n + 1) - b
の連立方程式を解いてみて,a, b が整数のものをピックアップするという作業が必要。あたりを付けて 1 つ,2 つ見つけて安心してはいけない。
これは結構大変である。数え落としのない,ブルートフォースで求めるのが安全確実。

for m = 19:1000, n = 6:m
    (m^2 + m - n^2 - n) == 312 &&
        println("m = $m, n = $n, (m - n) = $(m - n), (m + n + 1) = $(m + n + 1), $(m^2 + m - n^2 - n), N = $(Int(m*(m + 1)/2 - 18(18 + 1)/2))")
end

    m = 23, n = 15, (m - n) = 8, (m + n + 1) = 39, 312, N = 105
    m = 53, n = 50, (m - n) = 3, (m + n + 1) = 104, 312, N = 1260
    m = 156, n = 155, (m - n) = 1, (m + n + 1) = 312, 312, N = 12075

検算 N

(sum(19:23), sum(6:15)) |> println
(sum(19:53), sum(6:50)) |> println
(sum(19:156), sum(6:155)) |> println

    (105, 105)
    (1260, 1260)
    (12075, 12075)

 

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

算額(その1594)

2025年02月06日 | Julia

算額(その1594)

深川英俊,トニー・ロスマン:聖なる数学:算額,p. 114,森北出版株式会社,2010年4月22日.
キーワード:正12角形2個,対角線
#Julia, #SymPy, #算額, #和算, #数学

図に示すような正12角形の24本の対角線の長さの和が与えられたとき,内部にできる小さな正12角形の一辺の長さを求める術を述べよ。

外側の正12角形が内接する円の直径を R
内側の正12角形が内接する円の直径を r
とすると,
a = c*tan(15°), c = (a + b)*tan(15°)
対角線の長さの和を l とすると,
l = 12*2*(a + b)
より,以下の連立方程式を解き,a, b を求める。

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,
      l::positive
eq1 = (a + b)*tand(Sym(15)) - a/tand(Sym(15))
eq2 = 12*2(a + b) - l;

表には出てこないが,tan(15°) = 2 - √3 なので,
eq1 = -2*sqrt(3)*a - sqrt(3)*b + 2*b
である。

連立方程式を解いて,a, b を求める。

res = solve([eq1, eq2], (a, b));

# a
res[a] |> simplify |> println

    l*(7 - 4*sqrt(3))/24

l = 150 のとき,a = 0.897459621556135 である。

2*res[a](l => 150).evalf() |> println

    0.897459621556135

長精度演算すると以下のようになる。

2*150*(7 - 4*sqrt(big"3"))/24

    0.897459621556135323627682924706381...

# b
res[b] |> simplify |> println

    l*(-3 + 2*sqrt(3))/12

a = r*sin(15°) なので,
r = a/sin(15°)

c = R*sin(15°) なので,
R = c/sin(15°)

c = (a + b)*tan(15°) なので,
R = (a + b)*tan(15°)/sin(15°)

function draw(n, l, more=false)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    a = l*(7 - 4√3)/24
    b = l*(2√3 - 3)/12
    r = a/sind(15)
    R = (a + b)*tand(15)/sind(15)
    x = [R*cosd(z) for z = 0:360/n:360]
    y = [R*sind(z) for z = 0:360/n:360]
    c = R*sind(15)
    plot()
    for i = 1:n
        segment(x[i], y[i], x[mod(i - 2 + (n + 1)÷2, n) + 1], y[mod(i - 2 + (n + 1)÷2, n) + 1], :brown, lw=0.3)
    end
    circle(0, 0, R, :green)
    circle(0, 0, r, :red)
    circle(0, 0, c, :blue)
    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(r, 0, " r", :red, :left, :bottom, delta=delta/5)
        point(c, 0, "c ", :blue, :right, :bottom, delta=delta/5)
        lens!([0, 0.2], [-0.01, 0.2], inset = (1, bbox(0.7, 0.0, 0.3, 0.35)))
    end
end;

draw(12, 15, true)

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

算額(その1593)

2025年02月06日 | Julia

算額(その1593)

長野県 諏訪神社 文化2年(1805)
深川英俊,トニー・ロスマン:聖なる数学:算額,p. 114,森北出版株式会社,2010年4月22日.

キーワード:円1個,長方形,対角線
#Julia, #SymPy, #算額, #和算, #数学

長方形の中に円と対角線を引く。長方形の長辺 AB と短辺 BC が与えられたとき,円によって切り取られる対角線の長さ PQ を求めよ。

円の半径を r
a = AB, b = BC = 2r
P, Q の座標を (x1, y1), (x2, y2)
θ = ∠BDC
とおき,以下の連立方程式を解く。

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

using SymPy
@syms a::positive, b::positive, r::positive,
      x1::positive, x2::positive, PQ::positive
r = b/2
tanθ = b/a
y1 = (a - x1)*tanθ
y2 = (a - x2)*tanθ
eq1 = (x1 - (a - r))^2 + (y1 - r)^2 - r^2
eq2 = (x2 - (a - r))^2 + (y2 - r)^2 - r^2
eq3 = sqrt((x1 - x2)^2 + (y1 - y2)^2) - PQ;

res = solve([eq1, eq2, eq3], (PQ, x1, x2))[1];  # 1 of 2

# PQ
res[1] |> println

    sqrt(2)*sqrt(a)*b^(3/2)/sqrt(a^2 + b^2)

PQ の長さは sqrt((2a*b^3)/(a^2 + b^2)) である。

a = 185, b = 80 のとき,PQ = 68.2871764062511 である。

res[1](a => 185, b => 80).evalf() |> println

    68.2871764062511

function draw(a, b, more=false)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
    (x1, x2, PQ) = ((-sqrt(2)*a^(3/2)*b^(3/2) + a*(2*a^2 - a*b + b^2))/(2*(a^2 + b^2)), (sqrt(2)*a^(3/2)*b^(3/2) + a*(2*a^2 - a*b + b^2))/(2*(a^2 + b^2)), sqrt(2)*sqrt(a)*b^(3/2)/sqrt(a^2 + b^2))
    r = b/2
    tanθ = 2r/a
    θ = atand(tanθ)
    y1 = (a - x1)*tanθ
    y2 = (a - x2)*tanθ
    PQ2 = sqrt((x1 - x2)^2 + (y1 - y2)^2)
    println("PQ = $PQ = $PQ2")
    plot([0, a, a, 0, 0], [0, 0, 2r, 2r, 0], color=:magenta, lw=0.5)
    circle(a - r, r, r, :blue)
    segment(0, 2r, a, 0, :red)
    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 - r, r, "r,(a-r,r)", :blue, :center, delta=-delta/2)
        point(x1, y1, " Q:(x1,y1)", :blue, :left, :bottom, delta=delta/2)
        point(x2, y2, "P:(x2,y2)  ", :blue, :right, :vcenter)
        point(a, b, "A:(a,b)", :magenta, :right, :bottom, delta=delta/2)
        point(0, b, " B:(0,b)", :magenta, :left, :bottom, delta=delta/2)
        point(0, 0, " C:(0,0)", :magenta, :left, delta=-delta/2)
        point(a, 0, "D:(a,0)", :magenta, :right, delta=-delta/2)
        circle(a, 0, a/20, beginangle=180-θ, endangle=180)
        circle(a, 0, a/18, beginangle=180-θ, endangle=180)
        point(a - a/17, 0, "θ", :red, :right, :bottom, delta=delta/2, mark=false)
        ylims!(-5delta, b + 5delta)
    end
end;

draw(5, 3, true)

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

MacOS: ファイル・ディレクトリの変更日の時刻表示フォーマットの設定

2025年02月06日 | ブログラミング

MacOS Sequoia 15.x で,ディレクトリ表示中のファイル・ディレクトリの変更日の時刻表示が12時制になているのを24時制にする方法

ターミナルで 24時間制を強制する方法

1.ターミナルを開く(「アプリケーション」→「ユーティリティ」→「ターミナル」)

2.以下のコマンドを入力し、Enter を押す

defaults write NSGlobalDomain AppleICUForce24HourTime -bool true

3.Finder を再起動

killall Finder

または、ログアウト&ログインすることで設定を適用できます。

現在の設定を確認する

設定が適用されているか確認するには、以下のコマンドを実行してください。

defaults read NSGlobalDomain AppleICUForce24HourTime

これが 1 になっていれば 24時間制が有効になっています。

元に戻したい場合

12時間制に戻したい場合は、以下のコマンドを実行してください。

defaults write NSGlobalDomain AppleICUForce24HourTime -bool false
killall Finder
コメント
  • X
  • Facebookでシェアする
  • はてなブックマークに追加する
  • LINEでシェアする

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

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