楕円の共通内接線を求める式を教えてください

質問

大学の卒業研究で悩んでいます。
円柱の連続体を描画するプログラムを作成するにあたって、楕円同士の共通内接線を求めたいのです。
楕円単体の接線までは算出できました。でも、それを二つの楕円で共有させようとするとお手上げの状態です。
楕円の接線の式は、
円の中心を(x0,y0)、接線の座標を(x1,y1)、円の長径(x)をa、短径(y)をbとおいて、
X=(x1-x0)/a^2 、 Y=(y1-y0)/b^2 として
y=-(X/Y)x + (X/Y)x0+y0+(1/Y)
となりました。
己の計算力の低さが恨めしい今日この頃です。どうかどなたか、分かりやすく教えていただけないでしょうか？

回答No.4

tecchan22 · ベストアンサー

けっこう面倒くさそうですね。

まず、二つの楕円の式を（一方を中心が原点に来るように平行移動して）
x^2/a^2+y^2/b^2=1 , (x-x0)^2/c^2+(y-y0)^2/d^2=1
とおく。
さらにｘ座標を 1/c 倍、ｙ座標を 1/d 倍する変換（一次変換）
X=x/c ,Y=y/d　　　 ・・・(1)
を行うと、
(1)より x=cX ,y=dY を上の式に代入して、
{(c/a)^2}X^2+{(d/b)^2}Y^2=1 , (X-x0/c)^2+(Y-y0/d)^2=1
Ｘ，Ｙをｘ，ｙにもどし、c/a=A , d/b=B ,x0/c=s ,y0/d=t とおくと、
楕円 A^2x^2+B^2y^2=1 ，円 (x-s)^2+(y-t)^2=1　　　　　　・・・(2)
になる。
一次変換で接線は接線にうつるから、(2)で共通内接線を求めて、(1)の逆変換をしてやればよい。

(2)の楕円上の点を (cosθ/A,sinθ/B) とパラメーター表示する（０≦θ＜２π）と、
その点での接線の式は、A^2(cosθ/A)x + B^2(sinθ/B)y=1
整理して、Acosθx+Bsinθy=1　　　　　・・・(3)
これが(2)の円の接線となるためには、
(s,t) からの距離が１であればよいから、
|Ascosθ+Btsinθ-1|/√(Acosθ)^2+(Bsinθ)^2=1
分母を払って、|Ascosθ+Btsinθ-1|＝√{(Acosθ)^2+(Bsinθ)^2} 　　・・・(4)
両辺を二乗して、(Ascosθ+Btsinθ-1)^2=(Acosθ)^2+(Bsinθ)^2
これを展開して整理し、(sinθ)^2=1-(cosθ)^2・・(5)を代入して sinθ の一次式にし、残った sinθ を含む項たちと、そうでない項たちに左右に分けて再度二乗し、(5)を代入して cosθ だけの式にすると、
cosθの４次方程式が得られる。（この解の cosθ/A が、共通内接線及び共通外接線の、接点のｘ座標となる。だから４つ出てくる。）

さて、４次方程式の解の公式は四則とべき根だけで表されるから、それを用いて解を計算することは可能だが、多分現実的ではない。

サブルーチンとしてマセマティカにとばして計算させることなんて無理かな？（僕はプログラムには昏い。たとえ可能でも反則か？）

ならば、＃３さんの言われるように、数値解を近似で求めるのが良いか。

二つの楕円が離れた位置にあれば、通常異なる４解を－１～１の範囲に持つ（たまに重解になる）が、外側の２解α，βは、ニュートン法で１，－１から始めれば出るから、求めたα，βに対して４次方程式を (cosθ-α)(cosθ-β) で割れば、「ほとんど」割り切れて、cosθの２次方程式になり、めでたく４つの近似解が求まる。

さて、これらのうち、どれが共通「内」接線になるかだが、
(3)⇔ Acosθx+Bsinθy-1=0
より、f(x,y) = Acosθx+Bsinθy-1 とおけば、接線の式は f(x,y)=0 で、
これが（共通接線のとき）共通内接線であるためには、二つの楕円と円の中心 (0,0) と (s,t) が、ｆによって正領域と負領域に分かれればよい。
いま f(0,0)=-1<0 だから、f(s,t)>0であればよい。
つまり、Ascosθ+Btsinθ-1>0 であればいいので、
(4)から、、Ascosθ+Btsinθ-1＝√{(Acosθ)^2+(Bsinθ)^2}　・・★
を満たせばよい。
つまり、求めた cosθの近似解たち（ふつう４つ）に対して、sinθを±両方求め、
それら全て（ふつう八組）に対して、★を満たすものを求めれば、
めでたく cosθと sinθの組が二つ求まるという訳です。

あとはそれを(3)に代入し、逆変換を行い、（　(cosθ/a)x+(sinθ/b)y＝1 になりますね。）、最後にもとの位置に平行移動すれば、求まります。

ああ大変。
こんな方法は下手かも知れませんが、まあ一つの参考になれば。
長くなって済みません。

お礼

うわあ。すごいです……
やはり面倒だったんですね。途中の計算から頭がパンクしそうです。
でも、とりあえず気合入れて、解読させていただきます。
どうも、ありがとうございました。

回答No.5

kts2371148

＃３です。
二分法よりもいい方法を思いつきました。

２つの楕円を、楕円Ａ・楕円Ｂとします。
（傾いていてもＯＫです。）

まず、楕円Ａ上に、任意の点Ｐをとります。
同様に、楕円Ｂ上に、任意の点Ｑをとります。
（求める接線に近いと思われる点をとるとよいでしょう。）

直線ＰＱと楕円Ａの交点は２箇所あります。
接する場合は１箇所ですが、重解とみなして２箇所にします。
そこで、直線ＰＱと楕円Ａの交点をＰ１，Ｐ２とします。
同様に、直線ＰＱと楕円Ｂの交点をＱ１，Ｑ２とします。

楕円Ａ上におけるＰ１，Ｐ２の中間点を、新しく点Ｐとします。
楕円Ｂ上におけるＱ１，Ｑ２の中間点を、新しく点Ｑとします。

これを繰り返すと、点Ｐ，点Ｑは収束すると思われます。
その場合、直線ＰＱが接線になります。

これはかなりシンプルな方法なので、プログラムも容易だと思います。
また、概略だけを書きましたが、楕円Ａと楕円Ｂの位置関係によって
（接線の本数が０～４本のどれになるかによって）、工夫が必要になるでしょう。

上記の方法では中間点を新しい点にしていますが、
それ以外の点を新しい点にしたほうがいいかもしれません。
一つ、方法を考えてみました。

楕円Ａ上におけるＰ１，Ｐ２の中間点をＰ３とします。
楕円Ｂ上におけるＱ１，Ｑ２の中間点をＱ３とします。
そして、Ｐ１，Ｐ２，Ｐ３の３点を通る円弧Ａ’によって楕円Ａを近似し、
Ｑ１，Ｑ２，Ｑ３の３点を通る円弧Ｂ’によって楕円Ｂを近似します。
円Ａ’と円Ｂ’の共通接線を求めるのは容易です。
そして、その共通接線と同じ傾きで、
　楕円Ａに引いた接線が楕円Ｂと交わるなら、楕円Ａ上の接点を点Ｐ、楕円Ｂとの交点を点Ｑとする
　楕円Ｂに引いた接線が楕円Ａと交わるなら、楕円Ａとの交点を点Ｐ、楕円Ｂ上の接点を点Ｑとする
ようにした方が、収束が速いのではないかと思います。

あとは試行錯誤して、バランスのとれた方法を見つけるといいのではないかと思います。

お礼

ありがとうございます。
なんとかします。

回答No.3

kts2371148

よく考えると、今回の場合は厳密解を求める必要はないんですよね。
ということは、こんな方法はいかがでしょうか？

２つの楕円を、楕円Ａ・楕円Ｂとします。
そして、とりあえず楕円Ａと楕円Ｂの間を通る接線は考えないことにします。
また、とにかく突破口が欲しい状態ということですので、傾いた楕円は考えないことにします。

接線の式を、(cosθ)x + (sinθ)y = c とおきます。
そして、楕円上のすべての点は、(cosθ)x + (sinθ)y ≦ c の範囲内にあるものとします。
θの範囲は 0 ≦ θ < 2π ですが、とりあえず 0 ≦ θ ≦ π/2 の範囲を考えます。

楕円Ａ上のすべての点が (cos θ)x + (sin θ)y ≦ c を満たすときの
c の最小値はθの関数となりますので、これをを cA(θ) とおきます。
楕円Ｂ上のすべての点が (cos θ)x + (sin θ)y ≦ c を満たすときの
c の最小値はθの関数となりますので、これをを cB(θ) とおきます。
そして、もし cA(θ) = cB(θ) であれば、そのθにおいて共通接線が存在します。

さて、0 ≦ θ ≦ π/2 の範囲を考えていますので、
cA(0),cB(0),cA(π/2),cB(π/2) を求めます。

もし、cA(0) < cB(0) なのに cA(π/2) > cB(π/2) であれば、
どこか途中に共通接線があったことになります。
そこで、今度は cA(π/4) と cB(π/4)を比較します。
その結果から判断して、次は θ = π/8 か θ = 3π/8 のどちらかで比較します。
これを繰り返すと、どこで共通接線があったか、近似値を求めることができます。
これは、二分法と呼ばれる解き方です。
http://ja.wikipedia.org/wiki/%E4%BA%8C%E5%88%86%E6%B3%95

では、cA0 < cB0 かつ cA1 < cB1 の場合はどうでしょうか。
この場合は、共通接線がない場合もあれば、１本ある場合も、２本ある場合もあります。

0 ≦ θ ≦ π/2 において、cB(θ) - cA(θ) の値は、
　増加の後、減少する
　減少の後、増加する
　変化なし
のどれかです。もし、「減少の後、増加する」のであれば、
cB(θ) - cA(θ) の最小値が０または負になることもあり得ます。
そこで、二分法の考え方で、cB(θ) - cA(θ) の最小値がどこになるかを突き止めます。
もしその値が負であれば、cB(θ) - cA(θ) = 0 となるのは２箇所ありますので、
二分法を使って、その２箇所を求めます。

上記のような感じで、厳密解ではありませんが、近似値を求めることができます。
傾いた楕円の場合でも、範囲のとり方が少し変わるだけで、原理は同じです。
でも、これをプログラムするとなると大変そうですね。

お礼

うーん、ありがとうございます。
二分法ですかー。いちど兄にプログラムを組むときの豆知識として教えられましたけれど、確かに役立ちそうですね。
でも残念なことに、場合分けで上下にふるい分けをするという手法は、今回既に試してしまっているんです。確かに出来そうですが、もう少しスマートなやり方があったらいいなあと考えている状態です。
その上で、kts2371148さんのやり方のほうが堅実な気がします。ちょっと試してみようかな？　と思います。

回答No.2

info22

楕円の方程式を
(x-p)^2/a^2+(y-q)^2/b^2=1
(x-r)^2/c^2+(y-s)^2/d^2=1
、それぞれの楕円上の点(x0,y0),(x1,y1)とおけば
それぞれの楕円の接線は
(x-p)(x0-p)/a^2+(y-q)(y0-q)/b^2=1
(x-r)(x1-r)/c^2+(y-s)(y1-s)/d^2=1
ですが
この接線の式が一致する（係数が比例関係になる）条件から方程式を
立てます。接点の座標(x0,y0),(x1,y1)が楕円の方程式に代入した式の２つをあわせて4つの式が出来ますので、未知数が(x0,y0),(x1,y1)で4個ですから解けば良い事になります。

しかし、文字変数ばかりなので、色々な制約条件が見えにくいでしょう。
例えば2つの楕円が包含関係や交差する場合、一致する場合などを排除する条件も考えなくてはなりません。

(x0,y0),(x1,y1)の解から無効解（複素解も含む）を排除し、共通外接線の条件も排除し、残った共通内接線の条件を見出して、その条件付で、2本の共通内接線の方程式を決定する必要がなります。

最初から文字定数で一般化して解こうとすると、制約条件の具現化の見通しを立てることが困難です。

そのため、最初は、接点(x0,y0)と(x1,y1)以外の定数を全て文字変数ではなく、数値定数を使って問題を設定して、上記の方針に従って、方程式を立てて、接点を求めてみると良いでしょう。それでも無効解が出てきますが、具体的な座標である為、それらの点を結んでみれば、共通内接線以外が容易に排除できるでしょう。

途中、3次方程式が出てきますので、それだけでも大変です。
それが一般化して文字定数を使う場合、3次方程式の根が簡単に表せなくなるかと思います。
まず、定数を数値とした具体的な楕円の方程式から、確実に共通内接線を求める方法をマスターすることから手をつけて見て下さい。
十分それを習得したあと、楕円の位置関係を変えてみて試行し、
その次のステップとして、係数に文字変数を１つずつ追加しては、係数の制約条件を考慮しながら、一般化を進めて下さい。

頑張って下さい。

お礼

確かに。参考になります。
かなり面倒なプロセスに思えますが、今までやってきた計算ではどうしても出せませんでした。具体的な例から少しずつ一般式に近づけていくのですね。分かりました。
丁寧なアドバイス、ありがとうございました。やってみます。

回答No.1

tecchan22

二つの楕円の位置関係が、傾いているケースは考えなくて良いのですか？
（つまり、二つの楕円の長軸が、平行でも垂直でもない、という状況は考えなくて良いのですか？）

それが必要ないなら、あなたの求めた式を連立させればいいので、
計算は大変かも知れませんが、頑張って下さい！

※一方の中心を原点にすると計算が楽かもしれません。
※その上で、中心を原点に持ってきた楕円を一次変換して円にして計算すると、更に計算が楽になるかもしれません。
※ちなみに、
楕円 {(x-x0)^2/a^2}+{(y-y0)^2/b^2}=1 上の (x1,y1) における接線の式は、x の一つを x1 に、y の一つを y1 に置き換えた、
{(x1-x0)(x-x0)/a^2}+{(y1-y0)(y-y0)/b^2}=1
になります。
（あなたの求めたものと同じです。こうすると、覚えやすいでしょう？高校の教科書にのっていますが^^;）

頑張ってくださーい。

お礼

tecchan22さん、ありがとうございます。参考になります。
今のところ、傾いた円のことは考えていません。とにかく突破口が欲しい状態なので……
計算そのものはやっていたのですが、片方の円の中心を原点に置くという方法はすっかり忘れていました。恥ずかしい限りです。
もう少し計算を頑張ろうと思います。感謝します。