算額（その761）

算額（その761）

山形県鶴岡市羽黒町 出羽三山神社 文政6年(1823)
http://www.wasan.jp/yamagata/haguro.html

教材アーカイブス/算数・数学/文献/算額
http://www.s-soarer.jp/?plugin=attach&refer=%E6%95%99%E6%9D%90%E3%82%A2%E3%83%BC%E3%82%AB%E3%82%A4%E3%83%96%E3%82%B9%2F%E7%AE%97%E6%95%B0%E3%83%BB%E6%95%B0%E5%AD%A6%2F%E6%96%87%E7%8C%AE%2F%E7%AE%97%E9%A1%8D&openfile=DSCN1420.JPG

楕円内に2本の斜線と大円（半円）と甲円，乙円，丙円を入れる。乙円と丙円の直径がそれぞれ 3 寸と 2 寸（注）のとき，甲円の直径はいかほどか。

注：前者の URL では数字が判然とせず，検索の結果後者がヒットし，数字が確定できた。

楕円の長半径と短半径を a, b
甲円の半径と中心座標を r1, (0, b - r1)
乙円の半径と中心座標を r2, (0, r2 - b)
丙円の半径と中心座標を r3, (a - r3, 0)
半円の半径と中心座標を r0, (0, b - r0)
とおき，以下の連立方程式を解く。

include("julia-source.txt");

using SymPy
@syms a::positive, b::positive, r0::positive,
     r1::positive, r2::positive, r3::positive,
     x1::positive, x2::positive,
     y1::positive, y2::positive
eq1 = eq1 = r0^2/a^2 + (b - r0)^2/b^2 -1
eq2 = (a - r3)^2 + (b - r0)^2 - (r0 + r3)^2
eq3 = distance(x1, y1, -x2, -y2, 0, b - r1) - r1^2
eq4 = distance(x1, y1, -x2, -y2, 0, r2 - b) - r2^2
eq5 = distance(x1, y1, -x2, -y2, a - r3, 0) - r3^2
eq6 = x1^2/a^2 + y1^2/b^2 - 1
eq7 = x2^2/a^2 + y2^2/b^2 - 1
eq8 = 2r2 + r0 - 2b

using NLsolve

function nls(func, params...; ini = [0.0])
   if typeof(ini) <: Number
       r = nlsolve((vout, vin) -> vout[1] = func(vin[1], params..., [ini]), ftol=big"1e-40")
       v = r.zero[1]
   else
       r = nlsolve((vout, vin)->vout .= func(vin, params...), ini, ftol=big"1e-40")
       v = r.zero
   end
   return Float64.(v), r.f_converged
end;

function H(u)
   (a, b, r0, r1, x1, x2, y1, y2) = u
   return [
-1 + (b - r0)^2/b^2 + r0^2/a^2,  # eq1
(a - r3)^2 + (b - r0)^2 - (r0 + r3)^2,  # eq2
-r1^2 + (x1*(x1^2 + 2*x1*x2 + x2^2 + y1^2 + 2*y1*y2 + y2^2) - (x1 + x2)*(-b*y1 - b*y2 + r1*y1 + r1*y2 + x1^2 + x1*x2 + y1^2 + y1*y2))^2/(x1^2 + 2*x1*x2 + x2^2 + y1^2 + 2*y1*y2 + y2^2)^2 + (b - r1 - (b*y1^2 + 2*b*y1*y2 + b*y2^2 - r1*y1^2 - 2*r1*y1*y2 - r1*y2^2 - x1^2*y2 + x1*x2*y1 - x1*x2*y2 + x2^2*y1)/(x1^2 + 2*x1*x2 + x2^2 + y1^2 + 2*y1*y2 + y2^2))^2,  # eq3
-r2^2 + (x1*(x1^2 + 2*x1*x2 + x2^2 + y1^2 + 2*y1*y2 + y2^2) - (x1 + x2)*(b*y1 + b*y2 - r2*y1 - r2*y2 + x1^2 + x1*x2 + y1^2 + y1*y2))^2/(x1^2 + 2*x1*x2 + x2^2 + y1^2 + 2*y1*y2 + y2^2)^2 + (-b + r2 - (-b*y1^2 - 2*b*y1*y2 - b*y2^2 + r2*y1^2 + 2*r2*y1*y2 + r2*y2^2 - x1^2*y2 + x1*x2*y1 - x1*x2*y2 + x2^2*y1)/(x1^2 + 2*x1*x2 + x2^2 + y1^2 + 2*y1*y2 + y2^2))^2,  # eq4
-r3^2 + (x1 + x2)^2*(a*y1 + a*y2 - r3*y1 - r3*y2 - x1*y2 + x2*y1)^2/(x1^2 + 2*x1*x2 + x2^2 + y1^2 + 2*y1*y2 + y2^2)^2 + (a - r3 - ((a - r3)*(x1^2 + 2*x1*x2 + x2^2 + y1^2 + 2*y1*y2 + y2^2) - (y1 + y2)*(a*y1 + a*y2 - r3*y1 - r3*y2 - x1*y2 + x2*y1))/(x1^2 + 2*x1*x2 + x2^2 + y1^2 + 2*y1*y2 + y2^2))^2,  # eq5
-1 + y1^2/b^2 + x1^2/a^2,  # eq6
-1 + y2^2/b^2 + x2^2/a^2,  # eq7
-2*b + r0 + 2*r2,  # eq8
   ]
end;

(r2, r3) = (3, 2) .// 2
iniv = BigFloat[14, 7, 12, 6, 15, 6, 2, 7]
res = nls(H, ini=iniv)

   ([10.772733743669523, 6.13102420004887, 9.26204840009774, 4.383884823343673, 10.493949340307568, 5.4542418519895035, 1.3857702087399792, 5.287131918198976], true)

乙円と丙円の直径がそれぞれ 3 寸と 2 寸のとき，甲円の直径は 8.76777 寸となった。

「答」では「甲円の直径は十一寸零三分八厘六毛四一九有奇」となっており，計算結果と合わない。連立方程式に使った条件式のどれかが不適切なのかもしれない。

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

   甲円の直径 = 8.76777;  乙円の直径 = 3;  丙円の直径 = 2
   a = 10.7727;  b = 6.13102;  r0 = 9.26205;  r1 = 4.38388;  x1 = 10.4939;  x2 = 5.45424;  y1 = 1.38577;  y2 = 5.28713

function draw(more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   #(r2, r3) = (0.17, 0.12)
   (a, b, r0, r1, x1, x2, y1, y2) = res[1]
   @printf("甲円の直径 = %g;  乙円の直径 = %g;  丙円の直径 = %g\n", 2r1, 2r2, 2r3)
   @printf("a = %g;  b = %g;  r0 = %g;  r1 = %g;  x1 = %g;  x2 = %g;  y1 = %g;  y2 = %g\n", a, b, r0, r1, x1, x2, y1, y2)
   plot()
   ellipse(0, 0, a, b, color=:red)
   circle(0, b - r1, r1, :blue)
   circle(0, r2 - b, r2, :green)
   circle2(a - r3, 0, r3, :magenta)
   circle(0, b - r0, r0, :orange, beginangle=0, endangle=180)
   segment(-r0, 2r2 - b, r0, 2r2 - b, :orange)
   segment(x1, y1, -x2, -y2)
   segment(-x1, y1, x2, -y2)
   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 - r1, " 甲円:r1\n (0,b-r1)", :blue, :left, :vcenter)
       point(0, r2 - b, " 乙円:r2\n (0,r2-b)", :green, :left, :vcenter)
       point(r0, b - r0, "(r0,b-r0)", :orange, :right, :bottom, delta=delta/2)
       point(a - r3, 0, "丙円:r3\n(a-r3,0)", :black, :right, delta=-delta/2)
       point(a, 0, " a", :red, :left, :vcenter)
       point(0, b, " b", :red, :left, :bottom, delta=delta)
       point(x1, y1, "(x1,y1) ", :black, :right, :bottom, delta=delta)
       point(x2, -y2, " (x2,-y2)", :black, :left, :vcenter)
   end
end;

算額（その760）

算額（その760）

秋田県仙北市角館町西長野 熊野神社 弘化2年(1846)
http://www.wasan.jp/akita/kakunodatekumano2.html

直角三角形の中に正方形と，直径がそれぞれ 8 寸，6寸の甲円と乙円が入っている。正方形の一辺の長さはいかほどか。

直角三角形の直角を挟むに辺のうち，短い方を「鈎」，長い方を「股」
甲円の半径と中心座標を r1, (股 - a - r1, r1)
乙円の半径と中心座標を r2, (a - r2, a + r2)
正方形の一辺の長さを a
とおき，以下の連立方程式を解く。

include("julia-source.txt");

using SymPy
@syms r1::positive, r2::positive, a::positive,
     鈎::positive, 股::positive
eq1 = r1/a - r2/(鈎 - a)
eq2 = a/(股 - a) - 鈎/股
eq3 = dist(0, 0, 股, 鈎, 股 - a - r1, r1) - r1^2;
res = solve([eq1, eq2, eq3], (a, 鈎, 股))
res[2]

   (r1 + r2 + sqrt(r1^2 + r2^2), (r1 + r2)*(r1 + r2 + sqrt(r1^2 + r2^2))/r1, (r1 + r2)*(r1 + r2 + sqrt(r1^2 + r2^2))/r2)

2 組の解が得られるが，2 番目のものが適解である。
正方形の一辺の長さは 甲円と乙円の直径の二乗和の平方根を取り，甲円と乙円の直径を加える。
甲円と乙円の直径がそれぞれ 8 寸，6寸のとき，正方形の一辺の長さは 12 寸（1 尺 2 寸）である。

(8 + 6 + sqrt(8^2 + 6^2))/2

   12.0

ちなみに，鈎，股はそれぞれ 21, 28 寸である。

function draw(more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   (r1, r2) = (8, 6) .// 2
   (a, 鈎, 股) = (r1 + r2 + sqrt(r1^2 + r2^2), (r1 + r2)*(r1 + r2 + sqrt(r1^2 + r2^2))/r1, (r1 + r2)*(r1 + r2 + sqrt(r1^2 + r2^2))/r2)
   plot([0, 股, 股, 0], [0, 0, 鈎, 0], color=:black, lw=0.5)
   circle(股 - a - r1, r1, r1)
   circle(股 - r2, a + r2, r2, :blue)
   segment(股 - a, 0, 股 - a, a, :green)
   segment(股 - a, a, 股, a, :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(股 - a - r1, r1, "甲円:r1\n(股-a-r1,r1)", :red, :center, delta=-delta/2)
       point(股 - r2, a + r2, "乙円:r2\n(股-r2,a+r2)", :blue, :center, :bottom, delta=delta/2)
       point(股, 0, " 股", :black, :left, :bottom, delta=delta/2)
       point(股, 鈎, "(股,鈎)", :black, :right, :bottom, delta=delta)
       point(股, a, "(股,a)", :green, :right, delta=-delta)
       point(股 - a, 0, " 股-a", :green, :left, :bottom, delta=delta)
   end
end;

算額（その759）

算額（その759）

秋田県仙北市角館町西長野 熊野神社 弘化2年(1846)
http://www.wasan.jp/akita/kakunodatekumano2.html

直径が 4 寸の大円の中に，大円と同じ直径を持つ 2 個の弧と乙円 2 個が入っている。乙円の直径はいかほどか。

大円の半径と中心座標を r1, (0, 0), (0, r1), (0, -r1)
乙円の半径と中心座標を r2, (r1 - r2, 0)
として，以下の方程式を解く。

include("julia-source.txt");

using SymPy
@syms r1::positive, r2::positive
eq1 = (r1 - r2)^2 + r1^2 - (r1 + r2)^2
res = solve(eq1, r2)

   1-element Vector{Sym{PyCall.PyObject}}:
    r1/4

乙円の直径は，大円の直径の 1/4 である。
大円の直径が 4 寸のとき，乙円の直径は 1 寸である。

function draw(more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   r1 = 4//2
   r2 = r1/4
   @printf("乙円の直径 = %g;  大円の直径 = %g\n", 2r1, 2r2)
   plot()
   circle(0, 0, r1)
   circle(0, r1, r1, beginangle=210, endangle=330)
   circle(0, -r1, r1, beginangle=30, endangle=150)
   circle2(r1 - r2, 0, r2, :blue)
   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)
   end
end;

 

算額（その758）

算額（その758）

秋田県仙北市角館町西長野 熊野神社 弘化2年(1846)
http://www.wasan.jp/akita/kakunodatekumano2.html

半梯（半分の等脚台形の横倒し）内に，直角三角形と，天円，地円，人円が入っている。
鈎が 63 寸（注），地円と人円の直径がそれぞれ 28 寸，21 寸のとき，天円の直径はいかほどか。

注：算額の写真が不鮮明で，61 寸なのか 62 寸 なのか，はたまた 63 寸なのかが判然としない。63 寸としたときに，「答」と一致した。そうと分かれば，63 寸のようにも見える。

半梯の高さ，下底，上底がそれぞれ a, b, c
直角三角形の直角である頂点が x 軸と接する x 座標を d
天円の半径と中心座標を r1, (x1, y1)
地円の半径と中心座標を r2, (r2, r2)
人円の半径と中心座標を r3, (a - r3, r3)
として，以下の連立方程式を解く。

include("julia-source.txt");

using SymPy
@syms a::positive, c::positive, cd::positive, bd::positive,
     bc::positive, b::positive, d::positive, r1::positive,
     r2::positive, r3::positive, x1::positive, y1::positive
eq1 = c + (a - d) - cd - 2r3
eq2 = cd + bd - bc - 2r1
eq3 = d + b - bd - 2r2
eq4 = a^2 + (b - c)^2 - bc^2
eq5 = c^2 + (a - d)^2 - cd^2
eq6 = b^2 + d^2 - bd^2
eq7 = (bc + cd + bd)*r1 + d*b + (a - d)*c - (c + b)*a
res = solve([eq1, eq2, eq3, eq4, eq5, eq6, eq7], (r1, a, b, c, d, bc, bd));
res[3]

   (cd*(r2 + r3 - sqrt(r2^2 + r3^2))/(2*r3), cd*r2/(2*r3) + cd/2 + r2 - r2*sqrt(cd^2 - 4*cd*r3 - 4*r3^2)/(2*r3) + r3 + sqrt(cd^2 - 4*cd*r3 - 4*r3^2)/2, r2*(cd + 2*r3 + sqrt(cd^2 - 4*cd*r3 - 4*r3^2))/(2*r3), cd/2 + r3 - sqrt(cd^2 - 4*cd*r3 - 4*r3^2)/2, r2*(cd + 2*r3)/(2*r3) - r2*sqrt(cd^2 - 4*cd*r3 - 4*r3^2)/(2*r3), cd*sqrt(r2^2 + r3^2)/r3, cd*r2/r3)

4 組の解が得られるが，3 番目のものが適解である。
なお，4 番目の解は 3 番目の解と線対称ではないが，上底と下底の位置関係が算額のものと逆である。

「術」は「地円と人円それぞれ直径の二乗の和の平方根を，地円と人円の直径の和から引き，鈎を掛け人円の直径で割る」となっており，cd*(r2 + r3 - sqrt(r2^2 + r3^2))/(2*r3) と一致する。

次に，図を描くために天円の中心座標を求める。
4 組の解が得られるが，最初のものが適解である。

天円の直径は 42 寸である。

using SymPy
@syms a::positive, c::positive, cd::positive, bd::positive,
     bc::positive, b::positive, d::positive, r1::positive,
     r2::positive, r3::positive, x1::positive, y1::positive
(r1, a, b, c, d, bc, bd) = (21, 98 - 7*sqrt(Sym(2))/2, 14*sqrt(Sym(2)) + 56, 42 - 21*sqrt(Sym(2))/2, 56 - 14*sqrt(Sym(2)), 105, 84)
eq8 = distance(0, b, a, c, x1, y1) - r1^2
eq9 = distance(0, b, d, 0, x1, y1) - r1^2
res2 = solve([eq8, eq9], (x1, y1));
res2[1]

   (56 - 7*sqrt(2), 28)

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

   天円の直径 = 42;  地円の直径 = 28;  人円の直径 = 21;  鈎 = 63
   r1 = 21;  a = 93.0503;  b = 75.799;  c = 27.1508;  d = 36.201;  bc = 105;  bd = 84;  cd = 63

function draw(more=false)
   pyplot(size=(500, 500), grid=false, aspectratio=1, label="", fontfamily="IPAMincho")
   cd = 63
   r2 = 28//2
   r3 = 21//2
   (r1, a, b, c, d, bc, bd) = (cd*(r2 + r3 - sqrt(r2^2 + r3^2))/(2*r3), cd*r2/(2*r3) + cd/2 + r2 - r2*sqrt(cd^2 - 4*cd*r3 - 4*r3^2)/(2*r3) + r3 + sqrt(cd^2 - 4*cd*r3 - 4*r3^2)/2, r2*(cd + 2*r3 + sqrt(cd^2 - 4*cd*r3 - 4*r3^2))/(2*r3), cd/2 + r3 - sqrt(cd^2 - 4*cd*r3 - 4*r3^2)/2, r2*(cd + 2*r3)/(2*r3) - r2*sqrt(cd^2 - 4*cd*r3 - 4*r3^2)/(2*r3), cd*sqrt(r2^2 + r3^2)/r3, cd*r2/r3)
   (x1, y1) = (56 - 7*sqrt(2), 28)
   @printf("天円の直径 = %g;  地円の直径 = %g;  人円の直径 = %g;  鈎 = %g\n", 2r1, 2r2, 2r3, cd)
   @printf("r1 = %g;  a = %g;  b = %g;  c = %g;  d = %g;  bc = %g;  bd = %g;  cd = %g\n",
       r1, a, b, c, d, bc, bd, cd)
   plot([0, a, a, 0, 0], [0, 0, c, b, 0], color=:black, lw=0.5)
   circle(x1, y1, r1, :green)
   circle(r2, r2, r2)
   circle(a - r3, r3, r3, :blue)
   segment(0, b, d, 0)
   segment(d, 0, a, c)
   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(a, 0, " a", :black, :left, :bottom, delta=delta/2)
       point(d, 0, "d  ", :black, :right, :bottom, delta=delta/2)
       point(a, c, " (a,c)", :black, :left, :bottom, delta=delta/2)
       point(0, b, " b", :black, :left, :bottom, delta=delta/2)
       point(x1, y1, "天円:r1,(r1,r1)", :green, :center, delta=-delta/2)
       point(r2, r2, "地円:r2,(r2,r2)", :red, :center, delta=-delta/2)
       point(a - r3, r3, "人円:r3\n(a-r3,r3)", :blue, :center, delta=-delta/2)
   end
end;

算額（その757）

算額（その757）

山形県鶴岡市茨新田 鶴岡神明宮 文久2年(1862)
http://www.wasan.jp/yamagata/sinmei.html

盤の上に，日球，月球の2つの球が入っている円錐台がある。円錐台の麓に星球が複数個，互いに接しつつ円錐台の側面にも接して配置されている。
日球，月球の直径が与えられたとき，星球の数と直径を求めよ。

月球の半径と中心座標を r1, (0, 0, r1)
日球の半径と中心座標を r2, (0, 0, 2r1 + r2)
星球の半径と中心座標を r3, (0, x3, 0)
円錐台の底面と上面の半径を a, d
円錐台の側面の延長と z 軸の交点（円錐の頂点）を b
とおき，以下の連立方程式を解く。

まず，日球，月球の半径が与えられたとき，円錐台の底面と上面の半径と円錐の高さを求める。

include("julia-source.txt");

using SymPy
@syms a::positive, b::positive,
     r1::positive, r2::positive,
     c, d, h, n
@syms a, b, r1, r2, r3, x3, c, d, h
(r1, r2) = (9, 4).*10
c = sqrt(a^2 + b^2)
h = 2r1 + 2r2
eq1 = r1/(b - r1) - a/c
eq2 = r2/(b - 2r1 - r2) - a/c
eq3 = h/(a - d) - b/a;
res1 = solve([eq1, eq2, eq3], (a, b, d))

   1-element Vector{Tuple{Sym{PyCall.PyObject}, Sym{PyCall.PyObject}, Sym{PyCall.PyObject}}}:
    (135, 324, 80/3)

次に，星球の数と直径，そのうちの一つが円錐台の側面に外接するときの中心座標を求める。

@syms a::positive, b::positive,
     r1::positive, r2::positive,
     c, d, h, n
@syms a, b, r1, r2, r3, x3, c, d, h, n
(r1, r2) = (9, 4).*10
(a, b, d) = (135, 324, 80//3)
eq4 = distance(0, b, a, 0, x3, r3) - r3^2
eq5 = x3*sind(Sym(360)/2n) - r3
res2 = solve([eq1, eq2, eq3, eq4, eq5], (a, b, d, r3, x3))

   2-element Vector{NTuple{5, Sym{PyCall.PyObject}}}:
    (135, 324, 80/3, -405*sin(pi/n)/(2*sin(pi/n) - 3), -405/(2*sin(pi/n) - 3))
    (135, 324, 80/3, 270*sin(pi/n)/(3*sin(pi/n) + 2), 270/(3*sin(pi/n) + 2))

2 組の解が得られるが，最初のものが適解である。
n を指定すれば，星球の半径と中心座標が求まる。
たとえば n = 8 とすれば，以下のようになる。

res2[1][4](n => 8).evalf() |> println
res2[1][5](n => 8).evalf() |> println

   69.3567045353722
   181.237803023581

n = 8
(135, 324, 80/3, -405*sin(pi/n)/(2*sin(pi/n) - 3), -405/(2*sin(pi/n) - 3))

   (135, 324, 26.666666666666668, 69.3567045353722, 181.23780302358148)

n = 8 のとき，星球の半径は 69.3567045353722，中心座標は 181.237803023581 である。

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

   a = 135;  b = 324;  d = 26.6667;  r3 = 69.3567;  x3 = 181.238

上方からみた星球の様子は以下のようになる。青の円は星球の中心を通る円である。

n = 8
(x3, r3) = (-405/(2*sin(pi/n) - 3), -405*sin(pi/n)/(2*sin(pi/n) - 3))
println((x3, r3))
pyplot(size=(500, 500), grid=false, aspectratio=1, label="")
plot()
delta = (fontheight = (ylims()[2]- ylims()[1]) / 500 * 10 * 2) /3  # size[2] * fontsize * 2
circle(0, 0, x3, :blue)
rotate(x3, 0, r3, angle=360/n)
point(x3, 0, "星球:r3,(0,x3)", :red, :center, delta=-delta/2)

円錐台と日球，月球，星級と円錐台の内接・外接の状況は下図の様になる。

pyplot(size=(500, 500), grid=false, aspectratio=1, label="")
(r1, r2) = (9, 4).*10
h = 2r1 + 2r2
(a, b, d, r3, x3) = (135, 324, 80/3, -405*sin(pi/n)/(2*sin(pi/n) - 3), -405/(2*sin(pi/n) - 3))
@printf("a = %g;  b = %g;  d = %g;  r3 = %g;  x3 = %g\n", a, b, d, r3, x3)
plot([a, 0, -a, a], [0, b, 0, 0], color=:black, lw=0.5)
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)
segment(-d, h, d, h)
circle(0, r1, r1)
circle(0, 2r1 + r2, r2, :blue)
circle(x3, r3, r3, :green)
point(0, r1, "月球:r1,(0,0,r1)", :red, :center, delta=-delta)
point(0, 2r1 + r2, "日球:r2\n(0,0,2r1+r2)", :blue, :center, :bottom, delta=delta)
point(x3, r3, "星球:r3,(x3,0,r3)", :green, :center, delta=-delta)
point(a, 0, "a ", :black, :right, :bottom, delta=delta)
point(0, b, " b", :black, :left, :vcenter)
point(d, 2r1 + 2r2, " (d,2r1+2r2)", :black, :left, :vcenter)

三次元空間の表示は以下のようになる。