算額（その190）

算額（その190）

中村信弥「改訂増補 長野県の算額」
http://www.wasan.jp/zoho/zoho.html
県内の算額1（87）
長野県中野市田上 田上観音堂 文化6年(1809)
県内の算額1（64）
長野県下高井郡木島平村往郷 水穂神社 寛政12年(1800)
県内の算額1（46）
長野県長野市若穂 清水寺観音堂 寛政6年(1794)

第6問 鉤股の中に，全，甲，乙，丙，丁の正方形を入れる。丙，丁の一辺が 192 寸, 108 寸 のとき，全の一辺はいかほどか。

全，甲，乙，丙，丁 の正方形の一辺を a, b, c, d, e とする。c = 192, e = 108 で，構成される複数の直角三角形の鉤股の比を t =　AC / AB = (a + d + e - c)/(a + b + c - e) として，連立方程式を解く。

using SymPy

@syms a::positive, b::positive, d::positive;
c = 192
e = 108
t = (a + d + e - c)/(a + b + c - e)
eq1 = t - (b - c)/c
eq2 = t - e/(d - e)
eq3 = t - (a - c)/(b + c)

res = solve([eq1, eq2, eq3], (a, b, d))

   1-element Vector{Tuple{Sym, Sym, Sym}}:
    (588, 336, 252)

全，甲，乙の一辺の長さはそれぞれ 588 寸，336 寸，252 寸である。

ちなみに，鉤股弦は 3:4:5, t = 3/4

using Plots
using Printf

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 segment(x1, y1, x2, y2, color=:black; linestyle=:solid, linewidth=0.5)
   plot!([x1, x2], [y1, y2], color=color, linestyle=linestyle, linewidth=linewidth)
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 draw(more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   (a, b, d) = res[1]
   (c, e) = (192, 108)
   t = (a + d + e - c)/(a + b + c - e)
   println("t = $t")
   plot()
   rect(0, 0, a, a)
   rect(a, 0, a + b, b)
   rect(a + b, 0, a + b + c, c)
   rect(0, a, d, a + d)
   rect(0, a + d, e, a + d + e)
   plot!([0, a + b + c + c/t, 0, 0], [0, 0, a + d + e + e*t, 0], color=:black, lw=0.5)
   if more == true
       point(a, 0, " a", :black, :left, :top)
       point(a + b, 0, " a+b", :black, :left, :top)
       point(a + b + c, 0, " a+b+c", :black, :left, :top)
       point(0, a, "a ", :black, :right, :bottom)
       point(0, a + d, "a+d ", :black, :right, :bottom)
       point(0, a + d + e, "a+d+e ", :black, :right, :bottom)
       plot!([e, a + b + c, e, e], [c, c, a + d + e, c], color=:red, lw=1)
       point(e, c, " A", :black, :left, :bottom)
       point(a + b + c, c, " B", :black, :left, :bottom)
       point(e, a + d + e, " C", :black, :left, :bottom)
       hline!([0], color=:black, lw=0.5, xlims=(-175, 1400))
       vline!([0], color=:black, lw=0.5, ylims=(-70, 1050))
   else
      plot!(showaxis=false)
   end
end;

   t = 3/4

 

算額（その189）

算額（その189）

中村信弥「改訂増補 長野県の算額」
http://www.wasan.jp/zoho/zoho.html
県内の算額1（86）
長野県中野市田上 田上観音堂 文化6年(1809)

第5問 内外2個の菱形により決まる黒の部分の面積が 144 寸（本当は平方寸）
外側，内側の菱形の辺の長さが 17 寸，10 寸のとき共通対角線の長さはいかほどか。

内側の菱形の OB，OA の長さを a，b とする。BC の長さを c とする。

以下の連立方程式を解く。

include("julia-source.txt");

using SymPy

@syms a::positive, b::positive, c::positive
eq1 = a^2 + b^2 - 10^2  # ⊿ AOB についてピタゴラスの定理
eq2 = (a + c)^2 + b^2 - 17^2  # ⊿ AOC についてピタゴラスの定理
eq3 = 2b*(a + c) - 2a*b - 144  # 外側の菱形の面積から内側のひし形の面積を引いたものが黒積

solve([eq1, eq2, eq3], (a, b, c))

   1-element Vector{Tuple{Sym, Sym, Sym}}:
    (6, 8, 9)

共通対角線の長さは 2b = 16

using Plots
using Printf

function draw(more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   (a, b, c) = (6, 8, 9)
   plot([15, 0, -a - c, 0, a + c], [0, b, 0, -b, 0], linecolor=:black, linewidth=0.5, seriestype=:shape, fillcolor=:black)
   plot!([a, 0, -a, 0, a], [0, b, 0, -b, 0], linecolor=:black, linewidth=0.5, seriestype=:shape, fillcolor=:white)
   if more == true
       point(0, 0, " O", :green, :left, :bottom)
       point(a, 0, " B: a", :cyan, :left, :bottom)
       point(a + c, 0, " C: a+c", :red, :left, :bottom)
       point(0, b, " A: b", :green, :left, :bottom)
       hline!([0], color=:black, lw=0.5)
       vline!([0], color=:black, lw=0.5)
       segment(a + c, 0, a, 0, :gray, linewidth=1)
       segment(-a - c, 0, -a, 0, :gray, linewidth=1)
   else
      plot!(showaxis=false)
   end
end;

算額（その188）

算額（その188）

中村信弥「改訂増補 長野県の算額」
http://www.wasan.jp/zoho/zoho.html
県内の算額1（86）
長野県中野市田上 田上観音堂 文化6年(1809)

第4問 鉤股（直角三角形）内に大中小の3円を入れる。大円と中円の径を掛けると 72，大円と小円の径を掛けると 48 であるとき，それぞれの径を求めよ。
注：算額の図は「算額（その187）https://blog.goo.ne.jp/r-de-r/e/2053d2ea5ef9100ee13ecf80457132b9」と違い垂直線，水平線に接するという条件がない場合の解である。

大円，中円，小円の半径を r1,r2,r3 とする。中円，小円の中心座標を (x2, r2), (r3, y3) とする。
中円，小円の中心から弦への距離，鉤股弦の面積の関係式を解く。

using SymPy

function distance(x1, y1, x2, y2, x0, y0)
   p1, p2 = sympy.Point(x1, y1), sympy.Point(x2, y2)
   l = sympy.Line(p1, p2)
   l.distance(sympy.Point(x0, y0))^2  # 二乗距離を返す！！
end;

@syms X::positive, Y::positive, x2::pisitive, y3::positive,
     r1::positive, r2::positive, r3::positive, l::positive;

eq1 = 2r1 * 2r2 - 72
eq2 = 2r1 * 2r3 - 48
eq3 = (x2 - r1)^2 + (r1 - r2)^2 - (r1 + r2)^2
eq4 = (r1 - r3)^2 + (y3 - r1)^2 - (r1 + r3)^2
eq5 = distance(0, Y, X, 0, r3, y3) - r3^2
eq6 = distance(0, Y, X, 0, x2, r2) - r2^2
eq7 = X^2 + Y^2 - l^2
eq8 = (X + Y + l)*r1 - X * Y;
# res = solve([eq1, eq2, eq3, eq4, eq5, eq6, eq7, eq8])

solve() では解けないので，nlsolve() を使用する。

println(eq1, ",")
println(eq2, ",")
println(eq3, ",")
println(eq4, ",")
println(eq5, ",")
println(eq6, ",")
println(eq7, ",")
println(eq8, ",")

   4*r1*r2 - 72,
   4*r1*r3 - 48,
   (-r1 + x2)^2 + (r1 - r2)^2 - (r1 + r2)^2,
   (-r1 + y3)^2 + (r1 - r3)^2 - (r1 + r3)^2,
   -r3^2 + (-X*(X*r3 + Y^2 - Y*y3)/(X^2 + Y^2) + r3)^2 + (-Y*(X^2 - X*r3 + Y*y3)/(X^2 + Y^2) + y3)^2,
   -r2^2 + (-X*(X*x2 + Y^2 - Y*r2)/(X^2 + Y^2) + x2)^2 + (-Y*(X^2 - X*x2 + Y*r2)/(X^2 + Y^2) + r2)^2,
   X^2 + Y^2 - l^2,
   -X*Y + r1*(X + Y + l),

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=1e-14)
       v = r.zero[1]
   else
       r = nlsolve((vout, vin)->vout .= func(vin, params...), ini, ftol=1e-14)
       v = r.zero
   end
   return v, r.f_converged
end;

function H(u)
   (r1, r2, r3, X, Y, l, x2, y3) = u
   return [
4*r1*r2 - 72,
4*r1*r3 - 48,
(-r1 + x2)^2 + (r1 - r2)^2 - (r1 + r2)^2,
(-r1 + y3)^2 + (r1 - r3)^2 - (r1 + r3)^2,
-r3^2 + (-X*(X*r3 + Y^2 - Y*y3)/(X^2 + Y^2) + r3)^2 + (-Y*(X^2 - X*r3 + Y*y3)/(X^2 + Y^2) + y3)^2,
-r2^2 + (-X*(X*x2 + Y^2 - Y*r2)/(X^2 + Y^2) + x2)^2 + (-Y*(X^2 - X*x2 + Y*r2)/(X^2 + Y^2) + r2)^2,
X^2 + Y^2 - l^2,
-X*Y + r1*(X + Y + l),        
      ]
end;

iniv = [big"6.1", 4.2, 3.3, 24, 17, 31, 15, 13]
res = nls(H, ini=iniv)
println(res)

   (BigFloat[5.748772323345360879553386071530534002107771135435234241146154053030068382725282, 3.131103301291523318846022596310963240657188641774061152337086444494206214742888, 2.08740220086101554589734764786891092285620654321669230693296468944296425038439, 24.38365332206977707890963182337988451889403111451717974855438930320079233053448, 16.62684846651738071918895370581186617641229049594184966902009096464686858816747, 29.51295714189643603899180810391979664284976102613029457759886594916383321385045, 14.23405369758393117236351841678872247352580238769692473134187249862639602572867, 12.67697555362087005366317143755402346987899215067675675336938195307101120791379], true)

using Printf
@printf("大円径 = %6.3f; 中円径 = %6.3f; 小円径 = %6.3f\n", 2res[1][1], 2res[1][2], 2res[1][3])

   大円径 = 11.498; 中円径 =  6.262; 小円径 =  4.175

大円，中円，小円の径は元の単位で，11.498 寸，6.262 寸，4.175 寸である。

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, linecolor=:black, linewidth=0.5, seriestype=:shape, fillcolor=color)
   else
       plot!(ox .+ x, oy .+ y, color=color, linewidth=0.25)
   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 segment(x1, y1, x2, y2, color=:black; linestyle=:solid, linewidth=0.5)
   plot!([x1, x2], [y1, y2], color=color, linestyle=linestyle, linewidth=linewidth)
end;

function draw(more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   (r1, r2, r3, X, Y, l, x2, y3) = res[1]
   @printf("X = %6.3f; Y = %6.3f; l = %6.3f; x2 = %6.3f; y3 = %6.3f\n",
       X, Y, l, x2, y3)
   @printf("r1 = %6.3f; r2 = %6.3f; r3 = %6.3f\n", r1, r2, r3)
   plot([0, X, 0, 0], [0, 0, Y, 0], color=:black, lw=0.5)
   circle(r1, r1, r1, :red)
   circle(x2, r2, r2, :blue)
   circle(r3, y3, r3, :green)
   if more == true
       point(r1, r1, "大円(r1,r1,r1)")
       point(x2, r2, "中円(x2,r2,r2)")
       point(r3, y3, "小円(r3,y3,r3)")
       point(0, 0, " O", :green, :left, :bottom)
       point(X, 0, " X", :green, :left, :bottom)
       point(0, Y, " Y", :green, :left, :bottom)
       hline!([0], color=:black, lw=0.5)
       vline!([0], color=:black, lw=0.5)
   else
      plot!(showaxis=false)
   end
end;

   X = 24.384; Y = 16.627; l = 29.513; x2 = 14.234; y3 = 12.677
   r1 =  5.749; r2 =  3.131; r3 =  2.087

算額（その187）

算額（その187）

中村信弥「改訂増補 長野県の算額」
http://www.wasan.jp/zoho/zoho.html
県内の算額1（86）
長野県中野市田上 田上観音堂 文化6年(1809)

第4問 鉤股（直角三角形）内に大中小の3円を入れる。3円はそれぞれ鉤股弦に接し，鉤股弦内の垂直線，水平線にも接している。大円と中円の径を掛けると 72，大円と小円の径を掛けると 48 であるとき，それぞれの径を求めよ。
注：算額の図は垂直線，水平線に接するという条件がない図（「算額（その188）」https://blog.goo.ne.jp/r-de-r/e/674613630ed8ca8d49e2a557b185a4de 参照）になっているが，答えはその図と異なり上述の問いの場合の答えになっている。

大円，中円，小円の半径を r1,r2,r3 とする。AX = a, BX = b, CY = c, DY = d とする。
⊿AOD, ⊿AXB, ⊿CYD は相似で，相似比が r1:r2:r3 である。
また，AD^2 = l^2 = AO^2 + DO^2, 直角三角形の面積の関係について方程式を立てて解く。

using SymPy

@syms a::positive, b::positive, c::positive, d::positive,
     r1::positive, r2::positive, r3::positive, l::positive;

eq1 = 2r1 * 2r2 - 72
eq2 = 2r1 * 2r3 - 48
eq3 = d/b - c/a
eq4 = d/b - r3/r2
eq5 = d/(d + 2r1) - c/(a + 2r1)
eq6 = d/(d + 2r1) - r3/r1
eq7 = b/(d + 2r1) - r2/r1
eq8 = (d + 2r1)^2 + (a + 2r1)^2 - l^2
eq9 = (d + 2r1 + a + 2r1 + l)*r1 - (d +2r1) * (a + 2r1)
res = solve([eq1, eq2, eq3, eq4, eq5, eq6, eq7, eq8, eq9],
           (a, b, c, d, r1, r2, r3, l))

   1-element Vector{NTuple{8, Sym}}:
    (12, 9, 8, 6, 6, 3, 2, 30)

大円，中円，小円の半径は 6, 3, 2 で，元の単位でいえば 直径 12 寸, 6 寸, 4 寸 である。

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, linecolor=:black, linewidth=0.5, seriestype=:shape, fillcolor=color)
   else
       plot!(ox .+ x, oy .+ y, color=color, linewidth=0.25)
   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 segment(x1, y1, x2, y2, color=:black; linestyle=:solid, linewidth=0.5)
   plot!([x1, x2], [y1, y2], color=color, linestyle=linestyle, linewidth=linewidth)
end;

function draw(more)
    pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   (a, b, c, d, r1, r2, r3, l) = res[1]
   @printf("a = %6.3f; b = %6.3f; c = %6.3f; d = %6.3f\n", a, b, c, d)
   @printf("r1 = %6.3f; r2 = %6.3f; r3 = %6.3f\n", r1, r2, r3)
   plot([0, a + 2r1, 0, 0], [0, 0, d + 2r1, 0], color=:black, lw=0.5)
   circle(r1, r1, r1, :red)
   circle(2r1 + r2, r2, r2, :blue)
   circle(r3, 2r1 + r3, r3, :green)
   segment(0, 2r1, c, 2r1, :black, linewidth=0.25)
   segment(2r1, 0, 2r1, b, :black, linewidth=0.25)
   if more == true
       point(0, 0, " O")
       point(2r1, 0, " X")
       point(2r1+a, 0, " A")
       point(2r1, b, " B", :green, :left, :bottom)
       point(c, 2r1, " C", :green, :left, :bottom)
       point(0, 2r1, "Y ", :green, :right)
       point(0, 2r1+d, "D ", :green, :right)
       point(r1, r1, "大円(r1,r1,r1)")
       point(2r1+r2, r2, "中円(2r1+r2,r2,r2)")
       point(r3, 2r1+r3, "小円(r3, 2r1+r3,r3)")
       annotate!(2r1, 2r1+r3, text("AX = a\nBX = b\nCY = c\nDY = d\nAD = l", 10, :left, :green))
       hline!([0], color=:black, lw=0.5)
       vline!([0], color=:black, lw=0.5)
   else
      plot!(showaxis=false)
   end
end;

   a = 12.000; b =  9.000; c =  8.000; d =  6.000
   r1 =  6.000; r2 =  3.000; r3 =  2.000