福島県三春町 龍穏院 明治26年(1893)
http://www.wasan.jp/fukusima/ryuonin.html
街角の数学 Street Wasan ～落書き帳「○△□」～ 216. 二円の縁
街角の数学 Street Wasan ～落書き帳「○△□」～ 217. 二円の縁（その2）
http://streetwasan.web.fc2.com/math17.2.1.html
http://streetwasan.web.fc2.com/math17.2.2.html
キーワード：円3個，菱形
#Julia #SymPy #算額 #和算 #数学
大円 2 個が交わっており，中に小円と菱形が入っている。大円の直径が与えられたとき，黒積が最大となるときの菱形の対角線の長い方の長さはいかほどか。
大円の半径と中心座標を \(r_1,\ (r_1 - r_2,\ 0)\)
小円の半径と中心座標を \(r_2,\ (0,\ 0)\)
黒積は，灰色部分 = 扇形ABC - 三角形OAB - 四分円OCD の4倍である。
∠OAB = \(\theta\)，\(\theta = \arcsin(\sqrt{r_2(2r_1 - r_2)}/r_1)\)
黒積は小円の大きさ（半径）によって変化する。まず黒積が最大になるときの小円の半径を求める。
菱形の対角線の長い方の長さは \(4r_1 - 2r_2\) である。
include("julia-source.txt"); # julia-source.txt ソース
using SymPy
@syms r1, r2, y, θ, 円積率
円積率 = 0.79
円周率 = 4円積率
#円周率 = PI
y = sqrt(2r2*r1 - r2^2) # sqrt(r1^2 - (r1 - r2)^2)
θ = asin(y/r1)
S = 4(円周率*r1^2*(θ/2PI) - y*(r1 - r2)/2 - 円周率*r2^2/4)
S |> println
6.32*r1^2*asin(sqrt(2*r1*r2 - r2^2)/r1)/pi - 3.16*r2^2 - 2*(r1 - r2)*sqrt(2*r1*r2 - r2^2)
pyplot(size=(300, 150), grid=false, aspectratio=:none, label="")
plot(S(r1 => 1/2), xlims=(0, 0.5), xlabel="r2", ylabel="黒積")
vline!([0.28746844208856126])
g = diff(S, r2)
g |> println
6.32*r1*(r1 - r2)/(pi*sqrt(1 - (2*r1*r2 - r2^2)/r1^2)*sqrt(2*r1*r2 - r2^2)) - 6.32*r2 - 2*(r1 - r2)^2/sqrt(2*r1*r2 - r2^2) + 2*sqrt(2*r1*r2 - r2^2)
plot(g(r1 => 1/2), xlims=(0, 0.35), xlabel="r2", ylabel="黒積の導関数")
hline!([0])
solve(g(r1 => 1/2), r2)[1] |> println
0.287468442088561
S(r1 => 1/2, r2 => 0.287468442088561).evalf() |> println
0.115685731044897
大円の半径が \(r_1 = 1/2\) のときに，小円の半径が \(r_2 = 0.287468442088561\) のとき，黒積は最大値 \(0.115685731044897\) をとる。
菱長 = (4r1 - 2r2)(r1 => 1/2, r2 => 0.287468442088561) |> println
1.42506311582288
大円の半径が \(r_1 = 1/2\) のときに，小円の半径が \(r_2 = 0.287468442088561\) のとき，菱長は \(1.42506311582288\) である。
術は \(2r_1/(0.125/0.79^2 + 0.5)\) で，\(r_1 = 1/2\) のとき \(1.42798306829882\) となるが，違いがどこから来るのかわからない。
(2r1/(0.125/0.79^2 + 0.5))(r1 => 1/2) |> println
1.42798306829882
描画関数プログラムのソースを見る
function draw(r1, r2, more)
pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
y = sqrt(r1^2 - (r1 - r2)^2)
θ = asind(y/r1)
println("θ = ", θ)
θ1 = 180 - θ:0.1:180
xa = @. r1*cosd(θ1) + (r1 - r2)
ya = @. r1*sind(θ1)
θ2 = 180:-0.1:90
append!(xa, r2.*cosd.(θ2))
append!(ya, r2.*sind.(θ2))
append!(xa, xa[1])
append!(ya, ya[1])
plot([2r1 - r2, 0, r2 - 2r1, 0, 2r1 - r2], [0, y, 0, -y, 0], color=:green, lw=0.5)
plot!(xa, ya, color=:orange, lw=0.5, seriestype=:shape, fillcolor=:gray70)
circle2(r1 - r2, 0, r1)
circle(0, 0, r2, :blue)
segment(0, y, r1 - r2, 0, :gray70)
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(0, 0, "O", :blue, :center, delta=-delta)
point(0, r2, "D", :blue, :center, delta=-delta)
point(-r2, 0, "C ", :blue, :right, :vcenter)
point(0, sqrt(r2*(2r1 - r2)), " B:√(r2*(2r1-r2))", :red, :left, :vcenter)
point(r1 - r2, 0, " A:r1-r2 ", :red, :left, :vcenter)
point(2r1 - r2, 0, "2r1-r2 ", :green, :right, :vcenter)
point(r2, 0, " r2", :blue, :left, :vcenter)
end
end;
draw(1/2, 0.15, true)
以下のアイコンをクリックして応援してください
