裏 RjpWiki

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

算額(その700)

2024年02月16日 | ブログラミング

算額(その700)

埼玉の算額ほか
https://gunmawasan.web.fc2.com/files/saitama-sangaku-h24.html
愛宕神社の復元算額 明治13年(部分拡大図)(加須市)
https://gunmawasan.web.fc2.com/files/sangak-corner/atago-3s.jpg

正方形の中に四分円 2 個,大円,小円それぞれ 2 個ずつ入っている。小円の直径が 1 寸のとき,大円の直径はいかほどか。

正方形の一辺の長さを 2a
大円の半径と中心座標を r1, (r1, r1)
小円の半径と中心座標を r2, (a - r2, y2)
とおき,以下の連立方程式を解く

include("julia-source.txt");

using SymPy
@syms a::positive, r1::positive, r2::positive, y2::positive

eq1 = r2^2 + y2^2 - (2a - r2)^2
eq2 = (r1 + a)^2 + r1^2 - (2a - r1)^2
eq3 = (2a - r2)^2 + y2^2 - (2a + r2)^2
res = solve([eq1, eq2, eq3], (a, r1, y2))

   1-element Vector{Tuple{Sym, Sym, Sym}}:
    (3*r2, 3*r2*(-3 + 2*sqrt(3)), 2*sqrt(6)*r2)

大円の半径は,小円の半径の (6√3 - 9) 倍である。
小円の直径が 1 寸の場合,大円の直径は 6√3 - 9 = 1.392304845413264 寸である。

その他のパラメータは以下のとおりである。

a = 1.5;  r1 = 0.696152;  y2 = 2.44949

function draw(more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   r2 = 1//2
   (a, r1, y2) = (3*r2, 3*r2*(-3 + 2*sqrt(3)), 2*sqrt(6)*r2)
   @printf("大円の直径 = %g;  a = %g;  r1 = %g;  y2 = %g\n", 2r1, a, r1, y2)
   plot([a, a, -a, -a, a], [0, 2a, 2a, 0, 0], color=:magenta, lw=0.5)
   circle(a, 0, 2a, :blue, beginangle=90, endangle=180)
   circle(-a, 0, 2a, :blue, beginangle=0, endangle=90)
   circle(r1, r1, r1)
   circle(-r1, r1, r1)
   circle(a - r2, y2, r2, :green)
   circle(r2 - a, y2, r2, :green)
   if more
       delta = (fontheight = (ylims()[2]- ylims()[1]) / 500 * 10 * 2) /3  # size[2] * fontsize * 2
       hline!([0], color=:black, lw=0.5)
       vline!([0], color=:black, lw=0.5)
       point(a, 0, " a", :magenta, :left, :bottom, delta=delta/2)
       point(0, 2a, " 2a", :magenta, :left, :bottom, delta=delta/2)
       point(r1, r1, "大円:r1,(r1, r1)", :red, :center, delta=-delta/2)
       point(a - r2, y2, "小円:r2,(a-r2,y2)", :green, :center, delta=-delta/2)
   end
end;

 

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

算額(その699)

2024年02月16日 | Julia

算額(その699)

山口正義: やまぶき4,第58号,2018/12/06.
千葉県君津市鹿野山 鹿野山神野寺 万延二年(1861)

https://yamabukiwasan.sakura.ne.jp/ymbk58.pdf
山口正義:やまぶき,和算と歴史随想 からリンク
https://yamabukiwasan.sakura.ne.jp/page3.html

高さ 64 寸の正三角形内に,短径 6 寸の甲楕円,長径 48 寸の乙楕円が入っている。
甲楕円の長径,乙楕円の短径を求めよ。

正三角形の一辺の長さを 2a とする。
甲楕円の長半径,短半径,中心座標を a1, b1, (0, 2b2 + b1)
乙楕円の長半径,短半径,中心座標を a2, b2, (0, b2)
正三角形の斜辺と甲楕円,乙楕円との接点座標を (x1, y1), (x2, y2)
とおき,以下の連立方程式を解く。

include("julia-source.txt");

using SymPy
@syms a1::positive, b1::positive, x1::positive, y1::positive,
     a2::positive, b2::positive, x2::positive, y2::positive
@syms a::positive

s3 = sqrt(Sym(3))
(a1, a2) = (6, 48) .// 2
a = 64/s3
y1 = (a - x1)*s3
y2 = (a - x2)*s3
eq1 = x1^2/a1^2 + (y1 - (2b2 + b1))^2/b1^2 - 1
eq3 = x2^2/a2^2 + (y2 - b2)^2/b2^2 - 1
eq2 = -b1^2*x1/(a1^2*(y1 - (2b2 + b1))) + s3
eq4 = -b2^2*x2/(a2^2*(y2 - b2)) + s3
res = solve([eq1, eq2, eq3, eq4], (b1, x1, b2, x2))

   1-element Vector{NTuple{4, Sym}}:
    (13, 9*sqrt(3)/14, 37/2, 1152*sqrt(3)/91)

甲楕円の長径は 26 寸,乙楕円の短径は 37 寸である。

include("julia-source.txt");
function draw(more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   s3 = sqrt(3)
   h = 64
   a = h/s3
   (a1, a2) = (6, 48) .// 2
   (b1, x1, b2, x2) = (13, 9*sqrt(3)/14, 37/2, 1152*sqrt(3)/91)
   @printf("甲楕円の長径 = %g;  乙楕円の短径 = %g\n", 2b1, 2b2)
   y1 = (a - x1)*s3
   y2 = (a - x2)*s3
   plot([a, 0, -a, a], [0, h, 0, 0], color=:magenta, lw=0.5)
   ellipse(0, b2, a2, b2, color=:blue)
   ellipse(0, b1 + 2b2, a1, b1, color=:red)
   if more
       delta = (fontheight = (ylims()[2]- ylims()[1]) / 500 * 10 * 2) /3  # size[2] * fontsize * 2
       hline!([0], color=:black, lw=0.5)
       vline!([0], color=:black, lw=0.5)
       point(0, b1 + 2b2, " 甲楕円:a1,b1,(0,2b2+b1)", :black, :left, :vcenter)
       point(0, b2, " 乙楕円:a2,b2,(0,b2)", :black, :left, :vcenter)
       point(x1, y1, " (x1,y1)", :green, :left, :vcenter)
       point(x2, y2, " (x2,y2)", :green, :left, :vcenter)
   end
end;

 

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

算額(その698)

2024年02月16日 | Julia

算額(その698)

新潟県長岡市上岩井 根立寺 嘉永2年(1849)
http://www.wasan.jp/niigata/konryuji.html
涌田和芳,外川一仁: 三島根立寺の算額,長岡工業高等専門学校研究紀要,第53巻(2017)
https://kinpoku.nagaoka-ct.ac.jp/lib/kiyo/vol_53/53_17wakuta.pdf

直線上に楕円と天円,地円,人円がある。地円の直径が 1 寸のとき,極大となる天円の直径はいかほどか。

楕円の長半径,短半径,中心座標を a, b, (0, b)  注: この問題にいおいては,a < b である
天円の半径と中心座標を r1, (0, y1)
人円の半径と中心座標を r3, (0, 2r2 + r3)
地円の半径と中心座標を r2, (0, r2)
とおく。
算法助術の公式 85, 86, 86 を使用し,以下の連立方程式を解く。
天円は楕円の曲率円にあたるので,その半径は a^2/b である。
なお,楕円の左右の人円はこの問題を解く上では無関係である。

include("julia-source.txt");

using SymPy
@syms a::positive, b::positive, r1::positive, r2::positive, r3::positive
r2 = 1//2
eq1 = r1 - a^2/b
eq2 = r3 - (3r1 - 4r1^2/b)
eq3 = r2^2 - (a^2*(a^2 - r3^2))/(b^2 - a^2)
eq4 = r1 + r2 + r3 - b
res = solve([eq1, eq2, eq3, eq4], (r1, r3, a, b))

   1-element Vector{NTuple{4, Sym}}:
    (3/2, 5/2, 3*sqrt(3)/2, 9/2)

天円の半径は 2/3 となり,地円の 3 倍である。すなわち,天円の直径は 3 寸である。
その他のパラメータは以下のとおりである。

r1 = 1.5;  r3 = 2.5;  a = 2.59808;  b = 4.5
y1 = 7.5;  y3 = 3.5

include("julia-source.txt");
function draw(more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   r2 = 1//2
   (r1, r3, a, b) = (3/2, 5/2, 3*sqrt(3)/2, 9/2)
   y3 = 2r2 + r3
   y1 = 2r2 + 2r3 + r1
   @printf("r1 = %g;  r3 = %g;  a = %g;  b = %g\n", r1, r3, a, b)
   @printf("y1 = %g;  y3 = %g\n", y1, y3)
   plot()
   ellipse(0, b, a, b, color=:magenta)
   circle(0, y1, r1)
   circle(0, r2, r2, :blue)
   circle(0, y3, r3, :green)
   d = 2sqrt(r3^2 - r2^2)
   circle(d, r3, r3, :green)
   circle(-d, r3, r3, :green)
   if more
       delta = (fontheight = (ylims()[2]- ylims()[1]) / 500 * 10 * 2) /3  # size[2] * fontsize * 2
       hline!([0], color=:black, lw=0.5)
       vline!([0], color=:black, lw=0.5)
       point(0, b, " b", :black, :left, :vcenter)
       point(0, y3, " 人円:r3\n (0,y3)", :green, :left, :vcenter)
       point(d, r3, "(d,r3)", :green, :center, delta=-delta/2)
       point(-d, r3, "(-d,r3)", :green, :center, delta=-delta/2)
       point(0, y1, " 天円:r1\n (0,y1)", :black, :left, :vcenter)
       point(0, r2, "   地円:r2,(0,r2)", :black, :left, :vcenter)
   end
end;

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

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

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