算額(その457)
長野県長野市 久保寺観音堂 享和3年(1803)
中村信弥「改訂増補 長野県の算額」
http://www.wasan.jp/zoho/zoho.html
円に正三角形,正三角形に正方形が内接している。円の直径,正三角形の辺,正方形の辺のそれぞれの 3 乗の和を A,円の直径,正方形の辺の 3 乗の差を B とするとき,円の直径を A, B で表わせ。
include("julia-source.txt")
using SymPy
円の半径,正三角形の辺,正方形の辺をそれぞれ r, 2a, 2c とおく。
正方形の下辺と上辺がy軸と交差するy座標を b, d とする。
a = r√3/2, b = -r/2 から,c, d を求める。
@syms a, b, c, d, r
a = r*√Sym(3)/2
b = -r/2
eq1 = d - b - 2c
eq2 = 2c/(a - c) - sqrt(Sym(3))
solve([eq1, eq2], (c, d))
Dict{Any, Any} with 2 entries:
d => -3*sqrt(3)*r + 11*r/2
c => -3*sqrt(3)*r/2 + 3*r
r, a, c と A, B の関係式を記述する。
@syms r1::positive, r2::positive,
r3::positive, x3::positive,
A::positive, B::positive;
a = r*√Sym(3)/2
c = 3r*(2 - √Sym(3))/2
eq1 = (2r)^3 + (2a)^3 + (2c)^3 - A
eq2 = (2r)^3 - (2c)^3 - B;
eq1 + eq2 としてA, B, r の関係式を得る。
eq3 = eq1 + eq2
eq3 |> println
-A - B + 3*sqrt(3)*r^3 + 16*r^3
eq3 を r について解く(r を A, B で表す式を求める)。
solve(eq3, r)[1] |> factor |> println
(A + B)^(1/3)/(3*sqrt(3) + 16)^(1/3)
(中身)^(1/3) の形なので,3 乗根の中身を簡約化する。
@syms x, A::positive, B::positive;
apart((A + B)/(3sqrt(Sym(3)) + 16), x) |> simplify |>factor |> println
-(-16 + 3*sqrt(3))*(A + B)/229
求める関係式は (16 - 3√3)*(A + B)/229 の 3 乗根である。
求まるのは半径なので,2倍して直径にする関数 f を定義する。
f(A, B) = 2(A + B)^(1/3)/(3*sqrt(3) + 16)^(1/3);
検算してみる。
r = 1; a = r*√3/2; c = 3r*(2 - √3)/2 のとき,
A = 13.715575357311328; B = 7.480577065395304
である。
r = 1
a = r*√3/2
c = 3r*(2 - √3)/2
A = (2r)^3 + (2a)^3 + (2c)^3
B = (2r)^3 - (2c)^3
(A, B) |> println
(13.715575357311328, 7.480577065395304)
計算された A, B から f() によって円の直径を求める。
確かに 2.0 になる。
f(13.715575357311328, 7.480577065395304)
2.0
using Plots
function draw(more=false)
pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
r = 1
a = r*√3/2
b = -r/2
c = 3r*(2 - √3)/2
d = r*(11 - 6*√3)/2
plot([a, 0, -a, a], [-r/2, r, -r/2, -r/2], color=:blue, lw=0.5)
plot!([c, c, -c, -c, c], [-r/2, d, d, -r/2, -r/2], color=:green, lw=0.5)
circle(0, 0, r)
if more
delta = (fontheight = (ylims()[2]- ylims()[1]) / 500 * 10 * 2) / 3 # size[2] * fontsize * 2
point(r, 0, "r ", :red, :right, :bottom, delta=delta/2)
point(0, b, " b", :blue, :left, :bottom, delta=delta/2)
point(0, d, " d", :green, :left, :bottom, delta=delta/2)
point(c, 0, "c ", :green, :right, :bottom, delta=delta/2)
point(a, 0, "a ", :blue, :right, :bottom, delta=delta/2)
hline!([0], color=:gray, lw=0.5)
vline!([0], color=:gray, lw=0.5)
else
plot!(showaxis=false)
end
end;
※コメント投稿者のブログIDはブログ作成者のみに通知されます