二维导热问题的ADI-TDMA算法

基本原理

对于多维非稳态导热问题，其时间项的离散格式有两种：显式格式和隐式格式。在显式格式中，下一时间步的温度场由上一时间步的已知量求出，不需要求解代数方程组，但其时间步的选取受稳定性条件的限制，故通常只能取很小的时间步。隐式格式不受稳定性限制，但其需要对五对角矩阵（二维）或七对角矩阵（三维）进行求解，故每一时间步的计算量较大。

相比之下，交替方向隐式方法^[1]（alternating direction implicit, ADI）是在显式与隐式之间折中的一种方法。对于二维问题，先将$x$方向的温度值按隐式处理，$y$方向的温度值按显式处理，这样就将二维问题转化为一系列并列的一维问题，可以直接使用三对角阵算法（TDMA）进行计算；然后再倒过来，将$y$方向隐式，$x$方向显式处理。这种方法被称为Peaceman-Rachford ADI格式。可以证明，该方法对于二维问题是绝对稳定的，故没有稳定性限制，可以取较大的时间步长，同时每一时间步的计算量较小，可大幅加快计算速度。

下面介绍Peaceman-Rachford ADI格式下二维非稳态导热方程的求解过程。

基于有限体积法的全隐离散方程可表示为：
$$a_P^0(T_P-T_P^0)=a_E(T_E-T_P)-a_W(T_P-T_W)+a_N(T_N-T_P)-a_S(T_P-T_S)+S\Delta x\Delta y$$

(1)沿$x$方向求解：

先取$\Delta t/2$的时间步长，对x方向做隐式处理，求解温度场的中间量$U$，则全隐离散方程可改成：
a P 0 ( U P − T P 0 ) = a E ( U E − U P ) − a W ( U P − U W ) + a N ( T N 0 − T P 0 ) − a S ( T P 0 − T S 0 ) + S Δ x Δ y a_P^0(U_P-T_P^0)=a_E(U_E-U_P)-a_W(U_P-U_W)+a_N(T_N^0-T_P^0)-a_S(T_P^0-T_S^0)+S\Delta x\Delta y aP0​(UP​−TP0​)=aE​(UE​−UP​)−aW​(UP​−UW​)+aN​(TN0​−TP0​)−a