ガウシアンビーム光学の基礎 — ビームウエスト・広がり角・ABCD則

レーザーポインタを壁に当てると、小さな輝点が見えます。では、その壁を100m先に、1km先に、さらに月まで遠ざけたら、輝点はどこまで大きくなるでしょうか。

「レーザー光は平行光線だから広がらない」と習った方も多いと思いますが、これは正しくありません。波長 $\lambda = 1064\,\mathrm{nm}$、出射時の半径 $1\,\mathrm{mm}$ という、ごく標準的なレーザーを考えてみます。この記事で導く式を使うと、1km先でのビーム半径は約 $34\,\mathrm{cm}$、月面（38万km先）では約 $129\,\mathrm{km}$ になります。光は絶対に広がるのです。しかも厄介なことに、出口を細く絞るほど、遠方での広がりは大きくなります。この「絞ると広がる」というトレードオフの正体が、これから扱うガウシアンビーム光学です。

この知識が効いてくる場面を、具体的に2つ挙げておきます。

ひとつは宇宙光通信のリンク設計です。静止軌道の衛星から地上局へレーザーで通信するとき、送信望遠鏡の口径をいくつにすれば受信口径にどれだけの電力が届くのか。この計算は、まさにビーム半径 $w(z)$ の式そのものです。口径を大きくすると広がり角が小さくなり受信電力が増えますが、その分だけ指向精度の要求が厳しくなる — 設計トレードオフが式から直接読めます。

もうひとつはレーザー加工や光ピンセットの集光設計です。レンズでどこまで細く絞れるのか（スポット径）、絞った状態が前後どれくらいの距離で保たれるのか（焦点深度）。これらは幾何光学では「焦点は点」で終わってしまい、まったく答えが出ません。ガウシアンビームの枠組みで初めて、$2w_f \approx 4\lambda f/(\pi D)$ という設計式が手に入ります。

さらに、レーザー共振器の安定性解析、干渉計の光路設計、光ファイバとの結合効率計算 — 光を扱う工学のほぼ全域で、ガウシアンビームは共通言語になっています。

本記事の内容

• ヘルムホルツ方程式の近軸近似から、基本ガウシアンビーム解を省略なしで導出する
• ビーム半径 $w(z)$、曲率半径 $R(z)$、レイリー長 $z_R$、グーイ位相 $\zeta(z)$、$M^2$ 因子の意味を整理する
• 複素ビームパラメータ $q$ と ABCD 則を使い、レンズ集光のスポット径と焦点シフトを導く
• Pythonで $q$ 伝搬を実装し、ビームの広がり・強度分布・集光系の挙動を可視化する

前提知識

この記事を読む前に、以下の記事を読んでおくと理解が深まります。

• 波動光学入門（干渉と回折） — 光を波として扱う枠組み。回折の考え方が本記事の出発点です
• レンズの公式と結像 — 幾何光学の焦点。本記事ではこれが「近似でしかない」ことを見ます
• 光ファイバ通信の基礎 — ファイバの基本モードはほぼガウシアンで、本記事の結果が直接使えます

ガウシアンビームとは — 平面波でも球面波でもない「第三の解」

これから扱う対象の全体像です。ビームはある1点（ビームウエスト）で最も細くなり、そこから前後へ双曲線を描いて太くなっていきます。ウエスト付近では波面がほぼ平面なのに、離れるにつれて曲がり、遠方では中心がウエストにある球面波に近づく — この「場所によって性格が変わる」ところがガウシアンビームの核心です。横軸の単位に使ったレイリー長 $z_R$ と、遠方の広がりを決める半頂角 $\theta$ が、以降の主役になります。

波動光学の教科書に最初に出てくる解は、平面波と球面波です。平面波は無限に広い波面がまっすぐ進むもの、球面波は1点から放射状に広がるもの。どちらも数学的には美しいのですが、実際のレーザービームはそのどちらでもありません。

平面波は「横方向に無限に広がっている」ので、レーザーのような太さ数ミリの光束を表せません。球面波は「点から出て一様に全方向へ広がる」ので、前方に細く絞られたビームを表せません。実際のレーザーは、横方向に有限の太さを持ちながら、ほぼ真っ直ぐ進む — この中間的なふるまいをします。

3つの解の波面（等位相面）を並べたものです。平面波の波面は等間隔の平面がどこまでも続き、横方向の「太さ」という概念そのものを持ちません。球面波の波面は同心円で、点光源からあらゆる方向へ均等に広がってしまいます。右のガウシアンビームだけが、ウエスト付近では平面波のような平らな波面を、遠方では球面波のような曲がった波面を持ち、しかも横方向は有限の太さに収まっています。

イメージとしては、細長い砂時計を横に寝かせた形を思い浮かべてください。ビームはどこかで最も細くなり（くびれ）、そこから前後に向かってなだらかに太くなっていく。くびれのあたりでは波面はほぼ平面で、遠く離れると波面は球面に近づく。つまりガウシアンビームは、近くでは平面波のように、遠くでは球面波のようにふるまう解なのです。

そしてその横方向の強度分布は、名前のとおりガウス関数（正規分布と同じ釣鐘型）になります。中心が最も明るく、外側へ向かって滑らかに減衰する。この「ガウス型の横分布を持ち、砂時計型に伝搬する解」がガウシアンビームです。

なぜガウス型なのか、という問いには明確な答えがあります。後で見るように、これは近軸波動方程式の基本解（最低次モード）であり、レーザー共振器のなかで最も損失が小さいモードだからです。だからこそ、多くのレーザーは自然にこの形の光を出力します。「たまたまガウス型だった」のではなく、「方程式がそう要求している」わけです。

ここまでは絵の話でした。次はこの直感を、マクスウェル方程式から出発して数式として導きます。

出発点 — ヘルムホルツ方程式と近軸近似

真空中（または一様な媒質中）の電磁波は、マクスウェル方程式から波動方程式に帰着します。偏光の向きが場所によって変わらない状況では、電場をスカラー量 $E$ で近似して扱えます。

$$ \nabla^2 E – \frac{1}{c^2}\frac{\partial^2 E}{\partial t^2} = 0 $$

レーザー光は単一波長（角周波数 $\omega$）とみなせるので、時間依存性を $e^{j\omega t}$ と置いて分離します。

$$ E(\bm{r}, t) = \mathrm{Re}\left[ U(\bm{r})\, e^{j\omega t} \right] $$

これを波動方程式に代入すると、時間微分が $-\omega^2$ に置き換わり、空間だけの方程式が残ります。

$$ \begin{equation} \nabla^2 U + k^2 U = 0, \qquad k = \frac{\omega}{c} = \frac{2\pi}{\lambda} \end{equation} $$

これがヘルムホルツ方程式です。$k$ は波数で、1メートル進むあいだに位相が何ラジアン回るかを表します。

「ほぼ真っ直ぐ進む」を数式にする

ヘルムホルツ方程式は一般解が複雑ですが、レーザービームには強い制約があります。それは「$z$ 軸方向にほとんど真っ直ぐ進む」ということです。この制約を式に持ち込みましょう。

$z$ 方向に進む平面波は $e^{-jkz}$ です（時間依存を $e^{j\omega t}$ と取ったので、前進波は $e^{-jkz}$ になります）。レーザービームは「平面波に、ゆっくり変化する形の補正がかかったもの」と考えられます。そこで解を次の形に置きます。

$$ \begin{equation} U(x, y, z) = u(x, y, z)\, e^{-jkz} \end{equation} $$

ここで $u$ は包絡線（envelope）で、ビームの太さや波面の曲がりといった「ゆっくりした情報」を担います。速い位相回転はすべて $e^{-jkz}$ に押し込めた、という分担です。

これをヘルムホルツ方程式に代入します。$z$ による2階微分を丁寧に展開すると、

$$ \frac{\partial^2 U}{\partial z^2} = \left( \frac{\partial^2 u}{\partial z^2} – 2jk\frac{\partial u}{\partial z} – k^2 u \right) e^{-jkz} $$

横方向のラプラシアンを $\nabla_T^2 = \partial^2/\partial x^2 + \partial^2/\partial y^2$ と書くと、ヘルムホルツ方程式は

$$ \left( \nabla_T^2 u + \frac{\partial^2 u}{\partial z^2} – 2jk\frac{\partial u}{\partial z} – k^2 u + k^2 u \right) e^{-jkz} = 0 $$

$-k^2 u$ と $+k^2 u$ が相殺します。これは平面波の位相を先に取り出しておいた効果で、この相殺があるからこそ話が簡単になります。残りは

$$ \nabla_T^2 u + \frac{\partial^2 u}{\partial z^2} – 2jk\frac{\partial u}{\partial z} = 0 $$

緩慢変化包絡線近似（SVEA）

ここで最後の一手です。包絡線 $u$ は「$z$ 方向にゆっくり変化する」と仮定しました。ゆっくりとは、波長スケールで見て、1波長進むあいだの $u$ の変化がわずかであるという意味です。数式にすると

$$ \left| \frac{\partial^2 u}{\partial z^2} \right| \ll \left| 2k \frac{\partial u}{\partial z} \right| $$

となります。左辺は変化の変化（2階）、右辺は $k$ という大きな数がかかった1階の項。$k$ は可視光で $10^7\,\mathrm{m^{-1}}$ 程度という巨大な数なので、この不等式は「ビームの太さが波長よりずっと大きい限り」よく成り立ちます。この近似を緩慢変化包絡線近似（Slowly Varying Envelope Approximation, SVEA）と呼びます。

$\partial^2 u/\partial z^2$ を落とすと、目的の式が得られます。

$$ \begin{equation} \nabla_T^2 u – 2jk \frac{\partial u}{\partial z} = 0 \end{equation} $$

これが近軸ヘルムホルツ方程式（パラキシャル波動方程式）です。形をよく見てください。$z$ を時間、$\nabla_T^2$ を空間2階微分と読み替えると、これは量子力学のシュレディンガー方程式（自由粒子）とまったく同じ構造をしています。「光の横方向の広がりが $z$ とともに進行する」ことと、「波束が時間とともに広がる」ことが、数学的に同じ現象だというのは面白い対応です。

左図は $U = u\,e^{-jkz}$ という分解の意味を絵にしたものです。薄い青の細かい振動が $e^{-jkz}$ による速い位相回転、赤い曲線がそれを包む包絡線 $u$ で、$u$ はレイリー長のスケールでしかゆっくり変化しません（実際の波長は $z_R$ の 300 万分の 1 程度で、図よりはるかに細かい振動になります）。右図は、SVEA で捨てる項 $|\partial^2 u/\partial z^2|$ と残す項 $|2k\,\partial u/\partial z|$ を対数目盛で比べたもので、$\lambda = 1064\,\mathrm{nm}$、$w_0 = 1\,\mathrm{mm}$ では両者の比は最大でも $5.7\times10^{-8}$ です。捨てた項が8桁も小さいことが、この近似の安全性を裏づけています。

方程式が手に入りました。次はこれを実際に解き、砂時計型の解を取り出します。

基本ガウシアン解の導出

近軸ヘルムホルツ方程式を解くにあたり、先ほどの直感を思い出しましょう。求めたい解は「横方向にガウス型の振幅」と「なだらかに曲がった波面」を持つはずです。ガウス型の振幅は $\exp(-r^2/w^2)$、球面波に近い波面は $\exp(-jkr^2/2R)$ と書けます（球面波 $e^{-jkR_{\text{sph}}}$ を軸近傍で $R_{\text{sph}} \approx z + r^2/2z$ と展開したものです）。この2つは、ひとまとめにできます。

$$ \exp\left( -\frac{r^2}{w^2} \right) \exp\left( -j\frac{kr^2}{2R} \right) = \exp\left[ -j\frac{k r^2}{2}\left( \frac{1}{R} – j\frac{\lambda}{\pi w^2} \right) \right] $$

右辺の括弧の中身を確認しておきます。第2項は $-j\frac{kr^2}{2}\cdot\left(-j\frac{\lambda}{\pi w^2}\right) = -\frac{k\lambda r^2}{2\pi w^2}$ で、$k = 2\pi/\lambda$ を代入すると $k\lambda/(2\pi) = 1$ ですから、ちょうど $-r^2/w^2$ になります。振幅の減衰と波面の曲がりが、ひとつの複素数にまとめられたわけです。

そこで、この複素数を $1/q$ と名づけます。

$$ \begin{equation} \frac{1}{q(z)} = \frac{1}{R(z)} – j\frac{\lambda}{\pi w^2(z)} \end{equation} $$

これが複素ビームパラメータです。実部がビームの曲がり具合、虚部が太さを担っている、と覚えてください。

試行解を代入する

以上を踏まえ、次の形の解を仮定します。

$$ u(r, z) = A \exp\left[ -j\left( P(z) + \frac{k r^2}{2 q(z)} \right) \right] $$

$P(z)$ は $r$ に依存しない項で、全体の振幅と位相を担当します。$q(z)$ が先ほどの複素ビームパラメータです。$P$ と $q$ という2つの未知関数を、方程式に代入して決めていきます。

まず横方向のラプラシアンを計算します。円対称なので $\nabla_T^2 = \partial^2/\partial r^2 + (1/r)\partial/\partial r$ です。1階微分は

$$ \frac{\partial u}{\partial r} = -j\frac{kr}{q}\, u $$

これをもう一度 $r$ で微分すると、積の微分から2項出ます。

$$ \frac{\partial^2 u}{\partial r^2} = -j\frac{k}{q}\,u + \left( -j\frac{kr}{q} \right)^2 u = -j\frac{k}{q}u – \frac{k^2r^2}{q^2}u $$

$(1/r)\partial u/\partial r = -j(k/q)u$ を足し合わせると、

$$ \nabla_T^2 u = \left( -\frac{2jk}{q} – \frac{k^2 r^2}{q^2} \right) u $$

次に $z$ 微分です。$1/q$ の微分が $-q’/q^2$ になることに注意すると、

$$ \frac{\partial u}{\partial z} = -j\left( P’ – \frac{k r^2 q’}{2q^2} \right) u $$

この2つを近軸ヘルムホルツ方程式 $\nabla_T^2 u – 2jk\,\partial u/\partial z = 0$ に入れます。第2項は $(-2jk)\times(-j) = -2k$ となることに注意して、

$$ \left( -\frac{2jk}{q} – \frac{k^2r^2}{q^2} \right) u – 2k\left( P’ – \frac{kr^2 q’}{2q^2} \right) u = 0 $$

$u \neq 0$ なので括弧の中身がゼロです。ここで$r$ のべきごとに項を分けるのがポイントです。この式が任意の $r$ で成り立つには、$r^2$ の係数と $r^0$ の係数がそれぞれ独立にゼロでなければなりません。

$r^2$ の係数から:

$$ -\frac{k^2}{q^2} + \frac{k^2 q’}{q^2} = 0 \quad \Longrightarrow \quad \boxed{q'(z) = 1} $$

$r^0$ の係数から:

$$ -\frac{2jk}{q} – 2kP’ = 0 \quad \Longrightarrow \quad \boxed{P'(z) = -\frac{j}{q}} $$

驚くほど簡単な2本の常微分方程式に落ちました。とくに $q’ = 1$ は絶品で、複素ビームパラメータは自由空間では単に距離ぶんだけ足されていくことを意味します。

$$ q(z) = q_0 + z $$

初期条件を決める — ビームウエスト

積分定数 $q_0$ を決めましょう。ビームが最も細くなる位置（ビームウエスト）を $z = 0$ に取ります。そこでは波面が平面なので $R(0) = \infty$、太さは $w(0) = w_0$ です。定義式に入れると

$$ \frac{1}{q_0} = \frac{1}{\infty} – j\frac{\lambda}{\pi w_0^2} = -j\frac{\lambda}{\pi w_0^2} $$

したがって $q_0 = j\pi w_0^2/\lambda$ です。この虚部に現れる量には特別な名前がついています。

$$ \begin{equation} z_R \equiv \frac{\pi w_0^2}{\lambda} \end{equation} $$

これをレイリー長（Rayleigh range）と呼びます。すると $q_0 = j z_R$ となり、

$$ \begin{equation} q(z) = z + j z_R \end{equation} $$

という、拍子抜けするほど簡潔な形が得られました。ガウシアンビームの伝搬は、複素平面上で虚軸から出発して実軸方向へ平行移動するだけなのです。

左が $q(z) = z + jz_R$ の軌跡で、虚部が $z_R$ に固定されたまま実部だけが増える水平な直線です。ウエストでは純虚数 $jz_R$ から出発し、あとは進んだ距離ぶんだけ右へずれるだけ — これが $q’ = 1$ の図的な意味です。右はその逆数 $1/q$ の軌跡で、メビウス変換の性質から円になります。実部が $1/R$（波面の曲がり）、虚部が $-\lambda/(\pi w^2)$（太さ）を表すので、虚部が最も深い $z=0$ が最も細いウエスト、実部がゼロになるのもそこ（波面が平面）だと読み取れます。

$w(z)$ と $R(z)$ を取り出す

$q(z) = z + jz_R$ から、実際のビーム半径と曲率半径を読み出します。逆数を有理化します。

$$ \frac{1}{q(z)} = \frac{1}{z + jz_R} = \frac{z – jz_R}{z^2 + z_R^2} = \frac{z}{z^2+z_R^2} – j\frac{z_R}{z^2+z_R^2} $$

これを定義式 $1/q = 1/R – j\lambda/(\pi w^2)$ と実部・虚部で突き合わせます。実部から

$$ \frac{1}{R(z)} = \frac{z}{z^2+z_R^2} \quad \Longrightarrow \quad \boxed{R(z) = z\left[ 1 + \left(\frac{z_R}{z}\right)^2 \right]} $$

虚部から

$$ \frac{\lambda}{\pi w^2(z)} = \frac{z_R}{z^2+z_R^2} \quad \Longrightarrow \quad w^2(z) = \frac{\lambda (z^2+z_R^2)}{\pi z_R} $$

ここで $z_R = \pi w_0^2/\lambda$、すなわち $\lambda z_R/\pi = w_0^2$ という関係を使って整理すると、

$$ w^2(z) = \frac{\lambda z_R}{\pi}\left( 1 + \frac{z^2}{z_R^2} \right) = w_0^2\left[ 1 + \left( \frac{z}{z_R} \right)^2 \right] $$

$$ \begin{equation} \boxed{w(z) = w_0 \sqrt{1 + \left( \frac{z}{z_R} \right)^2}} \end{equation} $$

これがガウシアンビームの中心公式です。$w^2$ が $z^2$ の1次式になっているので、$w(z)$ のグラフは双曲線になります。冒頭で述べた「砂時計」の正体がこれです。

$P(z)$ を積分する — グーイ位相の出現

残る $P(z)$ を求めます。$P’ = -j/q = -j/(z+jz_R)$ を積分します。

$$ P(z) = -j \ln(z + jz_R) + \text{const} $$

ウエストで位相の基準を取る（$P(0) = 0$）と、

$$ P(z) = -j \ln\!\left( \frac{z + jz_R}{jz_R} \right) = -j \ln\!\left( 1 – j\frac{z}{z_R} \right) $$

対数の中身を極形式に直すのがコツです。$1 – jz/z_R$ の絶対値は $\sqrt{1+(z/z_R)^2} = w(z)/w_0$、偏角は $-\arctan(z/z_R)$ です。ここで

$$ \begin{equation} \zeta(z) \equiv \arctan\left( \frac{z}{z_R} \right) \end{equation} $$

と定義すると、$\ln(1-jz/z_R) = \ln(w/w_0) – j\zeta$ となります。代入すると

$$ P(z) = -j\left[ \ln\frac{w}{w_0} – j\zeta \right] = -j\ln\frac{w}{w_0} – \zeta $$

求めたいのは $\exp(-jP)$ なので、

$$ e^{-jP} = \exp\left( j \cdot j \ln\frac{w}{w_0} \right) \exp(j\zeta) = \frac{w_0}{w(z)} e^{j\zeta(z)} $$

$P$ の虚部が振幅の $w_0/w$ を、実部が余分な位相 $\zeta$ を生みました。振幅が $w_0/w$ に比例するのは、ビームが太くなるほど電場が薄まる — つまり全パワーが保存されることの現れです（後で確認します）。

完成した解

すべてを合わせると、基本ガウシアンビームの電場が得られます。

$$ \begin{equation} E(r,z) = E_0 \frac{w_0}{w(z)} \exp\left( -\frac{r^2}{w^2(z)} \right) \exp\left[ -j\left( kz + \frac{kr^2}{2R(z)} – \zeta(z) \right) \right] \end{equation} $$

4つの因子それぞれに意味があります。$w_0/w(z)$ は振幅の希釈、$\exp(-r^2/w^2)$ は横方向のガウス分布、$kz + kr^2/2R$ は伝搬と波面の曲がりによる位相、そして $-\zeta(z)$ がグーイ位相と呼ばれる余分な位相です。

導出は完了しました。次に、この式に登場する各パラメータが物理的に何を意味するのかを、一つずつ噛み砕いていきます。

各パラメータの意味を読み解く

数式が出そろったので、それぞれの記号が現場で何を意味するのかを整理します。

ビーム半径 $w(z)$ と $1/e^2$ の約束

$w(z)$ は「振幅が中心の $1/e$ になる半径」です。強度は振幅の2乗なので、$I \propto |E|^2 \propto \exp(-2r^2/w^2)$ となり、強度が中心の $1/e^2 \approx 13.5\%$ になる半径が $w$ です。レーザーのカタログで「ビーム径」と書かれていたら、ほぼ確実にこの $1/e^2$ 全幅 $2w$ を指しています。半径なのか直径なのか、$1/e$ なのか $1/e^2$ なのかで係数が2倍・$\sqrt{2}$ 倍と変わるので、実務ではここが事故の温床です。

半径 $a$ の丸アパーチャを通過するパワーの割合は、強度を積分して

$$ \frac{P(a)}{P_{\text{total}}} = 1 – \exp\left( -\frac{2a^2}{w^2} \right) $$

となります。$a = w$ で $86.5\%$、$a = 1.5w$ で $98.9\%$、$a = 2w$ で $99.97\%$。「レンズの口径はビーム半径の1.5倍以上取れ」という現場の経験則は、この数字から来ています。

左図で、振幅（緑）と強度（青）の裾の落ち方が違うことに注目してください。同じ $r = w$ の位置で、振幅は $1/e = 0.368$、強度は $1/e^2 = 0.135$ です。カタログの「ビーム径」がどちらの基準かを取り違えると $\sqrt{2}$ 倍の誤差が入る理由がここにあります。右図の透過率カーブは $a/w = 1.5$ を超えたあたりで急速に $100\%$ に飽和しており、口径をそれ以上大きくしても得られるパワーはほとんど増えません。逆に $a = w$ で切ると $13.5\%$ を捨てることになります。

レイリー長 $z_R$ — 焦点深度の物差し

$z = z_R$ を $w(z)$ の式に入れると $w(z_R) = \sqrt{2}\,w_0$、つまりビーム半径がウエストの $\sqrt{2}$ 倍（断面積で2倍、ピーク強度で半分）になる距離がレイリー長です。

$z \ll z_R$ ではビームはほとんど太くならず、実質的に平行光線とみなせます。$z \gg z_R$ では直線的に広がります。つまり $z_R$ は「この距離までなら光は絞られたまま」という焦点深度の物差しです。$\pm z_R$ の範囲（全長 $2z_R$）を共焦点パラメータ（confocal parameter）と呼び、レーザー加工では「焦点が効いている厚み」として使われます。

そして $z_R = \pi w_0^2/\lambda$ という式が示すのは、細く絞るほど焦点深度は2乗で短くなるという残酷な事実です。ウエスト半径を半分にすれば、焦点が保たれる距離は4分の1になります。

左図は同じ波長でウエスト半径だけを変えた3本の $w(z)$ です。$w_0 = 0.25\,\mathrm{mm}$ の赤い曲線は出発点こそ最も細いのに、$1\,\mathrm{m}$ も進まないうちに $w_0 = 1\,\mathrm{mm}$ の青を追い越してしまいます。丸印（$z = z_R$）の位置が右へ大きくずれていく様子が、焦点深度の差そのものです。右図はその関係を直接プロットしたもので、$w_0$ に対して放物線（2乗）になっています。半径を $1\,\mathrm{mm}$ から $0.5\,\mathrm{mm}$ に絞ると $z_R$ は $2.95\,\mathrm{m}$ から $73.8\,\mathrm{cm}$ へ、ちょうど4分の1に落ちます。

広がり角 $\theta$ — 遠方界の直線

$z \gg z_R$ での挙動を見ます。$w(z) = w_0\sqrt{1+(z/z_R)^2}$ の平方根の中で1が無視できるので、

$$ w(z) \approx w_0 \frac{z}{z_R} = w_0 \cdot \frac{z\lambda}{\pi w_0^2} = \frac{\lambda z}{\pi w_0} $$

$z$ に比例する — つまりビームは円錐状に広がります。その半頂角が広がり角です。

$$ \begin{equation} \boxed{\theta = \lim_{z\to\infty}\frac{w(z)}{z} = \frac{\lambda}{\pi w_0}} \end{equation} $$

分母に $w_0$ があることに注目してください。ウエストを細くするほど広がり角は大きくなる。これが冒頭で予告した「絞ると広がる」の正体です。回折の一般法則（開口が小さいほど回折角が大きい）が、ガウシアンビームでも同じ形で現れています。

$w(z)$ の双曲線は、この直線 $w = \theta z$ を漸近線に持ちます。近傍では平坦、遠方では直線 — 砂時計の形はこの2つの極限をなめらかにつないだものです。

左図は両対数目盛なので、$\theta = \lambda/(\pi w_0)$ の反比例関係が傾き $-1$ の直線として現れています。同じウエスト半径なら波長が短いほど広がり角が小さい（633 nm の破線が下）ことも読み取れます。右図は $\lambda = 1064\,\mathrm{nm}$、$w_0 = 1\,\mathrm{mm}$ のビームがどこまで広がるかを距離 $0.1\,\mathrm{m}$ から月までプロットしたものです。$z_R = 2.95\,\mathrm{m}$ を境に曲線が水平から傾き $1$ の直線へ折れ曲がり、そこから先は $w \approx \theta z$ にぴたりと乗ります。1 km で $33.9\,\mathrm{cm}$、月面で $129\,\mathrm{km}$ — 冒頭で予告した数字がこの直線の上にあります。

曲率半径 $R(z)$ — 平面波と球面波のあいだ

$R(z) = z + z_R^2/z$ のふるまいを追ってみます。

• $z \to 0$（ウエスト）: $R \to \infty$。波面は平面です。
• $z = z_R$: $R = 2z_R$ で最小値を取ります。ここが波面の曲がりが最もきついポイントです。
• $z \to \infty$: $R \to z$。中心がウエストにある球面波に漸近します。

$R$ が最小値を持つのは意外に感じるかもしれませんが、微分すればすぐわかります。$dR/dz = 1 – z_R^2/z^2 = 0$ より $z = z_R$ です。「ウエストでは平面、遠方では球面、その中間で最も曲がる」— ガウシアンビームが平面波と球面波の中間的存在であることを、$R(z)$ が端的に示しています。

左図の $R(z)$ が $z \to 0$ で発散し、$z = z_R$ で最小値 $2z_R$ を取り、その後は破線（$R = z$）に漸近していく様子が確認できます。曲率半径が「大きい＝平ら」なので、この谷が最も波面の曲がりが強い点です。右図は同じことを波面の形で描いたもので、$z=0$ の緑は直線（平面）、$z=z_R$ のオレンジが最も強く反り、$z=3z_R$ の青は半径の大きい球面に近づいて再び平らに見えはじめます。ウエストでも遠方でも「平ら」なのに理由が正反対（前者は曲がっていない、後者は半径が大きすぎる）という点が、この図の面白いところです。

グーイ位相 $\zeta(z)$ — 焦点で消える半波長

$\zeta(z) = \arctan(z/z_R)$ は、$z = -\infty$ から $+\infty$ まで伝搬するあいだに $-\pi/2$ から $+\pi/2$ まで、合計 $\pi$ だけ変化します。電場の位相は $-\zeta$ で入るので、平面波と比べて焦点を通過する際に $\pi$ だけ位相が遅れる（正確には進む向きの符号は規約次第ですが、$\pi$ ずれることは変わりません）。これがグーイ位相です。

小さな効果に見えますが、実害があります。レーザー共振器の共振条件は往復位相が $2\pi$ の整数倍という条件で決まるので、グーイ位相の分だけ共振周波数がずれます。共振器設計では必ず勘定に入れる項です。また、非線形光学の位相整合や、干渉計で焦点をまたぐ光路を扱うときにも効いてきます。

直感的には、ビームが横方向に閉じ込められている（$\Delta x$ が有限）ため、不確定性原理により横方向の運動量成分 $\Delta k_x$ を持ち、その分だけ縦方向の波数 $k_z = \sqrt{k^2 – k_T^2}$ がわずかに小さくなる — と説明されます。閉じ込めの代償として位相が遅れる、というわけです。

左図の $\zeta(z)$ は $-\pi/2$ と $+\pi/2$ を漸近線に持つ arctan 曲線で、全行程での変化量がちょうど $\pi$ になっていることが目で確認できます。変化の大半は $|z| < 2z_R$ の焦点近傍で起きており、遠方ではほとんど増えません。右図はその微分 $d\zeta/dz$ を波数 $k$ で割った量（実効波数の相対的な減り方）で、ピークでも $0.06\,\mathrm{ppm}$ 程度しかありません。1回通過あたりの効果は極めて小さいのに、共振器のように何万回も往復する系では無視できなくなる — この桁感がグーイ位相の扱いにくさです。

パワー保存の確認

強度を全断面で積分してみます。$I(r,z) = I_0 (w_0/w)^2 \exp(-2r^2/w^2)$ として、

$$ P = \int_0^\infty I(r,z)\, 2\pi r\, dr = I_0 \frac{w_0^2}{w^2} \cdot 2\pi \cdot \frac{w^2}{4} = \frac{\pi w_0^2 I_0}{2} $$

$z$ に依存しない定数になりました。振幅に $w_0/w$ を付けたのは、まさにこの保存を成り立たせるためです。逆にピーク強度は

$$ I(0, z) = \frac{2P}{\pi w^2(z)} $$

と、断面積に反比例して薄まります。ウエストで最大、レイリー長で半分です。

数値例で感覚をつかむ

式だけでは実感が湧かないので、実際の数字を入れてみましょう。Nd:YAGレーザー（$\lambda = 1064\,\mathrm{nm}$）で、ウエスト半径 $w_0 = 1\,\mathrm{mm}$ とします。

$$ z_R = \frac{\pi (1\times10^{-3})^2}{1.064\times10^{-6}} = 2.95\,\mathrm{m}, \qquad \theta = \frac{1.064\times10^{-6}}{\pi \times 1\times10^{-3}} = 0.339\,\mathrm{mrad} $$

広がり角は $0.339\,\mathrm{mrad}$、角度にして約 $0.019^\circ$。非常に小さいのですが、距離をかけると効いてきます。$z = 1\,\mathrm{km}$ で $w = 33.9\,\mathrm{cm}$、$z = 100\,\mathrm{km}$ で $33.9\,\mathrm{m}$、月面（$3.8\times10^8\,\mathrm{m}$）では $129\,\mathrm{km}$ です。アポロ計画の月レーザー測距で、地上の望遠鏡で拡大して送ってなお月面のスポットが数kmになるのは、この式の通りの帰結です。

一方 He-Ne レーザー（$\lambda = 633\,\mathrm{nm}$）で $w_0 = 0.5\,\mathrm{mm}$ なら、$z_R = 1.24\,\mathrm{m}$、$\theta = 0.403\,\mathrm{mrad}$。実験室のテーブル（数メートル）のスケールで、すでにビームが目に見えて太くなる領域に入っていることがわかります。

パラメータの意味が見えてきました。ここで一つ、深い性質に気づきます。$w_0$ と $\theta$ は独立に選べるでしょうか。次節でこれを詰めます。

ビームパラメータ積と $M^2$ 因子

$w_0$ を小さくすれば $\theta$ が大きくなる、という関係を見てきました。実は両者の積は波長だけで決まる定数になります。

$$ \begin{equation} w_0 \theta = w_0 \cdot \frac{\lambda}{\pi w_0} = \frac{\lambda}{\pi} \end{equation} $$

$w_0$ がきれいに消えました。この積をビームパラメータ積（Beam Parameter Product, BPP）と呼びます。$\lambda = 1064\,\mathrm{nm}$ なら $0.339\,\mathrm{mm\cdot mrad}$、$\lambda = 1550\,\mathrm{nm}$ なら $0.493\,\mathrm{mm\cdot mrad}$。

この式が語っているのは強い制約です。レンズをどう組み合わせても、BPP は変えられません。レンズはビームを細く絞れますが、その代わり必ず広がり角を増やします。ちょうど、粘土の総量を変えずに形だけ変えるようなものです。統計力学のリウヴィルの定理や、量子力学の不確定性関係と同じ構造で、光学系における「保存量」になっています。

ちなみに、不確定性関係との対応は形式的なアナロジーにとどまりません。ガウシアンビームは近軸波動方程式（＝自由粒子シュレディンガー方程式と同型）の最小不確定性波束に相当し、BPP が下限値 $\lambda/\pi$ を取るのはそのためです。

実際のレーザーは理想通りではない — $M^2$

理論上の下限が $\lambda/\pi$ ということは、実在のビームはこれより悪いということです。高次モードの混入、共振器の収差、増幅媒質の熱レンズ効果などにより、実際のレーザーの BPP は $\lambda/\pi$ より大きくなります。この劣化度合いを表すのが$M^2$ 因子（ビーム品質因子）です。

$$ \begin{equation} W_0 \Theta = M^2 \frac{\lambda}{\pi}, \qquad M^2 \ge 1 \end{equation} $$

大文字 $W_0, \Theta$ は実ビームのウエスト半径と広がり角です。$M^2 = 1$ が理想の基本ガウシアンビーム（$\mathrm{TEM}_{00}$）、$M^2 > 1$ が実ビーム。単一モードの He-Ne やファイバレーザーで $M^2 \approx 1.0\text{–}1.1$、高出力の固体レーザーで $1.2\text{–}2$、マルチモードのダイオードレーザーバーでは $20$ を超えることもあります。

$M^2$ の便利なところは、理論式をほぼそのまま流用できる点です。実ビームの伝搬は

$$ W(z) = W_0 \sqrt{1 + \left( \frac{M^2 \lambda z}{\pi W_0^2} \right)^2} = W_0\sqrt{1 + \left( \frac{z}{z_R’} \right)^2}, \qquad z_R’ = \frac{\pi W_0^2}{M^2 \lambda} $$

と書けます。つまり「波長を $M^2\lambda$ に置き換えた理想ビーム」として扱えばよい。レイリー長は $1/M^2$ 倍に短くなり、集光スポット径は $M^2$ 倍に太くなります。$M^2 = 2$ のレーザーは、同じレンズで絞ってもスポット面積が4倍になるということです。レーザー加工でビーム品質が価格に直結する理由がここにあります。

左図の各曲線は $W_0\Theta = M^2\lambda/\pi$ という双曲線で、レンズを入れても同じ曲線の上を滑るだけ（左下へ移動できない）というのが BPP 保存の意味です。$M^2$ が増えると曲線ごと右上へ平行に押し上げられ、どんなに工夫しても $M^2=1$ の曲線より内側には入れません。右図は同じウエスト半径 $0.2\,\mathrm{mm}$ から出発した3本の伝搬で、$M^2 = 4$ では $z_R’$ が $118\,\mathrm{mm}$ から $29.5\,\mathrm{mm}$ へ4分の1に縮み、同じ距離で4倍太くなります。ビーム品質がそのまま加工精度と価格に効く理由が読み取れます。

$M^2$ は実測でも求められます。ビームを一度レンズで絞り、その周りで $W(z)$ を複数点測って双曲線をフィットする（ISO 11146 の手順）と、$W_0$ と $\Theta$ が得られ、$M^2$ が算出できます。

さて、ここまでは自由空間の伝搬だけを扱ってきました。実際の光学系にはレンズやミラーが入ります。それらを通過したビームがどうなるかを、いちいち解き直すのは大変です。次節では、驚くほど簡単にこれを扱う道具 — ABCD 則 — を導入します。

複素ビームパラメータ $q$ と ABCD 則

幾何光学では、光線を「高さ $r$ と傾き $r’$」の2成分ベクトルで表し、光学素子を $2\times2$ 行列で表します。

$$ \begin{pmatrix} r_{\text{out}} \\ r’_{\text{out}} \end{pmatrix} = \begin{pmatrix} A & B \\ C & D \end{pmatrix} \begin{pmatrix} r_{\text{in}} \\ r’_{\text{in}} \end{pmatrix} $$

代表的な行列を挙げておきます。

複数素子を通るときは、行列を光が通る逆順に掛け合わせます。ここまでは幾何光学の話です。

ABCD 則

ガウシアンビームに対しては、次の驚くべき法則が成り立ちます。

$$ \begin{equation} \boxed{q_{\text{out}} = \frac{A q_{\text{in}} + B}{C q_{\text{in}} + D}} \end{equation} $$

これがABCD 則（Kogelnik の ABCD 法則）です。幾何光学の光線行列を、そのままガウシアンビームの複素パラメータに適用できてしまう。ビーム半径と曲率半径を別々に追いかける必要はなく、複素数 $q$ を1つ更新するだけで光学系全体を通せます。

なぜこれが成り立つのか、直感的な理由を示しておきます。曲率半径 $R$ の球面波を考えると、その光線は高さ $r$ に対して傾き $r’ = r/R$ を持ちます。この光線を行列に通すと

$$ r_{\text{out}} = A r + B\frac{r}{R}, \qquad r’_{\text{out}} = C r + D\frac{r}{R} $$

出射側の曲率半径は $R_{\text{out}} = r_{\text{out}}/r’_{\text{out}}$ ですから、

$$ R_{\text{out}} = \frac{A + B/R}{C + D/R} = \frac{AR + B}{CR + D} $$

球面波の曲率半径については、確かにこの形の変換則が成り立つわけです。ガウシアンビームは「曲率半径が複素数になった球面波」でしたから、$R \to q$ と置き換えればよい — というのが直感的な筋道です。厳密には、回折積分（Collins の一般化ホイヘンス・フレネル積分）を ABCD 行列で書き下し、ガウス型入力を代入して評価することで証明されます。

2つの基本素子で確認する

自由空間伝搬（$A=1, B=d, C=0, D=1$）:

$$ q_{\text{out}} = \frac{q + d}{1} = q + d $$

先に導いた $q’ = 1$、すなわち $q(z) = q_0 + z$ と完全に一致します。ABCD 則が既知の結果を再現しました。

薄肉レンズ（$A=1, B=0, C=-1/f, D=1$）:

$$ q_{\text{out}} = \frac{q}{-q/f + 1} \quad \Longrightarrow \quad \frac{1}{q_{\text{out}}} = \frac{1}{q} – \frac{1}{f} $$

逆数で書くとすっきりします。実部・虚部に分けると、

$$ \frac{1}{R_{\text{out}}} = \frac{1}{R_{\text{in}}} – \frac{1}{f}, \qquad w_{\text{out}} = w_{\text{in}} $$

物理的にもっともな結果です。薄いレンズは波面の曲率だけを変え、その場のビーム太さは変えない。レンズの直前と直後で光束の太さが変わらないのは当たり前で、変わるのは「これからどちらへ向かうか」だけ、というわけです。

2つの基本素子を、複素平面上の平行移動として描いたものです。左の自由空間伝搬は $q$ 平面で右へ $d$ だけずらす操作で、実部（ウエストからの距離 $z$）だけが増え、虚部（レイリー長 $z_R$）は伝搬しても変わりません。右の薄肉レンズは $1/q$ 平面で左へ $1/f$ だけずらす操作で、実部（$1/R$、波面の曲がり）だけが減り、虚部（$-\lambda/\pi w^2$、その場のビーム半径）は変わりません。後者の「虚部不変」が、まさに $w_{\text{out}} = w_{\text{in}}$ という結論そのものです。光学系の計算がこれほど単純な図形操作に還元されるのが、ABCD 則のありがたみです。

ABCD 則という強力な道具が手に入りました。これを使って、実務で最も知りたい量 — レンズでどこまで絞れるか — を導きます。

レンズによる集光 — スポット径と焦点シフト

平行に近いビーム（ウエストがちょうどレンズ位置にある、半径 $w$ のビーム）を焦点距離 $f$ のレンズで絞る、という最も基本的な状況を解きます。

レンズ直前では波面が平面（$R = \infty$）なので、

$$ q_{\text{in}} = j z_{R} , \qquad z_{R} = \frac{\pi w^2}{\lambda} $$

レンズを通すと $1/q_{\text{out}} = 1/q_{\text{in}} – 1/f = -j/z_R – 1/f$ です。通分すると

$$ \frac{1}{q_{\text{out}}} = -\frac{z_R + jf}{f z_R} \quad \Longrightarrow \quad q_{\text{out}} = -\frac{f z_R}{z_R + jf} $$

分母を実数化するため、分子・分母に $(z_R – jf)$ を掛けます。

$$ q_{\text{out}} = -\frac{f z_R (z_R – jf)}{z_R^2 + f^2} = -\frac{f z_R^2}{z_R^2+f^2} + j\frac{f^2 z_R}{z_R^2+f^2} $$

レンズから距離 $z$ 進んだ位置では $q(z) = q_{\text{out}} + z$ です。新しいウエストは実部がゼロになる位置（そこで $R = \infty$）なので、

$$ \begin{equation} \boxed{z_w = \frac{f z_R^2}{z_R^2 + f^2} = \frac{f}{1 + (f/z_R)^2}} \end{equation} $$

そしてそのときの虚部が新しいレイリー長 $z_R’$ です。

$$ z_R’ = \frac{f^2 z_R}{z_R^2+f^2} \quad \Longrightarrow \quad w_f^2 = \frac{\lambda z_R’}{\pi} = \frac{\lambda}{\pi}\cdot\frac{f^2 z_R}{z_R^2+f^2} $$

$z_R = \pi w^2/\lambda$ を代入すると $\lambda z_R/\pi = w^2$ なので、

$$ \begin{equation} \boxed{w_f = \frac{f w}{\sqrt{f^2 + z_R^2}}} \end{equation} $$

実用的な近似 — $2w_f \approx 4\lambda f/(\pi D)$

多くの実験では、入射ビームが太くレイリー長が非常に長い（$z_R \gg f$）状況です。例えば $\lambda = 633\,\mathrm{nm}$、$w = 2\,\mathrm{mm}$ なら $z_R = 19.9\,\mathrm{m}$ で、焦点距離 $100\,\mathrm{mm}$ のレンズに対して $200$ 倍も長い。このとき $\sqrt{f^2+z_R^2} \approx z_R$ として、

$$ w_f \approx \frac{f w}{z_R} = \frac{f w \lambda}{\pi w^2} = \frac{\lambda f}{\pi w} $$

入射ビームの $1/e^2$ 直径を $D = 2w$ と書き直すと、集光スポットの $1/e^2$ 直径は

$$ \begin{equation} \boxed{2w_f \approx \frac{4\lambda f}{\pi D} = \frac{4}{\pi}\lambda \cdot F_\# \approx 1.27\,\lambda\, F_\#} \end{equation} $$

ここで $F_\# = f/D$ は F ナンバーです。この形は覚えやすく、「$F_\#$ にほぼ波長を掛けたものが最小スポット径」という設計の目安になります。$F_\# = 2$ の光学系なら、波長 $1\,\mu\mathrm{m}$ の光は直径 $2.5\,\mu\mathrm{m}$ 程度までしか絞れません。回折限界の一般則（エアリーディスク径 $\approx 2.44\lambda F_\#$）と桁が揃っているのも安心材料です。

具体例で確認します。$\lambda = 633\,\mathrm{nm}$、入射 $1/e^2$ 直径 $D = 4\,\mathrm{mm}$、$f = 100\,\mathrm{mm}$ のとき、

$$ 2w_f = \frac{4 \times 633\times10^{-9} \times 0.1}{\pi \times 4\times10^{-3}} = 20.1\,\mu\mathrm{m} $$

厳密式 $w_f = fw/\sqrt{f^2+z_R^2}$ で計算しても $2w_f = 20.15\,\mu\mathrm{m}$ で、相対差は $0.0013\%$（相対誤差はおよそ $\frac{1}{2}(f/z_R)^2$）。この条件では近似式で十分です。

左図は $2w_f \approx 1.27\lambda F_\#$ が $F_\#$ に対する直線であることと、波長が長いほど傾きが急になることを示しています。灰色破線のエアリーディスク径 $2.44\lambda F_\#$ とは係数が2倍ほど違うだけで、同じ「回折限界」の家族であることが見て取れます。右図は $f = 100\,\mathrm{mm}$ 固定で入射半径を振ったもので、$w$ が大きいほどスポットが細くなる反比例の関係（赤破線の近似式）に、厳密式（青）がぴたりと重なります。ただし $w \approx 0.14\,\mathrm{mm}$ を境に青が反転して下がる — ちょうど $z_R = f$ となる点で、これより細い入射ビームでは $z_R < f$ となり近似式が破綻します（次節の焦点シフトの領域です）。

この式は「絞りたければ入射ビームを太くせよ」と言っています。直感に反するようですが、太いビームほど回折の影響が小さく、レンズがかける角度が大きくなるためです。実際、レーザー加工機ではビームエキスパンダで一度太くしてから集光レンズに入れます。

焦点シフト — 幾何光学の焦点にはピントが合わない

ここで $z_w$ の式に戻ります。

$$ z_w = \frac{f}{1+(f/z_R)^2} $$

幾何光学なら、平行光線は必ず焦点距離 $f$ の位置に集まるはずです。ところが分母が $1$ より大きいので、$z_w < f$ — 実際のウエストはレンズ寄りに来ます。これを焦点シフト（focal shift）と呼びます。

シフト量は $f/z_R$ の大小で決まります。

• $z_R \gg f$（太いビーム、短焦点）: $z_w \approx f$。幾何光学がよく当たります。先ほどの例では $z_w = 99.9975\,\mathrm{mm}$ で、ずれは $2.5\,\mu\mathrm{m}$ しかありません。
• $z_R \sim f$: 無視できないずれが出ます。$\lambda=1064\,\mathrm{nm}$、$w = 0.5\,\mathrm{mm}$（$z_R = 0.738\,\mathrm{m}$）を $f = 200\,\mathrm{mm}$ で集光すると $z_w = 186.3\,\mathrm{mm}$ で、$6.8\%$ 手前にずれます。
• $z_R \ll f$（細いビーム、長焦点）: $z_w \approx z_R^2/f \ll f$。$\lambda=633\,\mathrm{nm}$、$w = 0.1\,\mathrm{mm}$（$z_R = 49.6\,\mathrm{mm}$）を $f=50\,\mathrm{mm}$ で集光すると $z_w = 24.8\,\mathrm{mm}$。なんと焦点距離の半分の位置にウエストが来ます。

なぜこうなるのか。幾何光学の「平行光線」は $z_R = \infty$ に相当します。有限の $z_R$ を持つ実ビームは、レンズに入る時点ですでにわずかに発散しており、その発散ぶんだけ焦点が手前に引き寄せられる — と読めます。あるいは、レンズによる収束と回折による発散が競合し、両者が釣り合う位置が $f$ より手前になる、と考えてもよいでしょう。

これは実務で効きます。共焦点顕微鏡や光ピンセットのように焦点位置をミクロン精度で管理する系、あるいは焦点距離の長いビームエキスパンダの調整では、幾何光学の計算で位置を出すとピントがずれます。$q$ 伝搬で計算するのが正解です。

理論は出そろいました。次はこれらをPythonで実装し、目で確かめます。

Pythonでの実装

ここからは、$q$ パラメータの伝搬を実装して4つのことを可視化します。(1) 自由空間での $w(z)$ 曲線と双曲線・漸近線、(2) 強度分布の2Dマップ、(3) レンズ集光系のビーム径推移、(4) 焦点シフト現象です。

まずは共通の道具立てとして、$q$ パラメータを扱う関数群を用意します。

import numpy as np
import matplotlib
import matplotlib.pyplot as plt

# 日本語フォント設定
for cand in ["Hiragino Sans", "Yu Gothic", "Noto Sans CJK JP", "IPAexGothic", "Meiryo"]:
    if any(cand == f.name for f in matplotlib.font_manager.fontManager.ttflist):
        plt.rcParams["font.family"] = cand
        break
plt.rcParams["axes.unicode_minus"] = False


def rayleigh_range(w0, lam):
    """レイリー長 z_R = pi*w0^2/lambda"""
    return np.pi * w0**2 / lam


def q_at_waist(w0, lam):
    """ウエスト位置での複素ビームパラメータ q = j*z_R"""
    return 1j * rayleigh_range(w0, lam)


def w_R_from_q(q, lam):
    """複素パラメータ q からビーム半径 w と曲率半径 R を取り出す"""
    inv = 1.0 / q
    w = np.sqrt(-lam / (np.pi * inv.imag))   # 虚部がビーム半径を担う
    R = np.inf if inv.real == 0 else 1.0 / inv.real  # 実部が曲率を担う
    return w, R


def apply_abcd(M, q):
    """ABCD則 q' = (Aq+B)/(Cq+D)"""
    A, B, C, D = M[0, 0], M[0, 1], M[1, 0], M[1, 1]
    return (A * q + B) / (C * q + D)


def M_free(d):
    """自由空間伝搬の光線行列"""
    return np.array([[1.0, d], [0.0, 1.0]])


def M_lens(f):
    """薄肉レンズの光線行列"""
    return np.array([[1.0, 0.0], [-1.0 / f, 1.0]])

w_R_from_q が導出のキモをそのままコードにした部分です。$1/q$ の実部が曲率、虚部がビーム半径という対応を、数式どおりに実装しています。ABCD 則を使うと、レンズも自由空間も同じ1行（apply_abcd）で処理できるのが分かります。

(1) 自由空間での $w(z)$ と双曲線

まず、ウエストから前後に伝搬させて $w(z)$ を描きます。理論式と $q$ 伝搬の結果が一致することも同時に確認します。

lam = 1064e-9      # 波長 1064 nm (Nd:YAG)
w0 = 1.0e-3        # ウエスト半径 1 mm
zR = rayleigh_range(w0, lam)
theta = lam / (np.pi * w0)
print(f"レイリー長 z_R = {zR:.3f} m")
print(f"広がり角 theta = {theta*1e3:.4f} mrad")

z = np.linspace(-4 * zR, 4 * zR, 801)
w_theory = w0 * np.sqrt(1 + (z / zR) ** 2)                  # 理論式
w_q = np.array([w_R_from_q(q_at_waist(w0, lam) + zi, lam)[0] for zi in z])  # q伝搬
print(f"理論式とq伝搬の最大差: {np.max(np.abs(w_theory - w_q)):.2e} m")

fig, ax = plt.subplots(figsize=(10, 5))
ax.plot(z, w_theory * 1e3, color="royalblue", lw=2, label="ビーム半径 $w(z)$")
ax.plot(z, -w_theory * 1e3, color="royalblue", lw=2)
ax.plot(z, np.abs(theta * z) * 1e3, "--", color="crimson", lw=1.5, label="遠方界の漸近線 $w=\\theta z$")
ax.plot(z, -np.abs(theta * z) * 1e3, "--", color="crimson", lw=1.5)
ax.axvline(zR, color="gray", ls=":"); ax.axvline(-zR, color="gray", ls=":")
ax.axhline(w0 * 1e3, color="green", ls=":", lw=1)
ax.annotate("ウエスト $w_0$", xy=(0, w0 * 1e3), xytext=(0.6 * zR, 2.2),
            arrowprops=dict(arrowstyle="->"), fontsize=11)
ax.annotate("$z=z_R$ で $\\sqrt{2}w_0$", xy=(zR, np.sqrt(2) * w0 * 1e3),
            xytext=(1.3 * zR, 3.4), arrowprops=dict(arrowstyle="->"), fontsize=11)
ax.set_xlabel("伝搬距離 $z$ [m]"); ax.set_ylabel("ビーム半径 [mm]")
ax.set_title("ガウシアンビームの広がり（砂時計型の双曲線）")
ax.legend(loc="upper center"); ax.grid(alpha=0.3)
plt.tight_layout(); plt.show()

このグラフから3つのことが読み取れます。第一に、$w(z)$ は原点で平坦、遠方で直線に漸近する典型的な双曲線を描きます。「近くでは平面波、遠くでは球面波」という中間的性質が、形として見えているわけです。第二に、$z = z_R$（点線）で半径がちょうど $\sqrt{2}w_0 = 1.41\,\mathrm{mm}$ になっており、レイリー長の定義が確認できます。第三に、理論式と $q$ 伝搬の差は $10^{-18}\,\mathrm{m}$ オーダー（浮動小数点誤差）で完全に一致します。ABCD 則が正しく実装できている証拠です。

(2) 強度分布の2Dマップと横断面

次に、$(r, z)$ 平面での強度分布を描き、ビームの「中身」を見ます。

r = np.linspace(-3.5e-3, 3.5e-3, 400)
Z, R_grid = np.meshgrid(z, r)
W = w0 * np.sqrt(1 + (Z / zR) ** 2)
I = (w0 / W) ** 2 * np.exp(-2 * R_grid**2 / W**2)   # 規格化強度

fig, axes = plt.subplots(1, 2, figsize=(13, 4.8))
ax = axes[0]
pcm = ax.pcolormesh(z, r * 1e3, I, shading="auto", cmap="inferno")
ax.plot(z, w_theory * 1e3, "w--", lw=1.2, label="$\\pm w(z)$（強度 $1/e^2$ 境界）")
ax.plot(z, -w_theory * 1e3, "w--", lw=1.2)
ax.set_xlabel("伝搬距離 $z$ [m]"); ax.set_ylabel("横方向 $r$ [mm]")
ax.set_title("ガウシアンビームの強度分布 $I(r,z)$")
ax.legend(loc="upper right", fontsize=9)
fig.colorbar(pcm, ax=ax, label="規格化強度")

ax = axes[1]
for zi, c in zip([0, zR, 2 * zR, 3 * zR], ["k", "royalblue", "seagreen", "darkorange"]):
    wi = w0 * np.sqrt(1 + (zi / zR) ** 2)
    Ii = (w0 / wi) ** 2 * np.exp(-2 * r**2 / wi**2)
    ax.plot(r * 1e3, Ii, color=c, lw=2, label=f"$z={zi/zR:.0f}z_R$（ピーク {(w0/wi)**2:.2f}）")
ax.set_xlabel("横方向 $r$ [mm]"); ax.set_ylabel("規格化強度")
ax.set_title("横断強度プロファイルの変化（面積＝全パワーは一定）")
ax.legend(fontsize=9); ax.grid(alpha=0.3)
plt.tight_layout(); plt.show()

左のマップでは、ビームが中央でくびれ、前後に向かって明るさを失いながら太くなる様子がはっきり見えます。白破線の $\pm w(z)$ が、明るい領域の境界をきれいになぞっていることに注目してください。右の横断面では、$z = z_R$ でピーク強度がちょうど $0.50$ に落ちています。これは $I(0,z) \propto (w_0/w)^2 = 1/(1+(z/z_R)^2)$ という式そのままの結果で、$z = z_R$ で半分になるというレイリー長の別の顔です。各曲線の下の面積（回転体として積分した値）は等しく、パワーが保存されていることも読み取れます。

(3) レンズ集光系のビーム径推移

いよいよ ABCD 則の本領です。平行ビームをレンズで絞る系を、$q$ 伝搬で最初から最後まで追いかけます。

lam2 = 633e-9        # He-Ne 633 nm
w_in = 2.0e-3        # レンズ入射時のビーム半径 2 mm（D = 4 mm）
f_lens = 0.10        # 焦点距離 100 mm

q = q_at_waist(w_in, lam2)          # レンズ位置にウエストがある平行ビーム
zR_in = rayleigh_range(w_in, lam2)
q_after = apply_abcd(M_lens(f_lens), q)   # レンズ通過（1/q が -1/f だけシフト）

z_out = np.linspace(0, 0.16, 3201)
w_out = np.array([w_R_from_q(apply_abcd(M_free(zi), q_after), lam2)[0] for zi in z_out])

# 理論値との比較
w_f_exact = f_lens * w_in / np.sqrt(f_lens**2 + zR_in**2)
w_f_approx = lam2 * f_lens / (np.pi * w_in)
z_w_theory = f_lens / (1 + (f_lens / zR_in) ** 2)
print(f"入射ビームのレイリー長 : {zR_in:.2f} m")
print(f"スポット半径（厳密式） : {w_f_exact*1e6:.3f} um")
print(f"スポット半径（近似式） : {w_f_approx*1e6:.3f} um")
print(f"スポット直径 2w_f      : {2*w_f_exact*1e6:.2f} um")
print(f"ウエスト位置（理論）   : {z_w_theory*1e3:.4f} mm")
# q の実部がゼロになる位置がウエスト → z_w = -Re(q_after)
print(f"ウエスト位置（q から） : {-q_after.real*1e3:.4f} mm")
print(f"集光後のレイリー長     : {q_after.imag*1e6:.1f} um")

出力は次のようになります。

入射ビームのレイリー長 : 19.85 m
スポット半径（厳密式） : 10.074 um
スポット半径（近似式） : 10.075 um
スポット直径 2w_f      : 20.15 um
ウエスト位置（理論）   : 99.9975 mm
ウエスト位置（q から） : 99.9975 mm
集光後のレイリー長     : 503.7 um

厳密式と近似式 $\lambda f/(\pi w)$ が小数第3位まで一致しており、$z_R \gg f$ の条件下では近似式で十分実用になることが確認できます。ウエスト位置も、公式 $z_w = f/[1+(f/z_R)^2]$ と「$q$ の実部がゼロになる位置」という定義から求めた値がぴたり一致しました。注目すべきは集光後のレイリー長で、わずか $504\,\mu\mathrm{m}$ しかありません。半径 $10\,\mu\mathrm{m}$ まで絞ると、その状態が保たれるのは前後 $0.5\,\mathrm{mm}$ の範囲だけ — 「細く絞るほど焦点深度が2乗で縮む」という関係が数字で効いています。レーザー加工でZ軸のフォーカス制御が精密になる理由がこれです。

続いて、この集光の様子を図にします。

fig, axes = plt.subplots(1, 2, figsize=(13, 4.8))
ax = axes[0]
ax.plot(z_out * 1e3, w_out * 1e3, color="royalblue", lw=2)
ax.plot(z_out * 1e3, -w_out * 1e3, color="royalblue", lw=2)
ax.axvline(f_lens * 1e3, color="crimson", ls="--", lw=1.5, label="幾何光学の焦点 $z=f$")
ax.axvline(0, color="k", lw=3, alpha=0.5)
ax.text(1, 1.7, "レンズ", fontsize=11)
ax.set_xlabel("レンズからの距離 [mm]"); ax.set_ylabel("ビーム半径 [mm]")
ax.set_title(f"レンズ集光（$f$={f_lens*1e3:.0f} mm, 入射半径 {w_in*1e3:.0f} mm）")
ax.legend(); ax.grid(alpha=0.3)

# 焦点近傍を拡大
ax = axes[1]
mask = (z_out > 0.0975) & (z_out < 0.1025)
ax.plot((z_out[mask] - f_lens) * 1e3, w_out[mask] * 1e6, color="royalblue", lw=2)
ax.plot((z_out[mask] - f_lens) * 1e3, -w_out[mask] * 1e6, color="royalblue", lw=2)
ax.axhline(w_f_exact * 1e6, color="green", ls=":", label=f"$w_f$={w_f_exact*1e6:.1f} $\\mu$m")
ax.axhline(-w_f_exact * 1e6, color="green", ls=":")
zRf = rayleigh_range(w_f_exact, lam2)
ax.axvspan(-zRf * 1e3, zRf * 1e3, color="orange", alpha=0.2,
           label=f"$\\pm z_R'$={zRf*1e6:.0f} $\\mu$m（焦点深度）")
ax.set_xlabel("焦点からのずれ [mm]"); ax.set_ylabel("ビーム半径 [$\\mu$m]")
ax.set_title("焦点近傍の拡大 — スポット径と焦点深度")
ax.legend(fontsize=9); ax.grid(alpha=0.3)
plt.tight_layout(); plt.show()

左の全体図では、$2\,\mathrm{mm}$ の太いビームがレンズで折り曲げられ、$100\,\mathrm{mm}$ 先で一点に収束していく様子が見えます。幾何光学の焦点（赤破線）とほぼ重なっており、この条件では両者の差が肉眼で分からないほど小さいことがわかります。右の拡大図が本題で、焦点は「点」ではなく半径 $10\,\mu\mathrm{m}$ のくびれを持っています。オレンジの帯（$\pm z_R’ = \pm 504\,\mu\mathrm{m}$）の外では、ビーム半径が $\sqrt{2}$ 倍を超えて急速に太くなります。幾何光学では消えてしまう「焦点の構造」が、ガウシアンビームでは定量的に扱えることが一目でわかります。

(4) 焦点シフト — 幾何光学とのずれ

最後に、焦点シフトが顕著に現れる条件を作って可視化します。入射ビームを細くし（$w = 0.1\,\mathrm{mm}$）、$z_R$ と $f$ を同程度にします。

w_in2 = 0.1e-3        # 細い入射ビーム 100 um
f2 = 0.05             # f = 50 mm
zR_in2 = rayleigh_range(w_in2, lam2)
q2 = apply_abcd(M_lens(f2), q_at_waist(w_in2, lam2))
z2 = np.linspace(0, 0.08, 8001)
w2 = np.array([w_R_from_q(apply_abcd(M_free(zi), q2), lam2)[0] for zi in z2])
z_w2 = f2 / (1 + (f2 / zR_in2) ** 2)
print(f"入射レイリー長 z_R = {zR_in2*1e3:.2f} mm （f = {f2*1e3:.0f} mm と同程度）")
print(f"実際のウエスト位置 : {z2[w2.argmin()]*1e3:.3f} mm（理論 {z_w2*1e3:.3f} mm）")
print(f"幾何光学の焦点     : {f2*1e3:.1f} mm → シフト率 {(1-z_w2/f2)*100:.1f} %")

# 焦点シフト率を f/z_R の関数として描く
ratio = np.logspace(-2.5, 1.5, 300)        # f / z_R
shift = 1.0 / (1.0 + ratio**2)             # z_w / f

fig, axes = plt.subplots(1, 2, figsize=(13, 4.8))
ax = axes[0]
ax.plot(z2 * 1e3, w2 * 1e6, color="royalblue", lw=2, label="ビーム半径 $w(z)$")
ax.axvline(f2 * 1e3, color="crimson", ls="--", lw=1.5, label=f"幾何光学の焦点 {f2*1e3:.0f} mm")
ax.axvline(z_w2 * 1e3, color="seagreen", ls="-.", lw=1.5, label=f"実際のウエスト {z_w2*1e3:.1f} mm")
ax.set_xlabel("レンズからの距離 [mm]"); ax.set_ylabel("ビーム半径 [$\\mu$m]")
ax.set_title("焦点シフト：ウエストは幾何光学の焦点より手前に来る")
ax.legend(fontsize=9); ax.grid(alpha=0.3)

ax = axes[1]
ax.semilogx(ratio, shift, color="darkorange", lw=2)
ax.axhline(1.0, color="gray", ls=":", label="幾何光学の予測 $z_w=f$")
ax.scatter([f_lens / zR_in], [1 / (1 + (f_lens / zR_in) ** 2)], color="royalblue",
           zorder=5, s=60, label="(3)の条件：ほぼ一致")
ax.scatter([f2 / zR_in2], [1 / (1 + (f2 / zR_in2) ** 2)], color="seagreen",
           zorder=5, s=60, label="(4)の条件：50%手前")
ax.set_xlabel("$f / z_R$（レンズ焦点距離 ÷ 入射ビームのレイリー長）")
ax.set_ylabel("$z_w / f$（ウエスト位置の比）")
ax.set_title("焦点シフト率は $f/z_R$ だけで決まる")
ax.legend(fontsize=9); ax.grid(alpha=0.3, which="both")
plt.tight_layout(); plt.show()

出力は次のとおりです。

入射レイリー長 z_R = 49.63 mm （f = 50 mm と同程度）
実際のウエスト位置 : 24.810 mm（理論 24.814 mm）
幾何光学の焦点     : 50.0 mm → シフト率 50.4 %

左の図が示す事実は強烈です。幾何光学は「焦点は $50\,\mathrm{mm}$」と言うのに、実際にビームが最も細くなるのは $24.8\,\mathrm{mm}$ — 半分の位置です。ここにスクリーンを置いてピント合わせをすれば、教科書どおりに計算した位置は完全に外れます。右の図は、このずれが $f/z_R$ という1つの無次元量だけで決まることを示しています。$f/z_R \lesssim 0.1$（太いビーム・短焦点）なら幾何光学で十分ですが、$f/z_R \sim 1$ を超えると急速に破綻します。青点が (3) の条件、緑点が (4) の条件で、同じレンズの話でも入射ビームの太さ次第でここまで結果が変わるわけです。

実務的な教訓として、「入射ビームのレイリー長が焦点距離より十分長いか」を最初に確認する習慣をつけておくと、光学系の設計ミスを大幅に減らせます。判定は $z_R = \pi w^2/\lambda$ を電卓で叩くだけです。

まとめ

本記事では、ガウシアンビーム光学の基礎を、方程式の導出から実用設計式まで通して解説しました。

• 近軸ヘルムホルツ方程式 — ヘルムホルツ方程式に $U = u\,e^{-jkz}$ を代入し、緩慢変化包絡線近似で $\partial^2 u/\partial z^2$ を落とすと $\nabla_T^2 u – 2jk\,\partial u/\partial z = 0$ が得られます。これはシュレディンガー方程式と同型です
• 基本ガウシアン解 — 試行解を代入し $r$ のべきで分離すると $q’ = 1$、$P’ = -j/q$ という2本の簡単な式に落ちます。結果は $E \propto (w_0/w)\exp(-r^2/w^2)\exp[-j(kz + kr^2/2R – \zeta)]$
• 3つのパラメータ — ビーム半径 $w(z) = w_0\sqrt{1+(z/z_R)^2}$（双曲線）、曲率半径 $R(z) = z + z_R^2/z$（$z=z_R$ で最小 $2z_R$）、グーイ位相 $\zeta = \arctan(z/z_R)$（全行程で $\pi$）
• レイリー長と広がり角 — $z_R = \pi w_0^2/\lambda$ が焦点深度、$\theta = \lambda/(\pi w_0)$ が遠方界の広がり。細く絞るほど焦点深度は2乗で縮み、広がり角は増えます
• 不変量とビーム品質 — $w_0\theta = \lambda/\pi$ はレンズでは変えられない保存量。実ビームは $W_0\Theta = M^2\lambda/\pi$ で表し、$M^2$ を掛けた波長として扱えば理論式がそのまま使えます
• 複素パラメータと ABCD 則 — $1/q = 1/R – j\lambda/(\pi w^2)$ と定義すると、幾何光学の光線行列がそのまま $q’ = (Aq+B)/(Cq+D)$ として使えます。自由空間は $q+d$、レンズは $1/q – 1/f$
• 集光の設計式 — $2w_f \approx 4\lambda f/(\pi D) \approx 1.27\lambda F_\#$。ウエスト位置は $z_w = f/[1+(f/z_R)^2]$ で、$f/z_R$ が大きいほど幾何光学の焦点から手前にずれます

ガウシアンビームは、光を扱うあらゆる工学分野の共通言語です。ここで得た $w(z)$ と $q$ の道具は、レーザー共振器の安定性解析（往復 ABCD 行列の固有値問題）、光ファイバとの結合効率（モードフィールド径の重なり積分）、宇宙光通信のリンクバジェット（受信口径での電力捕捉率）へ、そのまま持ち込めます。

次のステップとして、以下の記事も参考にしてください。

• 波動光学入門（干渉と回折） — 回折の一般論。ガウシアンビームの広がりが回折の一形態であることが見えてきます
• レンズの組み合わせと光学系設計 — ABCD 行列の元になる幾何光学。ビームエキスパンダの設計へ
• 光ファイバ通信の基礎 — 単一モードファイバの基本モードはほぼガウシアン。結合効率の計算に本記事の式を使います
• 光衛星間リンク（光ISL）の原理 — 広がり角と指向精度のトレードオフが、そのままリンク設計になります
• 収差の種類と補正 — 実レンズの不完全さ。$M^2$ 劣化の原因のひとつです