裏 RjpWiki

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

算額(その550)

2023年12月10日 | Julia

算額(その550)

一〇一 大宮市高鼻町 氷川神社 明治31年(1898)
埼玉県立図書館:埼玉県史料集 第二集『埼玉の算額』,昭和44年,誠美堂印刷所,埼玉県与野市.

群馬の算額 143−10 榛名神社 明治33年

http://takasakiwasan.web.fc2.com/gunnsann/g143−10.html

 

一〇八 加須市騎西町 玉敷神社 大正4年(1915)
埼玉県立図書館:埼玉県史料集 第二集『埼玉の算額』,昭和44年,誠美堂印刷所,埼玉県与野市.

図のように外円をずらして描き,その隙間に甲円 5 個,乙円 2 個を入れる。乙円の直径が 1 寸のとき,甲円の直径はいかほどか。

外円の半径は甲円の半径の 3 倍,下の外円の中心は一番下の甲円の中心と同じになる。
外円の半径と中心座標を 3r1, (0, 0), (0, -2r1)
乙円の半径と中心座標を r2, (x2, y2)
とおき,以下の連立方程式を解く。なお,右端の甲円の中心座標は何の情報も与えない。

include("julia-source.txt");

using SymPy

@syms r1::positive, r2::positive, x2::positive, y2::positive
eq1 = x2^2 + (2r1 + y2)^2 - (3r1 + r2)^2
eq2 = x2^2 + y2^2 - (3r1 - r2)^2
eq3 = x2^2 + (2r1 - y2)^2 - (r1 + r2)^2;

この連立方程式から r2, x2, y2 を求める(式に r1 が変数として含まれる)と以下の解が得られる。

res = solve([eq1, eq2, eq3], (r2, x2, y2))

   1-element Vector{Tuple{Sym, Sym, Sym}}:
    (4*r1/5, 6*sqrt(2)*r1/5, 7*r1/5)

この連立方程式から r1, x2, y2 を求める(式に r2 が変数として含まれる)と解が得られない。なぜだろうか。

solve([eq1, eq2, eq3], (r1, x2, y2))

   Any[]

ともあれ,前の解から r2 = r1*4/5 が得られるので,r1 = r2*5/4 である。
乙円の直径が 1 寸のとき,甲円の直径は 1.25 寸である。

その他のパラメータは以下の通り。

   r1 = 0.625;  r2 = 0.5;  x2 = 1.06066;  y2 = 0.875
   甲円の直径 = 1.25;  乙円の直径 = 1

以下のプログラムは r2 を与えて r1, x2, y2 を求めて作図するものである。

function draw(r2, more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   r1 = 5r2/4
   r0 = 3r1
   (x2, y2) = r1 .* (6*sqrt(2)/5, 7/5)
   @printf("r1 = %g;  r2 = %g;  x2 = %g;  y2 = %g\n", r1, r2, x2, y2)
   @printf("甲円の直径 = %g;  乙円の直径 = %g\n", 2r1, 2r2)
   θ = atand(r1/sqrt(9r1^2 - r1^2))
   plot()
   circle(0, 0, r0)
   circle(0, -2r1, r0, beginangle=θ, endangle=180-θ)
   circle(0, 2r1, r1, :blue)
   circle(0, 0, r1, :blue)
   circle(0, -2r1, r1, :blue)
   circle(√3r1, -r1, r1, :blue)
   circle(-√3r1, -r1, r1, :blue)
   circle(x2, y2, r2, :magenta)
   circle(-x2, y2, r2, :magenta)
   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(x2, y2, "(x2,y2)", :magenta, :center, :top, delta=-delta/2)
       point(0, 2r1, " 2r1", :blue, :left, :vcenter)
       point(0, -2r1, " -2r1", :red, :left, :vcenter)
       point(√3r1, -r1, "(√3r1,-r1)", :blue, :center, :top, delta=-delta/2)
   end
end;

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

算額(その549)

2023年12月10日 | Julia

算額(その549)

群馬の算額 77-3 宮子神社 嘉永6年
http://takasakiwasan.web.fc2.com/gunnsann/g077-3.html

長方形内に半円,甲円,乙円および等円 2 個を入れる。長方形の短辺の長さが 3 寸のとき,等円の直径はいかほどか。

長方形の短辺,長辺を a, b とする。
半円の半径と中心座標を r0, (a/2, 0); r0 = a/2
甲円の半径と中心座標を r1, (r1, y1)
乙円の半径と中心座標を r2, (a - r2, y1)
等円の半径と中心座標を r3, (a - r3, y3), (a - r3, b - r3)
とおき,以下の連立方程式を解く。

include("julia-source.txt");

using SymPy

@syms a::positive, b, r0, r1, y1::positive, r2, r3, y3::positive

r0 = a//2
eq1 = (r0 - r1)^2 + y1^2 - (r0 + r1)^2
eq2 = (a - r3 - r0)^2 + y3^2 - (r0 + r3)^2
eq3 = 2r1 + 2r2 - a
eq4 = (a - r3 - r1)^2 + (y3 - y1)^2 - (r3 + r1)^2
eq5 = (a - r3 - r1)^2 + (b - r3 - y1)^2 - (r3 + r1)^2
eq6 = (r3 - r2)^2 + (y3 - y1)^2 - (r3 + r2)^2

res = solve([eq1, eq2, eq3, eq4, eq5, eq6], (b, r1, y1, r2, r3, y3))

   2-element Vector{NTuple{6, Sym}}:
    (-sqrt(3)*a/6 + a*(1 + 3*sqrt(3))/6, 3*a/8, sqrt(3)*a/2, a/8, a/6, sqrt(3)*a/3)
    (sqrt(3)*a/6 + a*(1 + 3*sqrt(3))/6, 3*a/8, sqrt(3)*a/2, a/8, a/6, sqrt(3)*a/3)

2 組の解が得られるが,2 番目の解が適解である。

b  = a*(1 + 4*sqrt(3))/6
r1 = a*3/8
y1 = a*sqrt(3)/2
r2 = a/8
r3 = a/6
y3 = a*sqrt(3)/3

等円の半径は長方形の短辺の 1/6 である。
長方形の短辺が 3 寸のとき,等円の半径は 1/2 寸,直径は 1 寸である。

その他のパラメータは以下の通り。

   b = 3.9641;  r1 = 1.125;  y1 = 2.59808;  r2 = 0.375;  r3 = 0.5;  y3 = 1.73205
   等円の直径 = 1;  甲円の直径 = 2.25;  乙円の直径 = 0.75

function draw(more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   a = 3
   r0 = a//2
   (b, r1, y1, r2, r3, y3) = a .* ((1 + 4*sqrt(3))/6, 3/8, sqrt(3)/2, 1/8, 1/6, sqrt(3)/3)
   @printf("b = %g;  r1 = %g;  y1 = %g;  r2 = %g;  r3 = %g;  y3 = %g\n", b, r1, y1, r2, r3, y3)
   @printf("等円の直径 = %g;  甲円の直径 = %g;  乙円の直径 = %g\n", 2r3, 2r1, 2r2)
   plot([0, a, a, 0, 0], [0, 0, b, b, 0], color=:black, lw=0.5)
   circle(r0, 0, r0, beginangle=0, endangle=180)
   circle(r1, y1, r1, :blue)
   circle(a - r3, y3, r3, :green)
   circle(a - r3, b - r3, r3, :green)
   circle(a - r2, y1, r2, :magenta)
   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(r0, 0, "半円:r0,(a/2,0)", :red, :center, :bottom, delta=delta)
       point(a, 0, "a ", :red, :right, :bottom, delta=delta/2)
       point(0, b, " b", :black, :left, :bottom, delta=delta/2)
       point(r1, y1, "甲円:r1,(r1,y1)", :blue, :center, :top, delta=-delta/2)
       point(a - r2, y1, "乙円:r2\n(r2,y1)", :magenta, :center, :top, delta=-delta/2)
       point(a - r3, b - r3, "等円:r3\n(a-r3,b-r3)", :green, :center, :top, delta=-delta/2)
       point(a - r3, y3, "等円:r3\n(a-r3,y3)", :green, :center, :top, delta=-delta/2)
   end
end;

 

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

算額(その548)

2023年12月10日 | Julia

算額(その548)

三十六 群馬県多野郡新町 稲荷神社 文政3年(1820)
群馬県和算研究会:群馬の算額,上武印刷株式会社,高崎市,1987年3月31日.

群馬の算額 36-3 稲荷神社 文政3年(1280)
http://takasakiwasan.web.fc2.com/gunnsann/g036-3.html

図のような等脚台形で,3 つの甲斜が 1 寸のとき,等脚台形の面積を最大にする乙斜の長さはいかほどか。

甲斜,乙斜を a, b とする。
台形の高さ h は,sqrt(a^2 - ((b - a)/2)^2) である。
面積は (a + b)*sqrt(4*a^2 - (a - b)^2)/4 である。

include("julia-source.txt");

using SymPy

@syms a, b
h = sqrt(a^2 - ((b - a)/2)^2)
S = (a + b)*h/2 |> simplify
S |> println

   (a + b)*sqrt(4*a^2 - (a - b)^2)/4

a = 1 のときに S のグラフを描くと以下のようになる。

pyplot(size=(400, 300), grid=true, aspectratio=:none, label="")
plot((1 + b)*sqrt(4*1^2 - (1 - b)^2)/4, xlims=(0,3), xlabel="b", ylabel="S")

b が 2 前後のときに S が最大値になることがわかる。

S を b について微分すると,接線の傾きを表す式 (2*a^2 + a*b - b^2)/(2*sqrt(3*a^2 + 2*a*b - b^2)) になる。

s = diff(S, b) |> simplify
s |> println

   (2*a^2 + a*b - b^2)/(2*sqrt(3*a^2 + 2*a*b - b^2))

接線の傾きが 0 ということは元の関数が極値(最大値または最小値)になるということである。
b = 2 のときに接線の傾きが 0 になる(b = 2 のときに最大値をとる)。

solve(s, b)[1] |> println

   2*a

S に a = 1, b = 2 を代入すると,3*sqrt(3)/4 ≒ 1.29903810567666 が最大値である。

S(a => 1, b => 2) |> println
S(a => 1, b => 2).evalf() |> println

   3*sqrt(3)/4
   1.29903810567666

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

算額(その547)

2023年12月10日 | Julia

算額(その547)

七十二 群馬県富岡市一ノ宮 貫前神社 嘉永2年(1849)
群馬県和算研究会:群馬の算額,上武印刷株式会社,高崎市,1987年3月31日.

群馬の算額 72−15 貫前神社 嘉永2年
http://takasakiwasan.web.fc2.com/gunnsann/g072-4.html

大円の中に弦 1 本,中円 2 個,小円 4 個が入っている。小円の直径が 7 寸のとき,大円の直径はいかほどか。

大円の半径と中心座標を R, (0, 0)
中円の半径と中心座標を 2r2, (r2, y1)
小円の半径と中心座標を r2, (2r2, y + r2)
とおき,以下の連立方程式を解く。

include("julia-source.txt");

using SymPy

@syms R::positive, y::negative, y1::positive, r2::positive

eq1 = r2^2 + y1^2 - (R - 2r2)^2
eq2 = (2r2 - r2)^2 + (y1 - y - r2)^2 - (2r2 + r2)^2
eq3 = (2r2)^2 + (y + r2)^2 - (R - r2)^2
res = solve([eq1, eq2, eq3], (R, y, y1))

   1-element Vector{Tuple{Sym, Sym, Sym}}:
    (25*r2/7, -r2*(7 + 8*sqrt(2))/7, 6*sqrt(2)*r2/7)

大円の半径は小円の半径の 25/7 倍である。
小円の直径が 7 寸のとき,大円の直径は 25 寸である。

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

r2 = 7/2
(25*r2/7, -r2*(7 + 8*sqrt(2))/7, 6*sqrt(2)*r2/7)

   (12.5, -9.15685424949238, 4.242640687119286)

function draw(more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   r2 = 7//2
   (R, y, y1) = (25*r2/7, -r2*(7 + 8*sqrt(2))/7, 6*sqrt(2)*r2/7)
   @printf("R = %g;  y = %g;  y1 = %g\n", R, y, y1)
   @printf("大円の直径 = %g\n", 2R)
   plot()
   circle(0, 0, R)
   circle(r2, y1, 2r2, :blue)
   circle(-r2, y1, 2r2, :blue)
   circle(0, y + r2, r2, :green)
   circle(2r2, y + r2, r2, :green)
   circle(-2r2, y + r2, r2, :green)
   circle(0, y1, r2, :green)
   segment(-sqrt(R^2 - y^2), y, sqrt(R^2 - y^2), y, :magenta)
   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, y1, " y1", :green, :left, :vcenter)
       point(0, y, " y", :magenta, :left, :top, delta=-delta/2)
       point(R, 0, "R ", :red, :right, :bottom, delta=delta/2)
       point(r2, y1, " 中円:2r2,(r2,y1)", :blue, :left, :vcenter)
       point(2r2, y + r2, "小円:r2\n(2r2,y+r2)", :green, :center, :top, delta=-delta)
   end
end;

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

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

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