算額（その151）

算額（その151）

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

第2問: 直角三角形の中に同じ大きさの8小の正方形を入れる。斜辺とそれら正方形に内接する円の直径を求めよ。

正方形の一辺の長さを a，円の半径を r とする。

⊿ABC の斜辺の長さ AC は 5a である。
⊿ABC の面積は 2 通りの方法で求めることができる。
(1) BC*AB/2 = 4a * 3a / 2 = 6a^2
(2) (AB*r + BC*r + CA*r)/2 = (3a + 4a + 5a)r/2 = 6a*r
(1), (2) は等しいので，r について解き 2 倍して直径を得る。

using SymPy
@syms a::positive, r::positive;

solve(6a^2 - 6a*r, r)[1] * 2

   2a

すなわち，円の直径は，正方形の一辺の長さの 2 倍である。

using Plots
using Printf

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, linewidth=0.5, seriestype=:shape, linecolor=color, fillcolor=color)
   else
       plot!(ox .+ x, oy .+ y, color=color, linewidth=0.5)
   end
end;

function rect(x1, y1, x2, y2, color=:pink)
   plot!([x1, x2, x2, x1, x1], [y1, y1, y2,  y2, y1], color=color, linewidth=0.5, seriestype=:shape, fillcolor=color)
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(a, more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   (x, y) = (19a/3, 19a/4)
   plot()
   circle(2a, 2a, a, :pink, fill=true)
   plot!([0, x, 0, 0], [0, 0, y, 0], color=:black, lw=0.5)
   for i = 0:a:4a
       rect(i, 0, i+a, a, :slategray1)
   end
   for j = 0:a:3a
       rect(0, j, a, j+a, :slategray1)
   end
   if more == true
       point(a, 4a, "  A", :black, :left, :bottom)
       point(a, a, "  B", :black, :left, :bottom)
       point(5a, a, "  C", :black, :left, :bottom)
       vline!([0], color=:black, lw=0.5)
       hline!([0], color=:black, lw=0.5)
       xticks!(0:5, ["0"; "a"; string.(2:5) .* "a"], tickfontsize=12)
       yticks!(0:4, ["0"; "a"; string.(2:4) .* "a"])
   else
       plot!(showaxis=false)
   end
end;

算額（その150）

算額（その150）

岐阜県 田代神社 天保12年(1841)8月
http://www.wasan.jp/gifu/tasiro.html

第5問: 正方形内に 2 個の半円弧を描き，楕円と甲乙丙の 10 個の円を入れる（乙円は半円）。楕円の短径が丙円の直径に等しいとき，楕円の長径wk求めよ。

一辺の長さが a = 1 の正方形を考える。大きい半円の半径は r0 = a/2，甲円の半径は r1 = a/4，乙円，丙円の半径を r2, r3 とする。丙円の中心の座標は (r2, r1) である。算額では「楕円の短径が丙円の直径に等しい」とあるが，「楕円の短径は丙円の半径に等しい」であろう（円の直径に対する用語として短軸・長軸と書いているときもあるが紛らわしい）。楕円の長径を r4 とする。また，大きい半円と楕円の接点座標を (x0, y0) とおき，連立方程式を解く。

using SymPy
@syms a::positive, r0::positive, r1::positive, r2::positive, r3::positive, r4::positive,
     x3::positive, y3::positive, x0::positive, y0::positive;
a = 1
r0 = a // 2  # 大きい半円の半径
r1 = a // 4  # 甲円の半径
eq1 = (r0 - r2)^2 + r1^2 - (r2 + r1)^2
eq2 = r2 + r3 - r1
eq5 = r4^2*(x0 - r3)/(r3^2*(2r1 - y0)) + (x0 - r0)/y0  # (x0, y0) における円の接線と楕円の接線の傾きが等しい
eq6 = (x0 - r0)^2 + y0^2 - r0^2 # (x0, y0) が円周上にある
eq7 = (x0 - r3)^2/r3^2 + (y0 - 2r1)^2/r4^2 - 1; # (x0, y0) が楕円周上にある

res = solve([eq1, eq2, eq5, eq6, eq7], (r2, r3, r4, x0, y0))

   1-element Vector{NTuple{5, Sym}}:
    (1/6, 1/12, sqrt(299/960 - 41*sqrt(41)/960), sqrt(41)/244 + 23/244, 9/244 + 11*sqrt(41)/244)

正方形の一辺の長さを a としたとき，楕円の長軸（長径の2倍）は

sqrt(299/960 - 41*sqrt(41)/960) * 2a

   0.38982852141066027

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 ellipse(ox, oy, ra, rb; φ=0, beginangle=0, endangle=360,
                    color=:black, lty=:solid, lwd=0.5, fcolor="")
"""
(ox, oy) を中心，ra, rb を半径とする楕円（楕円弧）。
fcolor を指定すると塗りつぶし。
"""
   θ = beginangle:0.1:endangle
   if φ == 0
       if fcolor == ""
           plot!(ra .* cosd.(θ) .+ ox, rb .* sind.(θ) .+ oy,
                 linecolor=color, linestyle=lty, linewidth=lwd)
       else
           plot!(ra .* cosd.(θ) .+ ox, rb .* sind.(θ) .+ oy,
                 linecolor=color, linestyle=lty, linewidth=lwd,
                 seriestype=:shape, fillcolor=fcolor)
       end
   else
       x = ra .* cosd.(θ)
       y = rb .* sind.(θ)
       cosine = cosd(φ)
       sine = sind(φ)
       if fcolor == ""
           plot!(cosine .* x .- sine .* y .+ ox,
                 sine .* x .+ cosine .* y .+ oy,
                 linecolor=color, linestyle=lty, linewidth=lwd)
       else
           plot!(cosine .* x .- sine .* y .+ ox,
                 sine .* x .+ cosine .* y .+ oy,
                 linecolor=color, linestyle=lty, linewidth=lwd,
                 seriestype=:shape, fillcolor=fcolor)
       end
   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, r3, r4, x0, y0) = res[1]
   a = 1
   r0 = a / 2
   r1 = a / 4
   @printf("甲円の半径 = %.5f;  乙円の半径 = %.5f;  丙円の半径 = %.5f;  楕円の長径 = %.5f;  接点の座標 = (%.5f, %.5f)\n",
       r1, r2, r3, r4, x0, y0)
   plot([0, a, a, 0, 0], [0, 0, a, a, 0], color=:black, lw=0.5)
   circle(r0, r1, r1)  # 甲
   circle(r0, a - r1, r1)  # 甲
   circle(r2, 0, r2, :orange, beginangle=0, endangle=180)  # 乙
   circle(r2, a, r2, :orange, beginangle=180, endangle=360)  # 乙
   circle(a - r2, 0, r2, :orange, beginangle=0, endangle=180)  # 乙
   circle(a - r2, a, r2, :orange, beginangle=180, endangle=360)  # 乙
   circle(r0, 0, r0, :blue, beginangle=0, endangle=180)  # 大きい半円
   circle(r0, a, r0, :blue, beginangle=180, endangle=360)  # 大きい半円
   circle(r2, r1, r3, :green)  # 丙
   circle(r2, a - r1, r3, :green)  # 丙
   circle(a - r2, r1, r3, :green)  # 丙
   circle(a - r2, a - r1, r3, :green)  # 丙
   ellipse(r3, r0, r3, r4)
   ellipse(a - r3, r0, r3, r4)
   if more == true
       point(r0, r1, "甲(r0,r1)", :red, :center)
       point(r2, 0, "乙(r2, 0)", :orange, :center, :bottom)
       point(r2, r1, "丙(r2,r1)", :green, :center)
       point(r3, r0, "(r3,r0)", :black, :center)
       point(x0, y0, "(x0,y0)")
       vline!([0], color=:black, lw=0.5)
       hline!([0], color=:black, lw=0.5)
   else
       plot!(showaxis=false)
   end
end;

   甲円の半径 = 0.25000;  乙円の半径 = 0.16667;  丙円の半径 = 0.08333;  楕円の長径 = 0.19491;  接点の座標 = (0.12050, 0.32555)