算額（その460）

算額（その460）

埼玉県鴻巣市 薬師堂 明治23年
山口正義(2016): やまぶき2, 第36号
https://yamabukiwasan.sakura.ne.jp/ymbk36.pdf

一〇一 大宮市高鼻町 氷川神社 明治31年(1898)
埼玉県立図書館：埼玉県史料集 第二集『埼玉の算額』，昭和44年，誠美堂印刷所，埼玉県与野市.

外円の中に正方形，大円，中円，小円，弦がある。小円の直径が 1 寸のとき，中円の直径はいかほどか。

外円の半径と中心座標を r0, (0, 0)
大円の半径と中心座標を r1, (x1, 2r2 - r0 + r1)
中円の半径と中心座標を r2, (0, 3r2 - r0), (0, r2 - r0)
小円の半径と中心座標を r3, (x3, 2r2 - r0 - r3)
として以下の連立方程式を解く。

include("julia-source.txt")

using SymPy

@syms r0::positive, r1::positive, x1::positive, r2::positive, r3::positive, x3::positive;

eq1 = x1^2 + (r1 - r2)^2 - (r1 + r2)^2
eq2 = x1^2 + (2r2 - r0 + r1)^2 - (r0 - r1)^2
eq3 = x3^2 + (r2 - r3)^2 - (r2 + r3)^2
eq4 = x3^2 + (2r2 - r0 - r3)^2 - (r0 - r3)^2
eq5 = distance(0, 4r2 - r0, r0 - 2r2, 2r2, x1, 2r2 - r0 + r1) - r1^2
res = solve([eq1, eq2, eq3, eq4, eq5], (r0, r1, x1, r2, x3))

   2-element Vector{NTuple{5, Sym}}:
    (25*r3/4, 5*r3/2, 5*sqrt(2)*r3/2, 5*r3/4, sqrt(5)*r3)
    (r3*(59 - 30*sqrt(2))/4, -4*sqrt(2)*r3 + 13*r3/2, -5*r3 + 3*sqrt(2)*r3/2, r3*(2*sqrt(2) + 7)/4, r3*sqrt(2*sqrt(2) + 7))

2 通りの解が得られるが，最初のものが適解である。

   r0 = 3.125;  r1 = 1.25;  x1 = 1.76777;  r2 = 0.625;  x3 = 1.11803
   中円の直径 = 1.25

中円の半径は 5*r3/4 で，小円の半径 r3 の 5/4 倍である。
r3 = 1/2寸のとき，小円の半径は 5/8寸，直径は 5/4寸 = 1寸2分5厘である。

using Plots

function draw(more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   r3 = 1/2
   (r0, r1, x1, r2, x3) = (25*r3/4, 5*r3/2, 5*sqrt(2)*r3/2, 5*r3/4, sqrt(5)*r3)
   @printf("r0 = %g;  r1 = %g;  x1 = %g;  r2 = %g;  x3 = %g\n", r0, r1, x1, r2, x3)
   @printf("中円の直径 = %g\n", 2r2)
   plot()
   circle(0, 0, r0, :black)
   circle(x1, 2r2 + r1 - r0, r1)
   circle(-x1, 2r2 + r1 - r0, r1)
   circle(0, r2 - r0, r2, :orange)
   circle(0, 3r2 - r0, r2, :orange)
   circle(x3, 2r2 - r0 - r3, r3, :green)
   circle(-x3, 2r2 - r0 - r3, r3, :green)
   segment(-sqrt(r0^2 - (2r2 - r0)^2), 2r2 - r0, sqrt(r0^2 - (2r2 - r0)^2), 2r2 - r0)
   plot!([0, r0 - 2r2, 0, 2r2 - r0, 0], [4r2 - r0, 2r2, r0, 2r2, 4r2 - r0], color=:blue, lw=0.5)
   if more
       delta = (fontheight = (ylims()[2]- ylims()[1]) / 500 * 10 * 2) /  3  # size[2] * fontsize * 2
       point(x1, 2r2 - r0 + r1, "大円:r1,(x1,2r2-r0+r1)", :red, :center, delta=-delta)
       point(0, 3r2 - r0, " 中円:r2,(0,3r2-r0)", :black, :left, :vcenter)
       point(0, r2 - r0, " r2-r0", :orange, :left, :vcenter)
       point(x3, 2r2 - r0 - r3, "小円:r3,(x3,2r2-r0-r3)", :black, :center, :bottom, delta=delta)
       point(0, 2r2 - r0, " 2r2-r0", :black, :left, :bottom, delta=delta/2)
       point(0, 4r2 - r0, " 4r2-r0", :black, :left, :vcenter)
       point(0, 2r2, " 2r2", :blue, :left, :vcenter)
       point(0, r0, " r0", :blue, :left, :bottom, delta=delta/2)
       point(r0, 0, "r0 ", :black, :right, :bottom, delta=delta/2)
       hline!([0], color=:gray, lw=0.5)
       vline!([0], color=:gray, lw=0.5)
   else
      plot!(showaxis=false)
   end
end;

算額（その459）

算額（その459）

長野県長野市 美和神社 文化10年(1813)

中村信弥(1999)：算額への招待
http://www.wasan.jp/syotai/syotai.html

上下 2 個の長方形の上段には大円，下段には小円を並べたところ，上下の長方形の横の長さの差が 5 寸になった。大円と小円の個数を求めよ。

算額では大円径，小円径を 91, 75 として解を求めているが，一般解を求めればよい。
どちらの長方形が大きいかの制約条件はない。

術は「剰一術（盈一術），現代解法（中村）は「法による除算」を用いているが，以下に示すコンピュータプログラムで全探索してもたかが知れている。

function prog(r1=91, r2=75; n=100)
   for i in 1:n
       for j in 1:n
           if abs(i*r1 - j*r2) == 5
               println("大円の個数 = $i, 小円の個数 = $j, $(i*r1) - $(j*r2)")
           end
       end
   end
end;

prog()

   大円の個数 = 5, 小円の個数 = 6, 455 - 450
   大円の個数 = 70, 小円の個数 = 85, 6370 - 6375
   大円の個数 = 80, 小円の個数 = 97, 7280 - 7275

円の個数の上限を 100 個までにして探索して複数の答えを得る。
どちらの長方形が大きいかの制約条件はないので，大円 5 個，小円 6 個が最小解である。
術，中村は大円 70 個，小円 85 個としているが，小円を入れる長方形の方が大きい最小の解である。

大円径，小円径が 67, 47 のとき，以下のような解が求まる。

prog(67, 47)

   大円の個数 = 12, 小円の個数 = 17, 804 - 799
   大円の個数 = 35, 小円の個数 = 50, 2345 - 2350
   大円の個数 = 59, 小円の個数 = 84, 3953 - 3948

大円径，小円径が 1009, 997 のときは探索範囲を 1000 個まで広げて，以下のような解が求まる。

prog(1009, 997, n=1000)

   大円の個数 = 415, 小円の個数 = 420, 418735 - 418740
   大円の個数 = 582, 小円の個数 = 589, 587238 - 587233

include("julia-source.txt")

using Plots

function draw(r1, r2, n1, n2, more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   (r1, r2) = (r1, r2) ./ 2
   (w1, w2) = (n1*2r1, n2*r2)
   plot([0, n1*2r1, n1*2r1, 0, 0], [0, 0, 2r1, 2r1, 0], color=:black, lw=0.5)
   for i in 1:n1
       circlef((2i - 1)r1, r1, r1, :cyan)
   end
   plot!([0, n2*2r2, n2*2r2, 0, 0], 2r1 .+ [0, 0, 2r2, 2r2, 0], color=:black, lw=0.5)
   for i in 1:n2
       circlef((2i - 1)r2, r2 + 2r1, r2, :red)
   end
   if more
       delta = (fontheight = (ylims()[2]- ylims()[1]) / 500 * 10 * 2) /  3  # size[2] * fontsize * 2
       annotate!(r1, r1, text(@sprintf("%i寸 × %i個 = %i寸", 2r1, n1, n1*2r1), :left, :vcenter, 10))
       annotate!(r1, 2r1 + r2, text(@sprintf("%i寸 × %i個 = %i寸", 2r2, n2, n2*2r2), :left, :vcenter, 10))
       hline!([0], color=:gray, lw=0.5)
       vline!([0], color=:gray, lw=0.5)
   else
      plot!(showaxis=false)
   end
end;