円板の固有振動数

下のURLで出している円板の固有振動数の式のKmnの出し方を教えてください
http://www.math.ryukoku.ac.jp/~iida/lecture/gr/g …

No.3 ベストアンサー 回答者： endlessriver 回答日時： 2014/06/05 19:46

掲題のサイトの議論は以下のような疑問があり、ここでの質問に意味がないように

思います。

１．境界条件
　詳しくないのですが、R(a)=0はありません。周囲はフリーです。
　中心でとめる場合も固定していないようなので境界条件は|R(r)|<∞

　しかない。したがって、R'(a)=0 だけになります（これも詳しくないが）。

２．まず、方程式はパラメータの正負によって次の２つの式がでます。
　R''+R'/x+(1-n^2/x^2)R=0 (1)
　R''+R'/x-(1+n^2/x^2)R=0 (2)
　R''+R'/x-(n^2/x^2)R=0 (3)

３．(1)の解は　　R=AJn+BYn
　　(2)の解は　　R=AIn+BKn
　　(3)の解は　　R=Ar^n+Br^(-n)

　掲題のサイトでは違うパラメータの方程式の解(1)(2)をごっちゃにしています。

４．r=0で有限という条件から
　　(1)の解は　　R=AJn
　　(2)の解は　　R=AIn
　　(3)の解は　　R=Ar^n

　となり、R'(a)=0を使うと(2)(3)はA=0　となります(In'>0,r^(n-1)>0)。結局、解は
　　R=AJn

　だけとなる。

５．以上の解でよければ求めるknmは（R'(a)=0を使うと）
Jn'(knm)=0

　です。これは、＃１のようにして求めることもでききますが、円形導波管の
　電磁界の教科書に載っています。

参考URL：http://tpweb2.phys.konan-u.ac.jp/~susa/keisan_bu …

この回答へのお礼

回答ありがとうございました

お礼日時：2014/06/05 22:15

No.2 回答者： endlessriver 回答日時： 2014/06/05 17:22

＃１です。

下記で求めたものはkaなので得られた値をaで割ったものが　kmnになります。

No.1 回答者： endlessriver 回答日時： 2014/06/05 17:16

数値計算するしかありません。

フリーソフトMaximaの結果を示します。
基のサイトはJm, kmn としていますが、以下ではJn, knmと考えています。

まず、g(n,x)=Jn'(x)/Jn(x)-In'(x)/In(x)
としてg(n,x)=0の解をニュートン法で求めます。

関数を変形したのは、（解が求めやすいと思って？）素直なグラフにするためです。
ベッセル関数の形は次のコマンドを実行すれば表示します。

　wxplot2d([bessel_j(1,x),bessel_j(2,x)], [x,0,50]);
　wxplot2d([bessel_i(1,x),bessel_i(2,x)], [x,0,5]);

１．まず、JnとInの微分は公式を使って、関数bj(n,x)とbi(n,x)を定義する。
　bj(n,x):=(bessel_j(n-1,x)-bessel_j(n+1,x))/2;
　bi(n,x):=(bessel_i(n+1,x)+bessel_i(n-1,x))/2;

　次に計算関数として
　g(n,x):=bj(n,x)/bessel_j(n,x)-bi(n,x)/bessel_i(n,x);

　を定義する。つぎのコマンドで、グラフを描いて大体の零点をもとめ、
　ニュートン法の開始値とする。

　[wxplot2d([g(0,x),g(1,x),g(2,x)],[x,0,20],[y,-10,10],
　　[legend, "g(0,x)", "g(1,x)","g(2,x)"])]$

２．ニュートン法のパッケージをロードする。
　load(mnewton);

３．変数に値をセット
　[n:0, x1:3, x2:7, x3:10]$

　ニュートン法の計算
　[mnewton(g(n,x),x,x1),mnewton(g(n,x),x,x2),mnewton(g(n,x),x,x3)];

　結果は次のように表示されます。
　(%o39) [[[x=3.196220616582541]],[[x=6.306437047688424]],[[x=9.439499137876405]]]

　すなわち、k01=3.196220616582541
k02=6.306437047688424
k03=9.439499137876405

　となります。

　同様に、次を実行すれば各値が求まります。

　[n:1, x1:4, x2:7.5, x3:11]$
　[mnewton(g(n,x),x,x1),mnewton(g(n,x),x,x2),mnewton(g(n,x),x,x3)];
　(%o41) [[[x=4.610899879049056]],[[x=7.799273800811232]],[[x=10.9580671919195]]]

　[n:2, x1:6, x2:9, x3:12.5]$
　[mnewton(g(n,x),x,x1),mnewton(g(n,x),x,x2),mnewton(g(n,x),x,x3)];
　(%o43) [[[x=5.905678235420522]],[[x=9.196882599635321]],[[x=12.40222096686439]]]