裏 RjpWiki

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

算額(その1185)

2024年08月05日 | Julia

算額(その1185)

(10) 京都市中京区三条大宮西二筋目下ル 武信稲荷神社 嘉永6年(1853)
近畿数学史学会:近畿の算額「数学の絵馬を訪ねて」,平成4年5月16日 初版第一刷,大阪教育図書株式会社,大阪市.
キーワード:円4個,正方形,等脚台形

正方形内に,水平線・斜線2本(等脚台形),日円 1 個,月円 1 個,星円 2 個をいれる。日円,月円の直径が 1 寸,4 寸のとき,星円の直径はいかほどか。

正方形の一辺の長さを 2a; 2a = 2r1 + 2r2
日円の半径と中心座標を r1, (0, 2r2 + r1)
月円の半径と中心座標を r2, (0, r2)
星円の半径と中心座標を r3, (a - r3, 2r2 - r3)
とおき,以下の連立方程式を解く。

include("julia-source.txt");

using SymPy
@syms a::positive, b::positive, r1::positive, r2::positive, r3::positive
a = r1 + r2
eq1 = (a - b) + 2r2 - sqrt((a - b)^2 + 4r2^2) - 2r3
eq2 = dist2(a, 0, b, 2r2, 0, r2, r2)
res = solve([eq1, eq2], (r3, b))[1]

   ((r1^2/2 + 2*r1*r2 + r2^2 - sqrt(4*r2^2*(r1 + r2)^2 + (r2^2 - (r1 + r2)^2)^2)/2)/(r1 + r2), r2^2/(r1 + r2))

res[1]

「術」は,「日円と月円の直径の積を,日円と月円の直径の和で割る」と実に簡潔である。
(1*4)/(1 + 4) = 4/5 = 0.8

sqrt(...) は,r1^2 + 2*r1*r2 + 2*r2^2 と簡約化できるので,それを代入して,simplify すれば,術と同じ式になる。

sqrt(4*r2^2*(r1 + r2)^2 + (r2^2 - (r1 + r2)^2)^2) |> factor |> println

   r1^2 + 2*r1*r2 + 2*r2^2

(r1^2/2 + 2*r1*r2 + r2^2 - (r1^2 + 2*r1*r2 + 2*r2^2)/2)/(r1 + r2)

日円,月円の直径が 1 寸,4 寸のとき,星円の直径は 1*4/(1 + 4) = 0.8 である。

function draw(r1, r2, more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   (r3, b) = ((r1^2 + 4*r1*r2 - r1*sqrt((4*r2^2*(r1 + r2)^2 + (r2^2 - (r1 + r2)^2)^2)/(r1 + r2)^2) + 2*r2^2 - r2*sqrt((4*r2^2*(r1 + r2)^2 + (r2^2 - (r1 + r2)^2)^2)/(r1 + r2)^2))/(2*(r1 + r2)), r2^2/(r1 + r2))
   a = r1 + r2
   @printf("日円,月円の直径が %g, %g のとき,星円の直径は %g である。\n", 2r1, 2r2, 2r3)
   @printf("r1 = %g;  r2 = %g;  r3 = %g;  2a = %g, b = %g\n", r1, r2, r3, 2a, b)
   plot([a, a, -a, -a, a], [0, 2a, 2a, 0, 0], color=:green, lw=0.5)
   circle(0, 2a - r1, r1)
   circle(0, r2, r2, :blue)
   circle2(a - r3, 2a - 2r1 - r3, r3, :magenta)
   segment(-a, 2a - 2r1, a, 2a - 2r1, :orange)
   segment(a, 0, b, 2a - 2r1, :orange)
   segment(-a, 0, -b, 2a - 2r1, :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(0, 2a - r1, "日円:r1\n(0,2a-r1)", :red, :center, delta=-delta/2)
       point(0, r2, "月円:r1,(0,r2)", :blue, :center, delta=-delta/2)
       point(a - r3, 2r2 - r3, "星円:r3,(a-r3,2r2-r3) ", :magenta, :right, :vcenter)
       point(b, 2r2, "(b,2r2)", :black, :center, :bottom, delta=delta/2)
   end
end;

draw(1/2, 4/2, true)

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

算額(その1184)

2024年08月05日 | Julia

算額(その1184)

(9) 安井金比羅宮
近畿数学史学会:近畿の算額「数学の絵馬を訪ねて」,平成4年5月16日 初版第一刷,大阪教育図書株式会社,大阪市.
キーワード:階差数列

あらたまの春し来ぬれば,おのづから空もうららに少女子が手毬つく音を数うれば,
はじめ,ひふみよいつむななやつここのつと聞こえたり,
二度目は二十三にして,三度は四十五ぞと云う,
四度は七十五なるぞと,
次第に増して終わりに成りて三千四百ふたそ五つまで続けつきたる手振りこそ幾たびならん聞かまほしけれ

意訳すれば,手毬をついた。最初は 9 回,2回目は 23 回,3回目は 45 回,4回目は 75 回ついた。だんだん増えていき,終いには 3425 回もついた!
3425 回ついたのは何回目か?

n 回目の手毬をついた回数 a[n] の第 2 階差が定数 8 なので,a[n] は n の二次式で表すことができる。
a[n] = α + n*β + n^2*γ とおき,α,β,γ を決定する。

include("julia-source.txt");

using SymPy
@syms α, β, γ, n
eq1 = α + β + γ - 9
eq2 = α + 2β + 4γ - 23
eq3 = α + 3β + 9γ - 45
# eq4 = α + 4β + 16γ - 75

res = solve([eq1, eq2, eq3], (α, β, γ))
res |> println

   Dict{Sym{PyCall.PyObject}, Sym{PyCall.PyObject}}(α => 3, γ => 4, β => 2)

eq = 3 + 2n + 4n^2 - 3425
res = solve(eq, n)[2]  # 2 of 2
res |> println

   29

n = 29
3 + 2n + 4n^2

   3425

一般化すると,

@syms n,  m
eq02 = 3 + 2n + 4n^2 - m
res2 = solve(eq02, n)[2]  # 2 of 2

「術」は,これを,「其数をよつにかえしてたいらかにひらきつきぬはすててこそしれ」と書いている。
「その数を 4 で割って,平方根を求め,整数部を取れ」ということだ。

sqrt(3425/4)

   29.261749776799064

n = 1:29
a = @. 3 + 2n + 4n^2;
plot(n, a, label="")
scatter!(n, a, label="")

 

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

算額(その1183)

2024年08月05日 | Julia

算額(その1183)

(8) 京都府宮津市天橋立文殊 智恩寺文殊堂 天保8年(1837)
近畿数学史学会:近畿の算額「数学の絵馬を訪ねて」,平成4年5月16日 初版第一刷,大阪教育図書株式会社,大阪市.
キーワード:円5個,正五角形

正五角形の中に,正五角形の中心を通り辺に接する円を 5 個描く。円が交わってできる中央部の円弧で囲まれた細長い部分の面積(図の黄色で示した部分)を求めよ。

正五角形の一辺の長さを a とする。
正五角形の中心から辺までの長さを h とする(h = BD)。
求める面積は,「扇形 CAB から三角形 ABC の面積を引いたものの 2*5 倍」である。

include("julia-source.txt");

using SymPy
@syms a, h
a = 10
eq = h*tand(Sym(36)) - (a//2)
ans_h = solve(eq, h)[1]
h2 = ans_h/2
h2 |> println

   a/(4*sqrt(5 - 2*sqrt(5)))

三角形 ABC の面積

S1 = h2*cosd(Sym(18))*h2*sind(Sym(18))
S1 |> println

   a^2*(-1/4 + sqrt(5)/4)*sqrt(sqrt(5)/8 + 5/8)/(16*(5 - 2*sqrt(5)))

扇形 CAB の面積 = (半径 h2 = BC の円の面積の 1/10)

∠ACB = 36°

S2 = PI*h2^2*(36//360)
S2 |> println

   pi*a^2/(160*(5 - 2*sqrt(5)))

S = 10(S2 - S1)
S |> println

   -5*a^2*(-1/4 + sqrt(5)/4)*sqrt(sqrt(5)/8 + 5/8)/(8*(5 - 2*sqrt(5))) + pi*a^2/(16*(5 - 2*sqrt(5)))

正五角形の一辺の長さが 10 のとき,求める面積は 2.39960452467531 である。

S(a => 10).evalf() |> println

   2.39960452467531

点 A の 座標は以下のとおり。

@syms x0::positive, y0::negative
h2 = 1
eq1 = x0^2 + (y0 + h2)^2 - h2^2
eq2 = (x0 - h2*cosd(54))^2 + (y0 - h2*sind(54))^2 - h2^2
res = solve([eq1, eq2], (x0, y0))

   1-element Vector{Tuple{Sym{PyCall.PyObject}, Sym{PyCall.PyObject}}}:
    (0.587785252292474, -0.190983005625053)

function fill(h2)
   θ = range(54, 90, 10)
   x = @. h2*cosd(θ)
   y = @. h2*sind(θ)
   cx = (x[1] + x[end])/2
   cy = (y[1] + y[end])/2
   for i in 1:length(x)
       append!(x, 2cx - x[i])
       append!(y, 2cy - y[i])
   end
   append!(x, x[1])
   append!(y, y[1])
   y = @. y - h2
   for i = 1:5
       plot!(x, y, seriestype=:shape, linecolor=:yellow, fillcolor=:yellow, color=:yellow)
       (x, y) = transform(x, y, deg=72)
   end
end;

function draw(a, more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   r = (a/2)/sind(36)
   h = (a/2)/tand(36)
   h2 = h/2
   θ = 18:72:378
   x = @. r*cosd(θ)
   y = @. r*sind(θ)
   plot(x, y, color=:blue, lw=0.5)
   circle(0, 0, r, :green)
   fill(h2)
   rotate(0, -h2, h2, angle=72)
   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.588h2, -0.191h2, "A", :black, :left, :bottom, delta=delta/2)
       point(0, 0, "B  ", :black, :right, :bottom, delta=delta/2)
       point(0, -h2, "C", :black, :center, delta=-delta/2)
       plot!([0, 0, 0.588h2, 0], [-h2, 0, -0.192h2, -h2], color=:black, lw=0.5)
       point(0, -h, "D", :blue, :center, delta=-delta/2)
       point(x[4], y[4], "E ", :blue, :right, delta=-delta/2)
   end
end;

draw(1, true)

 

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

算額(その1182)

2024年08月05日 | Julia

算額(その1182)

(7) 京都府宮津市天橋立文殊 智恩寺文殊堂 文政12年(1829)
近畿数学史学会:近畿の算額「数学の絵馬を訪ねて」,平成4年5月16日 初版第一刷,大阪教育図書株式会社,大阪市.
キーワード:円13個

台風のため破損して肝心な部分が読み取れないとのことである。
「問」の「如圖有□□重圓 只云赤青黄黒之四色(以下 8 行判読できず」,「答」の「圓径二尺五寸(かなりの部分不明)」,「術」の「廻リヲ圓径一尺(1 行不明)」という断片から,重なり合った 7 個の大円と 6 個の小円があるとき,小円の直径が 1 尺のとき,大円の直径はいかほどか」というような内容と思うが,「愚術曰四色接合之和積ヲ...」から始まる長い文章があるので,円の重なった部分を 4 色色分けされており,そのいずれかの色で塗られた面積を求めるというような複雑な内容であったとも思われる。

取り合えず,図だけを再現する。



大円の半径と中心座標を r1, (r1*cos(π/6), r1*sin(π/6))
小円の半径と中心座標を r2, (0, 2r1 - r2)
とおき,以下の方程式を解く。

include("julia-source.txt");

using SymPy
@syms r1, r2
(y1, x1) = r1.*sincos(PI/6)
eq1 = x1^2 + (2r1 - r2 -y1)^2 - (r2 + r1)^2
res = solve(eq1, r2)[1]
res |> println

   2*r1/5

小円は大円の 2/5 倍である。
小円が 1 尺のとき,大円は 2 尺 5 寸である。

function draw(r1, more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   r2 = 2r1/5
   plot()
   circle(0, 0, r1)
   rotate(0, r1, r1, angle=60)
   rotate(0, 2r1 - r2, r2, :blue, angle=60)
   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)
       (y1, x1) = r1.*sincos(pi/6)
       point(x1, y1, "(x1,y1) ", :red, :right, :bottom, delta=delta/2)
       point(0, 2r1 - r2, "2r1-r2", :blue, :center, delta=-delta/2)
   end
end;

draw(1, true)

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

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

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