算額（その421）

算額（その421）

長野県伊那市高遠町藤沢 諏訪大明神社 寶歴元年(1751)

中村信弥(2001)：幻の算額
http://www.wasan.jp/maborosi/maborosi.html

問題　直角三角形があり，その面積を S，(鈎 + 股 + 弦)/(股 - 鈎) = A として，各辺の長さを求めよ。

算額および中村も，S, A をもとにして鈎の 4 次方程式を解けばよいとしているが，ピタゴラスの定理を第 3 の条件式とすれば，SymPy で連立方程式を解くことで解が（一度に）求まる。

include("julia-source.txt");

using SymPy

@syms 鈎::positive, 股::positive, 弦::positive, S::positive, A::positive;

eq1 = 股*鈎/2 - S
eq2 = (鈎 + 股 + 弦)/(股 - 鈎) - A
eq3 = 鈎^2 + 股^2 - 弦^2

res = solve([eq1, eq2, eq3], (鈎, 股, 弦))

   4-element Vector{Tuple{Sym, Sym, Sym}}:
    (sqrt(S)*(-A*sqrt((A^5 - 5*A^3 - 2*A^2*sqrt(2*A^2 + 1) + 4*A + 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)) - A*sqrt((2*A^7 - 9*A^5 - 4*A^4*sqrt(2*A^2 + 1) + 3*A^3 + 14*A^2*sqrt(2*A^2 + 1) + 4*A + 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)) - sqrt((A^5 - 5*A^3 - 2*A^2*sqrt(2*A^2 + 1) + 4*A + 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)) + sqrt((2*A^7 - 9*A^5 - 4*A^4*sqrt(2*A^2 + 1) + 3*A^3 + 14*A^2*sqrt(2*A^2 + 1) + 4*A + 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)))/(sqrt(A)*(A^2 + 1)), sqrt(S)*(A*sqrt((A^5 - 5*A^3 - 2*A^2*sqrt(2*A^2 + 1) + 4*A + 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)) - A*sqrt((2*A^7 - 9*A^5 - 4*A^4*sqrt(2*A^2 + 1) + 3*A^3 + 14*A^2*sqrt(2*A^2 + 1) + 4*A + 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)) - sqrt((A^5 - 5*A^3 - 2*A^2*sqrt(2*A^2 + 1) + 4*A + 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)) - sqrt((2*A^7 - 9*A^5 - 4*A^4*sqrt(2*A^2 + 1) + 3*A^3 + 14*A^2*sqrt(2*A^2 + 1) + 4*A + 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)))/(sqrt(A)*(A^2 + 1)), 2*sqrt(S)*sqrt((A^3 - A - 2*sqrt(2*A^2 + 1))/(A^2 - 4))/sqrt(A))
    (sqrt(S)*(A*sqrt((A^5 - 5*A^3 - 2*A^2*sqrt(2*A^2 + 1) + 4*A + 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)) + A*sqrt((2*A^7 - 9*A^5 - 4*A^4*sqrt(2*A^2 + 1) + 3*A^3 + 14*A^2*sqrt(2*A^2 + 1) + 4*A + 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)) + sqrt((A^5 - 5*A^3 - 2*A^2*sqrt(2*A^2 + 1) + 4*A + 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)) - sqrt((2*A^7 - 9*A^5 - 4*A^4*sqrt(2*A^2 + 1) + 3*A^3 + 14*A^2*sqrt(2*A^2 + 1) + 4*A + 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)))/(sqrt(A)*(A^2 + 1)), sqrt(S)*(-A*sqrt((A^5 - 5*A^3 - 2*A^2*sqrt(2*A^2 + 1) + 4*A + 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)) + A*sqrt((2*A^7 - 9*A^5 - 4*A^4*sqrt(2*A^2 + 1) + 3*A^3 + 14*A^2*sqrt(2*A^2 + 1) + 4*A + 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)) + sqrt((A^5 - 5*A^3 - 2*A^2*sqrt(2*A^2 + 1) + 4*A + 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)) + sqrt((2*A^7 - 9*A^5 - 4*A^4*sqrt(2*A^2 + 1) + 3*A^3 + 14*A^2*sqrt(2*A^2 + 1) + 4*A + 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)))/(sqrt(A)*(A^2 + 1)), -2*sqrt(S)*sqrt(-(-A^3 + A + 2*sqrt(2*A^2 + 1))/(A^2 - 4))/sqrt(A))
    (sqrt(S)*(-A*sqrt((A^5 - 5*A^3 + 2*A^2*sqrt(2*A^2 + 1) + 4*A - 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)) + A*sqrt((2*A^7 - 9*A^5 + 4*A^4*sqrt(2*A^2 + 1) + 3*A^3 - 14*A^2*sqrt(2*A^2 + 1) + 4*A - 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)) - sqrt((A^5 - 5*A^3 + 2*A^2*sqrt(2*A^2 + 1) + 4*A - 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)) - sqrt((2*A^7 - 9*A^5 + 4*A^4*sqrt(2*A^2 + 1) + 3*A^3 - 14*A^2*sqrt(2*A^2 + 1) + 4*A - 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)))/(sqrt(A)*(A^2 + 1)), sqrt(S)*(A*sqrt((A^5 - 5*A^3 + 2*A^2*sqrt(2*A^2 + 1) + 4*A - 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)) + A*sqrt((2*A^7 - 9*A^5 + 4*A^4*sqrt(2*A^2 + 1) + 3*A^3 - 14*A^2*sqrt(2*A^2 + 1) + 4*A - 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)) - sqrt((A^5 - 5*A^3 + 2*A^2*sqrt(2*A^2 + 1) + 4*A - 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)) + sqrt((2*A^7 - 9*A^5 + 4*A^4*sqrt(2*A^2 + 1) + 3*A^3 - 14*A^2*sqrt(2*A^2 + 1) + 4*A - 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)))/(sqrt(A)*(A^2 + 1)), 2*sqrt(S)*sqrt((A^3 - A + 2*sqrt(2*A^2 + 1))/(A^2 - 4))/sqrt(A))
    (sqrt(S)*(A*sqrt((A^5 - 5*A^3 + 2*A^2*sqrt(2*A^2 + 1) + 4*A - 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)) - A*sqrt((2*A^7 - 9*A^5 + 4*A^4*sqrt(2*A^2 + 1) + 3*A^3 - 14*A^2*sqrt(2*A^2 + 1) + 4*A - 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)) + sqrt((A^5 - 5*A^3 + 2*A^2*sqrt(2*A^2 + 1) + 4*A - 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)) + sqrt((2*A^7 - 9*A^5 + 4*A^4*sqrt(2*A^2 + 1) + 3*A^3 - 14*A^2*sqrt(2*A^2 + 1) + 4*A - 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)))/(sqrt(A)*(A^2 + 1)), sqrt(S)*(-A*sqrt((A^5 - 5*A^3 + 2*A^2*sqrt(2*A^2 + 1) + 4*A - 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)) - A*sqrt((2*A^7 - 9*A^5 + 4*A^4*sqrt(2*A^2 + 1) + 3*A^3 - 14*A^2*sqrt(2*A^2 + 1) + 4*A - 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)) + sqrt((A^5 - 5*A^3 + 2*A^2*sqrt(2*A^2 + 1) + 4*A - 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)) - sqrt((2*A^7 - 9*A^5 + 4*A^4*sqrt(2*A^2 + 1) + 3*A^3 - 14*A^2*sqrt(2*A^2 + 1) + 4*A - 8*sqrt(2*A^2 + 1))/(A^4 - 8*A^2 + 16)))/(sqrt(A)*(A^2 + 1)), -2*sqrt(S)*sqrt((A^3 - A + 2*sqrt(2*A^2 + 1))/(A^2 - 4))/sqrt(A))

3番めの組が適解である。

res[3][1] |> simplify |> display

res[3][2] |> simplify |> display

res[3][3] |> simplify |> display

解が適正であることを，例で見てみよう。
鈎, 股, 弦がそれぞれ，3.9, 5.2, 6.5 のとき，S = 10.14, A = 12 である。

(鈎, 股, 弦) = (3.9, 5.2, 6.5)
println("S = $(股*鈎/2),  A = $((鈎 + 股 + 弦)/(股 - 鈎))")

   S = 10.14,  A = 11.999999999999996

res[3][1], res[3][2], res[3][3] に S = 10.14, A = 12 を代入すると，鈎，股，弦が求まる。

println([res[3][i](S => 10.14, A => 12).evalf() for i in 1:3])

   Sym[3.90000000000000, 5.20000000000000, 6.50000000000000]