裏 RjpWiki

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

算額(その170)

2023年03月18日 | Julia

算額(その170)

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

第29問: 直行するに直線に接し互いに外接する 3 個の等円(赤)を作り,青円と黒円を置く。青円の直径を知って黒円の直径を求めよ。

赤円,青円,黒円の直径をそれぞれ r1, r2, r3 として,更に,青円,黒円の中心の y 座標を y2, y3 として連立方程式を解く。なお,2 個の青円の中心が x 軸になるように座標を定めないと,連立方程式の係数が複雑になり解けないので注意。

using SymPy

@syms r1::positive, r2::positive, r3::positive, y2::positive, y3::positive;

eq1 = (r1 - r3)^2 + y3^2 - (r1 + r3)^2
eq3 = (2r1 - r3)^2 + (sqrt(Sym(3))r1 - y3)^2 - (r1 + r3)^2
eq4 = (r2 - r3)^2 + (y2 - y3)^2 - (r2 + r3)^2
eq5 = (2r1 - r2)^2 + (y2 - sqrt(Sym(3))r1)^2 - (r1 + r2)^2;

res = solve([eq1, eq3, eq4, eq5], (r1, r3, y2, y3))

   1-element Vector{NTuple{4, Sym}}:
    (4*sqrt(2)*r2/9 + r2, -2*sqrt(2)*r2 + 11*r2/3, 4*sqrt(6)*r2/9 + 16*sqrt(3)*r2/9, -49*sqrt(3)*r2*(116/441 - 580*sqrt(2)/441)/58)

黒円の半径は -2*sqrt(2)*r2 + 11*r2/3 である。
つまり,「黒円の半径 = (11 - 6√2)/3 × 青円の半径」である。

-2*sqrt(Sym(2))*r2 + 11*r2/3 |> simplify |> println

   r2*(11 - 6*sqrt(2))/3

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 = 1  # 青円の半径を 1 として,左右反転した図を描く
   (r1, r3, y2, y3) = (4*sqrt(2)*r2/9 + r2, -2*sqrt(2)*r2 + 11*r2/3, 4*sqrt(6)*r2/9 + 16*sqrt(3)*r2/9, -49*sqrt(3)*r2*(116/441 - 580*sqrt(2)/441)/58)
   plot()
   circle(r1, 0, r1, :indianred1, fill=true)
   circle(3r1, 0, r1, :indianred1, fill=true)
   circle(2r1, sqrt(3)r1, r1, :indianred1, fill=true)
   circle(r2, y2, r2, :steelblue1, fill=true)
   circle(r3, y3, r3, :snow4, fill=true)
   hline!([0], color=:black, lw=0.5)
   vline!([0], color=:black, lw=0.5)
   if more == true
       point(r2, y2, "(r2,y2)", :black)
       point(r3, y3, "(r3,y3)", :black)
       point(r1, 0, "r1", :black)
       point(2r1, sqrt(3)r1, "(2r1,√3r1)", :black)
   end
end;

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

算額(その169)

2023年03月18日 | Julia

算額(その169)

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

第28問: 互いに中心を通る 2 等円(黃)に 3 個の等円(青)を入れ,更に赤円 2 個,黒円 4 個を入れる。青円径を知って黒円径を求めよ。

算額(その85を一段階複雑にしたものである。

青円,赤円,黒円の半径を r1, r2, r3,黒円の中心座標を (x3, y3) とする。黃円の半径は 2r1 である。未知数 5 個に対して方程式は 4 本なので,r2, r3, x3, y3 が r1 を含む解が求まる。

using SymPy

@syms r1::positive, r2::positive, r3::positive, x3::positive, y3::positive;

eq1 = x3^2 + y3^2 - (r1 + r3)^2  # 青円と黒円が外接
eq2 = x3^2 + (r1 + r2 - y3)^2 - (r2 + r3)^2  # 赤円と黒円が外接
eq3 = r1^2 + (r1 + r2)^2 - (2r1 - r2)^2  # 赤円が黃円に内接
eq4 = (r1 + x3)^2 + y3^2 - (2r1 - r3)^2  # 黒円が黃円に内接
res = solve([eq1, eq2, eq3, eq4], (r2, r3, x3, y3))

   1-element Vector{NTuple{4, Sym}}:
    (r1/3, 2*r1/11, 5*r1/11, 12*r1/11)

黒円の半径は青円の半径の 2/11 である。
赤円の半径は青円の半径の 1/3 である。
黒円の中心座標は青円の半径の (5/11, 12/11) である。

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(R, more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   r1 = 1  # 青円の半径を 1 として図を描く
   (r2, r3, x3, y3) = (r1/3, 2*r1/11, 5*r1/11, 12*r1/11)
   plot()
   circle(r1, 0, 2r1, :yellow1, fill=true)
   circle(-r1, 0, 2r1, :yellow1, fill=true)
   circle(r1, 0, 2r1, :snow4)
   circle(-r1, 0, 2r1, :snow4)
   circle(0, 0, r1, :steelblue1, fill=true)
   circle(2r1, 0, r1, :steelblue1, fill=true)
   circle(-2r1, 0, r1, :steelblue1, fill=true)
   circle(0, r1 + r2, r2, :indianred1, fill=true)
   circle(0, -r1 - r2, r2, :indianred1, fill=true)
   circle(x3, y3, r3, :snow4, fill=true)
   circle(-x3, y3, r3, :snow4, fill=true)
   circle(x3, -y3, r3, :snow4, fill=true)
   circle(-x3, -y3, r3, :snow4, fill=true)
   if more == true
       point(0, r1, " r1", :black)
       point(0, r1+r2, " r1+r2", :black, :left, :bottom)
       point(x3, y3, "(x3,y3,r3)", :black, :left, :bottom)
       hline!([0], color=:black, lw=0.5)
       vline!([0], color=:black, lw=0.5)
   else
       plot!(showaxis=false)
   end
end;

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

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

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