2.4　橢圓方程和穩定狀態問題

考慮橢圓式（2.3.18），利用中央差分方案進行求導，結果得到一個二階的精度方案為

將式（2.3.18）降低成一個代數方程組為

通過替換I=j+（i-1）J；I=1，（I，J），該方程組可以轉換為一個單一矩陣方程為

其中A是N×N的系數矩陣，X是一個長度為N的向量，它包含未知變量φ從i=1，j=1到i=I，j=J所有N=IJ個格網點處的值，B是相同長度的向量，它包含式（2.4.2）右側的應力項。矩陣A是大而稀疏的，主要是有邊緣的三對角矩陣。

許多CFD問題最終都降低為矩陣方程的解，例如式（2.4.3）。因此一種有效的“求逆”方法顯得非常重要，尤其是矩陣非常大的時候（對于一個包含100×100格網的問題來說，矩陣尺寸為104×104），然而這種矩陣通常是稀疏矩陣，它的求解與拉普拉斯方程求解一樣，因此可以對其進行簡化。

對代數式（2.4.3）求解的兩種主要方法是直接法和迭代法。直接法例如Gauss消元法和Gauss-Jordan消元法是計算密集型的，矩陣大而密集時，舍入誤差也大，它們很少被用于大的三維CFD問題中。在矩陣稀疏且為帶狀矩陣的情形中，廣泛使用的方法是特別版的Gauss消元法。但是大多數情況下，所選擇的方法通常是迭代法，這種方法從一個初始估計開始，使用有效的加速手段計算到真實解。

2.4.1　直接解法

對式（2.4.3）直接求解的一種可用方法是傳統的Gauss消元法或Gauss-Jordan消元法，這種消元過程使用系統標準化和系統減法將矩陣的所有子元素和超對角元素消除?？梢詫ψ訉窃剡M行高斯前向消元然后使用回帶法，Press等（1992）對這些內容進行了討論。消元過程要求每一步除以一個元素，稱為樞軸點。為了使該方法盡可能精確，必須選擇最大可能元素作為樞軸點，這通常要求在每一步中對行和列重排。將行進行交換后使最大元素轉換到期望位置并不會改變最終的解，而將列交換則會使解變得雜亂。行交換稱為部分軸元法，它執行起來較為簡單，而列交換稱為全部軸元法，它需要額外的計算，因為必須對雜亂的解重新組合。通常使用的方法是將原始方程標準化并使其最大系數為1的情況下選擇可能是最大的元素作為樞軸點，這種方法稱為隱式旋轉。然而，這些方法要求方程右側的向量B為已知。對全部軸元法的Gauss-Jordan消元法程序可以參考Press等（1992）。

2.4.1.1　LU分解法

LU分解法為直接法中的首選，由于它所得三角集的解的簡單性質：替換。不論矩陣A是什么結構，無論它是稀疏的（大量元素為0）或者密集的（極少數或者沒有0元素）或者是帶狀的（只有對角線上或其附近的元素為非零值），LU分解都要求矩陣A被分解成下三角陣L和上三角陣U：A=L·U，那么A·X=（L·U）X=L（U·X）=B，因此可以首先用簡單的前向替換對L·Y=B進行求解，使得從L的結構可以找到向量Y，然后用Y求解U·X=Y，使可以用簡單的后向替換從U的結構中得到X，也就是求解U·X=L-1B，之后的任務就是有效地進行LU分解。Press等（1992）列出了這樣的算法，令

要得到系數αij和βij，寫出N+N2方程如下

通過采用特殊順序重組很容易對這些方程進行求解，使克勞特算法（Press等，1992）為

對每一個j=1，2，…，N ：

計算下一個j

繞軸旋轉就是將一個矩陣的行和列重新排列，以選出最大可能的矩陣元素對矩陣的其他元素進行標準化，這對克勞特算法很重要。然而，只需對行進行交換的部分軸元法就已經足夠，因為在矩陣求逆的過程中，對矩陣的行進行交換不會對解產生影響，它執行起來比全部軸元法簡單，全部軸元法會使解變得混亂［對于該問題的討論可以參考Press等（1992）］。如果L和U的元素事先已知，可以用前向替換對向量Y進行求解如下

然后可以用后向替換求未知的向量X

Press等（1992）和Fletcher（1998a）都列出了進行LU分解的FORTRAN源代碼（分別為LUDCMP和FACT），還列出了前向替換和后向替換（分別為LUBKSB和SOLVE）。Press等使用更為高效的克勞特算法進行LU分解。這些實現都是對密集矩陣而言的，它們對于稀疏矩陣來說在存儲空間和計算時間上都存在浪費問題，因為零元素也會被存儲和計算。有專門為稀疏矩陣而設計的方法，它使得算術運算量和存儲需求最小化。Press等（1992）為三對角矩陣提出了一種特別簡單的最優化方法（TRIDAG），該方法是一種將分解和替換相結合的方法（該方法大部分情況下都不必使用繞軸旋轉）。

LU分解的一個優勢是，它不要求方程右側的向量B為已知，因此對矩陣進行分解后，替換代碼可以利用它來有效地對不同B向量進行重復求解，這也意味著對于很大的方程組來說，由舍入誤差積累引起的嚴重精度喪失可以被修正，迭代改進對機器精度進行存儲（Press等，1992）。令Xc為修正解，Ww為誤差解，δX為誤差。因為A·Xc=B，乘以A和Xw后得不到B，而是得到B+δB。因此，只通過替換就可以對A·δX=δB=（A·Xw-B）進行求解，從而得到修正的δX，而不需要再次對矩陣A進行分解，Press等（1992）為此提出了一個算法。

2.4.1.2　 Thomas算法

對單一未知變量的系統使用三點有限差分方法或者線性內插的有限元素法，結果得到具有三對角結構的矩陣，用分步法（本章2.5節）對多維問題進行求解也會產生三對角系統。線松弛迭代法也要求三對角矩陣的解（本章2.4.2節），因此三對角矩陣的有效算法對CFD很重要，一個特別有效的（但由于其具有遞歸性，矢量化和并行化比較困難）算法是Thomas算法，該算法由兩個步驟組成：①是前向搜索，用矩陣主對角線上的單位元素將矩陣轉化為上三角形式（在這個過程中向量B被修改為B′）②后向搜索，用簡單回帶法得到未知向量X，優化的Gauss消元法必須使用遞歸關系進行求解。

令aixi-1+bixi+cixi+1=di為三對角矩陣方程。進行擴展后得

矩陣A是非常稀疏的，只有主對角、次對角和超對角線上有值填充，不需要對整個N×N矩陣進行存儲，只需要存儲ai、bi和ci三個長度向量。一般的LU分解需要進行計算的次數是O（N）次，而不是O（N3/3）次。

Thomas算法使用的回歸關系式為

其中系數αi和βi在前向搜索的過程中確定。將式（2.4.9）中的i改為i+1并將其替換為第i個方程，aixi-1+bixi+cixi+1=di，進行重新排列后得到

從而得到

這些是遞歸關系，αi和βi的值（可以由第一個方程b1x1+c1x2=d1進行確定）為

由此完成了前向搜索，遞歸關系式（2.4.9）中的所有系數都是已知的。從元素N進行的后向搜索提供了所要求的解，將式（2.4.9）帶入到最后一個方程bNxN-1+bNxN=dN，得到

從式（2.4.14）中得到xN后，當i=N-1，N-2，…，1時，可以用式（2.4.9）確定xi，因為式（2.4.11）到式（2.4.12）中的系數都是已知的。

和所有的Gauss消元法一樣，矩陣的病態條件會引起嚴重的舍入誤差，為了避免這種情況發生，必須使|bi|＞|ai|+|ci|。對Thomas算法進行編碼非常簡單，Fletcher（1988a）列出了一個事例代碼。

對包含單一未知變量的系統使用高階精度方案時，會生成一個矩陣方程，該矩陣的帶寬比三對角矩陣的帶寬大。然而，Thomas算法可以被泛化從而處理這些仍然具有稀疏性的矩陣。Fletcher（1988a）對五對角矩陣的泛化Thomas算法進行T描述和討論，如果未知變量多于一個，方程不止一個，那么求解過程通常會產生一個分塊三對角矩陣，是M×M的子矩陣形成了三對角矩陣的對角線、次對角線和超對角線元素，而不是由標量形成。也可以將Thomas算法泛化到這種情形（Fletcher，1988a），它要求進行O（NM3）次運算。

在簡單均勻矩形格網上對橢圓方程進行求解，矩陣元素是常數，在這種情形中，直接解法采用一個稱為循環約化的過程，然后用傅里葉級數表示，這種方法是很高效的，因為它使用了現代的快速傅里葉轉換（FFT）方法，不過這不是本書的討論范圍，這些方法的討論讀者可以參考Fletcher（1988a）和Press等（1992）。

2.4.2　迭代解法和松弛法

對于CFD和海洋模擬問題中的巨大稀疏矩陣特征來說，最常用的方法是迭代方法或松弛法。它將式（2.4.3）的矩陣A分裂成兩部分，其中一部分很容易求逆（E），另一部分為殘余部分（F）：A=E+F。那么方程變為

松弛法包括對如下方式進行連續迭代

相當于

迭代從初始估計X1。Rn是第n個迭代步驟的殘差向量，當解收斂時它的范數趨于0。M=E-1F稱為迭代矩陣，收斂速度取決于該矩陣的性質，尤其是它的譜半徑σ，它是矩陣的最大特征值λi的大小。特征值是代數方程的根，通過對行列式進行歸零處理得到

其中I是單位矩陣。然后對實數對稱矩陣，所有的特征值都是實數。

譜半徑的大小M決定收斂特性。如果σ＜1，則可以保證收斂性，這可以用誤差en=Xt-Xn進行證明，該誤差由方程en+1=Men=Mn·e0決定，因此在初始估計值e0處誤差接近0，如果Mn接近0，則解是收斂的，如果譜半徑小于1，則這一點可以保證，該條件通常轉變為對矩陣A為對角占優矩陣的要求。譜半徑的值越小，收斂速度越快。

2.4.2.1　Jacobi和Gauss-Seidel法

為了達到說明的目的，對拉普拉斯方程進行考慮，簡單起見，假設x方向和y方向上的格網大小相同，那么式（2.4.2）變為

經典的點Jacobi（或Richardson）方法由迭代方案組成如下：

如果對應的矩陣A可分解為一個嚴格的下三角矩陣L（只有對角線下方的元素非零）、一個對角矩陣D（只有對角線上的元素非零）和一個嚴格的上三角矩陣U（只有對角線上方的元素非零），那么Jacobi方法等價于

或者等價于

對于Dirichlet邊界條件，迭代矩陣M=D-1（L+U）的譜半徑通常情況下格網大小不相等，則有

其中Δx和Δy是格網大小。通常情況下該值雖小于1但很接近于1，因此，實際上Jacobi方法的收斂速度較低。當I，J＞＞1時，假設Δx=Δy，σJ～1-（π2/2N），其中N=IJ為格網點的數量。收斂速度定義為將誤差降低一個數量級所需的迭代次數，為-（logσJ）-1，因此它約為0.5N。因為每次迭代需要5N次運算，要使誤差減少一個數量級需要2.5N2次運算。

加權Jacobi方法是式（2.4.20）的簡單變體，即

其中可以選擇ω使收斂在一定程度上得到加速。

點Gauss-Seidel（or Liebmann）法是Jacobi方法的一個簡單變形如下

等價于

或

迭代矩陣M=（L+D）-1U的譜半徑等于Jacobi譜半徑的平方，因此，其收斂速度要快兩倍，但是仍然較慢而不實用。Jacobi方法和Gauss-Seidel方法，收斂性都可以保證。

松弛法可以是上述的點松弛法，其中在每次迭代循環中對一個格網點進行更新。利用Gauss-Seidel松弛法作為例子，在式（2.4.25）的左側唯一未知的點是（i，j），因此這是一種點迭代方法。假如認為格網點分布于一個棋盤上，很容易看出紅色（偶數）點處的解只取決于黑色（奇數）點，反之亦然。因此，紅色點集和黑色點集的交替更新是可能的，在每次松弛循環中的每次半搜索中對偶數點或奇數點進行更新，這一紅黑點松弛法（或奇偶排序法）是一種最為常用的方法。

松弛法也可是線性松弛，其中每次對整行或整列的點進行松弛。對式（2.4.25）稍作修改得到

對每一個j處的所有i都可以獲得一個解，也就是對j列的每一行所得三對角矩陣逐行求解，同樣地，也可以進行逐列求解。這就是線Gauss-Seidel松弛法，它比點Gauss-Seidel法速度快兩倍，收斂速度得到很大改善，而這是通過在每一步中引入點的解的信息，而不是相鄰點的信息。

2.4.2.2　逐次超松弛（SOR）法

如果認識到這些松弛方法的本質是對之前的解添加一個校正這一事實的話，就很容易理解逐次超松弛（SOR）法的工作原理。SOR一直是進行橢圓求解的主要方法，直到最近才被現代多格網和共軛梯度法所取代。當使用點Gauss-Seidel方法時SOR的優勢便顯露出來，每步迭代中不是對校正項，而是對下方的式子進行

這等價于

或者

這是Gauss-Seidel值和之前的值的加權平均。α是松弛參數，收斂性要求0＜α＜2。對于超松弛法，要使它比Gauss-Seidel方法收斂快，則需滿足α＞1。松弛參數的最優選擇是此時譜半徑為（Press等，1992）。需要注意，這僅是漸進有效的，只有在O（J）次迭代后得到。在初始階段，誤差逐漸增加。對于Jacobi和Gauss-Seidel方法采用的Dirichlet例，當Δx=Δy，I=J=N1/2時，將誤差降低一個數量級所需要的迭代次數約為0.33N1/2。顯然，該參數的選擇很關鍵，因為在該最優值周圍的一個狹窄窗口內收斂速度達到最好，離開該窗口后立即變差。不過，對一般問題而言，很難事先決定最優值，因此進行一些實驗是很必要的。點松弛和線松弛，Jacobi和Gauss-Seidel松弛都可以是超松弛的，但對點Jacobi方法來說，松弛參數的最優值是1，也就是常規Jacobi方法已經是最優的。

每次迭代中α可以選擇一個不同的值，因為上述最優值在開始時不是最優的。當α保持不變時，SOR方法稱為靜態超松弛法，而當α在每次迭代中不斷改變時，稱為非靜態超松弛法。在迭代過程中使用一個特定的α序列使收斂速度得到極大增加，這種方法稱為Chebyshev加速方法，Press等（1992）對該方法進行了描述，他在每次半搜索中改變α以保證誤差指標單調減小。α值序列從1開始，以最優值結束：SOR的輕微改變形式稱為對稱SOR，式（2.4.17）中的矩陣E為（αL+D）D-1（U+D）/（2α-α2），它允許SOR可以被Chebyshev加速法進行加速（Fletcher，1988a）。SOR也可以是點SOR或線SOR，在最優松弛參數相同的情況下后者比前者快倍（Hirsch，1988）。

一個重要的問題是緩和，在該過程中可以在高性能計算機上對松弛方案進行向量化和并行化處理，比如常規的線Jacobi松弛法的向量化比線Ganss-Seidel SOR方法的向量化簡單，因此，盡管它效率稍差，但也可能綜合效果會更好。就像點松弛方法中的紅黑方案那樣，可以在連續搜索過程中在交替的線上執行Jacobi線松弛法，這稱為斑馬線松弛法，可以對它進行逐行或逐列執行，斑馬線松弛法非常適合于現代計算機。

2.4.3　預條件共軛E梯度法

當使用SOR時，渴望使用一些方法使收斂得到加速。如前所述，收斂的速度取決于迭代矩陣的特征值λ，每個特征值表示一個因子，它使一次迭代中誤差的特定本征模減小。顯然，模的最大特征值決定收斂的速度。比值|λmax|/|λmin|是矩陣的條件數，表示特征值譜的帶寬，如果它比1大得多，則特征值接近λmax占主導并降低收斂速度，對收斂進行改進的一種方法是用一個預處理矩陣乘以預條件矩陣，從而預條件矩陣的特征值譜的條件數較小，而在λ為λmax時最小。

另一種收斂加速方法是共軛梯度法，它要對如下的函數進行最小化

相當于令其梯度等于0

它的求解與式（2.4.3）一樣。因此，方法簡化為對n維空間中的函數進行最小化，這種方法稱為共軛梯度法，它只適用于對稱的正定矩陣（對所有的z來說，zTAz＞0），這是對橢圓方程求解的一般情形。這種方法由以下步驟組成（Fletcher，1988a）

其中Pn組成一個搜索方向（開始時等于初始半徑R1），式（2.4.34）中的σn+1要保證這些搜索方向是相互正交的，而λn要保證式（2.4.32）在Pn方向上最小。在每一步中，基于A的特征函數，該方法通常使半徑R的級數展開式的特征向量貢獻消除。在沒有舍入誤差的情況下，該方法得到一個計算次數不超過N的魯棒算法，如果A的一些特征值聚集在一起，它的計算次數更少，不過，舍入誤差通常會降低這種方法的有效性。

如果對方程進行修改，使矩陣A可以轉換為一個近似單位矩陣，就可以避免上述問題，這稱為預處理，預條件共軛梯度法是從偏微分方程（PDEs）得到的有限差分方程進行求解的極有效方法（Fletcher，1998a；Press等，1992）。預處理減少了特征值的傳播，因此收斂速度更快，這種方法等價于求解下列方程，而不是解式（2.4.3）

其中

選擇矩陣W使WTW=E，選擇不同的E會使W產生不同的方案，那么預條件共軛梯度方案歸結為

其中S和代表偽殘差。與之前一樣，用σn+1的值表示收斂。

2.4.4　多格網法

使用現代多格網法對解加速的速度可能比使用SOR方法還快?？梢哉J為，松弛法本質上是迭代地使不正確的初始估計所引起的誤差場得到衰減，那么收斂的速度就取決于松弛方法使誤差衰減的確切程度，松弛法在衰減誤差場的高波數分量方面非常有效，通常只需要幾次迭代就能使波數分量衰減到接近格網分辨率波數（2Δx）-1。不過，低波數的衰減非常緩慢，因此在幾次迭代以后，總誤差的進一步減小變得緩慢。多格網方法是通過對一個事實的巧妙探索而實現對松弛法的收斂進行加速的，這個事實就是誤差的高頻分量在松弛過程中被迅速消除，某個分辨率格網下的低頻誤差在另一較低分辨率格網下會表現為高頻誤差。取代在一個單一格網上進行求解，它們在分辨率為各格網分辨率乘積的多格網上成功地進行了求解。

可以使用一個簡單的一維邊界值問題對多格網方法進行說明（Briggs，1987）為

利用中央二階精度方案，對應的有限差分方程變為

式（2.4.39）顯然是內部格網點j=2，…，n-1處的解的一個（n-2）×（n-2）三對角系統，可以用算法對其進行有效求解，例如Thomas算法。更一般地說，多維橢圓問題的一個有限差分近似中得到的線性方程系統的形式如式（2.4.3），如果yj表示一個近似解，那么誤差ej=xj-yj由殘差方程A·E=R=B-A·Y控制，可以使B=0，然后得到精確解為X=0，誤差就是E=-Y。

考慮y的初始估計由m=1到n的不同波數構成

用加權Jacobi方法進行100次迭代的應用表明，波數越高，衰減速度越快。波數分量n衰減得特別快，而波數1則衰減緩慢。而且，盡管前幾次迭代中高波數分量的抑制使總誤差大大減小，但誤差的進一步降低卻很緩慢，高頻分量被迅速消除，而低頻分量依然保留下來，因此，在前幾次迭代后，收斂速度急劇降低。在整個迭代方案中，誤差的低頻分量消除緩慢使得收斂緩慢。松弛方法是有效的平滑器，然而，盡管它們能對初始誤差分量進行有效抑制，卻使平滑分量本質上保持不變。收斂的速度不會受到加權因子的太大影響，實際上反而會受格網空間h的任何減小的影響。因此，嘗試重新定義格網以獲得更精確的解是不可能的，反而只會起到反效果。

多格網方法的指導原則是，在松弛過程中，振蕩分量（m＞＞n/2）大幅度衰減，而平滑分量（m＜＜n/2）則不會（Brandt，1977；Briggs，1987；Hackbusch，1985；McCormick，1987）。假設誤差的振蕩分量在某一特定格網中已經被利用幾次松弛搜索進行消除。因為平滑分量相對保持不變，進一步迭代不會帶來什么優勢，但是對密格網來說已經平滑的分量對粗格網來說卻是振蕩的，在粗格網上進行松弛會使這些誤差得到一定程度的消除，至少會變得更快。因此，用粗格網可能會更好，在該格網中的點上進行松弛，同時也進行松弛搜索。

這個過程可以從獲得期望解的格網開始，在進行幾次迭代后，當收斂變慢時，殘差被轉換到一個粗格網中。由于密格網上誤差的低頻分量在粗格網上看起來似乎是高頻分量，在對這些誤差分量進行衰減時，在粗格網上進行幾次迭代比在密格網上進行更多次的連續迭代要有效得多，然后反過來將這一改進的殘差估計（對殘差方程進行松弛處理）轉換到密格網中，從而獲得更好的誤差估計并使最終得到的解更精確，這一在密格網和粗格網上進行連續求解的過程可以保持收斂性，這種方法稱為粗格網校正，該過程也可以持續到格網達到最粗糙的情況（2×2），或者也可在某些中間點處終止。

另一種方法是從粗格網開始進行松弛搜索，并將其解作為更密格網上的一個初始估計，接著又再作為下一個更密格網的初始估計，直到得到滿意的解為止，重復整個循環直到收斂，這被稱為嵌套迭代。兩種方法都要求在兩個格網間對信息進行正確且高效轉換的方法，即從粗糙到密集和從密集到粗糙，前者稱為延伸，后者稱為限制。通常對延伸來說最合適的是簡單的線性內插法，而對于限制來說，則較合適的是簡單注入法，這種方法中粗格網點假設對應的密格網點會更好（盡管周圍的密集格網點的加權平均也會很好）。一個簡單的雙重粗格網校正方案如下（一直重復直到收斂）：

（1）用初始估計Yh在密格網Ωh上對AhXh=Bh進行p次松弛。

（2）在粗格網Ω2h上計算殘差

（3）用初始估計E2h=0在粗格網Ω2h上對A2hE2h=R2h進行p次松弛。

（4）修正密格網分辨率

（5）用新的初始估計Yh在密格網Ωh上對AhXh=B h進行p次松弛。

其中上標h表示原始格網，2h表示格網尺寸比原始格網大兩倍的粗格網。這個雙重格網方法是非常高效的，漸進收斂速度與格網大小有關。工作量與格網點總數N=nd大小相當，n為d問題維之一的格網點數，對最優的SOR來說為N3N2，對共軛梯度方法來說為N5/4

2.4.4.1　V循環

為了對粗格網（或分辨率更小的格網）平滑特性的優勢進行探索，采用的格網可能大于2個：

用初始估計Yh在密格網Ωh上對AhXh=Bh進行p次松弛。

在粗格網Ω2h上計算殘差

用初始估計E2h=0在粗格網Ω2h上對A2hE2h=R2h進行p次松弛。

在粗格網Ω4h上計算殘差

用初始估計E4h=0在粗格網Ω4h上對A4h E4h=R4h進行p次松弛。

在粗格網Ω8h上計算殘差

…

在粗格網ΩLh上對ALhELh=RLh進行松弛。

…

修正密格網分辨率

用初始估計Y4h在密格網Ω4h上對A4h X4h=B4h進行q次松弛。

修正密格網分辨率

用新的初始估計Y2h在密格網Ω2h上對A2h E2h=R2h進行q次松弛。

修正密格網分辨率

用新的初始估計Yh在密格網Ωh上對AhXh=B h進行q次松弛。

其中，上標4h表示格網大小為原始格網的四倍的粗格網，以此類推。這種方法稱為V循環多格網方法，它非常高效，多格網方法有很多可能的變形，比如W循環，與效率最好的基于FFT的直接解法的效率O（ndlogn）相比，該方法的效率為O（nd）。然而，必須注意的是，SOR不能用在多格網方法中，因為對多格網方法很關鍵的高頻衰減能力比常規松弛方法（例如Gauss-Seidel法）差很多，通常最適合多格網方法的是黑紅點松弛。更多細節和軟件方面的信息，讀者可以參考Fletcher（1988a）和Press等（1992）的著作。注意，對復雜形狀域（比如海盆）來說，基于FF T的方法和多格網方法的應用都很復雜。

需要注意，雖然松弛可應用到線性系統的原始方程或殘差方程中，但如果系統是非線性的，則必須將它用到原始方程中，更多細節讀者可以參考Brandt（1977）、Hackbush（1985）、Briggs（1987）、Fletcher（1988a）和Press（1992）等。

2.4.5　偽瞬態方法

穩態邊界值問題迭代方法通常等價于推進一個對應的時間依賴問題到它的漸進狀態。要證明這一點，可以將時間導數添加到拉普拉斯式（2.3.18）以得到一個雙曲線初始值問題

該方程可以用一個任意初始狀態求解到它的穩定狀態，與原始邊界值問題一樣受相同邊界條件約束，因此漸進穩定狀態對應于原始邊界值問題的解。使用顯式FTCS（時間前向，空間中央）方案（并假設Δx=Δy），得到

式中：n表示當前時間，Δt是時間步長，要滿足穩定性，（本章2.5節），當它等于0.25時，除了n表示時間級而不是迭代之外，所得到的方程與式（2.4.20）Jacobi迭代方案完全一樣。

這種人工添加一個時間依賴項從而將邊界值問題轉換為初始值問題并求漸進穩態解的方法稱為偽瞬態方法，它是一種克服求解復雜穩態問題的有效方法，復雜穩態問題中方程的性質有可能會在部分域中從雙曲線到橢圓發生改變，該方法在20世紀60年代早期被成功應用于當時幾乎不可處理的鈍再入體周圍高馬赫數流的穩態問題，這種流在駐點附近加速回到超音速狀態之前立即減速到次音速狀態。因此，這里的域包含嵌入的橢圓區域，而激震前沿首先是未知的。不同域中跨越邊界求解的復雜性必須作為解的一部分，這種復雜性可以避免，方法是通過將穩態問題轉換為一個初始值問題并從估計的初始狀態向前推進到漸進穩態。求解漸進狀態的等價非穩定問題有一個很大優勢，即方程是雙曲線的，可以暫時使用多維時間依賴問題的強大分割技術，比如交替方向隱式技術（ADI）。