裏 RjpWiki

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

算額(その1305)

2024年09月20日 | Julia

算額(その1305)

百十六 群馬県吾妻郡吾妻町金井 一宮神社 明治5年(1872)
群馬県和算研究会:群馬の算額,上武印刷株式会社,高崎市,1987年3月31日.
キーワード:円1個,四分円,斜線

四分円と,四分円の端を通る傾きが -1 の斜線と,直線に囲まれた領域に円を容れる。四分円の半径が 0.5 寸のとき,円の直径はいかほどか。

注:斜線と四分円の交点の y 座標を「高」と呼んでいる。四分円の半径と同じ長さである。また,斜線を単に「斜」と呼んでいるが,算額の図からもわかるように,傾き -1 の斜線である。

四分円の半径と中心座標を r1, (-r1, r1)
円の半径と中心座標を r2, (x2, r2)
とおき,以下の連立方程式を解く。

include("julia-source.txt");

using SymPy
@syms r1::positive, r2::positive, x2::positive
eq1 = dist2(0, r1, r1, 0, x2, r2, r2)
eq2 = (x2 + r1)^2 + (r1 - r2)^2 - (r1 + r2)^2
res = solve([eq1, eq2], (r2, x2))[3]  # 3 of 3

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

res[1] |> simplify |> x -> x/r1 |> expand |> x -> r1*x |> println

   r1*(6 - 4*sqrt(2))

円の半径は四分円の半径の (6 - 4√2) 倍である。
四分円の半径が 0.5 寸のとき,円の直径は 2*0.5*(6 - 4√2) = 0.3431457505076194 寸である。

2 * 0.5*(6 - 4√2)

   0.3431457505076194

function draw(r1, more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   (r2, x2) = ((r1*(3 - 2*sqrt(2)) + r1)^2/(4*r1), r1*(3 - 2*sqrt(2)))
   plot()
   circle(-r1, r1, r1, beginangle=270, endangle=360)
   segment(0, r1, r1, 0, :green)
   circle(x2, r2, r2, :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(0, r1, "r1", :red, :center, :bottom, delta=delta)
       point(-r1, 0, "-r1", :red, :center, delta=-delta)
       point(r1, 0, "r1", :green, :center, delta=-delta)
       point(-r1, r1, "四分円の中心", :red, :left, delta=-delta)
       point(x2, r2, "r2,(x2,r2)", :blue, :center, delta=-delta)
       ylims!(-6delta, r1 + 4delta)
   end
end;

draw(5, true)

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

算額(その1304)

2024年09月20日 | Julia

算額(その1304)

百十五 群馬県富岡市一宮 貫前神社 明治4年(1871)
群馬県和算研究会:群馬の算額,上武印刷株式会社,高崎市,1987年3月31日.
キーワード:円1個,正三角形,正五角形

正三角形の中に円,長方形,五角形を容れる。円の直径が 10 寸のとき,正五角形の一辺の長さはいかほどか。

1. 円の半径を r とおく。IE = r
2. 円を内包する正三角形の一辺の長さを 2a とおく。EF = a = √3r
3. 正五角形を内包する円の半径を R とおく。AE = AB = AD = R
4. EF = AB*cosd(18°) = R*cosd(18°) = a = √3r
5. R*cosd(18°) = √3r より,R = √3r/cosd(18°)
6. 正五角形の一辺の長さ = 2CD = 2*R*sind(36°)
7. 2*R*sind(36°) = 2 * (√3r/cosd(18°)) * sind(36°) 
8 r = 10/2 のとき s = 2 * (√3(10/2)/cosd(18)) * sind(36) = 10.704662693192699
9. 図を描くために必要な正三角形の一辺の長さを 2b とすると,相似関係から CH = b = a + (yo + R)/√3
以上により,円の直径が 10 寸のとき,正五角形の一辺の長さは 10.704662693192699 寸である。

2 * (√3(10/2)/cosd(18)) * sind(36) 

   10.704662693192699

s は三角関数を含むが,cosd(18), sind(36) はルートを含む式で表せるので,以下のように簡略化される。
s = 5(√15 - √3) = 10.704662693192699

include("julia-source.txt");

using SymPy
@syms d, r
s = 2r * (√Sym(3)/cosd(Sym(18))) * sind(Sym(36)) |> simplify
s/r |> x -> apart(x, d) |> simplify |> x -> x*r |> println

   r*(-sqrt(3) + sqrt(15))

円の直径が 10 寸のとき,正五角形の一辺の長さは 5(√15 - √3) = 10.704662693192699。

function draw(r, more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   a = √3r
   R = a/cosd(18)
   xo = R*sind(36)
   yo = R*cosd(36)
   plot()
   circle(0, yo + R + r, r, :blue)
   plot!([a, 0, -a, a], (yo + R) .+ [0, √3a, 0, 0], color=:green, lw=0.5)
   θ = 18:72:378
   x = R.*cosd.(θ)
   y = R.*sind.(θ)
   circle(0, yo, R, :pink)
   plot!(x, yo .+ y, color=:red, lw=0.5)
   plot!([a, a, -a, -a, a], [0, yo + R, yo + R, 0, 0], color=:orange, lw=0.5)
   b = a + (yo + R)/√3
   plot!([b, 0, -b, b], [0, yo + R + √3a, 0, 0], color=:black, 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(0, yo, " A", :green, :left, :vcenter)
       point(x[1], yo + y[1], " B", :green, :left, :vcenter)
       point(0, 0, "C", :green, :center, delta=-delta)
       point(R*sind(36), 0, "D", :green, :center, delta=-delta)
       point(0, yo + R, "E", :green, :center, delta=-delta)
       point(R*cosd(18), yo + R, " F", :green, :left, :vcenter)
       point(a, 0, "G", :green, :left, delta=-delta)
       point(b, 0, "H", :green, :left, delta=-delta)
       point(0, yo + R + r, " I", :green, :left, :vcenter)
       point(0, yo + R + √3a, " J", :green, :left, :vcenter)
       point(0, yo + R + r, "")
   end
end;

draw(10/2, true)

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

算額(その1303)

2024年09月20日 | Julia

算額(その1303)

百十五 群馬県富岡市一宮 貫前神社 明治4年(1871)
群馬県和算研究会:群馬の算額,上武印刷株式会社,高崎市,1987年3月31日.
キーワード:円2個,外円,斜線,弦2本

外円の中に等斜を 2 本引き,全円を容れる。外円の直径が 10 寸,斜の長さが 8 寸のとき,全円の直径はいかほどか。

外円の半径と中心座標を R, (0, 0)
全円の半径と中心座標を r, (0, R - r)
斜と外円の交点座標を (0, -R), (x, sqrt(R^2 - x^2)
斜の長さを 「斜」
とおき,以下の連立方程式を解く。

include("julia-source.txt");

using SymPy
@syms R::positive, r::positive, x::positive,
     斜::positive;
eq1 = (x^2 + (R - sqrt(R^2 - x^2))^2) - 斜^2
eq2 = dist2(0, -R, x, sqrt(R^2 - x^2), 0, R - r, r)
res = solve([eq1, eq2], (r, x))[1]  # 1 of 2

   (2*R*斜*(2*R*sqrt(128*R^8 - 128*R^6*斜^2 + 64*R^6*sqrt(4*R^4 - 4*R^2*斜^2 + 斜^4) + 48*R^4*斜^4 - 32*R^4*斜^2*sqrt(4*R^4 - 4*R^2*斜^2 + 斜^4) - 6*R^2*斜^6 + 8*R^2*斜^4*sqrt(4*R^4 - 4*R^2*斜^2 + 斜^4) - 斜^6*sqrt(4*R^4 - 4*R^2*斜^2 + 斜^4)) - 斜*(2*R - 斜)*(2*R + 斜)*(2*R^2 + sqrt(4*R^4 - 4*R^2*斜^2 + 斜^4)))/(32*R^6 - 24*R^4*斜^2 + 16*R^4*sqrt(4*R^4 - 4*R^2*斜^2 + 斜^4) + 6*R^2*斜^4 - 4*R^2*斜^2*sqrt(4*R^4 - 4*R^2*斜^2 + 斜^4) + 斜^4*sqrt(4*R^4 - 4*R^2*斜^2 + 斜^4)), 斜*sqrt(4*R^2 - 斜^2)/(2*R))

全円の半径を表す式は長くなり,SymPy では簡約化できないが,式は正しい。
外円の直径が 10 寸,斜が 8 寸のとき,全円の直径は 15/2 = 7.5 寸である。

2res[1](R => 10//2, 斜 => 8) |> println

   15/2

function draw(R, 斜, more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   (r, x) = (2*R*斜*(2*R*sqrt(128*R^8 - 128*R^6*斜^2 + 64*R^6*sqrt(4*R^4 - 4*R^2*斜^2 + 斜^4) + 48*R^4*斜^4 - 32*R^4*斜^2*sqrt(4*R^4 - 4*R^2*斜^2 + 斜^4) - 6*R^2*斜^6 + 8*R^2*斜^4*sqrt(4*R^4 - 4*R^2*斜^2 + 斜^4) - 斜^6*sqrt(4*R^4 - 4*R^2*斜^2 + 斜^4)) - 斜*(2*R - 斜)*(2*R + 斜)*(2*R^2 + sqrt(4*R^4 - 4*R^2*斜^2 + 斜^4)))/(32*R^6 - 24*R^4*斜^2 + 16*R^4*sqrt(4*R^4 - 4*R^2*斜^2 + 斜^4) + 6*R^2*斜^4 - 4*R^2*斜^2*sqrt(4*R^4 - 4*R^2*斜^2 + 斜^4) + 斜^4*sqrt(4*R^4 - 4*R^2*斜^2 + 斜^4)), 斜*sqrt(4*R^2 - 斜^2)/(2*R))
   @printf("外円の直径が %g,等斜の長さが %g のとき,全円の直径は %g である。\n", 2R, 斜, 2r)
   y = sqrt(R^2 - x^2)
   plot([-x, 0, x], [y, -R, y], color=:green, lw=0.5)
   circle(0, 0, R, :blue)
   circle(0, R - r, 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, R - r, "全円:r,(0,R-r)", :red, :center, delta=-delta/2)
       point(x, sqrt(R^2 - x^2), "(x,sqrt(R^2-x^2)) ", :green, :right, :vcenter)
       point(0, R, "R", :blue, :center, :bottom, delta=delta/2)
   end
end;

draw(10/2, 8, true)

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

算額(その1302)

2024年09月20日 | Julia

算額(その1302)

百八 群馬県邑楽郡板倉町板倉 雷電神社 慶応3年(1867)
群馬県和算研究会:群馬の算額,上武印刷株式会社,高崎市,1987年3月31日.
キーワード:円3個,外円,二等辺三角形,弦3本,

外円の中に水平な弦と長さの同じ斜めの弦 2 本(倒立した二等辺三角形)と等円 2 個を容れる。外円の直径が 5 寸,斜線(二等辺三角形の斜辺)の長さが 3 寸のとき,等円の直径はいかほどか。

外円の半径と中心座標を R, (0, 0)
等円の半径と中心座標を r, (r, r - sqrt(R^2 - x^2))
水平な弦と外円の交点座標を (x -sqrt(R^2 - x^2)
斜の長さを「斜」
とおき,以下の連立方程式を解く。

include("julia-source.txt");

using SymPy
@syms R::positive, r::positive, x::positive,
     斜::positive;
eq1 = sqrt(x^2 + (R - sqrt(R^2 - x^2))^2) - 斜
eq2 = r^2 + (r - sqrt(R^2 - x^2))^2 - (R - r)^2
res = solve([eq1, eq2], (r, x))[2]  # 2 of 2

   (-R + sqrt(2*R^2 - sqrt(4*R^4 - 4*R^2*斜^2 + 斜^4)) + sqrt(4*R^4 - 4*R^2*斜^2 + 斜^4)/(2*R), 斜*sqrt(4*R^2 - 斜^2)/(2*R))

等円の半径は -R + sqrt(2*R^2 - sqrt(4*R^4 - 4*R^2*斜^2 + 斜^4)) + sqrt(4*R^4 - 4*R^2*斜^2 + 斜^4)/(2*R) である。
SymPy では sqrt(4*R^4 - 4*R^2*斜^2 + 斜^4) をうまく処理できないので,その部分を手作業で簡約化した後で simplify すると,斜*(2R - 斜)/(2R) になる。

問では「外円之径五寸□分」と 1 文字欠損のように見えるが,直径がちょうど 5 寸のとき(すなわち欠損文字は零か?),等円の直径が 2 寸 4 分になる。

2res[1](R => 5.0/2, 斜 => 3) |> println

   2.40000000000000

2(斜*(2R - 斜)/(2R))(R => 5.0/2, 斜 => 3) |> println

   2.40000000000000

function draw(R, 斜, more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   r = 斜*(2R - 斜)/(2R)
   x = 斜*sqrt(4R^2 - 斜^2)/(2R)
   @printf("外円の直径が %g,等斜の長さが %g のとき,等円の直径は %g である。\n", 2R, 斜, 2r)
   @printf("R = %g;  斜 = %g;  r = %g;  x = %g\n", R, 斜, r, x)
   plot([x, -x, 0, x], [-sqrt(R^2 - x^2), -sqrt(R^2 - x^2), -R, -sqrt(R^2 - x^2)], color=:green, lw=0.5)
   circle(0, 0, R, :blue)
   circle2(r, r - sqrt(R^2 - x^2), 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(r, r - sqrt(R^2 - x^2), "等円:r,(r,r-sqrt(R^2-x^2))", :red, :center, delta=-delta/2)
       point(x, -sqrt(R^2 - x^2), "(x,-sqrt(R^2-x^2))   ", :green, :right, delta=-delta/2)
       point(0, R, "R", :blue, :center, :bottom, delta=delta/2)
   end
end;

draw(5/2, 3, true)

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

算額(その1301)

2024年09月20日 | Julia

算額(その1301)

百八 群馬県邑楽郡板倉町板倉 雷電神社 慶応3年(1867)
群馬県和算研究会:群馬の算額,上武印刷株式会社,高崎市,1987年3月31日.
キーワード:円2個,直角三角形,正方形

直角三角形の中に正方形と中鈎(直角の頂点から斜辺への垂線)を入れる。区画された領域に,甲円と乙円を 1 つずつ容れる。甲円と乙円の直径がそれぞれ 8 寸,15 寸のとき,正方形の一辺の長さはいかほどか。

正方形の一辺の長さを a とする。
直角三角形の鈎と股(直角を挟む 2 辺のうちの短い方と長い方)を,「鈎」,「股」とする。
股は問題を解く上では必要ないが,図を描くために求める。

また,図に示すように,正方形と中鈎の延長線を補助線として,交点座標を (股 - a, a + x)
甲円と乙円を含む直角三角形は相似で,相似比が a:(鈎 - a) である。
甲円の半径と中心座標を r1, (股 - r1, a + r1)
乙円の半径と中心座標を r2, (股 - a + r2, r2)
とおき,以下の連立方程式を解く。

include("julia-source.txt");

using SymPy
@syms a::positive, x::positive,
     鈎::positive, 股::positive,
     r1::positive, r2::positive;
eq1 = a/(鈎 - a) - (a + x)/a
eq2 = (鈎 - a) + a - sqrt((鈎 - a)^2 + a^2) - 2r1
eq3 = a + (a + x) - sqrt(a^2 + ( a + x)^2) - 2r2
eq4 = (鈎 - a)/a - 鈎/股;

SymPy の能力的に,一度には解けないので逐次的に解いていく。

まず,eq1, eq2 を連立させて,鈎,股を求める。

res = solve([eq1, eq4], (鈎, 股))[1]

   (a*(2*a + x)/(a + x), 2*a + x)

鈎,股を eq2 に代入し,x を求める。

eq12 = eq2(鈎 => res[1], 股 => res[2]) |> simplify |> numerator |> expand
ans_x = solve(eq12, x)[1]
ans_x |> println

   a*(a^2 - 4*a*r1 + 2*r1^2)/(2*r1*(a - r1))

鈎,股,x を eq3 に代入し,a を求める。

eq13 = eq3(鈎 => res[1], 股 => res[2], x => ans_x)
ans_a = solve(eq13, a)[3]  # 3 of 3
ans_a |> println

   r1 + r2 + sqrt(r1^2 + r2^2)

正方形の一辺の長さは,r1 + r2 + sqrt(r1^2 + r2^2) である。
甲円と乙円の直径がそれぞれ 8 寸,15 寸のとき,正方形の一辺の長さは 20 寸である。

ans_a(r1 => 8/2, r2 => 15/2).evalf() |> println

   20.0000000000000

a が求まれば x が決まり,更に鈎,股も決まる。
a = r1 + r2 + sqrt(r1^2 + r2^2)
x = a*(a^2 - 4a*r1 + 2r1^2)/(2r1*(a - r1))
鈎 = a*(2a + x)/(a + x)
股 = 2a + x

すべてのパラメータは以下のとおりである。

   r1 = 4;  r2 = 7.5;  a = 20;  x = 17.5;  鈎 = 30.6667;  股 = 57.5

function draw(r1, r2, more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   a = r1 + r2 + sqrt(r1^2 + r2^2)
   x = a*(a^2 - 4a*r1 + 2r1^2)/(2r1*(a - r1))
   鈎 = a*(2a + x)/(a + x)
   股 = 2a + x
   @printf("甲円と乙円の直径がそれぞれ %g 寸,%g 寸のとき,正方形の一辺の長さは %g 寸である。\n", 2r1, 2r2,a)
   @printf("r1 = %g;  r2 = %g;  a = %g;  x = %g;  鈎 = %g;  股 = %g\n", r1, r2, a, x, 鈎, 股)
   plot([0, 股, 股, 0], [0, 0, 鈎, 0], color=:green, lw=0.5)
   rect(股 - a, 0, 股, a, :red)
   circle(股 - r1, a + r1, r1, :blue)
   circle(股 - a + r2, r2, r2, :magenta)
   segment(股 - a, a + x, 股, 0)
   segment(股 - a, a + x, 股 - a, a)
   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(股, 鈎, "(股,鈎) ", :green, :right, :bottom, delta=delta/2)
       point(股 - a, 0, "股-a ", :red, :right, :bottom, delta=delta)
       point(股 - a, a + x, "(股-a,a+x) ", :black, :right, delta=-delta/2)
       point(股 - r1, a + r1, "甲円:r1\n(股-r1,a+r1)", :black, :center, :bottom, delta=delta/2, deltax=4delta)
       point(股 - a + r2, r2, "乙円:r2\n(股-a+r2,r2)", :magenta, :center, delta=-delta/2)
       xlims!(-2delta, 股 + 5delta)
   end
end;

draw(8/2, 15/2, true)

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

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

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