算額（その788）
享和3年癸亥五月 丸山良玄門人 豫州松山 大西佐兵衛義全
藤田貞資(1807)：続神壁算法
http://www.wasan.jp/jinpeki/zokujinpekisanpo.pdf
平田浩一：愛媛の和算と算額，日本数学会秋季総合分科会市民講演会，2013/9/23.
https://www.mathsoc.jp/publication/tushin/1804/1804hirata.pdf
弧環減球に中球 2 個と小球 1 個を除いた体積を求めよ。
弧環減球は円を両側から 2 円でくり抜き，それを回転させて得られる回転体である。図の太い黒で示した図形を y 軸を中心として回転する。上の図は，それを三次元表示したものである。
曲線が切り替わる点 (x0, y0) を求める。
include("julia-source.txt");
# julia-source.txt ソース https://blog.goo.ne.jp/r-de-r/e/ad3a427b84bb416c4f5b73089ae813cf
using SymPy
@syms R, r1, r2, r3, x0, y0
(r2, r3) = (3, 2) .// 2
R = 2r2 + r3
eq1 = (r1 + r3)^2 + (r3 + r2)^2 - (r1 + r2)^2
eq2 = (x0 - r1 - r3)^2 + y0^2 - r1^2
eq3 = x0^2 + y0^2 - R^2
res = solve((eq1, eq2, eq3), (r1, x0, y0))
2-element Vector{Tuple{Sym{PyCall.PyObject}, Sym{PyCall.PyObject}, Sym{PyCall.PyObject}}}:
(5, 9/4, -5*sqrt(7)/4)
(5, 9/4, 5*sqrt(7)/4)
[0, 5√7/4] で回転体の体積 res1 を求める。
@syms x
res1= integrate(PI*(-sqrt(25 - x^2) + 6)^2, (x, 0, res[2][3]))
res1 |> println
res1.evalf() |> println
pi*(-150*asin(sqrt(7)/4) + 8365*sqrt(7)/192)
21.5487629153420
[5√7/4, 4]で回転体の体積 res2 を求める。
res2 = integrate(PI*sqrt(16 - x^2)^2, (x, res[2][3], 4))
res2 |> println
res2.evalf() |> println
-2965*sqrt(7)*pi/192 + 128*pi/3
5.68345793167528
3 個の球の体積の和 res3 を求める。
res3 = 4PI/3*(2*(3//2)^3 + (2//2)^3)
res3 |> println
res3.evalf() |> println
31*pi/3
32.4631240870945
求める立体の体積 res を求める。
res = 2res1 + 2res2 - res3
res |> println
res |> simplify |> println
res.evalf() |> println
-2965*sqrt(7)*pi/96 + 2*pi*(-150*asin(sqrt(7)/4) + 8365*sqrt(7)/192) + 75*pi
75*pi*(-16*asin(sqrt(7)/4) + 4 + 3*sqrt(7))/4
22.0013176069401
求める体積は，75π(4 + 3√7 - 16asin(√7/4))/4 = 22.001317606940137 である。
75π*(4 + 3√7 - 16asin(√7/4))/4
22.001317606940137
function draw(more)
pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
(r2, r3) = (3, 2) .// 2
r1 = 5
R = 2r2 + r3
println("R = $R")
(r1, x0, y0) = (5, 9/4, 5*sqrt(7)/4)
θ1 = atand(y0, x0)
θ2 = atand(y0, r1 + r3 - x0)
plot()
circle(0, 0, R, :palevioletred1)
circle(0, 0, R, :black, beginangle=θ1, endangle=90, lw=2)
circle(0, 0, R, :black, beginangle=270, endangle=360-θ1, lw=2)
circle22(0, r3 + r2, r2, :magenta)
circle(0, 0, r3, :orange)
circle(r1 + r3, 0, r1, :blue)
circle(r1 + r3, 0, r1, :black, beginangle=180-θ2, endangle=180+θ2, lw=2)
segment(0, r3 + 2r2, 0, -r3 - 2r2, :black, lw=2)
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(r1 + r3, 0, "r1+r3")
point(0, r3, " r3", :black, :left, :bottom, delta=delta)
point(0, r3+ r2, " r3+r2", :magenta, :left, :vcenter)
point(x0, y0, " (x0,y0)", :black, :left, :vcenter, mark=false)
end
end;
