算額（その1231）

算額（その1231）

京都府京都市上京区北野馬喰町 北野天満宮 明治12年(1879)
近畿数学史学会：近畿の算額「数学の絵馬を訪ねて」，平成4年5月16日 初版第一刷，大阪教育図書株式会社，大阪市.
キーワード：円4個，二等辺三角形

二等辺三角形の中に，甲円 2 個，乙円と丙円を 1 個ずつ容れる。上斜（斜辺）が 1014 寸，下斜（底辺）が 1428 寸のとき，甲円の直径はいかほどか。

斜辺，底辺を c, 2a とおく。高さは b = sqrt(c^2 - a^2) となる。
甲円の半径と中心座標を r1, (x1, r1)
乙円の半径と中心座標を r2, (0, 2r3 + r2)
丙円の半径と中心座標を r3, (0, r3)
とおき，以下の連立方程式を解く。

include("julia-source.txt");

using SymPy
@syms a::positive, b::positive, c::positive,
     r1::positive, x1::positive,
     r2::positive, r3::positive
@syms a, b, c, r1, x1, r2, r3
# b = sqrt(c^2 - a^2) 式には代入しない
eq1 = dist2(a, 0, 0, b, 0, 2r3 + r2, r2)
eq2 = dist2(a, 0, 0, b, x1, r1, r1)
eq3 = x1^2 + (2r3 + r2 - r1)^2 - (r1 + r2)^2
eq4 = x1^2 + (r1 - r3)^2 - (r1 + r3)^2;
# res = solve([eq1, eq2, eq3], (r1, r2, r3))

using NLsolve

function nls(func, params...; ini = [0.0])
   if typeof(ini) <: Number
       r = nlsolve((vout, vin) -> vout[1] = func(vin[1], params..., [ini]), ftol=big"1e-40")
       v = r.zero[1]
   else
       r = nlsolve((vout, vin)->vout .= func(vin, params...), ini, ftol=big"1e-40")
       v = r.zero
   end
   return Float64.(v), r.f_converged
end;

function H(u)
   (r1, x1, r2, r3) = u
   return [
       a^2*b^2 - 2*a^2*b*r2 - 4*a^2*b*r3 + 4*a^2*r2*r3 + 4*a^2*r3^2 - b^2*r2^2,  # eq1
       b*(a^2*b - 2*a^2*r1 - 2*a*b*x1 + 2*a*r1*x1 - b*r1^2 + b*x1^2),  # eq2
       x1^2 - (r1 + r2)^2 + (-r1 + r2 + 2*r3)^2,  # eq3
       x1^2 + (r1 - r3)^2 - (r1 + r3)^2,  # eq4
   ]
end;

a = 1428//2
c = 1014
b = sqrt(c^2 - a^2)

iniv = BigFloat[180, 290, 203, 115]
res = nls(H, ini=iniv)

   ([178.5, 285.6, 203.09333333333333, 114.24], true)

甲円，乙円，丙円の直径はそれぞれ 357 寸，406.18666666666667寸，228.48 寸である。

function draw(a, c, more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   #b = sqrt(c^2 - a^2)
   (r1, x1, r2, r3) = res[1]
   @printf("上斜が %g，下斜が %g のとき，甲円の直径は %g である。\n", 2a, c, 2r1)
   plot([a, 0, -a, a], [0, b, 0, 0], color=:green, lw=0.5)
   circle2(x1, r1, r1)
   circle(0, 2r3 + r2, r2, :blue)
   circle(0, r3, r3, :magenta)
   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)
       point(0, b, " b", :green, :left, :bottom, delta=delta)
       point(a, 0, " a", :green, :left, :bottom, delta=delta)
       point(x1, r1, "甲円:r1\n(x1,r1)", :red, :center, delta=-delta)
       point(0, 2r3 + r2, "乙円:r2\n(0,2r3+r2)", :blue, :center, delta=-delta)
       point(0, r3, "丙円:r3\n(0,r3)", :magenta, :center, delta=-delta)
   end
end;

draw(1428/2, 1014, true)

 

算額（その1230）

算額（その1230）

福島県小野町 東堂山満福寺観音堂 天保8年(1837)
https://www.nippon.com/ja/japan-topics/c12803/

キーワード：円4個，直角三角形

外円の中に，直角三角形と，その 3 辺の中点を中心とする 3 個の円を容れる。鈎，股が 3 寸，4 寸のとき，外円の直径はいかほどか。

外円の半径と中心座標を R, (x, y)
鈎，股，弦を a, b, c
それぞれの中点を中心とする円の半径を r1, r2, r3
とおき，以下の連立方程式を解く。

include("julia-source.txt");

using SymPy
@syms a::positive, b::positive, c::positive, R::positive, r1::positive, r2::positive, r3::positive, x::positive, y::positive
@syms a, b, c, R, r1, r2, r3, x, y
#c = sqrt(a^2 + b^2)
r1 = a/2
r2 = b/2
r3 = c/2
eq1 = x^2 + (r1 - y)^2 - (R - r1)^2
eq2 = (r2 - x)^2 + y^2 - (R - r2)^2
eq3 = (r2 - x)^2 + (r1 - y)^2 - (R - r3)^2
res = solve([eq1, eq2, eq3], (R, x, y))[1]  # 1 of 2

   ((a^5 - a^4*c - a^3*b^2 - a^3*c^2 - a^2*b^3 + a^2*c^3 - sqrt(2)*a^2*sqrt(a*b^3*(-a^2 + 2*a*b - b^2 + c^2)) - sqrt(2)*a*b*sqrt(a*b^3*(-a^2 + 2*a*b - b^2 + c^2)) + sqrt(2)*a*c*sqrt(a*b^3*(-a^2 + 2*a*b - b^2 + c^2)) + b^5 - b^4*c - b^3*c^2 + b^2*c^3)/(4*(a^4 - 2*a^3*c - a^2*b^2 + a^2*c^2 + b^4 - 2*b^3*c + b^2*c^2)), (a^4*b^2 + a^3*b^3 - 3*a^3*b^2*c + sqrt(2)*a^3*sqrt(a*b^3*(-a^2 + 2*a*b - b^2 + c^2)) - 3*a^2*b^3*c + 3*a^2*b^2*c^2 + sqrt(2)*a^2*b*sqrt(a*b^3*(-a^2 + 2*a*b - b^2 + c^2)) - 2*sqrt(2)*a^2*c*sqrt(a*b^3*(-a^2 + 2*a*b - b^2 + c^2)) - a*b^5 + a*b^4*c + a*b^3*c^2 - a*b^2*c^3 - sqrt(2)*a*b*c*sqrt(a*b^3*(-a^2 + 2*a*b - b^2 + c^2)) + sqrt(2)*a*c^2*sqrt(a*b^3*(-a^2 + 2*a*b - b^2 + c^2)) + b^6 - b^5*c - b^4*c^2 + b^3*c^3)/(4*b*(a^4 - 2*a^3*c - a^2*b^2 + a^2*c^2 + b^4 - 2*b^3*c + b^2*c^2)), a*(a^4 - a^3*b - a^3*c + a^2*b*c - a^2*c^2 + a*b^3 - 3*a*b^2*c + a*b*c^2 + a*c^3 + b^4 - 3*b^3*c + 3*b^2*c^2 - b*c^3)/(4*(a^4 - 2*a^3*c - a^2*b^2 + a^2*c^2 + b^4 - 2*b^3*c + b^2*c^2)) + sqrt(2)*sqrt(-a*b^3*(a - b - c)*(a - b + c))*(b - c)*(a + b - c)/(4*(a^4 - 2*a^3*c - a^2*b^2 + a^2*c^2 + b^4 - 2*b^3*c + b^2*c^2)))

R, x, y を a, b, c を含む式で表すと，SymPy では簡約化できない複雑な長い式（省略）になる。
鈎，股として（弦は別途計算），具体的な数値を代入すればよい。

鈎，股が 3 寸，4 寸のとき，外円の直径は 2*72/23 = 6 + 6/23 である。

res[1](a => 3, b => 4, c => 5) |> println 
res[2](a => 3, b => 4, c => 5) |> println 
res[3](a => 3, b => 4, c => 5) |> println 

   72/23
   36/23
   24/23

function draw(a, b, more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   c = sqrt(a^2 + b^2)
   (R, x, y) = ((a^5 - a^4*c - a^3*b^2 - a^3*c^2 - a^2*b^3 + a^2*c^3 - sqrt(2)*a^2*sqrt(a*b^3*(-a^2 + 2*a*b - b^2 + c^2)) - sqrt(2)*a*b*sqrt(a*b^3*(-a^2 + 2*a*b - b^2 + c^2)) + sqrt(2)*a*c*sqrt(a*b^3*(-a^2 + 2*a*b - b^2 + c^2)) + b^5 - b^4*c - b^3*c^2 + b^2*c^3)/(4*(a^4 - 2*a^3*c - a^2*b^2 + a^2*c^2 + b^4 - 2*b^3*c + b^2*c^2)), (a^4*b^2 + a^3*b^3 - 3*a^3*b^2*c + sqrt(2)*a^3*sqrt(a*b^3*(-a^2 + 2*a*b - b^2 + c^2)) - 3*a^2*b^3*c + 3*a^2*b^2*c^2 + sqrt(2)*a^2*b*sqrt(a*b^3*(-a^2 + 2*a*b - b^2 + c^2)) - 2*sqrt(2)*a^2*c*sqrt(a*b^3*(-a^2 + 2*a*b - b^2 + c^2)) - a*b^5 + a*b^4*c + a*b^3*c^2 - a*b^2*c^3 - sqrt(2)*a*b*c*sqrt(a*b^3*(-a^2 + 2*a*b - b^2 + c^2)) + sqrt(2)*a*c^2*sqrt(a*b^3*(-a^2 + 2*a*b - b^2 + c^2)) + b^6 - b^5*c - b^4*c^2 + b^3*c^3)/(4*b*(a^4 - 2*a^3*c - a^2*b^2 + a^2*c^2 + b^4 - 2*b^3*c + b^2*c^2)), a*(a^4 - a^3*b - a^3*c + a^2*b*c - a^2*c^2 + a*b^3 - 3*a*b^2*c + a*b*c^2 + a*c^3 + b^4 - 3*b^3*c + 3*b^2*c^2 - b*c^3)/(4*(a^4 - 2*a^3*c - a^2*b^2 + a^2*c^2 + b^4 - 2*b^3*c + b^2*c^2)) + sqrt(2)*sqrt(-a*b^3*(a - b - c)*(a - b + c))*(b - c)*(a + b - c)/(4*(a^4 - 2*a^3*c - a^2*b^2 + a^2*c^2 + b^4 - 2*b^3*c + b^2*c^2)))
   @printf("a = %g;  b = %g;  c = %g;  R = %g;  x = %g;  y = %g\n", a, b, c, R, x, y)
   plot([0, b, 0, 0], [0, 0, a, 0], color=:green, lw=0.5)
   circle(x, y, R)
   circle(0, a/2, a/2, :blue)
   circle(b/2, 0, b/2, :magenta)
   circle(b/2, a/2, c/2, :orange)
   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)
       point(x, y, "外円:R,(x,y)", :red, :center, delta=-delta/2)
       point(0, a/2, "小円:r1,(0,a/2)", :blue, :center, delta=-delta/2)
       point(b/2, 0, "中円:r2,(b/2,0)", :magenta, :center, delta=-delta/2)
       point(b/2, a/2, "大円:r3,(b/2,a/2)", :orange, :center, delta=-delta/2)
       point(0, a, "a ", :green, :right, :bottom, delta=delta)
       point(b, 0, " b", :green, :left, delta=-delta)
   end
end;

draw(3, 4, true)

算額（その1229）

算額（その1229）

(9) 兵庫県姫路市井ノ口宮山 荒川神社 明治3年(1870)
近畿数学史学会：近畿の算額「数学の絵馬を訪ねて」，平成4年5月16日 初版第一刷，大阪教育図書株式会社，大阪市.
キーワード：円2個，直角三角形，正方形

「問」は簡単で，「図のように大円の直径が 9 寸のとき，小円の直径はいかほどか」

これだけでは小円が内接する直角三角形大きさが定まらないが，特に断らない限り「直角三角形の辺（鈎，股，弦）の比は 3:4:5 である」ということであろう。

鈎 = 大円の直径
股 = 4鈎/3
弦 = 5鈎/3
そして，直角三角形に内接する円の直径は，「鈎 + 股 - 弦」である。
簡約化すると小円の直径は大円の直径の 2/3 である。
よって，大円の直径が 9 寸のとき，小円の直径は 6 寸である。

using SymPy
@syms 大円の直径, 鈎, 股, 弦
鈎 = 大円の直径
小円の直径 = (鈎 + 4鈎/3 - 5鈎/3)
小円の直径 |> println

   2*大円の直径/3

図を描くために直角三角形の頂点と，小円の中心座標を求める。

include("julia-source.txt");

using SymPy

@syms r1::positive, x::positive, y::negative,
     r2::positive, x2::positive, y2::positive,
     鈎::positive
鈎 = 2r1
eq1 =dist2(鈎/√Sym(2), 鈎/√Sym(2), x, y, x2, y2, r2)
eq2 =dist2(0, 0, x, y, x2, y2, r2)
eq3 =dist2(0, 0, 鈎/√Sym(2), 鈎/√Sym(2), x2, y2, r2)
eq4 = (x^2 + y^2) - (5*鈎/3)^2
eq5 = ((x - 鈎/√Sym(2))^2 + (y - 鈎/√Sym(2))^2) - (4*鈎/3)^2
res = solve([eq1, eq2, eq3, eq4, eq5], (x, y, r2, x2, y2))[1]  # 1 of 2

   (7*sqrt(2)*r1/3, -sqrt(2)*r1/3, 2*r1/3, sqrt(2)*r1, sqrt(2)*r1/3)

直角三角形の頂点座標 (7√2r1/3, -√2r1/3)
小円の中心座標 (√2r1, √2r1/3)

function draw(r1, more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   鈎 = 2r1
   (x, y, r2, x2, y2) = (7√2r1/3, -√2r1/3, 2r1/3, √2r1, √2r1/3)
   plot([0, 鈎/√2, 0, -鈎/√2, 0], [0, 鈎/√2, 2鈎/√2, 鈎/√2, 0], color=:green, lw=0.5)
   plot!([-鈎/√2, -x, 0, x, 鈎/√2], [鈎/√2, y, 0, y, 鈎/√2], color=:green, lw=0.5)
   circle(0, 鈎/√2, r1, :blue)
   circle2(x2, y2, r2)
   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)
       point(0, 鈎/√2, "大円:r1,(0,鈎/√2)", :blue, :center, delta=-delta)
       point(x2, y2, "小円:r2\n(x2,y2)", :red, :center, delta=-delta)
       point(x, y, "(x,y)", :green, :left, :bottom, delta=delta)
       point(鈎/√2, 鈎/√2, "(鈎/√2,鈎/√2)", :green, :left, :bottom, delta=delta)
       point(0, 2鈎/√2, "(0,2鈎/√2)", :green, :left, :bottom, delta=delta)
   end
end;

draw(9/2, true)