算額（その230）

算額（その230）

(12) 京都市伏見区御香宮門前町 御香宮神社（ごこうのみやじんじゃ）文久3年(1863)
近畿数学史学会：近畿の算額「数学の絵馬を訪ねて」，平成4年5月16日 初版第一刷，大阪教育図書株式会社，大阪市.

絵馬堂が意外におもしろい！（3）　－ 御香宮神社 －
http://flattravel.blog.fc2.com/blog-entry-16.html

キーワード：円3個，長方形

大円の一部を欠き取る長方形内に，天円と地円を入れる。長方形の長辺と短辺がそれぞれ 4寸，3寸で，大円の径が 13寸のとき，天円と地円の径はいかほどか。

長方形の左下角の座標を (x, y) とおく。
大円の半径，中心座標を r0, (0, 0)
天円の半径，中心座標を r1, (x1, y1) # (x + r1, y + 3 - r1)
地円の半径，中心座標を r2, (x2, y2) # (x + 4 - r2, y + r2)
とおき，以下の方程式を解く。
注: x, x1, x2 は ::negative である。

using SymPy

@syms r0::positive, x::negative, y::positive,
   r1::positive, x1::negative, y1::positive,
   r2::positive, x2::negative, y2::positive;
r0 = 13//2
(x1, y1) = (x+r1, y + 3 - r1)
(x2, y2) = (x + 4 - r2, y + r2)
eq1 = x1^2 + y1^2 - (r0 + r1)^2
eq2 = x2^2 + y2^2 - (r0 - r2)^2
eq3 = x^2 + y^2 - r0^2
eq4 = (x + 4)^2 + (y + 3)^2 - r0^2

res = solve([eq1, eq2, eq3, eq4], (r1, r2, x, y))

   2-element Vector{NTuple{4, Sym}}:
    (4/5, 6/5, -28/5, 33/10)
    (36, 6/5, -28/5, 33/10)

最初の解が妥当なものである。すなわち天円の径は 8/5 = 1.6 つまり 1寸6分，地円の径は 12/5 = 2.4 つまり 2寸4分である。

   r1 = 0.8000000;  r2 = 1.2000000;  x = -5.6000000;  y = 3.3000000
   天円径 = 1.6000000,  地円径 = 2.4000000

using Plots
using Printf

function circle(ox, oy, r, color=:red; beginangle=0, endangle=360)
   θ = beginangle:0.1:endangle
   x = r.*cosd.(θ)
   y = r.*sind.(θ)
   plot!(ox .+ x, oy .+ y, color=color, linewidth=0.5)
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")
   (r1, r2, x, y) = res[1]
   r0 = 13//2
   (x1, y1) = (x+r1, y + 3 - r1)
   (x2, y2) = (x + 4 - r2, y + r2)
   @printf("r1 = %.7f;  r2 = %.7f;  x = %.7f;  y = %.7f\n", r1, r2, x, y)
   @printf("天円径 = %.7f,  地円径 = %.7f\n", 2r1, 2r2)
   plot([x, x+4, x+4, x, x], [y, y, y+3, y+3, y], color=:blue, lw=0.5)
   circle(0, 0, r0)
   circle(x1, y1, r1, :green)
   circle(x2, y2, r2, :gold)
   if more == true
       point(x1, y1, "(x1,y1)", :green, :center, :top)
       point(x1, y1, "天円", :green, :center, :bottom)
       point(x2, y2, "(x2,y2)", :gold, :center, :top)
       point(x2, y2, "地円", :gold, :center, :bottom)
       point(x, y, "(x,y)", :blue)
       point(x+4, y+3, "(x+4,y+3)", :blue, :right, :bottom)
       point(0, r0, " r0", :red)
       hline!([0], color=:black, lw=0.5)
       vline!([0], color=:black, lw=0.5)
   else
       plot!(showaxis=false, xlims=(-6, -0.5), ylims=(3.1, 6.35))
   end
end;