数式処理ソフト DERIVE(デライブ) de ドライブ

74.石の安定性(球欠・楕円体球欠)(2)

目次

 1. 第38回「石の安定性(球欠・楕円体球欠)(1)」のおさらい
 2. カバリエリの原理
 3. カバリエリの原理の応用(楕円、楕円体)
 4. カバリエリの原理の応用(楕円弓形)
 5. カバリエリの原理の応用(楕円体球欠)とパップス-ギルダンの定理
 6. 楕円体球欠(斜めに切断)の体積
 7. 楕円体球欠(斜めに切断)の重心

1. 第38回「石の安定性(球欠・楕円体球欠)(1)」のおさらい

「2015年の2月の今月のご挨拶「宿題を片付けよう」の第1弾として、第38回で取り上げた球欠や楕円体球欠の体積や重心について考えてみよう」

「ミャンマーのゴールデンロックがきっかけだったわね。
それにしても、第38回「石の安定性(球欠・楕円体球欠)(1)」は、2008年4月のことだから、ずいぶん前のことね」

「ともちゃんも、暑いところ、ご苦労さん。
 まったく。光陰矢のごとしじゃな。
 最初に、記号を整理しておこうかの。
 
 上図が「球欠」で、球の半径をr、球欠の高さをh、重心の高さをZ、と記載している。
 h=2rの場合は、完全な球体となっておるのじゃな」

「そうだったわね。ちなみに、図のaは、球欠の切り欠いた底面の半径よ。
 このとき、Vを球欠の体積とすると、V=πh2(3r- h)/3、成り立つ。
 また、重心の高さ、Z=3(h- 2r)2/(4(3r- h)) 、という結果だったわね」

「次に「楕円体球欠」じゃ。

 上図が「楕円体球欠」じゃ。
 第38回では、回転楕円体を取り上げている。
 高さ方向の長半径c、水平方向の短半径b、重心の高さZ、楕円体球欠の高さをh、とする。
 このとき、体積をVとすると、V=πb2h2(3c- h)/(3c2)、
 また、重心位置Zは、 Z=3(h- 2c)2/(4(3c- h)) 、となった。
 今回は、楕円体球欠の体積や重心の別解について、触れることからはじめてみよう」

「というわけで、次節へGO!」
目次へ戻る

2. カバリエリの原理

「宿題となっているのは、楕円体球欠を斜めに切断した場合の体積や重心などを求めることじゃ。
 しかし、まずは、前回のおさらいをしておこうと思って、あらためて、考えてみたのじゃ」

「それが、このカバリエリの原理のことね。
 カバリエリ(カヴァリエリ)は、Bonaventura Cavalieri(1571~1630)という人の名前よ。
 矢野健太郎先生の「数学史」(科学新興社:1967年8月1日第1刷)によれば、カバリエリさんは、1635年に、次のような原理を示したの。
 

 それが、上図。
 平面上の図形MとNを、直線で切断したときに、その切断長の割合が常にm:nとなっているとき、図形の面積の比率も、m:nであるというもの。
 そもそも「原理」、って、「定理」よね?」

「ま、そうじゃ。
 原理も定理もほぼ似たようなものじゃな。
 原理の方が根本的なもののように感じられるが、このような名称は、多分に発見された時代の歴史からくるので、みだりに変えると混乱する。
 ここでも「原理」と書いておこう。
 カバリエリの原理は、立体図形についても、直線を平面に置き換えて、長さの比率を断面積の比率とすれば、同様に成り立つ。
 では、平面と立体について、ともちゃん、証明してごらん」

「わかりやした。
 とはいっても、比率って奴は、ちょっと苦手なので、座標軸を書いた方が分かりやすい。
 まずは、平面の場合のカバリエリの原理から説明するよ。
 
 上のように座標軸を導入すれば、図形Mの面積、S(M)=∫m(y)dy、
 Nの面積、S(N)=∫n(y)dy、書けるので、ε=n(y)/m(y)、とすれば、S(N)=∫n(y)dy=∫ε×m(y)dy、
 このとき、ε が定数である(これは比率が一定という意味よ)場合は、S(N)=ε∫m(y)dy=ε×S(M)、となってカバリエリの原理は、説明できる。
 えーと。先回りして言っておくけど、図形Nは、Mと比較して、y軸方向の伸縮は、ないとしている。
 そうでないと、m(y)なりn(y)がゼロとなってしまう区間が生じるので、そもそも、比率が一定ではなくなるでしょ。

 また、図形が、上の図のように、直線で切断したときに複数の部分(赤いところ)に分かれているときも、m(y)やn(y)を複数部分の長さの和を表すものと考えれば、いいのよ。
 次に立体の場合ね。
 平面の場合の長さを断面積に、面積を体積に置き換え、2つの立体の断面積をm(z)、n(z)で表し、
 V(M)=∫m(z)dz、
 V(N)=∫n(z)dz=∫ε×m(z)dz=ε×∫m(z)dz、ここで、ε=n(z)/m(z)=定数であることを使っている。
 よって、V(N)/V(M)=ε、であり、立体図形の場合のカバリエリの原理が証明された」

「では、次節で、楕円と回転楕円体の場合について、適用してみよう」 
目次へ戻る

3. カバリエリの原理の応用(楕円、楕円体)

「まず、楕円の場合ね。
 
 上の図では、分かりやすいように、円(黒)と楕円(青)を重ねて置いてあるわ。
 今、x軸に平行で、位置yの直線と両図形の交点のx座標値(の正の値)を円:m(y)、楕円:n(y)で表す。
 楕円の長軸の長さをc、短軸の長さをbと書いて、図では、c>bとしている。なお、円の場合の半径は、cとしている。
 具体的に確認すると、円の方程式は、x^2+y^2=c^2、楕円は、(x/b)^2+(y/c)^2=1、なので、
 m(y)=√(c^2-y^2)=c×√(1-(y/c)^2)、
、n(y)=b×√(1-(y/c)^2)=ε×c√(1-(y/c)^2)、ただし、ε=b/c、
 よって、n(y)/m(y)=ε=b/c、であり、これは、定数であることが分かる。
 従って、カバリエリの原理より、面積S(楕円)=面積S(円)×ε=πc^2×ε=πbc、となる。
 これは、既知の楕円の面積の公式(π×長半径×短半径)と一致する。

 次に回転楕円体の場合ね。
 

 上図では、球(黒:半径 c)、回転楕円体(青:長軸はz方向の長さ c、短軸はx及びy方向の長さ b)を示している。
 xy平面に水平な平面(位置z=z)と球との交点は、一つの円となり、その面積をm(z)、
 回転楕円体との交点も一つの円となるので、その面積をn(z)と置く。
 球の方程式:x^2+y^2+z^2=c^2、
 回転楕円体の方程式:(x/b)^2+(y/b)^2+(z/c)^2=1、から、
 m(z)=π(x^2+y^2)=πc^2×(1-(z/c)^2)、
 n(z)=πb^2×(x^2+y^2)=πb^2×(1-(z/c)^2)
 よって、n(z)/m(z)=(b/c)^2、となり、これは、定数なので、εと置く。
 すると、球の体積=(4π/3)c^3、から、立体図形の場合のカバリエリの原理を基にして、回転楕円体の体積は、
 V(回転楕円体)=V(球)×ε=(4π/3)c^3×(b/c)^2=(4π/3)b^2×c、なって、既知の公式と一致することが分かる」

「ま、そういうことじゃ。
 蛇足じゃが、楕円体は、回転楕円体でなくてもよい。そのときは、ε=ab/c^2、となる。
 さて、次節では、平面上の弓形や楕円弓形に対して応用してみることにしよう」
目次へ戻る

4. カバリエリの原理の応用(楕円弓形)

「第38回では、平面図形について、取り上げなかったので、まずは、そこから確認してみることにしよう。
 たとえば、円の一部を欠いたもの(ここでは、「弓形」と呼ぶ)、楕円の一部を欠いたもの(ここでは、「楕円弓形」と呼ぶ)などじゃ。
 そこで、下図を参照して欲しいんじゃが、弓形と楕円弓形じゃ。

 黒い曲線は、弓形、Mで、半径cの円をy=kの位置で欠いた弓形じゃ。(-c=<k<=c)
 また、青い曲線は、楕円弓形、Nで、長径は、c、短径は、bの楕円を同じく、y=kの直線で欠いたもの。
 それぞれについて、面積と重心を計算してご覧。
 結果は、h=c-kを弓形の高さとして整理することにしよう。
 0=<h=<2c、じゃ」

「まず、弓形の方だけど、y=kとの交点のx座標(正のもの。以下同様)は、√(c^2-k^2)、
 k=<y<=c、である、y=yの直線との交点 m(y)=√(c^2-y^2)なので、面積は、2倍の係数を忘れずに、
 S(M)=2×∫(y=k~c)(√(c^2-y^2)dy=- c^2ASIN((c- h)/c)- (2(c- h)√(h(2c- h))-πc^2)/2、
 楕円弓形は、カバリエリの原理に従うと、S(N)=(b/c)×S(M)。
 一応検算してみると、楕円:(x^2)^2+(y/c)^2=1、から、y=kとの交点は、ε×√(c^2-k^2)、
 S(N)=2×∫(y=k~c)(ε×√(c^2-y^2)dy=2ε×∫(y=k~c)(√(c^2-y^2)dy
 となるので、S(M)=2×∫(y=k~c)(√(c^2-y^2)dy、と比較すれば、S(N)/S(M)=εとなることは明らかでしょ。
 参考までにグラフを描いてみた。
 c=1、b=0.5、


 上図は、ε=1/2、の場合の弓形と楕円弓形の面積を縦軸で、横軸は、各弓形の高さhを示す。
 hの有効範囲は、h=0~2まで、なので注意して頂戴。

 次に重心について考える。
 対称性から、各弓形の重心は、y軸上にあるので、重心のy座標をG(M)、G(N)で表す。
 定義から、G(M)=2∫(y=k~c)y×(√(c^2-y^2)dy/S(M)
 また、G(N)=2ε×∫(y=k~c)(√(c^2-y^2)dy/S(N)であるが、S(N)=ε×S(M)を使えば、
結局、G(M)=G(N)であることが分かる。
 つまり、弓形と楕円弓形の重心の位置は、同一であることが示された。
 具体的な表式は、G(M)=G(N)=4(h(2c- h))^(3/2)/(3(2c^2ACOS((c- h)/c)- 2(c- h)√(h(2c- h))))、
 こいつは、少し、意外な結果」

「重心については、確かにな。
 
 しかし、G(M)=∫y×S(M)dy/S(M)、G(N)=∫y×S(N)dy/S(N)、なので、
 S(N)=ε×S(M)を使えば、G(N)=∫y×ε×S(M)dy/(εS(M))=∫y×S(M)dy=G(M)、
 よって、平面上の図形でカバリエリの原理が使える場合には、MとNの重心は、同一であることが一般に成り立つことが分かる。
 次節では、(回転)楕円体球欠について、取り上げてみよう」 
目次へ戻る

5. カバリエリの原理の応用(楕円体球欠)とパップス-ギルダンの定理

「立体の場合じゃ。
 第38回で、球欠と(回転)楕円体球欠の体積と重心を求めた。
 楕円体球欠については、長軸に垂直な平面で切断した場合を取り上げた。
 だいぶ、時間が経過しているので、まずは、おさらいからはじめることにしよう。

 上図では、球 M(黒:半径 c)、回転楕円体 N(青:z軸方向の長軸の長さ c、x及びy方向の短軸の長さ b)を示している。
 そして、z=kの水平な面でそれぞれを切断し、その上部を、球欠及び楕円体球欠と考える。
 kの変域は、-c~cまでとする。平面の場合と同様に、球欠、楕円体球欠の高さ h=c-k、と置こう。
 z=zの平面と球との交線は、円であり、その半径は、√(c^2-z^2)から断面積 m(z)=π×(c^2-z^2)、
 従って、球欠の体積 V(M)=∫(z=k~c)m(z)dz=π∫(z=k~c)(c^2-z^2)dz、

 一方、(回転)楕円体球欠では、z=zの水平面との交線は、円であり、その半径は、b×√(1-(z/c)^2)から、
 断面積n(z)=π×(b/c)^2×(c^2-z^2)、
 よって、楕円体球欠の体積 V(N)=∫(z=k~c)n(z)dz=π(b/c)^2∫(z=k~c)(c^2-z^2)dz、
 すなわち、楕円体球欠の体積=球欠の体積×ε、であることが分かる。ここで、ε=(b/c)^2は断面積の比率である。
 これは、立体の場合のカバリエリの原理から導かれる式と一致しておるじゃろう」

「念のため、具体的に計算すると、前回に出てきた式、V(N)=πh^2×ε×(3c - h)/3、と同一となる。※」
 ※第38回は、ε=(b/c)^2を使っていないので、上の式でεを(b/c)^2とすれば、同一となる。

「実はな、今回の記事を書くために、第38回の記事を読み返したところ、球欠の体積と楕円体球欠の式が定数倍しか違わないことに気がついたのじゃよ。
 それで、ようやく、ぴんときたんじゃ」

「それで、カバリエリの原理なのね」

「うん。昔、何かで読んだ記憶があったのじゃが、その人の名前を忘れてしまっていてな。
なにか樽の体積がどうのこうのという話題と一緒に記載されていたような記憶があったのじゃ。
 そこから調べたら、すぐに出てきた※」
※ 「数学史」(矢野健太郎:科学新興社:1967年8月1日初版) 103ページ。

「その樽がどうのこうのというのは、何だったの?」

「パップス-ギルダンの定理というものじゃ。
 お酒の樽(西洋のことじゃから、ワインとかウイスキー樽のようなものを思い浮かべて欲しい)のような回転対称性を持つ立体の表面積と体積を求める定理じゃ。
 
 上図のようなイメージじゃな。
 パップス-ギルダンの定理は、樽の体積=2π×z軸からの重心までの長さ G×断面の面積 S、
 また、側面積=2π×側面の長さ×重心までの長さ G、(表面積は、樽のふたとそこの部分の円の面積を加算すればよい)
 という2つからなる。
 ちなみに、パップスは、Pappus, 西暦340年頃のギリシャ-ローマ時代の数学者。
 ギルダンは、パップスの得たこの定理を再発見したことで知られている。Paul Guldin,1577-1643 オーストリアの数学者。
 円筒座標系(r,z,φ)を考える。
 z軸から断面の先端までの長さをf(z)とすれば、体積については、
 

 上の図のように、樽の体積要素は、微少な3角柱の体積で、(1/2)f(z)×f(z)×dφ×dz、
 従って、体積=2π×(1/2)∫f(z)^2×dz、ここで、zの変域は、z=k~cまでとする。
 一方、z軸から断面の重心までの距離をGとすれば、重心の定義により、G=∫∫(0~f(z))r×dr×dz/S、
 rについては、積分できて、L=(1/2)∫f(z)^2×dz/S、よって、(1/2)∫f(z)^2×dz=G×S、
 これから、体積=2π×(1/2)∫(z))r^2×dz=2π×G×S、
 以上じゃ。
 側面積は、ともちゃんが考えてごらん。
 ヒントは、厚さがものすごく薄い膜のような立体の側面積を求めることを体積の式からやってみるといい」

「厚さがΔ(極めて小さい定数)の膜からなる図形を考えて、側面の長さをLとすると、断面積 S=Δ×L、
 一方、体積=2π×S×G=2π×Δ×G×L、であるが、この薄い図形の体積=側面積×Δとも考えられるので、
 側面積=体積/Δ=2π×S×G、と求められる」

「ま、そういうことじゃ。
 今回の楕円体球欠は、回転楕円体なので、パップス-ギルダンの定理で体積を求めることもできる。
 L=(1/2)∫f(z)^2×dz/S、なんじゃが、f(z)は、楕円の方程式、(y/b)^2+(z/c)^2=1、のyにあたるので、
y^2=f(z)^2=b^2(1-(z/c)^2)=(b/c)^2×(c^2-z^2)、
 これから、L=(1/2)∫(z=k~c)((b/c)^2×(c^2-z^2))dz/S、
 よって、楕円体球欠の体積は、2π×L×S=π×∫(z=k~c)((b/c)^2×(c^2-z^2))dz、となって、
この節の前半で出てきた式と確かに、一致する。
 また、寄り道をしてしまったな。
 次節では、本題の楕円体の長軸に斜めに切断した場合の楕円体球欠の体積を求めよう」
目次へ戻る

6. 楕円体球欠(斜めに切断)の体積

「下図は、回転楕円体 (x/b)^2+(y/b)^2+(z/c)^2=1、をx軸の周りを反時計回りにθだけ回転した図形じゃ。


 元の楕円体上の点(x,y,z)は、傾けた新しい図形上の点(x',y',z')に写される。
 x'=x、
 y'=y×cosθ-z×sinθ
 z'=y×sinθ+z×cosθ
 これは、逆には、
 x=x'、
 y=y'×cosθ+z×sinθ、
 z=-y'×sinθ+z×cosθ、
 なので、(x/b)^2+(y/b)^2+(z/c)^2=1、に代入して、’をとると、
 ((y^2 + εz^2)COS(θ)^2 + 2yz(1- ε)SIN(θ)COS(θ) + (εy^2 + z^2)SIN(θ)^2 + x^2)/(c^2ε) = 1、
 ここで、ε≡b^2/c^2、である。
 この新しい楕円体を「斜楕円体」と呼ぶとすると、斜楕円体とz=zの水平面が交わる場合は、その交線もまた楕円じゃ。
 では、その楕円の式を導いてご覧」

「斜楕円体の方程式が、
 ((y^2 + εz^2)COS(θ)^2 + 2yz(1- ε)SIN(θ)COS(θ) + (εy^2 + z^2)SIN(θ)^2 + x^2)/(c^2ε) = 1、なので、
 これをxとyを変数とする楕円の標準形に変形する。
 この作業は、DERIVEでは、なかなか難しいんで、手作業に頼ると、次式を得ることができるの。
 x^2/(√(α×δ))^2+(y+β/α)^2/(√(δ))^2=1、
 ここで、α≡(cos(θ))^2+ε×(sin(θ))^2、
 β≡(1-ε)sin(θ)cos(θ)
 γ≡ε(cos(θ))^2+×(sin(θ))^2、
 δ≡(εc^2/α)+z^2×(β^2/α^2)-γz^2/α、
 従って、楕円の面積は、S=π√(αδ)×√δ=πδ√α、
 S=((z^2(ε^2 -2ε + 1)COS(θ)^4 + COS(θ)^2(z^2(ε -1)^2SIN(θ)^2 -z^2(ε^2 -ε + 1) + c^2ε) + ε(c^2ε -z^2)SIN(θ)^2)/(COS(θ)^2 + εSIN(θ)^2)^(3/2))、となる」

「複雑な式じゃが、cos(θ)をsin(θ)で表すよう、DERIVEのオプションメニューの「モード設定」で変更すると、
 S=ε(c^2(ε - 1)SIN(θ)^2 -z^2 + c^2)/((ε - 1)SIN(θ)^2 + 1)^(3/2)、と簡単化される。
 θ=0の場合、S=πε(c^2-z^2)、これは、半径=√(ε(c^2-z^2))の円である。zの変域も±cであり、正しい。
 θ=π/2の場合、S=π(εc^2-z^2)/√ε、これは、zの変域が、±c√εであり正しく、z=0と置いて、S=πc^2×√ε=π×c×c√ε、となるので、長径=c、短径=c√εの楕円であるので、これも合っているようじゃな。
 c=1、θ=π/4、ε=1/4の場合のグラフを描いてみた。
 縦軸がz、横軸がy、原点を通る斜線は、z=zの切断面の楕円の中心の軌跡じゃ」


「なるほどね。
 これで、ようやく、斜楕円体球欠の体積Vが次式で求められるわ。
 V=∫(z=k~zの最大値)S dz、ただし、kは、-zの最大値~zの最大値までのパラメータ、
 zの最大値/最小値は、δ=(εc^2/α)+z^2×(β^2/α^2)-γz^2/α=0から、z=±c×√((cos(θ))^2+ε(sin(θ))^2)、
 よって、斜楕円体球欠の体積は、
 V=ε(2c^3((ε -1)SIN(θ)^2 + 1)^(3/2) + 3c^2k(1 -ε)SIN(θ)^2 -k(3c^2 -k^2))/(3((ε -1)SIN(θ)^2 + 1)^(3/2))、
 斜楕円体球欠の高さ h=zの最大値-k、を使うと、
 V=h^2ε(3c√((ε -1)SIN(θ)^2 + 1) -h)/(3((ε -1)SIN(θ)^2 + 1)^(3/2))、となる。
 特別な場合、h=2×zの最大値、の時は、V=4πc^3ε/3、と定数であり、楕円体の体積と一致することも分かる」

「かなり、簡単な形となったな。
 h=2×zの最大値の時というのは、斜楕円体の体積全体なので、これは、回転していない場合の体積と一致するのは、当然であるから、検算の一つになる。
 Vの式をDERIVE画面では、
 
 となっている。この方が見やすいじゃろう。
 ここで、hは、斜楕円体球欠の高さ、cは、元の楕円体の縦軸方向の長さ、ε=短径^2/長径^2、θは、x軸の周りの反時計回りの回転角度。
 θ=0の時は、V=h^2ε(3c- h)/3、となり、第38回及び(5)の式と一致することも分かる」 
目次へ戻る

7. 楕円体球欠(斜めに切断)の重心

「前節で斜楕円体球欠の体積が求められた。
そこで、斜楕円体球欠の重心位置を計算して欲しい」

「対称性から考えて、重心は、x=0面にある。重心の位置を(GY、GZ)で表す。
 GZ=∫(z=zの最大値-h~zの最大値)z×S dz/V、ただし、Sは、前節で求めた水平面との交線が作る楕円の面積、Vは、斜楕円体球欠の体積、
 GZ=3(4ch√(COS(θ)^2 + εSIN(θ)^2) + 4c^2(1 - ε)SIN(θ)^2 - 4c^2 - h^2)/(4(h - 3c√(COS(θ)^2 + εSIN(θ)^2)))
 θ=0の場合は、GZ=3(2c-h)^2)/(4(3c - h))、となって、第38回の結果と一致する。
 GYを求めるために、水平面との断面の楕円の中心は、y=-βz/α=(z(ε-1)SIN(θ)COS(θ)/((ε-1)SIN(θ)^2 + 1))、であることを用いる。
 この断面の重心のy座標は、中心であるので、
 GY=∫(z=zの最大値-h~zの最大値)y×S dz/V、から、
 GY=h^2ε(1-ε)SIN(θ)COS(θ)(4ch√((ε-1)SIN(θ)^2 + 1) + 4c^2(1-ε)SIN(θ)^2-4c^2-h^2)/(4((ε-1)SIN(θ)^2 + 1)^(5/2))
 θ=0及びπ/2の場合は、GY=0となる」

「この後、斜楕円体球欠の静的安定性を調べる予定じゃったが、だいぶ長くなってきたので、省略しよう。
 いや、ともちゃんも、お疲れさんじゃった」
目次へ戻る

最終更新日 2015/8/17 一部追記 2017/4/7