裏 RjpWiki

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

算額(その166)

2023年03月16日 | Julia

算額(その166)

岐阜県大垣市西外側町 大垣八幡神社 天保年間
http://ryugen3.sakura.ne.jp/toukou3/wasankibousya.PDF

第17問: 長方形内に半円(萌黄)と円弧を描き,その中に赤円と黒円を入れる。赤円と黒円の直径を知って長方形の長い方の一辺の長さを求めよ。

図のように記号を定め,連立方程式を解く。

using SymPy

@syms  r0::positive, r1::positive, r2::positive, r3::positive, x2::positive, x::positive;

r2 = 24
r3 = 17
eq1 = x2^2 + (r0 - 2r1 + r2)^2 - (r0 - r2)^2
eq2 = (x - r3)^2 + (r0 - r3)^2 - (r0 + r3)^2
eq3 = x2^2 + (r1 - r2)^2 - (r1 + r2)^2
eq4 = (r0 - 2r1)^2 + x^2 - r0^2;

res = solve([eq1, eq2, eq3, eq4], (r0, r1, x2, x))

   2-element Vector{NTuple{4, Sym}}:
    (289*sqrt(102)/49 + 11849/196, -17*sqrt(186001 + 19528*sqrt(102))/392 + 289*sqrt(102)/98 + 11849/392, sqrt(-204*sqrt(186001 + 19528*sqrt(102))/49 + 13872*sqrt(102)/49 + 142188/49), 34*sqrt(102)/7 + 408/7)
    (289*sqrt(102)/49 + 11849/196, 17*sqrt(186001 + 19528*sqrt(102))/392 + 289*sqrt(102)/98 + 11849/392, sqrt(204*sqrt(186001 + 19528*sqrt(102))/49 + 13872*sqrt(102)/49 + 142188/49), 34*sqrt(102)/7 + 408/7)

赤円,黒円の径が 24, 17 のとき,長方形の長い方の一辺の長さは (408 + 34√102)/7 = 107.34045255775867 である。
術文では,「黒径 / (1 - sqrt(黒径 / 赤径))」 であるが,SymPy では 「(sqrt(赤径)*黒径^1.5 + 赤径*黒径)/(赤径 - 黒径)」までしか簡単にならない。

(408 + 34√102)/7, (sqrt(r2)*r3^1.5 + r2*r3)/(r2 - r3), r3 / (1 - sqrt(r3 / r2))

   (107.34045255775867, 107.34045255775865, 107.34045255775871)

using Plots

function circle(ox, oy, r, color=:red; beginangle=0, endangle=360, fill=false)
   θ = beginangle:0.1:endangle
   x = r.*cosd.(θ)
   y = r.*sind.(θ)
   if fill
       plot!(ox .+ x, oy .+ y, linecolor=color, linewidth=0.5, seriestype=:shape, fillcolor=color)
   else
       plot!(ox .+ x, oy .+ y, color=color, linewidth=0.5)
   end
end;

function point(x, y, string="", color=:green, position=:left, vertical=:top; mark=true)
   mark && scatter!([x], [y], color=color, markerstrokewidth=0)
   annotate!(x, y, text(string, 10, position, color, vertical))
end;

function draw(more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   (r2, r3) = (24, 17)
   (r0, r1, x2, x) = res[1]
   r0 < 2r1 && ((r0, r1, x2, x) = res[2])
   println("r0 = $r0;  r1 = $r1;  x2 = $x2\nx = $x = $(x.evalf())")
   plot()
   circle(0, 0, r0, beginangle=0, endangle=90)
   circle(0, r0 - r1, r1, :olivedrab1, beginangle=270, endangle=450, fill=true)
   circle(x2, r0 - 2r1 + r2, r2, :indianred1, fill=true)
   circle(x - r3, r0 - r3, r3, :snow4, fill=true)
   plot!([0, x, x, 0, 0], [r0 - 2r1, r0 - 2r1, r0, r0, r0 - 2r1], color=:black, lw=0.5)
   if more == true
       point(0, r0, " r0", :black, :left)
       point(0, r0 - r1, " r0-r1", :black, :left)
       point(0, r0 - 2r1, " r0-2r1", :black, :left)
       point(x2, r0 - 2r1 + r2, "(x2,r0-2r1+r2)", :black, :top)
       point(x - r3, r0 - r3, "(x-r3,r0-r3)", :black, :top)
       point(x, r0 - 2r1, " x", :black)
       hline!([0], color=:black, lw=0.5)
       vline!([0], color=:black, lw=0.5)
   else
       plot!(showaxis=false)
       plot!(xlims=(0,112), ylims=(50, 125))
   end
end;

draw(false)

   r0 = 289*sqrt(102)/49 + 11849/196;  r1 = -17*sqrt(186001 + 19528*sqrt(102))/392 + 289*sqrt(102)/98 + 11849/392;  x2 = sqrt(-204*sqrt(186001 + 19528*sqrt(102))/49 + 13872*sqrt(102)/49 + 142188/49)
   x = 34*sqrt(102)/7 + 408/7 = 107.340452557759

 


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

コメントを投稿

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

Julia」カテゴリの最新記事