算額（その472）

算額（その472）

宮城県角田市横倉 愛宕神社 明治15年(1882)1月
http://www.wasan.jp/miyagi/yokokuraatago.html

徳竹亜紀子，谷垣美保: 2021年度の算額調査，仙台高等専門学校名取キャンパス 研究紀要，第 58 号, p.7-28, 2022.
https://www.sendai-nct.ac.jp/natori-library/wp/wp-content/uploads/2022/03/kiyo2022-2.pdf

団扇の中に，正三角形，大円 2 個，小円 2 個が入っている。大円の直径が 23 寸 4 分のとき，小円の直径はいかほどか。

扇を構成する外円の半径と中心座標を r0, (0, 0)
大円の半径と中心座標を r1, (r1, r1 - r0/2)
小円の半径と中心座標を r2, (x2, -r2 - r0/2)
団扇の下部の円弧の半径と中心座標を r0/2, (0, -r0)
外円と円弧の交点座標を (x, y)
とおき，以下の連立方程式を解く。

include("julia-source.txt")

using SymPy
@syms r0::positive, r1::positive, r2::positive, x2::positive;
r1 = 234//20
eq1 = x2^2 + (r0/2 - r2)^2 - (r0/2 + r2)^2
eq2 = r1^2 + (r1 - r0/2)^2 - (r0 - r1)^2
eq3 = x2^2 + (-r0/2 - r2)^2 - (r0 - r2)^2
res = solve([eq1, eq2, eq3], (r0, r2, x2))

   1-element Vector{Tuple{Sym, Sym, Sym}}:
    (117/5, 351/100, 117*sqrt(30)/50)

外円の直径は大円の直径の 2 倍である（方程式によらなくてもわかることではあるが）。
小円の直径は大円の直径の 3/10 倍である。23.4*3/10 = 7.02 = 7寸0分2厘である。

術では「大円の直径を3倍して「1位ずらす」（10で割る）」ことで，7 寸としている（普通は 7 寸あまりありと記載される）。

図を描くために円弧の描き始めと描き終わりの角度を求めるために，外円と円弧の交点座標 (x, y) を求める。

using SymPy
@syms r0::positive, r1::positive, x, y;
r1 = 234//20
r0 = 117//5
eq11 = x^2 + y^2 - r0^2
eq12 = x^2 + (y + r0)^2 - (r0/2)^2
rs2 = solve([eq11, eq12], (x, y))

   2-element Vector{Tuple{Sym, Sym}}:
    (-117*sqrt(15)/40, -819/40)
    (117*sqrt(15)/40, -819/40)

using Plots

function draw(more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   r1 = 234//20
   (r0, r2, x2) = (2*r1, 3*r1/10, sqrt(30)*r1/5)
   @printf("小円の直径 = %g;  r0 = %g;  r1 = %g;  r2= %g;  x2= %g\n", 2r2, r0, r1, r2, x2)
   plot([√3r0/2, 0, -√3r0/2, √3r0/2], [-r0/2, r0, -r0/2, -r0/2], color=:black, lw=0.5)
   circle(0, 0, r0, :black)
   circle(r1, r1 - r0/2, r1, :red)
   circle(-r1, r1 - r0/2, r1, :red)
   circle(x2, -r2 - r0/2, r2, :blue)
   circle(-x2, -r2 - r0/2, r2, :blue)
   (x, y) = (117*sqrt(15)/40, -819/40)
   θ = round(Int, atand(Float64((y + r0)/x)))
   circle(0, -r0, r0/2, :green, beginangle=θ, endangle=180 - θ)
   if more
       delta = (fontheight = (ylims()[2]- ylims()[1]) / 500 * 10 * 2) /  3  # size[2] * fontsize * 2
       hline!([0], color=:black, lw=0.5)
       vline!([0], color=:black, lw=0.5)
       point(r1, r1 - r0/2, "大円:r1,(r1,r1-r0/2)", :red, :center, :top, delta=-delta)
       point(x2, -r2 - r0/2, "小円:r2,(x2,-r2-r0/2)", :blue, :center, :top, delta=-delta)
       point(0, -r0/2, " -r0/2", :black, :left, :bottom, delta=delta/3)
       point(0, -r0, " -r0", :green, :left, :bottom, delta=delta/3)
       point(√3r0/2, -r0/2, "(√3r0/2,-r0/2) ", :black, :right, :bottom, delta=delta/3)
       point(x, y, "(x,y)")
   end
end;

算額（その471）

算額（その471）

宮城県角田市横倉 愛宕神社 明治15年(1882)1月
http://www.wasan.jp/miyagi/yokokuraatago.html

徳竹亜紀子，谷垣美保: 2021年度の算額調査，仙台高等専門学校名取キャンパス 研究紀要，第 58 号, p.7-28, 2022.
https://www.sendai-nct.ac.jp/natori-library/wp/wp-content/uploads/2022/03/kiyo2022-2.pdf

二等辺三角形（圭）内に全円と菱形，および等円 2 個が入っている。全円の直径が 14 寸のとき，等円の直径を求めよ。

二等辺三角形の底辺の長さを 2a, 高さを b とする。また，菱形の右端の頂点座標を(c, r1) とする。
全円の半径と中心座標を r1，(0, r1)
等円の半径と中心座標を r2, (x2, 0)
とし，以下の連立方程式を解く。

include("julia-source.txt")

using SymPy
@syms a::positive, b::positive, c::positive,
     r1::positive, r2::positive, x2::positive;
r1 = 14//2
eq1 = x2^2 + (r1 - r2)^2 - (r1 + r2)^2
eq2 = (a - x2)r1 - r2*a
eq3 = (b - r1)a - b*c
eq4 = ((b - r1)a/r1)^2 - (b^2 + a^2)
eq5 = r1/(a - c) - b/a
eq5 = r1*a - (a - c)b
#eq5 = (2a +2sqrt(a^2 + b^2))r1/2 - a*b
#eq5 = distance(a, 0, 0, b, x2, r2) - r2^2
res = solve([eq1, eq2, eq3, eq4], (a, c, r2, x2))

   2-element Vector{NTuple{4, Sym}}:
    (sqrt(3*b - 2*sqrt(2)*sqrt(b^2 - 21*b + 98) - 28)*(7*b + 7*sqrt(2*b^2 - 42*b + 196) - 98)/(sqrt(b)*(b - 14)), sqrt(3*b - 2*sqrt(2)*sqrt(b^2 - 21*b + 98) - 28)*(7*b^2 + 7*b*sqrt(2*b^2 - 42*b + 196) - 147*b - 49*sqrt(2*b^2 - 42*b + 196) + 686)/(b^(3/2)*(b - 14)), (21*b - 14*sqrt(2*b^2 - 42*b + 196) - 196)/b, sqrt(588*b - 392*sqrt(2)*sqrt(b^2 - 21*b + 98) - 5488)/sqrt(b))
    (sqrt(3*b + 2*sqrt(2)*sqrt(b^2 - 21*b + 98) - 28)*(7*b - 7*sqrt(2*b^2 - 42*b + 196) - 98)/(sqrt(b)*(b - 14)), sqrt(3*b + 2*sqrt(2)*sqrt(b^2 - 21*b + 98) - 28)*(7*b^2 - 7*b*sqrt(2*b^2 - 42*b + 196) - 147*b + 49*sqrt(2*b^2 - 42*b + 196) + 686)/(b^(3/2)*(b - 14)), (21*b + 14*sqrt(2*b^2 - 42*b + 196) - 196)/b, sqrt(588*b + 392*sqrt(2)*sqrt(b^2 - 21*b + 98) - 5488)/sqrt(b))

2 組の解が得られるが，最初のものが適解である。しかし，解はすべて b が含まれている。つまり，b は任意の値を取りうるということである。

問の中の条件として「全円の直径が 14 寸のとき」とあるだけで b に関する条件がない。
図において，二等辺三角形の斜辺は全円の接線になるので，b が小さくなれば a は大きくなり，それに連れて等円の半径も大きくなる。「逆も真なり」である。

res[1][3] |> println

   (21*b - 14*sqrt(2*b^2 - 42*b + 196) - 196)/b

術では，「等円の直径は，全円の直径を 2 倍して 7 で割る」としている。等円の半径が 2 になるときの b は 784/31 = 25.2903225806452 である。

solve(res[1][3] - 2)[1].evalf() |> println

   25.2903225806452

res[1][3](b => 784/31) |> println

   2.00000000000000

等円の直径 = 4;  a = 10.4766;  b = 25.2903;  c = 7.57686;  r2= 2;  x2= 7.48331

たとえば b = 17 のときは，等円の直径は 7 になる。
等円の直径 = 7;  a = 19.799;  b = 16;  c = 11.1369;  r2= 3.5;  x2= 9.89949
これを術もどきでいうと「等円の直径は，全円の直径を 2 倍して 4 で割る」ことになってしまう。

using Plots

function draw(more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   r1 = 14//2
   b = 784/31  # 25.2903225806452
   (a, c, r2, x2) = (sqrt(3*b - 2*sqrt(2)*sqrt(b^2 - 21*b + 98) - 28)*(7*b + 7*sqrt(2*b^2 - 42*b + 196) - 98)/(sqrt(b)*(b - 14)), sqrt(3*b - 2*sqrt(2)*sqrt(b^2 - 21*b + 98) - 28)*(7*b^2 + 7*b*sqrt(2*b^2 - 42*b + 196) - 147*b - 49*sqrt(2*b^2 - 42*b + 196) + 686)/(b^(3/2)*(b - 14)), (21*b - 14*sqrt(2*b^2 - 42*b + 196) - 196)/b, sqrt(588*b - 392*sqrt(2)*sqrt(b^2 - 21*b + 98) - 5488)/sqrt(b))
   @printf("等円の直径 = %g;  a = %g;  b = %g;  c = %g;  r2= %g;  x2= %g\n", 2r2, a, b, c, r2, x2)
   plot([a, 0, -a, 0], [0, b, 0, 0], color=:black, lw=0.5)
   plot!([c, 0, -c, 0, c], [r1, 2r1, r1, 0, r1], color=:blue, lw=0.5)
   circle(0, r1, r1, :red)
   circle(x2, r2, r2, :magenta)
   if more
       delta = (fontheight = (ylims()[2]- ylims()[1]) / 500 * 10 * 2) /  3  # size[2] * fontsize * 2
       hline!([0], color=:black, lw=0.5)
       vline!([0], color=:black, lw=0.5)
       point(a, 0, "a", :black, :left, :bottom, delta=delta/2)
       point(0, b, " b", :black, :left, :vcenter)
       point(c, r1, " (c,r1)", :blue, :left, :vcenter)
       point(0, r1, " 全円:r1,(0,r1)", :red, :left, :vcenter)
       point(x2, r2, "等円:r2,(x2,r2)", :magenta, :center, :top, delta=-delta)
       point(0, 2r1, " 2r1", :red, :left, :bottom, delta=delta/2)
       
   end
end;