Excel VBA　数学実験室

≫ [Amazon書籍] 話してウケる!　不思議がわかる!（理系のネタ全書）

ベッセル関数の定義

　ベッセル関数 は三角関数のように周期的に振動する関数です。応用範囲はとても広く、古典力学の振動・波動に始まって、電磁気学における電磁波、さらには量子力学の波動方程式の解にまで登場します。エクセルでは BESSELJ() という形で組み込まれています。整数次ベッセル関数は一般的に次のような級数によって定義されます：


　n は全ての整数をとることができますが、Excel における BESSELJ() では負整数を入れるとエラーが出ます。しかし J _n と J _−n の間には、


という簡単な関係式があるので、グラフを描くことに支障はありません。ベッセル関数は n を実数にまで拡大定義することができますが、さすがにそこまでは Excel はサポートしていません。また、Excel のヘルプなどでは上式分母の (n + k)! のところが　Γ(n + k + 1) となっていたりしますが、n, k ともに整数である場合はあまり意味がなく、階乗記号で表す方が表記が簡単になります。
　このブログでは当面の間、n を整数の範囲にとどめておきます。時間があるときに実数で定義されたベッセル関数を自分で作りますので、それまで少しお待ちください。

　ベッセル関数は n を変化させると J ₀(x), J ₁(x), J ₂(x), …… というように異なる関数となります。さっそくエクセルで J ₀(x) と J ₁(x) のグラフを描いてみましょう。

　

　三角関数と似ていますが、x の増加と共に振幅は減少していきます。また、n が大きくなると関数の立ち上がりは非常に遅くなります。試しに n = 9, 10 の場合を見てみましょう。

　

　x = 5 付近でようやく立ち上がるという、のんびりとした関数です。最初の波の振幅も小さくなっています。

　このブログは何か変わった関数を探すという趣旨なので、ベッセル関数を紹介して終わりというわけにはいきません。 J ₁(x) において変数 x を xsinx としてみると・・・・・・

　

　特徴のあるグラフが出現しましたね。面白い法則が隠されていることにお気づきですか？　なだらかな山と山の間に、鋭いピークを持つ山がありますが、その鋭い山のほうは 0, 1, 2, 3, ・・・・・・というように数が増えていっています。もっと広範囲の x でも確認してみましたが、この法則に間違いはないようです。うーむ。不思議ですねえ。

ベッセル関数を足し合わせます

　sinx + sin2x のように三角関数同士を加えても、加法定理を使って sin と cos で表される三角関数になりますね。同じようにベッセル関数同士を足しても、やはり周期関数になります。試みに J ₀(x) と J ₉(x) を足し合わせてみましょう。

　

　上図に J ₀(x) と J ₉(x) がそれぞれ描かれていて、下図が足し合わせた結果です。グラフを良く見てイメージを掴んでください。

　最初にお話したように、番号 n が大きくなるほど関数の立ち上がりは遅くなります。 J ₉(x) は x = 5 付近でようやく立ち上がりますから、 J ₀ + J ₉ は x < 5 の範囲でほぼ J ₀ と同じように振る舞います。 x = 10 付近では山と谷が打ち消し合って、 0 に近い値をとりますね。ベッセル関数は、同符号であれば n が異なっていても、 x が大きくなると最終的には全てのグラフは重なってしまいます。 J ₀(x) と J ₉(x) も 25 < x あたりからはほぼ同じ関数として振る舞うので、両者を足し合わせると山と山、谷と谷が強め合って、およそ 2 倍の振幅で振動します。そのあたりのこともグラフを見て確認しておいてください。

ベッセル関数を掛け合わせます

　三角関数同士を掛けるとやはり周期関数になるように、ベッセル関数の積も周期関数を生み出します。 J ₀(x) と J ₉(x) を掛け合わせた場合のグラフを見てみましょう。

　

　1 以下の振幅同士を掛けているので、振幅はとても小さくなります（中図）。そのままでは見にくいので、縦軸スケールを拡大した図を一番下に載せています。 J ₉(x) が立ち上がるのは x = 5 あたりですから、x < 5 の領域では振幅はほぼ 0 です。 x = 10 の付近でそれぞれの関数は山と谷になりますから、掛けると負号となって、大きな谷を作ります。 x が大きな範囲では両関数のグラフはほとんど重なり合うので、谷のところは山となり、ほぼ全域で正値をとる関数となります（おおよそ y > 0 の領域で振動します）。
　選択する番号 n によって形は様々ですが、「ベッセル関数同士の掛け算は小さな波を作り出す」ということをイメージに定着させてください。

ベッセル概周期関数

　「三角関数の重ね合わせ sinax + cosbx において、 a と b の比が有理数であれば一定周期が保障されますが、無理数であると周期が乱れて概周期関数（おおよその周期関数）になる」ということを以前にお話しました。
　同じことがベッセル関数でも成り立ちます。第 42 回では J ₀(x) + J ₉(x) を考えましたが、変数はともに x であったために、足し合わせても単純な周期関数が現れるだけでした。しかし、たとえば次のように J ₁(x) + J ₃(√2x) のような関数を作ってみると・・・・・・

　

　こういう微妙な波形が生じます。「周期があるようなー、ないようなー」というグラフですね。本来のベッセル関数の振幅は 1 / √x に比例して小さくなってゆくのですが、そういう規則性も崩れてしまっています。振幅が小さくなったと思ったら、次はまた大きくなったりしています。係数の無理数が大きくなると関数はどんどんぎざぎざな形になります。 J ₁(x) + J ₃(√119x) としてみると・・・・・・

　

　なんだか得体の知れないグラフになりましたね。

ベッセル関数と三角関数

　ベッセル関数と三角関数は、よく似た周期関数です。級数展開の形も似ていますね。しかし似ているといっても、別種の関数であるので、足したり掛けたりしてみると、思わぬ波形が生じたりします。たとえば J ₁(x) + cos(2x) という関数を作ってみると・・・・・・

　

　振幅の大きさが交互に変わる周期関数です。面白い規則性があります。原点の山だけを除外して、他の波を図のように a 群と b 群に分けてみます。 a 群は少しずつ振幅を減少させ、逆に b 群は僅かずつ増加していますね。やがて a 群と b 群は同じ振幅となって混じりあい、 x が大きいところでは一定の波となります。その理由は J ₁(x) は振幅を減少させていく関数であるために、遠方ではその波が消えてしまって cos(2x) だけが残るからです。

　次は J ₁(x) と cos(2x) の積を作ってみます：

　

　異様な形のグラフが描かれていますね。鋭いピークが 2 つセットになっていて、正負の値を交互に繰り返します。この理由を解析するのは・・・・・・私にはちょっと難しすぎます。誰かわかる人がいたら、コメントで教えてください。

≫ Python ベッセル関数
　
　ベッセル関数 J_k(x) を k = 0 から n まで足し合わせた関数を作ってみました。


　B_n(x) を計算する VBA ユーザー定義関数（BESSELSUM関数）のコードは記事の一番下に載せておくので、自分のエクセルで試してみたい人はコピーしてください。n = 1 から 3 までのグラフは次のようになります。

　

　このように、n = 3 までは、あまり変化しません。

　

　しかし n = 4 以降は不規則な変化をみせるようになります。
　ところが n がある程度大きくなると次のような形に落ち着きます。

　　

　上の図は n = 28 としていますが、以降はほとんど形は変わりません（ k の大きなベッセル関数は微かな振動しかしないため、寄与がほとんどなくなります）。たぶんこれが最終形態と思われます（ n = 100 まで確認）。


　ワークシートで BESSELSUM関数を使用するときは


としてください。
　
　次は二変数ベッセル関数です。まず最初に距離 r を変数とする z = J₀(r) のグラフを見てみます：

　

　原点付近でひと息に振幅が落ちてしまいます。ベッセル関数の振幅は近似的に √x に逆比例するので、原点から遠ざかるほど振幅の減少率は緩やかになります。次は変数を x + y にしてみます。

　

　範囲は x, y > 0 としてあります。ベッセル関数の周期はおよそ 2*pi（ pi：円周率）ですから、


を満たす点で振幅が極大となります。言い換えると上式で表される直線上で波高が揃うのです。次は √x + y という変数をもつベッセル関数です：

　

　波が極大となる点を結ぶ線が曲がっていますね。今度は曲線


の上で振幅が極大となるので、このようなグラフとなるわけです。最後はベッセル関数と三角関数の重ね合わせて


という関数を作ってグラフを描いてみます：

　

　とても複雑な関数です。波高が揃っていないうえに、位相が変わるたびに波の形そのものを変えてゆきます。周期関数同士を重ね合わせると、思いもよらない形のグラフを生み出します。
　
　tanx は π/2 + nπ（n = ±1, ±2, ... ）で発散する関数ですが、この変数に値域の大きさが |π/2| を超えないような関数を入れると全域で連続な関数となります。今回はその１例としてベッセル関数 J_n(x) を組み込んでみます。 |J_n(x)| ≦ 1です。

　y = tan[J₀(x)] のグラフを描いてみます。

　

　原点で最大値となり、|x| の増加とともに減衰する振動関数です。
　tanx, J₀(x) ともに偶関数なので、その合成関数もまた偶関数となります。

　値域の大きさが|π/2| を超えない範囲でベッセル関数の線形結合をとれば、連続性が失われません。例として


という関数のグラフを描いてみます。

　

　係数によって原点にピークをもつ J₀(x) の寄与を抑えて、相対的に立ち上がりの遅い J₁(x) や J₂(x) の寄与を増やしているので、最大値の位置は原点からずれます。

　J₁(x) と cosx/2 を組み込んで、


という関数のグラフを描いてみます。

　

　cosx/2 の効果があるので、振幅が 0 まで落ちることはなく、原点から離れたところでも一定の振動があります。
　
　ベッセル関数を極方程式に組み込んだ


という方程式を扱ってみましょう。ベッセル関数は三角関数によく似た周期関数です。しかし x が増加するとと（1 / √x に比例して）振幅を減衰させていき、x → ∞ では 0 となります。また n によって立ち上がりの速さと振幅が異なります。この関数が極方程式に組み込まれると、とても美しい軌跡を描きます。さっそく r = J₀(θ) (0 ≦ θ)のグラフを見てみましょう。

　

　見事ですね。点 (1, 0) からスタートして、原点に巻きつくようにして円を描き、少しずつその半径を小さくしていきます。「周回ごとに必ず原点を通過する」というのがこのグラフの特徴です。螺旋とは異なる形です。 J₃(x) は原点から始まる関数ですから、 r = J₃(θ) のグラフは下のようになります：

　

　正葉曲線と重ね合わせてみます：


という関数を作ってグラフを描いてみると・・・・・・

　

　r = J₀(θ) と r = sin(3θ) 双方の特徴が合わさったようなグラフになっていますね！　周回ごとに描かれる正葉曲線が小さくなっていきます。次にベッセル関数の中に三角関数を入れてみます。つまり


というグラフを描いてみますと ...

　

　ぐにゃぐにゃっとした星形未満のような軌跡ですね。ついでにもう１つ。


としてみると・・・・・

　

　突起部分が増えました。sin(nθ) の n を大きくしていくとギザギザがどんどん増えます。
　

新しい関数に出会ったときは

　足してみたり、掛けてみたり、変数を変えてみたり・・・・・・こういうことを繰り返すことで、高度な関数も次第に明確なイメージとなって頭に定着するようになります。理工系の学生さんはガンマ関数、デルタ関数、球面調和関数、ルジャンドル関数、・・・・・・というように次から次へと新しい特殊関数の使い方を覚えていかなくてはなりません。目の回るような忙しさかもしれませんが、「急がば回れ」と少し気持ちにゆとりをもって、関数のイメージを掴むようにつとめてください。これはそういうことを怠っていた私自身の反省でもあります。大学院へ進むと基礎的な勉強のための時間を確保するのは本当に難しいです。私と同じ轍を踏まないように、何事も若いうちにやっておいたほうがいいです。

ラグビーにはまってます

　ラグビーワールドカップ面白いですね！　日本代表選手たちの活躍には本当に感動しました。 3 勝もしたのにベスト 8 に進めなかったのは本当に残念です。準決勝 2 試合も見ました。昨晩はオーストラリアとアルゼンチンの試合でした。序盤から選手たちが出血も厭わずにぶつかり合う大熱戦で、その迫力に圧倒されました。正直、私はまだルールすらきちんと把握していない「にわかファン」にすぎないのですが、「世の中にはこんな面白いスポーツがあったのか！」と知って興奮しました。「大男たちが戦列を作り、ぶつかりあって陣地を確保してゆく」という光景は、大昔の戦争の名残をかなりとどめていますね。現代においては、他にこういうスポーツはあまり残っていない気がします。
　英語ブログもやっているので、英国の歴史資料などを読むことが多いのですが、昔のゲルマン戦士たちの勇ましく戦う光景が再現されているようで、ついつい観ているほうも血が騒ぎますね。決勝戦も楽しみです！