波動方程式から回折現象を導くその2

プログラミング教えてください.

回折像の理論

前回（■）は \begin{align} u(x_\mathrm{P},y_\mathrm{P},z_\mathrm{P}) = \frac{1}{4\pi} \iint_{\mathrm{S}_2} dS \, \frac{e^{ikr}}{r}\left(\frac{1-ikr}{r^2}z_\mathrm{P} -\frac{\partial}{\partial z}\right)u \tag{7} \end{align} なる式を導出した（ただし $r=\sqrt{(x-x_\mathrm{P})^2+(y-y_\mathrm{P})^2+z_\mathrm{P}^2}$ ）.
この $(7)$ 式に従って開口部で波 $u$ を積分すると, $z>0$ での点 $\mathrm{P}(x_\mathrm{P},y_\mathrm{P},z_\mathrm{P})$ での波の振幅が計算できる.
よく知られた回折の近似理論に, Fresnel回折（Wikipedia: フレネル回折）とFraunhofer回折（Wikipedia: フラウンホーファー回折）というものがあるので, 以下でこれらを実際に導出してみたいと思う.

$(7)$ 式の次のステップとして, [1]では球面波 \begin{align} &u(x,y,z)=e^{ikr}/r \tag{8-1} \\ &\left(r=\sqrt{(x-x_\mathrm{Q})^2+(y-y_\mathrm{Q})^2+(z-z_\mathrm{Q})^2}~(\ne0)\right) \end{align} （ただし点 $\mathrm{Q}(x_\mathrm{Q},y_\mathrm{Q},z_\mathrm{Q})$ は $z<0$ の点）が回折するとして, これを $(7)$ 式に代入して計算している.
そこでここでは平面 $z=0$ に垂直な平面波 \begin{align} &u(x,y,z)=e^{ikz} \tag{8-2} \end{align} が入射するとして計算してみる（垂直でなくても計算はあまり変わらない）.
$(8\text{-}2)$ 式を $(7)$ 式に代入して, \begin{align} u(x_\mathrm{P},y_\mathrm{P},z_\mathrm{P}) &= \frac{1}{4\pi} \iint_{\mathrm{S}_2} dS \, \frac{e^{ikr}}{r}\left(\frac{1-ikr}{r^2}z_\mathrm{P} -\frac{\partial}{\partial z}\right)e^{ikz} \\ &= \frac{1}{4\pi} \iint_{z\,=\,0} dxdy \, \chi(x,y) \, \frac{e^{ikr}}{r}\left(\frac{1-ikr}{r^2}z_\mathrm{P} -ik\right)e^{ik0} \\ &= \frac{1}{4\pi} \iint_{z\,=\,0} dxdy \, \chi(x,y) \, e^{ikr}\frac{z_\mathrm{P}-ikr(r+z_\mathrm{P})}{r^3}. \tag{9} \end{align} ただし1行目から2行目の変形で, 開口部で1をとりそれ以外で0をとる定義関数 $\chi(x,y)$ を導入した.

$(9)$ 式において, 回折させる波を可視光であるとする.
光の波数 $k$ は非常に大きいので, $(9)$ 式は \begin{align} u(x_\mathrm{P},y_\mathrm{P},z_\mathrm{P}) &\approx \frac{k}{4\pi i} \iint_{z\,=\,0} dxdy \, \chi(x,y) \left[\frac{r+z_\mathrm{P}}{r^2} e^{ikr} \right], \tag{9$^\prime$} \\ r&=\sqrt{(x-x_\mathrm{P})^2+(y-y_\mathrm{P})^2+z_\mathrm{P}^2} \end{align} と近似できる.
この積分は大変なので, まず係数 $(r+z_\mathrm{P})/r^2$ について \begin{align} r\approx z_\mathrm{P}+\frac{1}{2}\frac{(x-x_\mathrm{P})^2+(y-y_\mathrm{P})^2}{z_\mathrm{P}} ~\Rightarrow~\frac{r+z_\mathrm{P}}{r^2} \approx \frac{2}{z_\mathrm{P}} \end{align} と近似する（これは回折像を狭い領域でしか考えないという近似）.
すると $(9^\prime)$ 式は \begin{align} &u(x_\mathrm{P},y_\mathrm{P},z_\mathrm{P}) \\ &\approx \frac{k}{2\pi iz_\mathrm{P}} \iint_{z\,=\,0} dxdy \, \chi(x,y) \exp{ik\left(z_\mathrm{P}+\frac{1}{2}\frac{(x-x_\mathrm{P})^2+(y-y_\mathrm{P})^2}{z_\mathrm{P}}\right)} \\ &= \frac{k}{2\pi iz_\mathrm{P}} e^{ikz_\mathrm{P}} \iint_{z\,=\,0} dxdy \, \chi(x,y) \exp{ik\frac{(x-x_\mathrm{P})^2+(y-y_\mathrm{P})^2}{2z_\mathrm{P}}} \tag{10-1} \end{align} となる.
これがFresnel回折の式である.
Fraunhofer回折では $(10\text{-}1)$ 式で更に \begin{align} (x-x_\mathrm{P})^2+(y-y_\mathrm{P})^2 \approx -2xx_\mathrm{P}+x_\mathrm{P}^2-2yy_\mathrm{P}+y_\mathrm{P}^2 \end{align} と近似して（これはスクリーンが遠くにあるとする近似）, \begin{align} &u(x_\mathrm{P},y_\mathrm{P},z_\mathrm{P}) \\ &\approx \frac{k}{2\pi iz_\mathrm{P}} e^{ikz_\mathrm{P}} \iint_{z\,=\,0} dxdy \, \chi(x,y) \exp{ik\frac{-2xx_\mathrm{P}+x_\mathrm{P}^2-2yy_\mathrm{P}+y_\mathrm{P}^2}{2z_\mathrm{P}}} \\ &= \frac{k}{2\pi iz_\mathrm{P}} e^{ik\left(z_\mathrm{P}+\frac{x_\mathrm{P}^2+y_\mathrm{P}^2}{2z_\mathrm{P}}\right)} \iint_{z\,=\,0} dxdy \, \chi(x,y) \exp{-ik\frac{xx_\mathrm{P}+yy_\mathrm{P}}{z_\mathrm{P}}} \tag{10-2} \end{align} が得られる.
この $(10\text{-}2)$ 式は定義関数のFourier変換となっており, 文献[2]などで回折像の計算にFourier変換が出てくるのはこれが理由である.

回折像の計算

最後に具体的な回折像をいくつか計算したい.
光の波数を $k=10^7~\mathrm{m}^{-1}$ （波長 $628~\mathrm{nm}$ の赤色光に相当）とし, 開口部は1辺の幅が $1~\mathrm{mm}$ の正方形であるとする（座標の原点はこの正方形の中心である）.
この条件下で, $z_\mathrm{P}=0.1~\mathrm{m}$ の場合と $z_\mathrm{P}=1.0~\mathrm{m}$ の場合の $u(x_\mathrm{P},y_\mathrm{P},z_\mathrm{P})$ の大きさを, $(10\text{-}1)$ 式から計算した.
同様の計算結果は[3]でも見ることができる.


図1 $z_\mathrm{P}=0.1~\mathrm{m}$ のときの回折像	図2 $z_\mathrm{P}=1.0~\mathrm{m}$ のときの回折像

図3 図1の原点付近の拡大（ $z_\mathrm{P}=0.1~\mathrm{m}$ ）	図4 図2の原点付近の拡大（ $z_\mathrm{P}=1.0~\mathrm{m}$ ）