算額（その885）

算額（その885）

六四 加須市不動岡 総願寺 慶応2年(1866)
埼玉県立図書館：埼玉県史料集 第二集『埼玉の算額』，昭和44年，誠美堂印刷所，埼玉県与野市.

外円の中に，大円 2 個，小円 6 個を入れる。大円の直径が 2 寸のとき，小円の直径はいかほどか。

外円の半径と中心座標を R, (0, 0)
大円の半径と中心座標を r1, (0, r1); r1 = R/2
小円の半径と中心座標を r2, (x2, 2r2)
とおき，以下の連立方程式を解く。

include("julia-source.txt");

using SymPy
@syms R::positive, r1::positive, r2::positive, x2::positive
R = 2r1
eq1 = x2^2 + 4r2^2 - (R - r2)^2
eq2 = x2^2 + (r1 - 2r2)^2 - (r1 + r2)^2
res = solve([eq1, eq2], (r2, x2))

   1-element Vector{Tuple{Sym{PyCall.PyObject}, Sym{PyCall.PyObject}}}:
    (2*r1/5, 4*sqrt(3)*r1/5)

小円の半径は大円の半径の 2/5 倍である。
大円の直径が 2 寸のとき，小円の直径は 4/5 寸 = 8分 である。

その他のパラメータは以下のとおりである。
   r1 = 1;  r2 = 0.4;  x2 = 1.38564;  R = 2

function draw(more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   r1 = 2/2
   (r2, x2) = r1 .* (2/5, 4√3/5)
   R = 2r1
   @printf("大円の直径が %g のとき，小円の直径は %g である\n", 2r1, 2r2)
   @printf("r1 = %g;  r2 = %g;  x2 = %g;  R = %g\n", r1, r2, x2, R)
   plot()
   circle(0, 0, R, :blue)
   circle22(0, r1, r1)
   circle4(x2, 2r2, r2, :green)
   circle2(x2, 0, r2, :green)
   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, R, " R", :blue, :left, :bottom, delta=delta/2)
       point(R, 0, " R", :blue, :left, :bottom, delta=delta/2)
       point(0, r1, "大円:r1\n(0,r1)", :red, :center, delta=-delta/2)
       point(0, -r1, "大円:r1\n(0,-r1)", :red, :center, delta=-delta/2)
       point(0, 0, "", :blue)
       point(x2, 2r2, "小円:r2\n(x2,2r2)", :green, :center, delta=-delta)
   end
end;

算額（その884）

算額（その884）

六四 加須市不動岡 総願寺 慶応2年(1866)
埼玉県立図書館：埼玉県史料集 第二集『埼玉の算額』，昭和44年，誠美堂印刷所，埼玉県与野市.

外円内に二等辺三角形と，大円を 2 個，小円を 2 個入れる。小円の直径が 1 寸のとき，大円の直径はいかほどか。

外円の半径と中心座標を R, (0, 0)
大円の半径と中心座標を r1, (0, -r1), (0, r1); r1 = R/2
小円の半径と中心座標を r2, (x2, y2); y2 < 0
二等辺三角形の底辺と外円の接点座標を (x0, y0); y0 < 0
とおき，以下の連立方程式を解く。

include("julia-source.txt");

using SymPy
@syms R::positive, x0::positive, y0::negative, r1::positive,
     r2::positive, x2::positive, y2::negative
R = 2r1
y0 = y2 + r2
L = sqrt(x0^2 + (2r1 - y0)^2)
eq1 = x0^2 + y0^2 - R^2
eq2 = r1/3r1 - x0/L
eq3 = x2^2 + (r1 + y2)^2 - (r1 + r2)^2
eq4 = x2^2 + y2^2 - (R - r2)^2;
res = solve([eq1, eq2, eq3, eq4], (r1, x2, y2, x0))

   1-element Vector{NTuple{4, Sym{PyCall.PyObject}}}:
    (9*r2, 8*r2, -15*r2, 8*sqrt(2)*r2)

大円の半径は小円の半径の 9 倍である。
小円の直径が 1 寸のとき，大円の直径は 9 寸である。

その他のパラメータは以下のとおりである。
r2 = 0.5;  r1 = 4.5;  x2 = 4;  y2 = -7.5;  R = 9;  x0 = 5.65685;  y0 = -7

function draw(more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   r2 = 1/2
   (r1, x2, y2, x0) = r2 .* (9, 8, -15, 8√2)
   R = 2r1
   y0 = y2 + r2
   @printf("小円の直径が %g のとき，大円の直径は %g である\n", 2r2, 2r1)
   @printf("r2 = %g;  r1 = %g;  x2 = %g;  y2 = %g;  R = %g;  x0 = %g;  y0 = %g\n", r2, r1, x2, y2, R, x0, y0)
   plot()
   circle(0, 0, R, :blue)
   circle22(0, r1, r1)
   plot!([x0, 0, -x0, x0], [y0, R, y0, y0], color=:green, lw=0.5)
   circle2(x2, y2, r2, :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, R, " R", :blue, :left, :bottom, delta=delta/2)
       point(R, 0, " R", :blue, :left, :bottom, delta=delta/2)
       point(0, r1, "大円:r1\n(0,r1)", :red, :center, delta=-delta/2)
       point(0, -r1, "大円:r1\n(0,-r1)", :red, :center, delta=-delta/2)
       point(0, 0, "", :blue)
       point(x0, y0, " (x0,y0)", :blue, :left, :vcenter)
       point(x2, y2, "  小円:r2,(x2,y2)", :magenta, :left, delta=-delta)
   end
end;

算額（その883）

算額（その883）

六四 加須市不動岡 総願寺 慶応2年(1866)
埼玉県立図書館：埼玉県史料集 第二集『埼玉の算額』，昭和44年，誠美堂印刷所，埼玉県与野市.

半円の中に二等辺三角形と等円 3 個を入れる。等円の直径が 1 寸のとき，半円の直径はいかほどか。

二等辺三角形の底辺の長さを 2a
半円の半径と中心座標を r1, (0, 0)
等円の半径と中心座標を r2, (0, r2), (x2, r2)
とおき，以下の連立方程式を解く。

include("julia-source.txt");

using SymPy
@syms d, a::positive, r1::positive, r2::positive,
     x2::positive
eq1 = numerator(apart(dist(a, 0, 0, r1, 0, r2) - r2^2, d))
eq2 = numerator(apart(dist(a, 0, 0, r1, x2, r2) - r2^2, d))
eq3 = x2^2 + r2^2 - (r1 - r2)^2
res = solve([eq1, eq2, eq3], (a, r1, x2))

   1-element Vector{Tuple{Sym{PyCall.PyObject}, Sym{PyCall.PyObject}, Sym{PyCall.PyObject}}}:
    (r2*sqrt(sqrt(2) + 2 + sqrt(2*sqrt(2) + 3))/sqrt(1 + sqrt(2)), r2*(1 + sqrt(2*sqrt(2) + 3)), r2*sqrt(2 + 2*sqrt(2)))

res[1][2] |> sympy.sqrtdenest |> simplify |> println

   r2*(sqrt(2) + 2)

半円の半径は等円の半径の 1 + sqrt(2√2 + 3) = √2 + 2 倍である。
等円の直径が 1 寸のとき，半円の直径は 3.414213562373095 寸である。

算額では「術曰置二個開平方加二個乗等円径得半円径合問」は，「等円の直径の √2 + 2 倍」である。

その他のパラメータは以下のとおりである。

r2 = 0.5;  a = 0.776887;  r1 = 1.70711;  x2 = 1.09868

function draw(more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   r2 = 1/2
   (a, r1, x2) = (r2*sqrt(sqrt(2) + 2 + sqrt(2*sqrt(2) + 3))/sqrt(1 + sqrt(2)), r2*(1 + sqrt(2*sqrt(2) + 3)), r2*sqrt(2 + 2*sqrt(2)))
   @printf("等円の直径が %g のとき，半円の直径は %g である\n", 2r2, 2r1)
   @printf("r2 = %g;  a = %g;  r1 = %g;  x2 = %g\n", r2, a, r1, x2)
   plot(ylims=(-0.15, 1.85))
   circle(0, 0, r1, beginangle=0, endangle=180)
   circle(0, r2, r2, :blue)
   circle2(x2, r2, r2, :blue)
   segment(0, r1, a, 0, :green)
   segment(0, r1, -a, 0, :green)
   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(a, 0, "a", :green, :center, delta=-delta)
       point(0, r1, "r1", :red, :center, :bottom, delta=delta)
       point(r1, 0, "r1", :red, :center, delta=-delta)
       point(0, r2, "等円:r2,(0,r2)", :blue, :center, delta=-delta)
       point(x2, r2, "等円:r2,(x2,r2)", :blue, :center, delta=-delta)
       point(0, 0, "", :red)
   end
end;

算額（その882）

算額（その882）

六四 加須市不動岡 総願寺 慶応2年(1866)
埼玉県立図書館：埼玉県史料集 第二集『埼玉の算額』，昭和44年，誠美堂印刷所，埼玉県与野市.

長方形の中に対角線と斜線を引き，区画された領域に甲円，乙円 2 個ずつを入れる。乙円の直径が 1 寸のとき，長方形の長辺の長さはいかほどか。

長方形の長辺，短辺をそれぞれ a, b
斜線と長辺の交点座標を (c, 0)
甲円の半径と中心座標を r1, (a - r1, y1), (a/2, r1)
乙円の半径と中心座標を r2, (r2, b/2 + r1)
とおき，以下の連立方程式の数値解を求める。

include("julia-source.txt");

using SymPy
@syms d, a::positive, b::positive, c::positive,
     r1::positive, y1::positive, r2::positive
eq1 = numerator(apart(dist(a, 0, 0, b, a - r1, y1) - r1^2, d))
eq2 = numerator(apart(dist(a, 0, 0, b, a/2, r1) - r1^2, d))
eq3 = numerator(apart(dist(a, 0, 0, b, r2, b/2 + r2) - r2^2, d))
eq4 = numerator(apart(dist(a, b, c, 0, a - r1, y1) - r1^2, d))
eq5 = numerator(apart(dist(a, b, c, 0, a/2, r1) - r1^2, d));

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)
   (a, b, c, r1, y1) = u
   return [
       -a^2*r1^2 + a^2*y1^2 - 2*a*b*r1*y1,  # eq1
       a^2*b^2 - 4*a^2*b*r1 - 4*b^2*r1^2,  # eq2
       a^2*b^2 - 4*a^2*b*r2 - 4*a*b^2*r2 + 8*a*b*r2^2,  # eq3
       a^2*b^2 - 2*a^2*b*y1 - a^2*r1^2 + a^2*y1^2 - 2*a*b^2*c - 2*a*b^2*r1 + 4*a*b*c*y1 + 2*a*b*r1*y1 + 2*a*c*r1^2 - 2*a*c*y1^2 + b^2*c^2 + 2*b^2*c*r1 - 2*b*c^2*y1 - 2*b*c*r1*y1 - c^2*r1^2 + c^2*y1^2,  # eq4
       a^2*b^2 - 4*a^2*b*r1 - 4*a*b^2*c + 12*a*b*c*r1 + 4*b^2*c^2 - 4*b^2*r1^2 - 8*b*c^2*r1,  # eq5
   ]
end;

r2 = 1/2
iniv = BigFloat[53, 28, 30, 6, 11] ./5
res = nls(H, ini=iniv)

   ([5.23606797749979, 2.618033988749895, 2.8944271909999157, 0.6180339887498949, 1.0], true)

長方形の長辺の長さは 5.23606797749979 である。

算額では「術曰置五個開平方加三個得直長合問」は √5 + 3 = 5.23606797749979 なので，数値解がこれに一致していることがわかる。

その他のパラメータは以下のとおりである。
r2 = 0.5;  a = 5.23607;  b = 2.61803;  c = 2.89443;  r1 = 0.618034;  y1 = 1

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

算額（その881）

算額（その881）

六四 加須市不動岡 総願寺 慶応2年(1866)
埼玉県立図書館：埼玉県史料集 第二集『埼玉の算額』，昭和44年，誠美堂印刷所，埼玉県与野市.

直線の上に乙円 4 個が載り，その上に甲円，丙円，丁円が 2 個ずつ載っている。乙円の直径が 2 寸のとき，甲円の直径はいかほどか。

甲円の半径と中心座標を r1, (2r2, y1)
乙円の半径と中心座標を r2, (r2, r2)
丙円の半径と中心座標を r3, (0, y3)
丁円の半径と中心座標を r4, (r4, y1)
とおき，以下の連立方程式を解く。

include("julia-source.txt");

using SymPy
@syms d, r1::positive, r1::positive, y1::positive,
     r2::positive, r3::positive, y3::positive, r4::positive
eq1 = r2^2 + (y1 - r2)^2 - (r1 + r2)^2
eq2 = 4r2^2 + (y1 - y3)^2 - (r1 + r3)^2
eq3 = r2^2 + (y3 - r2)^2 - (r2 + r3)^2
eq4 = r4^2 + (y1 - y3)^2 - (r3 + r4)^2
eq5 = r1 + 2r4 - 2r2
solve([eq1, eq2, eq3, eq4, eq5], (r1, y1, r3, y3, r4))

   1-element Vector{NTuple{5, Sym{PyCall.PyObject}}}:
    (3*r2/2, r2*(2 + sqrt(21))/2, 7*r2/10, r2*(1 + 3*sqrt(21)/10), r2/4)

甲円の半径は乙円の半径の 3/2 倍である。
乙円の直径が 2 寸のとき，甲円の直径は 3 寸である。

その他のパラメータは以下のとおりである。
r1 = 1.5;  y1 = 3.29129;  r2 = 1;  r3 = 0.7;  y3 = 2.37477;  r4 = 0.25

function draw(more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   r2 = 2//2
   (r1, y1, r3, y3, r4) = r2 .* (3/2, (2 + √21)/2, 7/10, 1 + 3√21/10, 1/4)
   @printf("乙円の直径が %g のとき，甲円の直径は %g である\n", 2r2, 2r1)
   @printf("r1 = %g;  y1 = %g;  r2 = %g;  r3 = %g;  y3 = %g;  r4 = %g\n", r1, y1, r2, r3, y3, r4)
   plot()
   circle2(r2, r2, r2)
   circle2(3r2, r2, r2)
   circle2(2r2, y1, r1, :green)
   circle2(r4, y1, r4, :blue)
   circle(0, y3, r3, :magenta)
   circle(0, 2y1 -y3, 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(2r2, y1, "甲円:r1,(2r2,y1)", :green, :center, delta=-delta)
       point(r2, r2, "乙円:r2\n(r2,r2)", :red, :center, delta=-delta)
       point(3r2, r2, "乙円:r2\n(3r2,r2)", :red, :center, delta=-delta)
       point(0, y3, "丙円:r3\n(0,y3)", :magenta, :center, delta=-delta)
       point(r4, y1, "丁円:r4\n(r4,y1)", :blue, :center, :bottom, delta=5delta, deltax=-3delta)
   end
end;

算額（その880）

算額（その880）

六四 加須市不動岡 総願寺 慶応2年(1866)
埼玉県立図書館：埼玉県史料集 第二集『埼玉の算額』，昭和44年，誠美堂印刷所，埼玉県与野市.

二等辺三角形内に，甲円 1 個，乙円 4 個，丙円 2 個が入っている。乙円の直径が 1 寸のとき，甲円の直径はいかほどか。

二等辺三角形の底辺と高さを 2a，b
甲円の半径と中心座標を r1, (0, y1)
乙円の半径と中心座標を r2, (r2, r2), (3r2, r2)
丙円の半径と中心座標を r3, (x3, y3)
とおき，以下の連立方程式の数値解を求める。

include("julia-source.txt");

using SymPy
@syms d, a::positive, b::positive,
     r1::positive, y1::positive,
     r2::positive, r3::positive, y3::positive
eq1 = numerator(apart(dist(a, 0, 0, b, 0, y1) - r1^2, d))
eq2 = numerator(apart(dist(a, 0, 0, b, 3r2, r2) - r2^2, d))
eq3 = numerator(apart(dist(a, 0, 0, b, 2r2, y3) - r3^2, d))
eq4 = r2^2 + (y1 - r2)^2 - (r1 + r2)^2 |> expand
eq5 = 4r2^2 + (y1 - y3)^2 - (r1 + r3)^2 |> expand
eq6 = r2^2 + (y3 - r2)^2 - (r2 + r3)^2 |> expand;

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)
   (a, b, r1, y1, r3, y3) = u
   return [
       a^2*b^2 - 2*a^2*b*y1 - a^2*r1^2 + a^2*y1^2 - b^2*r1^2,  # eq1
       a^2*b^2 - 2*a^2*b*r2 - 6*a*b^2*r2 + 6*a*b*r2^2 + 8*b^2*r2^2,  # eq2
       a^2*b^2 - 2*a^2*b*y3 - a^2*r3^2 + a^2*y3^2 - 4*a*b^2*r2 + 4*a*b*r2*y3 + 4*b^2*r2^2 - b^2*r3^2,  # eq3
       -r1^2 - 2*r1*r2 + r2^2 - 2*r2*y1 + y1^2,  # eq4
       -r1^2 - 2*r1*r3 + 4*r2^2 - r3^2 + y1^2 - 2*y1*y3 + y3^2,  # eq5
       r2^2 - 2*r2*r3 - 2*r2*y3 - r3^2 + y3^2,  # eq6
   ]
end;

r2 = 1/2
iniv = BigFloat[3, 3, 1, 1.5, 1, 2]
res = nls(H, ini=iniv)

   ([2.726702475495804, 2.665648262572091, 0.73841681234051, 1.6329943517456873, 0.3542486889354094, 1.1926332525571277], true)

乙円の直径が 1 のとき，甲円の直径は 1.47683362468102 である。

算額では「術曰置二十八個開平方加八個除以球個得甲円径合問」で，数式で表すと (sqrt(28) + 8)/9 = 1.47683362468102 となり，数値解が正しいことがわかる。

その他のパラメータは以下のとおりである。

a = 2.7267;  b = 2.66565;  r1 = 0.738417;  y1 = 1.63299;  r3 = 0.354249;  y3 = 1.19263

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