1、推导磁场分量的更新方程
对磁场分量的时间导数应用中心差分公式,将时间中心点取为 $ n\Delta t $。以
$$
\frac{\partial H_x}{\partial t} = \frac{1}{\mu_x} \left( \frac{\partial E_y}{\partial z} - \frac{\partial E_z}{\partial y} - \sigma_{mx} H_x - M_{ix} \right)
$$
为例,基于场位置用有限差分近似可得:
$$
\frac{H_x^{n + 1/2}(i, j, k) - H_x^{n - 1/2}(i, j, k)}{\Delta t}
= \frac{1}{\mu_x(i, j, k)} \cdot \frac{E_y^n(i, j, k + 1) - E_y^n(i, j, k)}{\Delta z}
- \frac{1}{\mu_x(i, j, k)} \cdot \frac{E_z^n(i, j + 1, k) - E_z^n(i, j, k)}{\Delta y}
- \frac{\sigma_{mx}(i, j, k)}{\mu_x(i, j, k)} H_x^n(i, j, k)
- \frac{1}{\mu_x(i, j, k)} M_{ix}^n(i, j, k)
$$
经过一些处理,将未来项 $ H_x^{n + 1/2}(i, j, k) $ 移到左边,其他项移到右边,得到:
$$
H_x^{n + 1/2}(i, j, k) = \frac{2\mu_x(i, j, k) - \Delta t \sigma_{mx}(i, j, k)}{2\mu_x(i, j, k) + \Delta t \sigma_{mx}(i, j, k)} H_x^{n - 1/2}(i, j, k)
+ \frac{2\Delta t}{[2\mu_x(i, j, k) + \Delta t \sigma_{mx}(i, j, k)] \Delta z} [E_y^n(i, j, k + 1) - E_y^n(i, j, k)]
$$
对于其他磁场分量 $ H_y $ 和 $ H_z $ 也可按相同方法推导更新方程,如 $ H_y $ 的更新方程:
$$
H_y^{n + 1/2}(i, j, k) = H_y^{n - 1/2}(i, j, k)
+ \frac{2\Delta t}{\mu \Delta x \ln(\Delta x / a)} [E_z^n(i + 1, j, k) - E_z^n(i, j, k)]
- \frac{\Delta t}{\mu \Delta z} [E_x^n(i, j, k + 1) - E_x^n(i, j, k)]
$$
并可写成与通用更新方程相同形式:
$$
H_y^{n + 1/2}(i, j, k) = C_{hyh}(i, j, k) \cdot H_y^{n - 1/2}(i, j, k)
+ C_{hyez}(i, j, k) \cdot [E_z^n(i + 1, j, k) - E_z^n(i, j, k)]
+ C_{hyex}(i, j, k) \cdot [E_x^n(i, j, k + 1) - E_x^n(i, j, k)]
$$
其中:
- $ C_{hyh}(i, j, k) = 1 $
- $ C_{hyez}(i, j, k) = \frac{2\Delta t}{\mu_y(i, j, k) \Delta x \ln(\Delta x / a)} $
- $ C_{hyex}(i, j, k) = -\frac{\Delta t}{\mu_y(i, j, k) \Delta z} $