楕円と直線の交点と曲率半径の関係

楕円曲線式　x^2/a^2+y^2/b^2=1　（原点を中心として、x軸長2a、y軸長2bの楕円）
のaが特定されているとして

上記楕円曲線と
直線　y=xtanθ　（原点を通る、特定の傾きθの直線）
との交点における
楕円の曲率半径が　特定値r　であるときに

上記楕円曲線式のbを算出する式を
教えてください　

No.1 回答者： 178-tall 回答日時： 2012/09/29 16:11

パラメータ表示 (x, y) = (a*cosθ, b*sinθ) なら、

曲率半径 r = {(a^2 + b^2)cos^2(θ) - a^2}^(3/2)/(ab)
になるらしい。

風聞だけで自分じゃ確かめていないが、似たような形になるのだろう。
{a, θ} を与えて、b^2 の 3 次方程式、かな？

　　

No.2 回答者： alice_44 回答日時： 2012/09/29 19:19

楕円の曲率の導出は、↓の pp.21-22 にあるよ。

www007.upp.so-net.ne.jp/tmh_ogaw/mymath/TheoryOfCurves.pdf

尤も、質問の問題では、
r = 1/κ = { (a sin φ)^2 + (b cos φ)^2 }^(3/2) / (ab)
の式から
(a cos φ, b sin φ) = (L cos θ, L sin θ)
を使って φ を消去しなければならないから、θ の使い方には気をつけて。

φ と L を消去すると、
(r^2){ (b cos θ)^2 + (a sin θ)^2 }^3 = ab{ (a^2 sin θ)^2 + (b^2 cos θ)^2 }^3
これは、b については 7 次方程式みたいだけど、こんなもの解けるのだろうか？


#　どうでもいいことだが、
#　楕円は、二次曲線のひとつであって、楕円曲線ではない。

No.3 回答者： 178-tall 回答日時： 2012/09/29 21:11

これは「実務系」の問題とお見受けしました。

どの座標系を選択しても、特殊なケース以外では高次方程式を解かざるを得ないようです。
手馴れた座標にて逐次解法を駆使されるのがベスト、だと愚考します。

直交系での曲率算式もあちこちあるようで。
一例だけ、参考 URL に。

　　

参考URL：http://flab.ces.kyutech.ac.jp/~fujio/index.files …

No.4 回答者： 178-tall 回答日時： 2012/09/30 21:16

参考 URL の曲率算式を盗用し、直線 y = t*x と楕円との交点における楕円の曲率半径が r であるとして、楕円曲線式の b を勘定する算式でも出しますか…。

曲率算式 (これは、ひたすら盲信)
　r^2 = {(a^4yo^2 + b^4xo^2)^3}/(a^8b^8)　　　…(G.1)

直線 y = t*x と楕円との交点 (yo/t, yo)
　(yo/t)^2/a^2 + yo^2/b^2 = 1
　→ yo^2 = (tab)^2/(ta^2 + b^2)　　　…(*)
　　　　　　↓
　xo^2 = (yo/t)^2 = (ab)^2/(ta^2 + b^2)　　…(**)

(*), (**) を (G.1) へ代入した結果。
　r^2 = {(a^4t^2 + b^4)^3}/[(a^2b^2){(at)^2 + b^2}]

…委細は「眉唾」ながら、いずれこんな形。
多項式展開して零点探し、なのでしょう。

　　

No.5 回答者： alice_44 回答日時： 2012/10/01 09:49

計算の経路は異なっても、求まる曲率は同じはずだけれど…

No.2 と No.4 で b の方程式の次数が違うのは、なぜだろう？

まあ、いづれにせよ、b についての方程式を正しく立てて、
近似解法に持ち込むだけです。実務系の目的であれば、
数値計算の部分はパソコンに任せてしまえばいい。
特にプログラミングしなくても、既存のソフトが使えるでしょう。
精度によっては、エクセルのソルバーでも済みます。

それ以前の、数学部分の話として、(a cosθ, b sinθ) の
偏角が θ でないことは、ちゃんと理解しておいたほうがよいです。
No.2 に書いた (a cosφ, b sinφ) = (L cosθ, L sinθ)
という式の意味をよく考えてみてください。

No.6 回答者： 178-tall 回答日時： 2012/10/01 11:09

「b の値が判ればいい」というので、スプレッドシートにて「不動点収束」を試行。

その前に、当方の習慣病である算式訂正。
　r^2 = {(a^4t^2 + b^4)^3}/[(a^2b^2){(at)^2 + b^2}^3]　　　…(*** 改)

これを例えば、
　b^2 = {(a^4t^2 + b^4)^3}/[(ar)^2{(at)^2 + b^2}^3]
と変形。
(*** 改) の右辺は b^2 の有理式なので、R(b^2) とでもします。
ここで、b^2 = R(b^2) を満たす「不動点収束」へ。

b の初期値 bo を (*** 改) の両辺へ代入すると、よほどの強運に恵まれぬかぎり、
　bo^2 ≠ R(bo^2)
このあと、次の初期値に R(bo^2) を使うのが「不動点収束」のやり口。

スプレッドシートに R(bo^2) の算式をペタペタはりつければ、準備完了。
bo の初期値に長径 a と比べ小さめの値を入れると収束。
有効桁内で収束するのは、繰り返し 30 回ほどを要しますが、手間要らずなのが便利。

お試しを。

　　

No.7 回答者： 178-tall 回答日時： 2012/10/01 11:23

露骨なミス、お気づきとは思いますが、一応訂正。

その前に、当方の習慣病である算式訂正。
　r^2 = {(a^4t^2 + b^4)^3}/[(a^2b^2){(at)^2 + b^2}^3]

これを例えば、
　b^2 = {(a^4t^2 + b^4)^3}/[(ar)^2{(at)^2 + b^2}^3]　　　…(*** 改)
と変形。
(*** 改) の右辺は b^2 の有理式なので、R(b^2) とでもします。
ここで、b^2 = R(b^2) を満たす「不動点収束」へ。

No.8 回答者： 178-tall 回答日時： 2012/10/01 12:47

蛇足を少々。

>パラメータ表示 (x, y) = (a*cosθ, b*sinθ) なら、
>曲率半径 r = {(a^2 + b^2)cos^2(θ) - a^2}^(3/2)/(ab)　になるらしい。

こちらから r^2 を出してみると？

原題の y = t を勘案すると、
　b*sinθ = ta*cosθ
これから、
　sin^2(θ) = b^2/{(ta)^2 + b^2}
　cos^2(θ) = (ta)^2/{(ta)^2 + b^2}
これを上記引用式　r^2 = {(a^2 + b^2)cos^2(θ) - a^2}^3/(ab)^2　へ代入すると、
　r^2 = {(a^4t^2 + b^4)^3}/[(ab)^2{(at)^2 + b^2}^3]
となるらしい。

無断借用した二つの曲率半径算式は同値といえそう。

　　

No.9 回答者： 178-tall 回答日時： 2012/10/01 12:56

習慣病は不治の病也。

原題の y = tx を勘案すると、…

　　

No.10 回答者： 178-tall 回答日時： 2012/10/01 22:14

蛇足のついで。

(1 次) Newton の一試行例を。

前掲の式を
　P(B) = B*(ra)^2*{(at)^2 + B}^3 - (a^4t^2 + B^2)^3
と変形。(B = b^2)
P(B) の零点に近づけるための微係数は、
　P'(B) = (ra)^2*{(at)^2 + 4B}{(at)^2 + B}^2 - 6B(a^4t^2 + B^2)^2
を使用。

これも、スプレッドシートにて算式をペタペタはりつければ、準備完了。
「不動点収束」させたサンプルで、繰り返し 20 回ほどで有効桁内で収束。
思ったほど速くならないのは、曲率半径の変化が少ないからか。(a/b = 1.5 の例題)

もっと円に近づけば、ますますカッたるくなるのでしょう。