算額（その1446）
十一 金ケ崎町 個人宅 元文6年(1741)
山村善夫：現存 岩手の算額，昭和52年1月30日，熊谷印刷，盛岡市.
http://www.wasan.jp/yamamura/yamamura.html
キーワード：正方形3個，円1個，正三角形1個，直角三角形1個
#Julia, #SymPy, #算額, #和算
ごちゃごちゃしているので，分かりやすく書き下すと，以下の 7 式が述べられている。
(1) 外側の大正方形の面積から，中正方形，小正方形，円，正三角形，直角三角形の面積を除いた面積（外積）と大方面の4乗の和が 65580.185
(2) 大正方形，中正方形，小正方形の一辺の長さ（大方面，中方面，小方面），円の直径（円径），直角三角形の直角を挟む二辺の短い方（鈎），正三角形の一辺の長さ（三角面）の 6 つを加えると 51
(3) 円径 - 中方面 = 1
(4) 中方面 - 小方面 = 1
(5) 小方面 - 鈎 = 1
(6) 鈎 - 三角面 = 1
(7) 股 - 鈎 = 2
である。小方面（およびその他のパラメータ）はいかほどか。ただし，円積率として 0.79 ≒ π/4 = 0.7853981633974483，三角積率として 0.43 ≒ √3/4 = 0.4330127018922193 を使え。
注1：直径が5の円の面積は 5^2*0.79 = 19.75 と計算する。正確には π*2.5^2 = 19.634954084936208
注2：一辺が5の正三角形の面積は 5^2*0.43 = 10.75 と計算する。正確には 5^2*√3/4 = 10.825317547305483
問には 0.43 を使えと書いているが，0.433 を使えば 10.825 になり，精度は格段に上がるのになぜ？以下では 0.433 を使うと期待される答えが高精度で得られる。また，山村も 0.433 を使っている。
以上 (1) 〜 (7) を連立方程式として解く。
using SymPy
@syms 大方面, 中方面, 小方面, 円径, 三角面, 鈎, 股
eq1 = (大方面^2 - 円径^2*0.79 - 中方面^2 - 小方面^2 - 三角面^2*0.433 - 鈎*股/2) + 大方面^4 - 65580.185
eq2 = 大方面 + 中方面 + 小方面 + 円径 + 鈎 + 三角面 - 51
eq3 = (円径 - 中方面) - 1
eq4 = (中方面 - 小方面) - 1
eq5 = (小方面 - 鈎) - 1
eq6 = (鈎 - 三角面) - 1
eq7 = (股 - 鈎) - 2
res = solve([eq1, eq2, eq3, eq4, eq5, eq6, eq7], (大方面, 中方面, 小方面, 円径, 三角面, 鈎, 股))[1] # of 4
(16.0000000000000, 8.00000000000000, 7.00000000000000, 9.00000000000000, 5.00000000000000, 6.00000000000000, 8.00000000000000)
三角積率として 0.433 を使った。
その結果として，求められた 7 個のパラメータは全て（誤差のない）整数値になった。というか，結果が整数値になるように 65580.185 のように妙な数値を持ち込んだのであろう。
小正方形の一辺の長さ（小方面）は 7 である。
三角積率として 0.43 を採用すると得られる解は，
(15.9999954343728, 8.00000091312543, 7.00000091312543, 9.00000091312543, 5.00000091312543, 6.00000091312543, 8.00000091312543)
になる。
算額出題者は三角積率として 0.433 を使ったのではなかろうか?
山村の解説では，術は小方面を未知数として高次方程式を解いているようだ。
eq2 ~ eq7 を連立方程式として，大方面, 中方面, 円径, 三角面, 鈎, 股 を求めると，それぞれは小方面を含む式で表される。
res = solve([eq2, eq3, eq4, eq5, eq6, eq7], (大方面, 中方面, 円径, 三角面, 鈎, 股));
得られた解を eq1 に代入する。
eq1 = eq1(円径 => res[円径], 股 => res[股], 中方面 => res[中方面], 三角面 => res[三角面], 大方面 => res[大方面], 鈎 => res[鈎])
eq1 |> println
-小方面^2 + (51 - 5*小方面)^4 + (51 - 5*小方面)^2 - 0.433*(小方面 - 2)^2 - (小方面 - 1)*(小方面 + 1)/2 - (小方面 + 1)^2 - 0.79*(小方面 + 2)^2 - 65580.185
簡約化すると，小方面についての 4 次式が得られる。
eq1 = simplify(eq1)
eq1 |> println
625.0*小方面^4 - 25500.0*小方面^3 + 390171.277*小方面^2 - 2653533.428*小方面 + 6702216.423
方程式を解いて得られる解は，連立方程式を解いて得られるものと同じである。
solve(eq1, 小方面)[1] |> println
7.00000000000000
係数を -2 倍して，山村の式と比較してみよう。
-2eq1 |> println
-1250.0*小方面^4 + 51000.0*小方面^3 - 780342.554*小方面^2 + 5307066.856*小方面 - 13404432.846
山村（及び術？）は，x^2 の係数として -780342.555 としているがこれは 780342.554 である。
@syms x
eq0 = -13404432.846 + 5307066.856x - 780342.555x^2 + 51000x^3 - 1250x^4
eq0 |> println
-1250*x^4 + 51000*x^3 - 780342.555*x^2 + 5307066.856*x - 13404432.846
この違いは大きい。解を求めると，x = 7 となるべきところが 7.00000029828747 になる。
solve(eq0, x)[1] |> println
7.00000029828747
