裏 RjpWiki

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

算額(その1367)

2024年10月22日 | Julia

算額(その1367)

四十一 群馬県高崎市下小鳥町 幸宮神社 文政7年(1824)
群馬県和算研究会:群馬の算額,上武印刷株式会社,高崎市,1987年3月31日.
キーワード:直角三角形,鈎股弦
#Julia, #SymPy, #算額, #和算

直角三角形の面積が 756 歩,股寸と弦寸の和が 1 丈 4 尺 7 寸のとき,釣,股,弦はいかほどか。

直角三角形の三辺を「鈎」,「股」,「弦」
直角三角形の面積を K1
股と弦の和を K2
とおき,以下の連立方程式を解く。

using SymPy
@syms 鈎, 股, 弦, K1, K2
eq1 = 鈎^2 + 股^2 - 弦^2  # こちらで与えると解けない
eq1 = sqrt(鈎^2 + 股^2) - 弦
eq2 = 鈎*股/2 - K1
eq3 = 股 + 弦 - K2;
res = solve([eq1, eq2, eq3], (鈎, 股, 弦))[3]  # 3 of 3

   (6^(2/3)*(12*2^(1/3)*K2^2 + 18^(1/3)*(1 + sqrt(3)*I)^2*(18*K1*K2 + sqrt(3)*sqrt(K2^2*(108*K1^2 - K2^4)))^(2/3))/(36*(1 + sqrt(3)*I)*(18*K1*K2 + sqrt(3)*sqrt(K2^2*(108*K1^2 - K2^4)))^(1/3)), 12*K1*(1 + sqrt(3)*I)*(108*K1*K2 + 6*sqrt(3)*sqrt(K2^2*(108*K1^2 - K2^4)))^(1/3)/(12*2^(1/3)*K2^2 + 18^(1/3)*(1 + sqrt(3)*I)^2*(18*K1*K2 + sqrt(3)*sqrt(K2^2*(108*K1^2 - K2^4)))^(2/3)), (12*K1*(1 + sqrt(3)*I)*(108*K1*K2 + 6*sqrt(3)*sqrt(K2^2*(108*K1^2 - K2^4)))^(1/3) + K2*(-12*2^(1/3)*K2^2 - 18^(1/3)*(1 + sqrt(3)*I)^2*(18*K1*K2 + sqrt(3)*sqrt(K2^2*(108*K1^2 - K2^4)))^(2/3)))/(-12*2^(1/3)*K2^2 - 18^(1/3)*(1 + sqrt(3)*I)^2*(18*K1*K2 + sqrt(3)*sqrt(K2^2*(108*K1^2 - K2^4)))^(2/3)))

解は 3 組得られるが,3 番目のものが適解である。虚数解として得られるものがあるが,虚部は実質的に 0 である。

鈎,股,弦はそれぞれ,21 寸,72 寸,75 寸である。

res[1](K1 => 756, K2 => 147).evalf() |> println
res[2](K1 => 756, K2 => 147).evalf() |> println
res[3](K1 => 756, K2 => 147).evalf() |> println

   21.0000000000000
   72.0 - 8.67361737988404e-19*I
   75.0000000000000

 

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

算額(その1366)

2024年10月22日 | Julia

算額(その1366)

三十六 群馬県多野郡新町 稲荷神社 文政3年(1820)
群馬県和算研究会:群馬の算額,上武印刷株式会社,高崎市,1987年3月31日.
キーワード:円8個,等脚台形
#Julia, #SymPy, #算額, #和算

等脚台形の中に,大円 5 個,小円 3 個を容れる。台形の上底が与えられたとき,小円の直径はいかほどか。

注:算額(その1364)の設問の不備を修正し,枝葉を取り除いたものである。

等脚台形の下底,上底,高さを 2a, 2b, h
大円の半径と中心座標を r1, (2r1, r1), (r1, h - r1)
小円の半径と中心座標を r2, (0, 2r1 + r2)
甲円の半径と中心座標を r3, (x3, r3)
乙円の半径と中心座標を r4, (0, h - r4)
とおき,以下の連立方程式を解く。
方程式は算額(その1364)と同じで,与える条件と求めるパラメータが違うだけである。
甲円と乙円のパラメータも求めておく(図には描かない)。

include("julia-source.txt");

using SymPy
@syms a::positive, b::positive, h::positive,
     r1::positive, r2::positive,
     r3::positive, x3::positive, r4::positive;
eq1 = r1^2 + (h - 2r1)^2 - 4r1^2
eq2 = r1/(a - 2r1) - 1/√Sym(3)
eq3 = r1/(b - r1) - √Sym(3);
# res1 = solve([eq1, eq2, eq3], (r1, a, h));

eq4 = r1^2 + (h - r1 - (2r1 + r2))^2 - (r1 + r2)^2
eq5 = r1^2 + (r1 - r4)^2 - (r1 + r4)^2;
# res = solve([eq1, eq2, eq3, eq4, eq5], (r1, r2, r4, a, h))

eq6 = r3/(a - x3) - 1/√Sym(3)
eq7 = x3 - 2r1 - 2sqrt(r1*r3)
res = solve([eq1, eq2, eq3, eq4, eq5, eq6, eq7],
   (r1, r2, r3, x3, r4, a, h))[1]  # 2 of 2

   (b*(3 - sqrt(3))/2, -989*b*sqrt(151316*sqrt(3) + 262087) - 3*b/2 + sqrt(3)*b/2 + 571*sqrt(3)*b*sqrt(151316*sqrt(3) + 262087), b*(3 - sqrt(3))/6, 2*b, b*(3 - sqrt(3))/8, b*(sqrt(3) + 3)/2, b*(-2967*sqrt(151316*sqrt(3) + 262087) - 2*sqrt(3) + 6 + 1713*sqrt(453948*sqrt(3) + 786261))/2)

パラメータによっては二重根号を外して簡約化できる。

r1
res[1] |> println

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

r2
res[2] |> sympy.sqrtdenest |> simplify |> println

   b*(-5 + 3*sqrt(3))/2

r3
res[3] |> sympy.sqrtdenest |> simplify |> println

   b*(3 - sqrt(3))/6

x3
res[4] |> sympy.sqrtdenest |> simplify |> println

   2*b

r4
res[5] |> println

   b*(3 - sqrt(3))/8

a
res[6] |> println

   b*(sqrt(3) + 3)/2

h
res[7] |> sympy.sqrtdenest |> simplify |> println

   b*(sqrt(3) + 3)/2

小円の半径は b*(3 - sqrt(3))/6 なので,上底が 14 寸のとき小円の直径は 14*(3 - sqrt(3))/6 = 2.95854811567262 寸である。
このとき,下底は 32*(3 - sqrt(3))/3 = 13.524791385931977 である。

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

   a = 16;  r1 = 4.28719;  r2 = 0.66323;  r3 = 1.42906;  x3 = 13.5248;  r4 = 1.0718;  b = 6.7624;  h = 16

---

なお,小円の直径を知るだけならば,以下のような手順で小円の直径を(紙と鉛筆だけで)求めることができる。

上底 2b と大円の半径 r1 の関係式は以下であり,大円の半径は r1 = b*(3 - √3)/2 である。

@syms b, r1
eq01 = (1 + 1/√Sym(3))r1 - b
eq01 |> println

   -b + r1*(sqrt(3)/3 + 1)

r1 = solve(eq01, r1)[1] |> simplify
r1 |> println

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

一辺の長さが 2r1 の正三角形の頂点を中心とする3個の大円が作る中央の隙間にピッタリハマる円の半径は以下のようになる(図中の灰色で描いた小さな直角三角形を参照)。

r1*2/√Sym(3) - r1 |> simplify |> println

    b*(-5 + 3*sqrt(3))/2

上底が 14 寸のとき,小円の直径は 1.3730669589464242 寸である。

b = 14/2
2*b*(-5 + 3*sqrt(3))/2

   1.3730669589464242

function draw(b, more)
   pyplot(size=(800, 800), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   (r1, r2, r3, x3, r4, a, h) = (b*(3 - sqrt(3))/2, -989*b*sqrt(151316*sqrt(3) + 262087) - 3*b/2 + sqrt(3)*b/2 + 571*sqrt(3)*b*sqrt(151316*sqrt(3) + 262087), b*(3 - sqrt(3))/6, 2*b, b*(3 - sqrt(3))/8, b*(sqrt(3) + 3)/2, b*(-2967*sqrt(151316*sqrt(3) + 262087) - 2*sqrt(3) + 6 + 1713*sqrt(453948*sqrt(3) + 786261))/2)
   r1 = b*(3 - √3)/2
   r2 = b*(3√3 - 5)/2
   r3 = b*(3 - √3)/6
   x3 = 2b
   r4 = b*(3 - √3)/8
   a = b*(√3 + 3)/2
   h = b*(√3 + 3)/2
   @printf("上底が %g のとき,小円の直径は %g である(下底は %g)。\n", 2b, 2r2, 2a)
   @printf("b = %g;  r1 = %g;  r2 = %g;  r3 = %g;  x3 = %g;  r4 = %g;  a = %g;  h = %g\n", b, r1, r2, r3, x3, r4, a, h)
   plot([a, b, -b, -a, a], [0, h, h, 0, 0], color=:green, lw=0.5)
   circle(0, r1, r1)
   circle2(2r1, r1, r1)
   circle2(r1, h - r1, r1)

   circle(0, 2r1 + r2, r2, :blue)
   circle2(r1, h - 2r1 - r2, r2, :blue)

   #=
   circle(0, h - r4, r4, :green)
   circle2(r1, r4, r4, :green)

   circle2(x3, r3, r3, :orange)
   l = h - (2r1 + r2) - r4
   circle2(r1 + l*cosd(30), h - (2r1 + r2) + l*sind(30), r4, :green)
   =#
   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(2r1, r1, "大円:r1,(2r1,r1)", :red, :center, delta=-delta)
       point(0, r1, "大円:r1,(0,r1)", :red, :center, delta=-delta)
       point(r1, (1 + √3)r1, "大円:r1,(r1,(1+√3)r1)", :red, :center, delta=-delta)
       point(0, 2r1 + r2, "小円:r2,(0,2r1+r2)", :blue, :left, :vcenter, deltax=4delta)
       #=
       point(0, h - r4, "乙円:(0,h-r4)", :green, :left, :vcenter, deltax=7delta)
       point(x3, r3, "甲円:(x3,r3)", :orange, :right, :vcenter, deltax=-9delta)
       =#
       point(b, h, "(b,h)", :brown, :right, :bottom, delta=delta/2)
       point(a, 0, "a", :brown, :left, :bottom, delta=delta/2)
       plot!([r1, 2r1, r1, r1], [r1, r1, h - 2r1 - r2, r1], color=:gray70, lw=0.5)
   end  
end;

draw(14/2, true)

 

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

算額(その1365)

2024年10月22日 | Julia

算額(その1365)

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

算額(その23),算額(その995)算額(その396)を組み合わせた問題である。

図形と求めるものは,算額(その806)と同じである。
算額(その806)では与えられた条件を使って解を求めた。
今回の問題では,一般解を求める。

直角三角形の 3 辺を,「鈎」,「股」,「弦」
甲円の半径と中心座標を r1, (r1, r1)
乙円の半径と中心座標を r2, (x2, r2)
丁円の半径と中心座標を r3, (x3, r3)
とおき,以下の連立方程式を解く。
SymPy の性能上,一度に解を求めることができないので,まずは eq1, eq2, eq3, eq4 を解いて 「股」,r1, r2, x2 を求める。
得られた解を既知として,eq5, eq6 を解いて r3, x3 を求める。

include("julia-source.txt");

using SymPy
@syms 鈎::positive, 股::positive, 弦::positive,
     r1::positive, r2::positive, x2::positive,
     r3::positive, x3::positive;
eq1 = 鈎^2 + 股^2 - 弦^2
eq2 = 鈎 + 股 - 弦 - 2r1
eq3 = (x2 - r1) - 2sqrt(r1*r2)
eq4 = r2/(股 - x2) - r1/(股 - r1);
res1 = solve([eq1, eq2, eq3, eq4], (股, r1, r2, x2))[2]  # 2 of 2

   (sqrt(弦^2 - 鈎^2), -弦/2 + 鈎/2 + sqrt(弦^2 - 鈎^2)/2, (-弦*(4*sqrt(弦)*sqrt(-(弦 - 鈎)*(-2*弦^2 + 2*弦*sqrt(弦^2 - 鈎^2) + 鈎^2)) + 4*弦^2 - 3*弦*鈎 - 4*弦*sqrt(弦^2 - 鈎^2) - 鈎^2 + 3*鈎*sqrt(弦^2 - 鈎^2)) + 2*鈎*(-弦^2 + 弦*sqrt(弦^2 - 鈎^2) + 鈎^2) + sqrt(弦^2 - 鈎^2)*(4*sqrt(弦)*sqrt(-(弦 - 鈎)*(-2*弦^2 + 2*弦*sqrt(弦^2 - 鈎^2) + 鈎^2)) + 4*弦^2 - 3*弦*鈎 - 4*弦*sqrt(弦^2 - 鈎^2) - 鈎^2 + 3*鈎*sqrt(弦^2 - 鈎^2)))/(2*鈎^2), (4*sqrt(弦)*sqrt(-(弦 - 鈎)*(-2*弦^2 + 2*弦*sqrt(弦^2 - 鈎^2) + 鈎^2)) + 4*弦^2 - 3*弦*鈎 - 4*弦*sqrt(弦^2 - 鈎^2) - 鈎^2 + 3*鈎*sqrt(弦^2 - 鈎^2))/(2*鈎))

股 = sqrt(弦^2 - 鈎^2)
r1 = -弦/2 + 鈎/2 + sqrt(弦^2 - 鈎^2)/2
r2 = (-弦*(4*sqrt(弦)*sqrt(-(弦 - 鈎)*(-2*弦^2 + 2*弦*sqrt(弦^2 - 鈎^2) + 鈎^2)) + 4*弦^2 - 3*弦*鈎 - 4*弦*sqrt(弦^2 - 鈎^2) - 鈎^2 + 3*鈎*sqrt(弦^2 - 鈎^2)) + 2*鈎*(-弦^2 + 弦*sqrt(弦^2 - 鈎^2) + 鈎^2) + sqrt(弦^2 - 鈎^2)*(4*sqrt(弦)*sqrt(-(弦 - 鈎)*(-2*弦^2 + 2*弦*sqrt(弦^2 - 鈎^2) + 鈎^2)) + 4*弦^2 - 3*弦*鈎 - 4*弦*sqrt(弦^2 - 鈎^2) - 鈎^2 + 3*鈎*sqrt(弦^2 - 鈎^2)))/(2*鈎^2)
x2 = (4*sqrt(弦)*sqrt(-(弦 - 鈎)*(-2*弦^2 + 2*弦*sqrt(弦^2 - 鈎^2) + 鈎^2)) + 4*弦^2 - 3*弦*鈎 - 4*弦*sqrt(弦^2 - 鈎^2) - 鈎^2 + 3*鈎*sqrt(弦^2 - 鈎^2))/(2*鈎)
eq5 = (x3 - r1) - 2sqrt(r1*r3)
eq6 = (x2 - x3) - 2sqrt(r2*r3)
res2 = solve([eq5, eq6], (r3, x3))[1]

   ((2*sqrt(弦)*sqrt(2*弦^3 - 2*弦^2*鈎 - 2*弦^2*sqrt(弦^2 - 鈎^2) - 弦*鈎^2 + 2*弦*鈎*sqrt(弦^2 - 鈎^2) + 鈎^3) + 2*弦^2 - 弦*鈎 - 2*弦*sqrt(弦^2 - 鈎^2) - 鈎^2 + 鈎*sqrt(弦^2 - 鈎^2))^2/(2*(鈎*sqrt(-弦 + 鈎 + sqrt(弦^2 - 鈎^2)) + sqrt(-4*弦^(3/2)*sqrt(2*弦^3 - 2*弦^2*鈎 - 2*弦^2*sqrt(弦^2 - 鈎^2) - 弦*鈎^2 + 2*弦*鈎*sqrt(弦^2 - 鈎^2) + 鈎^3) + 4*sqrt(弦)*sqrt(弦^2 - 鈎^2)*sqrt(2*弦^3 - 2*弦^2*鈎 - 2*弦^2*sqrt(弦^2 - 鈎^2) - 弦*鈎^2 + 2*弦*鈎*sqrt(弦^2 - 鈎^2) + 鈎^3) - 8*弦^3 + 4*弦^2*鈎 + 8*弦^2*sqrt(弦^2 - 鈎^2) + 5*弦*鈎^2 - 4*弦*鈎*sqrt(弦^2 - 鈎^2) - 鈎^3 - 鈎^2*sqrt(弦^2 - 鈎^2)))^2), -弦/2 + 鈎/2 + sqrt((2*sqrt(弦)*sqrt(2*弦^3 - 2*弦^2*鈎 - 2*弦^2*sqrt(弦^2 - 鈎^2) - 弦*鈎^2 + 2*弦*鈎*sqrt(弦^2 - 鈎^2) + 鈎^3) + 2*弦^2 - 弦*鈎 - 2*弦*sqrt(弦^2 - 鈎^2) - 鈎^2 + 鈎*sqrt(弦^2 - 鈎^2))^2/(鈎*sqrt(-弦 + 鈎 + sqrt(弦^2 - 鈎^2)) + sqrt(-4*弦^(3/2)*sqrt(2*弦^3 - 2*弦^2*鈎 - 2*弦^2*sqrt(弦^2 - 鈎^2) - 弦*鈎^2 + 2*弦*鈎*sqrt(弦^2 - 鈎^2) + 鈎^3) + 4*sqrt(弦)*sqrt(弦^2 - 鈎^2)*sqrt(2*弦^3 - 2*弦^2*鈎 - 2*弦^2*sqrt(弦^2 - 鈎^2) - 弦*鈎^2 + 2*弦*鈎*sqrt(弦^2 - 鈎^2) + 鈎^3) - 8*弦^3 + 4*弦^2*鈎 + 8*弦^2*sqrt(弦^2 - 鈎^2) + 5*弦*鈎^2 - 4*弦*鈎*sqrt(弦^2 - 鈎^2) - 鈎^3 - 鈎^2*sqrt(弦^2 - 鈎^2)))^2)*sqrt(-弦 + 鈎 + sqrt(弦^2 - 鈎^2)) + sqrt(弦^2 - 鈎^2)/2)

丁円の半径は,「鈎」と「弦」を含む長い式になる。

-弦/2 + 鈎/2 + sqrt((2*sqrt(弦)*sqrt(2*弦^3 - 2*弦^2*鈎 - 2*弦^2*sqrt(弦^2 - 鈎^2) - 弦*鈎^2 + 2*弦*鈎*sqrt(弦^2 - 鈎^2) + 鈎^3) + 2*弦^2 - 弦*鈎 - 2*弦*sqrt(弦^2 - 鈎^2) - 鈎^2 + 鈎*sqrt(弦^2 - 鈎^2))^2/(鈎*sqrt(-弦 + 鈎 + sqrt(弦^2 - 鈎^2)) + sqrt(-4*弦^(3/2)*sqrt(2*弦^3 - 2*弦^2*鈎 - 2*弦^2*sqrt(弦^2 - 鈎^2) - 弦*鈎^2 + 2*弦*鈎*sqrt(弦^2 - 鈎^2) + 鈎^3) + 4*sqrt(弦)*sqrt(弦^2 - 鈎^2)*sqrt(2*弦^3 - 2*弦^2*鈎 - 2*弦^2*sqrt(弦^2 - 鈎^2) - 弦*鈎^2 + 2*弦*鈎*sqrt(弦^2 - 鈎^2) + 鈎^3) - 8*弦^3 + 4*弦^2*鈎 + 8*弦^2*sqrt(弦^2 - 鈎^2) + 5*弦*鈎^2 - 4*弦*鈎*sqrt(弦^2 - 鈎^2) - 鈎^3 - 鈎^2*sqrt(弦^2 - 鈎^2)))^2)*sqrt(-弦 + 鈎 + sqrt(弦^2 - 鈎^2)) + sqrt(弦^2 - 鈎^2)/2

鈎 = 3 寸, 弦 = 6 寸を代入して,r3 = 0.20699488352959394 寸を得る。

---

鈎,弦として,与えられた数値そのままを用いて連立方程式を解けば,数値解が求まる。

using SymPy
@syms 鈎::positive, 股::positive, 弦::positive,
     r1::positive, r2::positive, x2::positive,
     r3::positive, x3::positive;
(鈎, 弦) = (3, 6)
eq1 = 鈎^2 + 股^2 - 弦^2
eq2 = 鈎 + 股 - 弦 - 2r1
eq3 = (x2 - r1) - 2sqrt(r1*r2)
eq4 = r2/(股 - x2) - r1/(股 - r1);
eq5 = (x3 - r1) - 2sqrt(r1*r3)
eq6 = (x2 - x3) - 2sqrt(r2*r3)
res1 = solve([eq1, eq2, eq3, eq4, eq5, eq6], (股, r1, r2, x2, r3, x3))[1]

   (3*sqrt(3), -3/2 + 3*sqrt(3)/2, -117/2 - 12*sqrt(14 - 8*sqrt(3)) + 6*sqrt(42 - 24*sqrt(3)) + 69*sqrt(3)/2, -15*sqrt(3)/2 + 6*sqrt(14 - 8*sqrt(3)) + 27/2, 3*(-3*sqrt(6) + 4*sqrt(7 - 4*sqrt(3)) + 5*sqrt(2))^2/(4*(sqrt(-39 - 8*sqrt(2)*sqrt(7 - 4*sqrt(3)) + 4*sqrt(6)*sqrt(7 - 4*sqrt(3)) + 23*sqrt(3)) + sqrt(-1 + sqrt(3)))^2), (3*sqrt(-2 + 2*sqrt(3))*(-3*sqrt(6) + 4*sqrt(7 - 4*sqrt(3)) + 5*sqrt(2))/2 + 3*(-1 + sqrt(3))*Abs(sqrt(-39 - 8*sqrt(2)*sqrt(7 - 4*sqrt(3)) + 4*sqrt(6)*sqrt(7 - 4*sqrt(3)) + 23*sqrt(3)) + sqrt(-1 + sqrt(3)))/2)/Abs(sqrt(-39 - 8*sqrt(2)*sqrt(7 - 4*sqrt(3)) + 4*sqrt(6)*sqrt(7 - 4*sqrt(3)) + 23*sqrt(3)) + sqrt(-1 + sqrt(3))))

Abs = abs
(3*sqrt(3), -3/2 + 3*sqrt(3)/2, -117/2 - 12*sqrt(14 - 8*sqrt(3)) + 6*sqrt(42 - 24*sqrt(3)) + 69*sqrt(3)/2, -15*sqrt(3)/2 + 6*sqrt(14 - 8*sqrt(3)) + 27/2, 3*(-3*sqrt(6) + 4*sqrt(7 - 4*sqrt(3)) + 5*sqrt(2))^2/(4*(sqrt(-39 - 8*sqrt(2)*sqrt(7 - 4*sqrt(3)) + 4*sqrt(6)*sqrt(7 - 4*sqrt(3)) + 23*sqrt(3)) + sqrt(-1 + sqrt(3)))^2), (3*sqrt(-2 + 2*sqrt(3))*(-3*sqrt(6) + 4*sqrt(7 - 4*sqrt(3)) + 5*sqrt(2))/2 + 3*(-1 + sqrt(3))*Abs(sqrt(-39 - 8*sqrt(2)*sqrt(7 - 4*sqrt(3)) + 4*sqrt(6)*sqrt(7 - 4*sqrt(3)) + 23*sqrt(3)) + sqrt(-1 + sqrt(3)))/2)/Abs(sqrt(-39 - 8*sqrt(2)*sqrt(7 - 4*sqrt(3)) + 4*sqrt(6)*sqrt(7 - 4*sqrt(3)) + 23*sqrt(3)) + sqrt(-1 + sqrt(3))))

   (5.196152422706632, 1.098076211353316, 0.6465370682525275, 2.7832432350115006, 0.2069948835295971, 2.0515879469406673)

以下は図を描くプログラムである。

function draw(鈎, 弦, more)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   股 = sqrt(弦^2 - 鈎^2)
   (r1, r2, x2) = (-弦/2 + 鈎/2 + 股/2, (-弦*(4*sqrt(弦)*sqrt(-(弦 - 鈎)*(-2*弦^2 + 2*弦*股 + 鈎^2)) + 4*弦^2 - 3*弦*鈎 - 4*弦*股 - 鈎^2 + 3*鈎*股) + 2*鈎*(-弦^2 + 弦*股 + 鈎^2) + 股*(4*sqrt(弦)*sqrt(-(弦 - 鈎)*(-2*弦^2 + 2*弦*股 + 鈎^2)) + 4*弦^2 - 3*弦*鈎 - 4*弦*股 - 鈎^2 + 3*鈎*股))/(2*鈎^2), (4*sqrt(弦)*sqrt(-(弦 - 鈎)*(-2*弦^2 + 2*弦*股 + 鈎^2)) + 4*弦^2 - 3*弦*鈎 - 4*弦*股 - 鈎^2 + 3*鈎*股)/(2*鈎))
   (r3, x3) = ((2*sqrt(弦)*sqrt(2*弦^3 - 2*弦^2*鈎 - 2*弦^2*股 - 弦*鈎^2 + 2*弦*鈎*股 + 鈎^3) + 2*弦^2 - 弦*鈎 - 2*弦*股 - 鈎^2 + 鈎*股)^2/(2*(鈎*sqrt(-弦 + 鈎 + 股) + sqrt(-4*弦^(3/2)*sqrt(2*弦^3 - 2*弦^2*鈎 - 2*弦^2*股 - 弦*鈎^2 + 2*弦*鈎*股 + 鈎^3) + 4*sqrt(弦)*股*sqrt(2*弦^3 - 2*弦^2*鈎 - 2*弦^2*股 - 弦*鈎^2 + 2*弦*鈎*股 + 鈎^3) - 8*弦^3 + 4*弦^2*鈎 + 8*弦^2*股 + 5*弦*鈎^2 - 4*弦*鈎*股 - 鈎^3 - 鈎^2*股))^2), -弦/2 + 鈎/2 + sqrt((2*sqrt(弦)*sqrt(2*弦^3 - 2*弦^2*鈎 - 2*弦^2*股 - 弦*鈎^2 + 2*弦*鈎*股 + 鈎^3) + 2*弦^2 - 弦*鈎 - 2*弦*股 - 鈎^2 + 鈎*股)^2/(鈎*sqrt(-弦 + 鈎 + 股) + sqrt(-4*弦^(3/2)*sqrt(2*弦^3 - 2*弦^2*鈎 - 2*弦^2*股 - 弦*鈎^2 + 2*弦*鈎*股 + 鈎^3) + 4*sqrt(弦)*股*sqrt(2*弦^3 - 2*弦^2*鈎 - 2*弦^2*股 - 弦*鈎^2 + 2*弦*鈎*股 + 鈎^3) - 8*弦^3 + 4*弦^2*鈎 + 8*弦^2*股 + 5*弦*鈎^2 - 4*弦*鈎*股 - 鈎^3 - 鈎^2*股))^2)*sqrt(-弦 + 鈎 + 股) + 股/2)
   @printf("鈎 が %g,弦 が %g のとき,丁円の直径は %g である。\n", 鈎, 弦, 2r3)
   @printf("鈎 = %g;  弦 = %g;  股 = %g;  r1 = %g;  r2 = %g;  x2 = %g;  r3 = %g;  x3 = %g\n", 鈎, 弦, 股, r1, r2, x2, r3, x3)
   println("r3 = $r3")
   plot([0, 股, 0, 0], [0, 0, 鈎, 0], color=:magenta, lw=0.5)
   circle(r1, r1, r1)
   circle(x2, r2, r2, :blue)
   circle(x3, r3, r3, :green)
   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, r1, "甲円:r1,(r1,r1)", :red, :center, delta=-delta/2)
       point(x2, r2, "乙円:r2,(x2,r2)", :blue, :center, delta=-delta/2)
       point(x3, r3, "丁円:r3,(x3,r3)", :green, :center, delta=-5delta)
       point(0, 鈎, "鈎 ", :magenta, :right, :vcenter)
       point(股, 0, " 股", :magenta, :center, delta=-delta)
       dimension_line(delta, 鈎 + 2delta, 股 + delta, 2delta, "弦", :gray50, delta=2delta, deltax=2delta)
       plot!(xlims=(-8delta, 股 + 3delta), ylims=(-6delta, 鈎 + 4delta))
   end  
end;

draw(3,6, true))

 

 

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

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

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