大連理工大學工業裝備結構分析優化與CAE軟件全國重點實驗室有限體積法基本原理夏廣慶,張軍軍,陳沖,

韓亞杰,鹿暢9.1非穩態控制方程在實際中，非穩態問題通常無法通過單獨求解連續性方程、動量方程或者能量方程解決，大多數時候需要聯立這些方程構成能夠描述指定問題的、復雜方程組，隨后通過各個參數之間的耦合計算得到最終期望的結果。常見不可壓縮流體非穩態輸運過程中，主要求解的方程有連續性方程、動量方程。

9.2控制方程離散——動量方程離散前述章節中介紹了穩態方程離散方法，對于非穩態方程而言其離散方法與前述方程離散過程相同，區別在于增加了非穩態項離散。采用有限體積法對動量方程離散，首先對方程在控制體積上進行積分，包括空間積分和時間積分：

9.2控制方程離散——動量方程離散（非穩態與對流項）對于非穩態項進行離散，可以采用Euler隱式格式進行離散，即：

對于對流項，由高斯散度定理可以得到：其中，????

=??·U??為面??上的對流通量，為一標量。對于2D情況有：9.2控制方程離散——動量方程離散（非穩態與對流項）

將上式代入前式中，可以得到：9.2控制方程離散——動量方程離散（擴散項）

上式中面上梯度項?U

的方向從單元體心??指向臨近網格單元??，即該矢量的方向和網格??與網格??體心連線共線，如果體心連線和面法向矢量不垂直，則需要引入非正交修正，這里網格為正規網格，無需此步驟。對于r面，單位面矢量(1,0):9.2控制方程離散——動量方程離散（擴散項）對于l

面，單位面矢量(-1,0):對于t面，單位面矢量(0,1):對于b面，單位面矢量(0,-1):所以可以得到擴散項的離散結果為：9.2控制方程離散——動量方程離散（壓力項及小結）對壓力項采用散度定理可以得到：對各微分項采用中心差分法進行離散，得到：將上述各項進行整理，采用如前所述非結構化網格求和描述方法，可以得到動量方程的離散形式為：9.2控制方程離散——動量方程離散（小結）如果采用正規網格，對流項、擴散項、壓力項均采用中心差分形式，將前述所得代入上式可以得到：可以看出，所得系數在??時刻均為已知，進一步展開可以得到：9.2控制方程離散——壓力修正法（連續性方程離散）前面推導得到了速度的離散格式，但是在實際計算中，速度與壓力之間是耦合的，僅靠動量方程無法求解壓力項，由此需要用到連續方程計算壓力。其基本思想是：首先根據給定的速度場（或者由前述推導得到的離散方程計算得到的速度場）修正壓力，使壓力滿足連續性方程。在此過程中能夠得到壓力的修正值??’，通常的做法有兩種：1）直接求解??’的方程；2）求解壓力的泊松方程。OpenFOAM采用第二種方法。對連續性方程進行積分，隨后運用散度定理可以得到：如果方程收斂，則有如下離散格式成立：9.2控制方程離散——壓力修正法（連續性方程離散）將上式進行整理，變換，可以得到：將上式代入連續性方程中：由高斯定理，將上式逆向推導，得到壓力p

的泊松方程：注意：方程兩端都為未知量，無法進行求解。由于引入了一個新的場量????????，CFD將一個物理問題轉化為了代數方程組的求解問題，而????????包含了新計算出的速度，因此會變。9.2控制方程離散——壓力修正法（壓力修正方程）

但在該式中包含前一步壓力，因此，該離散格式存在滯后性，實際收斂后速度應該滿足下式：

9.2控制方程離散——壓力修正法（壓力修正方程）將前式做如下變換：

上式中左側已知，右側未知，能夠求解，在得到壓力后，代入式壓力的泊松方程可以得到速度的修正值，最終得到速度的迭代值。9.2控制方程離散——壓力修正法（壓力泊松方程）如前所述壓力泊松方程為：對于左端HbyA

項：由此????????定義式先計算出????????|??和????????|??，再計算面上的值，再得到壓力泊松方程左端項：9.2控制方程離散——壓力修正法（壓力泊松方程）如前所述壓力泊松方程為：

得到壓力的更新方程：9.2控制方程離散——壓力修正法（壓力泊松方程）

2）由于邊界的通量?????和?????????????是固定的，所以不需要插值。這里由原本的連續性方程得到邊界的????????和????，然后將其代入壓力泊松方程中，可以得到邊界的壓力值。如對左上角網格而言，假設有：U??

=(0,0)，U??

=(1,0)，則可以的得到：9.2控制方程離散——邊界單元離散（左邊界不含角點）9.2控制方程離散——邊界單元離散（右邊界不含角點）9.2控制方程離散——邊界單元離散（上邊界不含角點）9.2控制方程離散——邊界單元離散（下邊界不含角點）9.2控制方程離散——邊界單元離散（左上角Dirichlet）9.2控制方程離散——邊界單元離散（右上角Dirichlet）9.2控制方程離散——邊界單元離散（左下角Dirichlet）9.2控制方程離散——邊界單元離散（右下角Dirichlet）9.3非穩態問題求解算法——SIMPLE算法非穩態問題中各物理量隨時間變化，物理量之間存在耦合關系，它們之間相互依賴，需要通過一定的算法保證問題求解的正確性。針對本問題常用的求解方法有：壓力耦合方程組的半隱式方法（Semi-ImplicitMethodforPressureLinkedEquations，SIMPLE），壓力的隱式算子分割算法（PressureImplicitwithSplittingofOperators，PISO）。SIMPLE算法開始

是否收斂？

SIMPLE算法結束否是

9.3非穩態問題求解算法——PISO算法PISO算法最初被用來處理瞬態問題，其方程的構建過程與SIMPLE算法基本一致，因此PISO算法的基本思想與SIMPLE算法相似，PISO算法與SIMPLE算法第1、2步相同，不同之處在與第3步：在該步驟中，設置修正次數nCorr，在每一次修正中，需要進行以下步驟：按照SIMPLE算法的第3步計算????????和????

按照SIMPLE算法的第4步計算???按照SIMPLE算法的第5步計算U??

如果nCorr=1則該算法與SIMPLE算法相同，如果nCorr>1則需要進行多次修正，直到速度場和壓力場收斂，可以看到PISO只是在SIMPLE的基礎上增加了一個內循環，即修正循環。9.4Example1——非結構化網格2D熱傳導本節以經典頂蓋驅動方腔流為例對2D非穩態問題進行算例演示，該算例模型與OpenFOAM中算例模型一致。算例文件夾路徑為：OpenFOAM-10/tutori