算額（その1091）

算額（その1091）

四十三 岩手県一関市真滝 熊野白山滝神社 弘化3年(1846)
山村善夫：現存 岩手の算額，昭和52年1月30日，熊谷印刷，盛岡市. http://www.wasan.jp/yamamura/yamamura.html
キーワード：円1個，四分円2個，正方形

正方形の中に四分円を 2 個，円を 1 個容れる。円の直径が与えられたときに黒積（灰色の部分の 2 倍）を求めよ。

正方形の一辺の長さを R，円の半径を r とおき，以下の方程式を解く。

include("julia-source.txt")

using SymPy

@syms R::positive, r::positive
eq1 = (R/2)^2 + r^2 - (R - r)^2
res = solve(eq1, R)[1]
res[1] |> println

   8*r/3

正方形の一辺の長さは円の半径の 8/3 倍である。

黒積の半分（図の灰色の部分）は，
扇形 OAR の面積（半径　R の円の面積の 1/6 = πR^2/6）から，
⊿AOC の面積（√3R/2 * R/2 / 2 = √3R^2/8）と，
円の面積の半分（πr^2/2）を引いたものの，
2倍である。

r = 1/2 のとき，R = 4/3 となり，黒積は 0.30648601314366886 である。

r = 1/2
R = 4/3
2(pi*R^2/6 - R^2*√3/8 - pi*r^2/2)

   0.30648601314366886

術は，円径^2*1.5593 とあるが，おかしい?

(2r)^2*1.5593

   1.5593

function draw(more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   r = 1/2
   R = 8r/3
   黒積 = 2(pi*R^2/6 - R^2*√3/8 - pi*r^2/2)
   @printf("円の直径が %g(正方形の一辺の長さが %g) のとき，黒積は %g である。\n", 2r, R, 黒積)
   θ = 0:0.1:60
   x = @. R*cosd(θ)
   y = @. R*sind(θ)
   append!(x, [R/2, R])
   append!(y, [0, 0])
   plot([0, R, R, 0, 0], [0, 0, R, R, 0], color=:blue, lw=0.5)
   plot!(x, y, seriestype=:shape, fillcolor=:gray80)
   circlef(R/2, r, r, :white, beginangle=-90, endangle=90)
   circle(0, 0, R, beginangle=0, endangle=90)
   circle(R, 0, R, beginangle=90, endangle=180)
   circle(R/2, r, r, :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)
       segment(R/2, 0, R/2, √3R/2)
       segment(0, 0, R/2, √3R/2)
       point(R/2, √3R/2, "A", :black, :center, :bottom, delta=delta/2)
       point(R, 0, " R", :black, :left, :bottom, delta=delta/2)
       point(0, 0, "O ", :black, :right, :bottom, delta=delta/2)
       circle(0, 0, 0.05R,  :black, beginangle=0, endangle=60)
       circle(0, 0, 0.055R,  :black, beginangle=0, endangle=60)
       point(0.13R, 0.01R, "60°", :black, :right, :bottom, delta=delta/2, mark=false)
       point(R/2, r, " B", :black, :left, :vcenter)
       point(R/2, 0, " C", :black, :left, :bottom, delta=delta)
       plot!(xlims=(-3delta, R + 3delta))
   end
end;