Prerequisite 薛定谔方程
，类氢原子的定态波函数
本文使用原子单位制．我们希望求解氢原子在电场中的薛定谔方程
\begin{equation}
-\frac{1}{2m} \boldsymbol{\nabla}^2 \Psi + V( \boldsymbol{\mathbf{r}} , t)\Psi = \mathrm{i} \frac{\partial}{\partial{t}} \Psi
\end{equation}
虽然最直观的方法是使用直角坐标，但计算效率很低．实际中一般使用球坐标系，用球谐函数展开波函数．如果哈密顿算符是轴对称的（例如线偏振电场），那么只需要 $m = 0$ 的球谐函数，即勒让德多项式．
\begin{equation}
\Psi( \boldsymbol{\mathbf{r}} , t) = \frac{1}{r}\sum_{l,m} \psi_{l,m}(r) Y_{l,m}( \boldsymbol{\mathbf{r}} )
\end{equation}
其中 $\psi_{l,m}(r)$ 是约化径向波函数（scaled radial wave function），式中每一项叫做一个
分波（partial wave）．
在薛定谔方程中，一种简单的做法是使用以下势能表示库仑势能以及电磁场对电子的作用
\begin{equation}
V( \boldsymbol{\mathbf{r}} , t) = -\frac{Z}{r} - q \boldsymbol{\mathbf{E}} (t) \boldsymbol\cdot \boldsymbol{\mathbf{r}}
\end{equation}
在原子单位中电子电荷 $q = -1$，该势能使用了长度规范
，此外我们还可以使用速度规范速度规范
，也是等效的．在长度规范中，我们近似认为原子的长度远小于电磁波的波长，所以电磁波不随位置变化．
线偏振光
若我们取电场极化方向为 $ \hat{\boldsymbol{\mathbf{z}}} $，则角动量 $L_z$ 是一个守恒量．假设初始波函数关于 $ \hat{\boldsymbol{\mathbf{z}}} $ 轴对称，那么在波函数的整个演化（propagation）过程中，我们只需要 $m=0$ 的球谐函数展开波函数，即
\begin{equation}
\Psi( \boldsymbol{\mathbf{r}} , t) = \frac{1}{r}\sum_{l'} \psi_{l'}(r) Y_{l', 0}( \hat{\boldsymbol{\mathbf{r}}} )
\end{equation}
另外薛定谔方程中 $V( \boldsymbol{\mathbf{r}} , t) = \boldsymbol{\mathbf{E}} (t) \boldsymbol\cdot \boldsymbol{\mathbf{r}} = E(t) z$，进而可以用球谐函数表示（
eq. 9 ）
\begin{equation}
V( \boldsymbol{\mathbf{r}} , t) = \sqrt{\frac{4\pi}{3}} E(t) r \cdot Y_{1,0}( \hat{\boldsymbol{\mathbf{r}}} )
\end{equation}
以上两式代入薛定谔方程，再把每一项与 $Y_{l,0}( \hat{\boldsymbol{\mathbf{r}}} )$ 做内积（放入积分 $\int Y_{1,0}^*( \hat{\boldsymbol{\mathbf{r}}} ) \,\mathrm{d}{\Omega} $ 中）可得一系列耦合的方程
\begin{equation}
H_0 \psi_{l}(r) + \sum_{l' = 0}^{\infty} E(t)rF_{l, l'} \psi_{l'}(r) = \mathrm{i} \frac{\partial \psi_{l}}{\partial t} \quad (l=0,1,\dots)
\end{equation}
其中无场哈密顿算符为
\begin{equation}
H_0 = -\frac{1}{2m} \frac{\mathrm{d}^{2}}{\mathrm{d}{r}^{2}} -\frac{Z}{r} + \frac{l(l+1)}{2mr^2}
\end{equation}
矩阵 $ \boldsymbol{\mathbf{F}} $ 就是跃迁偶极子矩阵的角向积分（
eq. 3 ）
\begin{equation}
\begin{aligned}
F_{l,l'} &= \sqrt{\frac{4\pi}{3}} \left\langle Y_{l,0} \middle| Y_{1,0} \middle| Y_{l',0} \right\rangle \\
&= \sqrt{(2l+1)(2l'+1)} \begin{pmatrix}l & 1 & l'\\ 0 & 0 & 0\end{pmatrix} ^2\\
&= \sqrt{\frac{2l+1}{2l'+1}} \begin{bmatrix}l & 1 & l'\\ 0 & 0 & 0\end{bmatrix} ^2
\end{aligned}
\end{equation}
可见，当没有外场的时候每一个项（即每一个分波）都可以独立演化，而电场将不同的分波耦合起来．根据氢原子的选择定则 $\Delta l = \pm 1$（
eq. 5 ），矩阵 $ \boldsymbol{\mathbf{F}} $ 中除了两条副对角线上的元都为零．另外由
eq. 2 易得 $ \boldsymbol{\mathbf{F}} $ 是一个对称矩阵．
任意偏振光
将所有 $(l,m)$ 按某种顺序排列，例如
\begin{equation}
(0, 0),\ (1,-1),\ (1,0),\ (1,1),\ (2,-2), \dots
\end{equation}
并将他们编号为 $\lambda = 1,2, \dots$，那么可以把
eq. 2 记为
\begin{equation}
\Psi( \boldsymbol{\mathbf{r}} , t) = \frac{1}{r}\sum_\lambda \psi_\lambda(r) Y_\lambda( \boldsymbol{\mathbf{r}} )
\end{equation}
势能函数为 $V( \boldsymbol{\mathbf{r}} , t) = \boldsymbol{\mathbf{E}} (t) \boldsymbol\cdot \boldsymbol{\mathbf{r}} = E_x x + E_y y + E_z z$．
eq. 6 的径向耦合方程拓展为
\begin{equation}
\begin{aligned}
H_0 \psi_{\lambda}(r) + \sum_{\lambda'} r \left[E_x(t) F_{\lambda, \lambda'}^{x} + E_y(t) F_{\lambda, \lambda'}^{y} + E_z(t) F_{\lambda, \lambda'}^{z} \right] \psi_{\lambda'}(r) \\
= \mathrm{i} \frac{\partial \psi_{\lambda}}{\partial t} \quad (\lambda=0,1,\dots)
\end{aligned}
\end{equation}
三个耦合矩阵分别为
\begin{equation}
\begin{aligned}
F_{\lambda,\lambda'}^{x} = \sqrt{\frac{2\pi}{3}} E_x \left( \left\langle Y_{l,m} \middle| Y_{1,-1} \middle| Y_{l',m'} \right\rangle - \left\langle Y_{l,m} \middle| Y_{1,1} \middle| Y_{l',m'} \right\rangle \right)
\end{aligned}
\end{equation}
\begin{equation}
\begin{aligned}
F_{\lambda,\lambda'}^{y} = \mathrm{i} \sqrt{\frac{2\pi}{3}} E_y \left( \left\langle Y_{l,m} \middle| Y_{1,-1} \middle| Y_{l',m'} \right\rangle + \left\langle Y_{l,m} \middle| Y_{1,1} \middle| Y_{l',m'} \right\rangle \right)
\end{aligned}
\end{equation}
\begin{equation}
\begin{aligned}
F_{\lambda,\lambda'}^{z} = \sqrt{\frac{4\pi}{3}} E_z \left\langle Y_{l,m} \middle| Y_{1,0} \middle| Y_{l',m'} \right\rangle
\end{aligned}
\end{equation}
其中（
eq. 3 ）
\begin{equation}
\left\langle Y_{l,m} \middle| Y_{1,m_1} \middle| Y_{l',m'} \right\rangle = (-1)^m\sqrt{\frac{3(2l+1)(2l'+1)}{4\pi}} \begin{pmatrix}l & 1 & l'\\ 0 & 0 & 0\end{pmatrix} \begin{pmatrix}l & 1 & l'\\ -m & m_1 & m'\end{pmatrix}
\end{equation}
算符拆分
在实际的程序中，我们可以把演化子 $ \exp\left(- \mathrm{i} H \Delta t\right) $ 拆成两项．虽然这么做会引入一定的误差，但是却大大提高了效率
\begin{equation}
\exp\left(- \mathrm{i} H \Delta t\right) = \exp\left(- \mathrm{i} H_0 \frac{\Delta t}{2}\right) \exp\left(- \mathrm{i} F \Delta t\right) \exp\left(- \mathrm{i} H_0 \frac{\Delta t}{2}\right) + \mathcal{O}\left(\Delta t^3 \right)
\end{equation}
也就是说，在每个时间步长 $\Delta t$ 中，我们可以把波函数先根据
\begin{equation}
H_0 \psi_{l}(r) = \mathrm{i} \frac{\partial \psi_{l}}{\partial t}
\end{equation}
演化 $\Delta t/2$，其中 $t_{mid}$ 取这段时间的中点．再对每个 $ \boldsymbol{\mathbf{r}} $ 根据
\begin{equation}
\sum_{l'} E(t_{mid})rF_{l, l'} \psi_{l'}(r) = \mathrm{i} \frac{\partial \psi_{l}}{\partial t}
\end{equation}
演化 $\Delta t$．最后再根据
eq. 17 演化 $\Delta t/2$．
理论上也可以把演化子划分成两项 $ \exp\left(- \mathrm{i} F \Delta t\right) \exp\left(- \mathrm{i} H_0 {\Delta t}\right) $ 但据说eq. 16 更精确．至于相邻两步之间产生的 $ \exp\left(- \mathrm{i} H_0 {\Delta t}/{2}\right) \exp\left(- \mathrm{i} H_0 {\Delta t}/{2}\right) $ 是否可以合并为 $ \exp\left(- \mathrm{i} H_0 {\Delta t}\right) $，取决于所使用的算法这么做以后是否会引入更大误差．
网格和演化算法
Prerequisite Crank-Nicolson 算法（一维）
未完成：以下内容应该放在一维薛定谔方程里面讲解
可以使用二维数组储存波函数，每一列（或行）是一个分波的 $\psi_l(r)$．径向网格可以使用等间距网格，但 FEDVR 网格效率要更高．
演化可以并使用 Crank-Nicolson 算法 Crank-Nicolson 算法（一维）演化．但是 Lanczos 算法效率更高．