裏 RjpWiki

Julia ときどき R, Python によるコンピュータプログラム,コンピュータ・サイエンス,統計学

算額(その2022)

2024年08月16日 | Julia

算額(その2022)

(3) 大阪府茨木市大字総持寺 総持寺 弘化3年(1846)
近畿数学史学会:近畿の算額「数学の絵馬を訪ねて」,平成4年5月16日 初版第一刷,大阪教育図書株式会社,大阪市.
キーワード:円4個,四分円4個,楕円

楕円の中に円 1 個,等円 4 個,四分円 4 個を容れる。楕円の短径が 58 寸のとき,等円の直径はいかほどか。

注:楕円の長径に関しては何の定義もない(何でも良い)。四分円は円と同じ半径を持つ円の 1/4 である。

楕円の短半径を b,等円の半径を r とおく。

方程式を立てるまでもなく,等円の半径 r は短半径 b の (2 - √2)/2 倍である。
短径が 58 寸のとき,等円の直径は 16.98780669118024 寸である。

include("julia-source.txt");

using SymPy

@syms b
r = (b/√Sym(2) - b/2)*√Sym(2) |> simplify
r |> println
r(b => 58/2).evalf() |> println

   b*(2 - sqrt(2))/2
   8.49390334559012

function draw(b, more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   a = 60
   r = b*(2 - √2)/2
   plot()
   ellipse(0, 0, a, b)
   circle4((b - r)/√2, (b - r)/√2, r)
   circle(0, 0, b, :blue)
   circle(b, b, b, :green, beginangle=180, endangle=270)
   circle(-b, b, b, :green, beginangle=270, endangle=360)
   circle(-b, -b, b, :green, beginangle=0, endangle=90)
   circle(b, -b, b, :green, beginangle=90, endangle=180)
   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, b, "b", :green, :center, :bottom, delta=delta)
       point((b - r)/√2, (b - r)/√2, "等円:r,((b-r)/√2,(b-r)/√2)", :red, :left, :bottom, delta=delta,  deltax=-7delta)
       point(b, b, " 四分円:b,(b,b)", :green, :left, :vcenter)
   end
end;

draw(58/2, true)

コメント
  • X
  • Facebookでシェアする
  • はてなブックマークに追加する
  • LINEでシェアする

算額(その2021)

2024年08月16日 | Julia

算額(その2021)

(2) 大阪府豊中市服部元町 服部天神社 天保14年(1843)
近畿数学史学会:近畿の算額「数学の絵馬を訪ねて」,平成4年5月16日 初版第一刷,大阪教育図書株式会社,大阪市.
キーワード:円4個,四分円,長方形

長方形の中に四分円,中円,甲円,乙円,丙円を容れる。四分円,中円の直径が 9 寸,3 寸のとき,甲円,乙円,丙円の直径はいかほどか。

1. デカルトの円定理による方法

四分円,中円,甲円,乙円,丙円の半径を r1, r2, r3, r4, r5 とおき,曲率を k1, k2, k3, k4, k5 とおき,デカルトの円定理を適用する。

r1 = 4.5;  r2 = 1.5;  a = 6.69615;  r3 = 0.602886;  x3 = 3.29423;  r4 = 0.323085;  x4 = 2.41154;  r5 = 0.200962;  x5 = 1.90192

r1 = 9//2
r2 = 3//2
k3 = 1/r1 + 1/r2 + 2sqrt(1/(r1*r2))
r3 = 1/k3
k4 = 1/r1 + 1/r3 + 2sqrt(1/(r1*r3))
r4 = 1/k4
k5 = 1/r1 + 1/r4 + 2sqrt(1/(r1*r4))
r5 = 1/k5
2 .*(r3, r4, r5)

   (1.2057713659400522, 0.6461709275204175, 0.40192378864668404)

四分円,中円の直径が 9 寸,3 寸のとき,甲円,乙円,丙円の直径は 1.2057713659400522 寸, 0.6461709275204175 寸, 0.40192378864668404 寸である。

2. 2円の関係式に基づき,順次パラメータを求める場合

上に加えて,甲円,乙円,丙円の中心の x 座標を定義し,逐次的に以下の連立方程式を解く。

include("julia-source.txt");

using SymPy

@syms a::positive, r1::positive,
     r2::positive, r3::positive, x3::positive,
     r4::positive, x4::positive, r5::positive,
     x5::positive
(r1, r2) = (9, 3) .// 2
eq1 = (a - r2)^2 + (r1 - r2)^2 - (r1 + r2)^2
eq2 = x3^2 + (r1 - r3)^2 - (r1 + r3)^2
eq3 = x4^2 + (r1 - r4)^2 - (r1 + r4)^2
eq4 = x5^2 + (r1 - r5)^2 - (r1 + r5)^2
eq5 = (a - r2 - x3)^2 + (r2 - r3)^2 - (r2 + r3)^2
eq6 = (x3 - x4)^2 + (r3 - r4)^2 - (r3 + r4)^2
eq7 = (x4 - x5)^2 + (r4 - r5)^2 - (r4 + r5)^2;
res = solve([eq1, eq2, eq3, eq4, eq5, eq6, eq7],
   (a, r3, x3, r4, x4, r5, x5))[2]

   (3/2 + 3*sqrt(3), 9/2 - 9*sqrt(3)/4, -9/2 + 9*sqrt(3)/2, 63/2 - 18*sqrt(3), 18 - 9*sqrt(3), 3/2 - 3*sqrt(3)/4, 9/2 - 3*sqrt(3)/2)

円の半径は,前述のデカルトの円定理を用いて得られたものと同じである。

function draw(r1, r2, more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   (a, r3, x3, r4, x4, r5, x5) = (3/2 + 3*sqrt(3), 9/2 - 9*sqrt(3)/4, -9/2 + 9*sqrt(3)/2, 63/2 - 18*sqrt(3), 18 - 9*sqrt(3), 3/2 - 3*sqrt(3)/4, 9/2 - 3*sqrt(3)/2)
   @printf("r1 = %g;  r2 = %g;  a = %g;  r3 = %g;  x3 = %g;  r4 = %g;  x4 = %g;  r5 = %g;  x5 = %g\n", r1, r2, a, r3, x3, r4, x4, r5, x5)
   plot([0, a, a, 0, 0], [0, 0, r1, r1, 0], color=:green, lw=0.5)
   circle(0, r1, r1, beginangle=270, endangle=360)
   circle(a - r2, r2, r2, :blue)
   circle(x3, r3, r3, :orange)
   circle(x4, r4, r4, :purple)
   circle(x5, r5, r5, :magenta)
   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, r1, "四分円:r1,(0,r1)", :red, :left, :bottom, delta=delta/2)
       point(a - r2, r2, "中円:r2,(a-r2,r2)", :blue, :center, delta=-delta)
       point(x3, r3, "甲円:r3\n(x3,r3)", :orange, :center, :bottom, delta=delta)
       point(x4, r4, " 乙円:r4,(x4,r4)", :purple, :left, :bottom, delta=-delta/2)
       point(x5, r5, " 丙円:r5,(x5,r5)", :magenta, :left, :vcenter, delta=-delta)
       point(a, r1, "(a,r1)", :green, :right, :bottom, delta=delta)
   end
end;

draw(9/2, 3/2, true)

コメント
  • X
  • Facebookでシェアする
  • はてなブックマークに追加する
  • LINEでシェアする

算額(その2020)

2024年08月16日 | Julia

算額(その2020)

(2) 大阪府豊中市服部元町 服部天神社 天保14年(1843)
近畿数学史学会:近畿の算額「数学の絵馬を訪ねて」,平成4年5月16日 初版第一刷,大阪教育図書株式会社,大阪市.
キーワード:円8個,外円

外円の中に 7 個の円を容れる。外円と甲円の直径が 3 寸,1.8 寸のとき,乙円,丙円,丁円の直径を求めよ。

算額(その868)」は累円であるが,これはそのうちの最初に 3 個の円についてである。
パラメータの求め方はそちらを参照(なお,SymPy のバージョンの違いで,当時「算額(その868)」で得られた解が今では違ったものになっている)。
ここでは,累円それぞれの数値解を求めるやりかたでパラメータを求める。

include("julia-source.txt");

using SymPy

@syms R::positive, r1::positive,
     r2::positive, x2::positive, y2::real,
     r3::positive, x3::positive, y3::real,
     r4::positive, x4::positive, y4::positive
R = 3//2
r1 = 18//20
x2 = r2
eq1 = r2^2 + y2^2 - (R - r2)^2
eq4 = r2^2 + (R - r1 - y2)^2 - (r1 + r2)^2;

res1 = solve([eq1, eq4], (r2, y2))[1]

   (9/16, -3/4)

@syms R::positive, r1::positive,
     r2::positive, x2::positive, y2::real,
     r3::positive, x3::positive, y3::real
R = 3//2
r1 = 18//20
(r2, y2) = (9//16, 3//4)
x2 = r2
eq2 = x3^2 + y3^2 - (R - r3)^2
eq5 = x3^2 + (R - r1 - y3)^2 - (r1 + r3)^2
eq7 = (x3 - x2)^2 + (y3 - y2)^2 - (r3 + r2)^2;

res2 = solve([eq2, eq5, eq7], (r3, x3, y3))[2]

   (3/8, 9/8, 0)

@syms R::positive, r1::positive,
     r2::positive, x2::positive, y2::real,
     r3::positive, x3::positive, y3::real
R = 3//2
r1 = 18//20
(r2, y2) = (9//16, 3//4)
x2 = r2
(r3, x3, y3) = (3/8, 9/8, 0)
eq3 = x4^2 + y4^2 - (R - r4)^2
eq6 = x4^2 + (R - r1 - y4)^2 - (r1 + r4)^2
eq8 = (x4 - x3)^2 + (y4 - y3)^2 - (r4 + r3)^2;

res2 = solve([eq3, eq6, eq8], (r4, x4, y4))[1]

   (0.225000000000000, 1.12500000000000, 0.600000000000000)

外円と甲円の直径が 3 寸,1.8 寸のとき,乙円,丙円,丁円の直径は 9/8 寸,3/4 寸,9/20 寸である。

function draw(more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   (R, r1) = (3/2, 1.8/2)
   (r2, x2, y2) = (0.5624999999999999, 0.5624999999999999, -0.75)
   (r3, x3, y3) = (0.3749999999999999, 1.125, -3.3306690738754696e-16)
   (r4, x4, y4) = (0.225000000000000, 1.12500000000000, 0.600000000000000)
   plot()
   circlef(0, 0, R, :cornsilk2)
   circle2f(0, R - r1, r1, 1)
   circle2f(x2, y2, r2, 2)
   circle2f(x3, y3, r3, 3)
   circle2f(x4, y4, r4, 4)
   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, R, "R", :black, :center, :bottom, delta=delta)
       point(0, R - r1, "甲円:r1,(0,R-r1)", :white, :center, delta=-delta)
       point(x2, y2, "乙円:r2,(x2,y2)", :white, :center, delta=-delta)
       point(x3, y3, "丙円:r3\n(x3,y3)", :white, :center, :vcenter)
       point(x4, y4, "丁円:r4\n(x4,y4)", :white, :center, :vcenter)
   end
end;

draw(true)

コメント
  • X
  • Facebookでシェアする
  • はてなブックマークに追加する
  • LINEでシェアする

PVアクセスランキング にほんブログ村

PVアクセスランキング にほんブログ村