氢原子的跃迁偶极子矩阵元（束缚态到连续态）

                     

贡献者： addis

  • 本文处于草稿阶段。
预备知识　跃迁概率（含时微扰）

　　 本文使用原子单位制。

图
图 1：$ \left\lvert \left\langle k_z \middle| z \middle| 1,0,0 \right\rangle \right\rvert $，实线、虚线：平面波近似，点：无近似（代码：plot_hydrogen_trans_dipole.m）

1. 长度规范（平面波近似）

　　 归一化的平面波和归一化的氢原子基态为（用平面波近似末态库仑函数）

\begin{equation} \left\lvert \boldsymbol{\mathbf{k}} \right\rangle = (2\pi)^{-3/2} \mathrm{e} ^{ \mathrm{i} \boldsymbol{\mathbf{k}} \boldsymbol\cdot \boldsymbol{\mathbf{r}} }~, \qquad \left\lvert 0 \right\rangle = \frac{1}{\sqrt{\pi}} \mathrm{e} ^{-r}~. \end{equation}
长度规范下的跃迁偶极子，可以在极坐标系中积分（令 $ \hat{\boldsymbol{\mathbf{k}}} $ 方向为极轴，$ \hat{\boldsymbol{\mathbf{r}}} $ 与其夹角为 $\theta$）
\begin{equation} \left\langle \boldsymbol{\mathbf{k}} \middle| \boldsymbol{\mathbf{r}} \middle| 0 \right\rangle = \frac{ \hat{\boldsymbol{\mathbf{k}}} }{(2\pi)^{3/2}\sqrt{\pi}} \int_0^{+\infty} \int_0^\pi \mathrm{e} ^{-r} \mathrm{e} ^{- \mathrm{i} k r \cos\theta} r \cos\theta \cdot 2\pi r^2 \sin\theta \,\mathrm{d}{\theta} \,\mathrm{d}{r} ~ \end{equation}
换元，令 $u = \cos\theta$，得1
\begin{equation} \begin{aligned} \left\langle \boldsymbol{\mathbf{k}} \middle| \boldsymbol{\mathbf{r}} \middle| 0 \right\rangle &= \frac{1}{\sqrt 2 \pi} \hat{\boldsymbol{\mathbf{k}}} \int_0^{+\infty} r^3 \mathrm{e} ^{-r} \int_{-1}^1 \mathrm{e} ^{- \mathrm{i} k r u} u \,\mathrm{d}{u} \cdot \,\mathrm{d}{r} \\ &= \mathrm{i} \frac{\sqrt2 \hat{\boldsymbol{\mathbf{k}}} }{k\pi} \int_0^{+\infty} r^2 \mathrm{e} ^{-r} \left[ \cos\left(kr\right) - \frac{1}{kr} \sin\left(kr\right) \right] \,\mathrm{d}{r} \\ &= - \mathrm{i} \frac{8 \sqrt2}{\pi} \frac{ \boldsymbol{\mathbf{k}} }{(k^2+1)^3}~, \end{aligned} \end{equation}
注意这是一个纯虚数。Matlab 代码如下：

代码 1：hydrogen_trans_dipole_plane_approx_z.m
% hydrogen transition dipole, approximate Coulomb plane wave with plane wave
% since middle and right parts are symmetric,
%     assuming \bvec k in z+ direction
% <\bvec k|\bvec r|n,0,0> = [0, 0, <\bvec k|z|n,0,0>]
% output numerical integration and analytical result (eq_HyIon2_1)
function [dipole_z, dipole_analy_z] = ...
    hydrogen_trans_dipole_plane_approx_z(kz, Z, n)
% === params ===
rmax = 20; Nr = 200; Nth = 100;
% ==============
k = [0, 0, kz];
r = linspace(0, rmax, Nr); dr = rmax / (Nr-1);
th = linspace(0, pi, Nth); dth = pi / (Nth-1);
ph = 0;
[R, Th] = ndgrid(r, th);
zz = R .* cos(Th);
Psi_n = hydrogen_Psi(Z, n, 0, 0, r', th, ph);
% <\bvec k|\bvec r|0>
Psi_k = 1/(2*pi)^(3/2) * exp(1i*kz.*zz);
dipole_z = sum(sum(conj(Psi_k).*zz.*Psi_n .*R.^2.*sin(Th))) ...
    * dr * dth * 2*pi;
% eq_HyIon2_1
if n == 1
    dipole_analy_z = -1i*(8*sqrt(2))/pi * kz / (dot(k,k) + 1)^3;
end
end

2. 速度规范

　　 注意一阶微扰理论中的初态和末态波函数都是无微扰（无外场）情况下的，与规范无关。要计算 $ \left\langle \boldsymbol{\mathbf{k}} \middle| \boldsymbol\nabla \middle| 0 \right\rangle $，先看积分

\begin{equation} \int \exp\left(- \mathrm{i} \boldsymbol{\mathbf{k}} \boldsymbol\cdot \boldsymbol{\mathbf{r}} \right) \exp\left(-r\right) \,\mathrm{d}^{3}{r} = \frac{8\pi }{(k^2 + 1)^2}~. \end{equation}
使用算符 $ \boldsymbol\nabla $ 的反厄米性得
\begin{equation} \begin{aligned} &\int \exp\left(- \mathrm{i} \boldsymbol{\mathbf{k}} \boldsymbol\cdot \boldsymbol{\mathbf{r}} \right) \boldsymbol\nabla \exp\left(-r\right) \,\mathrm{d}^{3}{r} = -\int [ \boldsymbol\nabla \exp\left( \mathrm{i} \boldsymbol{\mathbf{k}} \boldsymbol\cdot \boldsymbol{\mathbf{r}} \right) ]^* \exp\left(-r\right) \,\mathrm{d}^{3}{r} \\ &= \mathrm{i} \boldsymbol{\mathbf{k}} \int \exp\left(- \mathrm{i} \boldsymbol{\mathbf{k}} \boldsymbol\cdot \boldsymbol{\mathbf{r}} \right) \exp\left(-r\right) \,\mathrm{d}^{3}{r} = \mathrm{i} \frac{8 \pi \boldsymbol{\mathbf{k}} }{(k^2 + 1)^2}~, \end{aligned} \end{equation}
乘以归一化系数得
\begin{equation} \left\langle \boldsymbol{\mathbf{k}} \middle| \boldsymbol\nabla \middle| 0 \right\rangle = \mathrm{i} \frac{2\sqrt{2}}{\pi}\frac{ \boldsymbol{\mathbf{k}} }{(k^2 + 1)^2}~. \end{equation}
该式代入式 11 （$q = -1$）得微分截面为
\begin{equation} \frac{\partial \sigma}{\partial \Omega} = \frac{32}{mc\omega} \frac{k( \boldsymbol{\mathbf{k}} \boldsymbol\cdot \hat{\boldsymbol{\mathbf{e}}} )^2}{(k^2 + 1)^4} = \frac{64}{mc} \frac{k( \boldsymbol{\mathbf{k}} \boldsymbol\cdot \hat{\boldsymbol{\mathbf{e}}} )^2}{(k^2 + 1)^5}~. \end{equation}
对于质子数为 $Z$ 类氢原子有
\begin{equation} \frac{\partial \sigma}{\partial \Omega} = \frac{32 Z^5}{mc\omega} \frac{k( \boldsymbol{\mathbf{k}} \boldsymbol\cdot \hat{\boldsymbol{\mathbf{e}}} )^2}{(k^2 + Z^2)^4}~. \end{equation}

3. 两种规范对比

　　 如果 $ \left\lvert \boldsymbol{\mathbf{k}} \right\rangle $ 是库仑函数（能量本征态）应该有（式 2 ）

\begin{equation} \left\langle \boldsymbol{\mathbf{k}} \middle| \boldsymbol{\mathbf{r}} \middle| 0 \right\rangle = -\frac{ \left\langle \boldsymbol{\mathbf{k}} \middle| \boldsymbol\nabla \middle| 0 \right\rangle }{m\omega_{k0}}~. \end{equation}
其中 $\omega_{k0} = k^2/2 + 1/2$，但实际上式 3 和式 6 满足
\begin{equation} \left\langle \boldsymbol{\mathbf{k}} \middle| \boldsymbol{\mathbf{r}} \middle| 0 \right\rangle = -2\frac{ \left\langle \boldsymbol{\mathbf{k}} \middle| \boldsymbol\nabla \middle| 0 \right\rangle }{m\omega_{k0}}~. \end{equation}
这说明在使用平面波近似库仑函数时，长度规范的 transition amplitude 恰好是速度规范的 2 倍，截面是四倍（待求证）。

　　 教材中推导微分截面一般使用速度规范，因为速度规范的结果与实验吻合更好。

4. 使用库仑平面波

长度规范

　　 理论上若把上面的平面波换成库仑平面波（库仑势能中的精确散射态），那么理论上用不同的规范结果是一样的。先看球面波投影（类比式 2 ）：

\begin{equation} \begin{aligned} &\quad \left\langle C_{l',m'}(k) \middle| r\cos\theta \middle| \psi_{n,l,m} \right\rangle \\ &= \delta_{m,m'}(\delta_{l+1,l'}\mathcal C_{l,m} + \delta_{l,l'+1}\mathcal{C}_{l',m'}) \sqrt{\frac{2}{\pi}}\int_0^\infty F_{l'}(\eta, kr)( \hat{\boldsymbol{\mathbf{r}}} ) R_{n,l}(r)r^2 \,\mathrm{d}{r} ~. \end{aligned} \end{equation}
库仑平面波投影（式 11 ）：
\begin{equation} \langle{\psi_{ \boldsymbol{\mathbf{k}} }^{(-)}}|{z}|{\psi_{n,l,m}}\rangle = \sum_{l'} \frac{ \mathrm{i} ^{-l'}}{k} \mathrm{e} ^{ \mathrm{i} \sigma_{l'}} Y_{l',m}( \hat{\boldsymbol{\mathbf{k}}} ) \int \left\langle C_{l',m}(k) \middle| z \middle| n,l,m \right\rangle ~. \end{equation}
其中 $\eta = -Z/k$。对 $l=0$ 有

　　 笔者没有见过该积分的解析解，长度规范的 Matlab 数值积分代码如下：

代码 2：hydrogen_trans_dipole_sph.m
% <C_{l1,m1}(k1)|z|n2,l2,m2>
% eq_HyIon2_3
function ret = hydrogen_trans_dipole_sph(Z, k1, l1, m1, ...
                                          n2, l2, m2, r_max)
if m1 ~= m2 || abs(l1 - l2) ~= 1
    ret = 0; return;
end
eta = -Z/k1;
f = @(r) coulombF_sym(l1, eta, k1*r) .* hydrogen_Rnl(Z, n2, l2, r) .* r.^2;
I_r = sqrt(2/pi)*integral(f, 0, r_max);
l_ = min(l1, l2);
I_ang = sqrt(((l_+1)^2-m2^2)/((2*l_+1)*(2*l_+3)));
ret = I_r * I_ang;
end
代码 3：hydrogen_trans_dipole_plane.m
% <\bvec k|z|n,l,m>
% eq_HyIon2_4
% 未验证
function ret = hydrogen_trans_dipole_plane(...
    Z, k1, th1, ph1, n2, l2, m2, r_max)
if n2 <= 0 || l2 < 0 || l2 >= n2 || abs(m2) > l2
    error('illegal n2,l2,m2');
end
ret = 0;
for l1 = [l2-1, l2+1]
    m1 = m2;
    if l1 < 0 || abs(m1) > l1
        continue;
    end
    mel_sph = hydrogen_trans_dipole_sph(...
        Z, k1, l1, m1, n2, l2, m2, r_max);
    ret = ret + (1i^(-l1))/k1 * exp(1i*coulomb_sigma(l1, -Z/k1))...
        * SphHarm(l1,m1,th1,ph1) * mel_sph;
end
end

速度规范

　　

未完成：速度规范数值积分，验证和长度规范一样

画图代码

代码 4：plot_hydrogen_trans_dipole.m
% plot hydrogen transition dipole
% compare different methods

% == params ==
Z = 1;
n = 1; l = 0; m = 0;
kmin = 0.01; kmax = 3; Nk = 51;
r_max = 100;
% ============

kz = linspace(kmin, kmax, Nk);
dipole_z0_appr = zeros(1, Nk);
dipole_z0_ana = dipole_z0;
dipole_z = dipole_z0;

% plane wave approx
for i = 1:Nk
    [dipole_z0_appr(i), dipole_z0_ana(i)]...
        = hydrogen_trans_dipole_plane_approx_z(kz(i), Z, n);
end

% Coulomb wave
parfor i = 1:Nk
    disp(i);
    dipole_z(i) = hydrogen_trans_dipole_plane(...
        Z, kz(i), 0, 0, n, l, m, r_max);
end

figure; plot(kz, abs(dipole_z0_appr)); hold on;
plot(kz, abs(dipole_z0_ana), '--');
plot(kz, abs(dipole_z), '.-k');

legend({'plane wave approx', ...
    'plane wave approx (analytical)', 'exact (Coulomb plane wave)'});
grid on; axis([0, 3, 0, 2.5]);

xlabel 'k [au]'; ylabel 'abs dipole';

1. ^ 最后一步可通过 Wolfram Alpha 或 Mathematica 获得。

                     

© 小时科技 保留一切权利