有限體積法基本原理 課件 chap06-2D輸運方程及有限體積法

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

韓亞杰,鹿暢6.12D輸運方程通用溫度輸運的對流擴散方程：忽略時間導數(shù)項和對流項，可以得到：在直角坐標系下將梯度算子與點積算子展開，可以得到：6.2控制方程離散控制方程中各項在有限體積內(nèi)積分：

6.2控制方程離散根據(jù)散度定理將控制方程中體積積分項進行變換：單元具有有限數(shù)量的面（??），對于四邊形單元，??=4，而對于三角形單元，??=3。因此，對整個表面的曲面積分可以用對單元的??個面積分來代替：在二階有限體積法中，流動量在整個面上都是線性變化的，因此，整個面的積分可以簡化為面心的值（常數(shù)）：6.2控制方程離散——內(nèi)部單元離散格式為進一步簡化方程，參考圖6–1中的2D四邊形單元，該單元有左面（??）、右面（??）、頂面（??）和底面（??）。單位法向量總是指向單元外，??方向從左到右為正，??方向從下到上為正。這種離散方法適用于網(wǎng)格中所有單元形狀和方向與圖6–1中四邊形單元（正規(guī)網(wǎng)格）相同，為進一步簡化，接下來將對內(nèi)部單元和邊界單元的離散格式分別考慮。6.2控制方程離散——內(nèi)部單元離散格式所有擴散項都將使用中心差法，使用圖6–2中的符號表示單元中心點之間的距離??，上式可以表示為：

按照單元中心溫度（????、????、????、????和????）對上式進行排序，可以得到：6.2控制方程離散——內(nèi)部單元離散格式所以對于內(nèi)部單元，2D擴散方程的有限體積離散格式及系數(shù)如下：使用求和符號對離散格式進行表示：6.2控制方程離散——邊界單元（Dirichlet）右、上和下表面與內(nèi)部單元相連，因此，通過這些面的熱擴散可采用與內(nèi)部單元相同的處理方法來進行計算：本章將在左面應用定值（Dirichlet）邊界條件。左邊界面上的溫度為??????????，到壁面的距離為到下一個單元中心點距離的一半。因此，通過左邊界面的熱擴散為：將上式代入邊界單元對應離散格式中：6.2控制方程離散——邊界單元（Dirichlet）為與內(nèi)部單元離散格式保持一致，將上式變換為：需要注意的是，使用求和符號時對于邊界面????

=0，并且????的求和不包括邊界面。6.2控制方程離散——邊界單元（Neumann）邊界施加固定通量（Neumann）時的邊界單元如圖6-4所示，該單元右、上和下表面與內(nèi)部單元相連，因此，通過這些面的熱擴散可采用與內(nèi)部單元相同的處理方法來進行計算：

6.2控制方程離散——邊界單元（Neumann）為與內(nèi)部單元離散格式保持一致，將上式變換為：需要注意的是，使用求和符號時對于邊界面????

=0，并且????的求和不包括邊界面。6.3控制方程各單元系數(shù)與線性系統(tǒng)根據(jù)前述推導，可以得到控制方程對應單元各系數(shù)(Dirichlet邊界與內(nèi)部單元)：Neumann邊界與內(nèi)部單元系數(shù)表：單元類型Interior0Left0Right0Bottom0Top0單元類型Interior0Left00Right00Bottom00Top006.3控制方程各單元系數(shù)與線性系統(tǒng)將前述單元系數(shù)用求和符號表示，如下表：單元面類型Interior0Dirichlet0Neumann00

6.2控制方程離散——多邊界界單元對于2D情況，一些邊界單元可能具有不止一個邊界面，且不同邊界面的邊界條件也可能不同，如圖6-5所示，計算域左上角具有2個邊界面，該單元的右面與下面與內(nèi)部單元相連，因此，通過這些面的熱擴散通量可采用與內(nèi)部單元相同的處理方法來進行計算：：如果邊界面的左側為固定溫度??????????,??，頂部為固定溫度??????????,t，那么也可以在邊界面上使用中心差分法。6.2控制方程離散——多邊界界單元為與內(nèi)部單元離散格式保持一致，將上式變換為：使用求和符號表示：6.3Example1——帶熱源的熱傳導根據(jù)問題模型，可以知道桿件中溫度滿足如下控制方程：6.3Example1——帶熱源的熱傳導Step-1：網(wǎng)格劃分Step-2：添加材料屬性Step-3：計算方程中各系數(shù)單元1010100-4071006021010100-2051005031010100-205100504100100-409100605

0101010-202100506

101010100100407

101010100100408

1001010-204100506.3Example1——帶熱源的熱傳導Step-3：計算方程中各系數(shù)單元90101010-2021005010

101010100100401110101010010040121001010-2041005013

010010-40510060141010010-2031005015

1010010-2031005016100010-40710060Step-3：計算方程中各系數(shù)(采用求和符號表示各系數(shù))單元面類型Interior10100Dirichlet(L)2002000Dirichlet(B)2003000Dirichlet(R)2004000Dirichlet(T)20050006.4Example1——帶熱源的熱傳導OpenFOAM求解器構造——求解器結構及各文件功能簡介(復習)files文件主要完成兩個工作：1）指定需要編譯的文件，這里的文件是代碼文件，即.C文件，而不包括.H文件；2）指定編譯類型及對應類型名稱，如編譯成庫類型，其關鍵字為LIB=···；編譯成可執(zhí)行程序類型，其關鍵字為EXE=···。options文件也主要完成兩個工作：1）指定編譯的頭文件所在目錄，如果所用頭文件在當前文件夾或為標準C++頭文件，無需指定，EXE_INC和關鍵字-I用于指定頭文件所在目錄；2）指定編譯當前程序用到的庫，EXE_LIBS和關鍵字-L/-l用于指定用到的庫。6.4Example1——帶熱源的熱傳導createFields.H——初始化基本物理場待求溫度場速度場溫度擴散系數(shù)熱源OpenFOAM求解器構造OpenFOAM求解器構造6.4Example1——帶熱源的熱傳導convectionAndDiffusion.C——求解器源文件網(wǎng)格讀取線性方程構造方程求解及輸出算例文件路徑正確性檢測時間類初始化網(wǎng)格類初始化基本物理場初始化構造線性方程：各項與算例問題中對應線性方程求解結果文件輸出6.4Example1——帶熱源的熱傳導OpenFOAM求解器構造編輯完成OpenFOAM各求解器后，在求解器下執(zhí)行wmake對求解器進行編譯，得到可執(zhí)行程序。求解器源文件編譯得到求解器名稱g++之后的指令參數(shù)依次指定了采用的C++標準版本、編譯器參數(shù)、機器位數(shù)、整數(shù)位數(shù)等信息，之后以“-I”開頭的表示使用到的頭文件所在的目錄，以“-L”開頭的表示指定的庫，以“-1”開頭的表示鏈接到指定的庫，“-oapp”表示生成的程序名為app。6.4Example1——帶熱源的熱傳導OpenFOAM

求解器對應算例構造——帶熱源熱傳導算例構造算例文件較為方便方法修改已有的、與問題模型相似的算例，這里以修改cavity為例進行介紹，對應算例文件路徑：openfoam10/tutorials/incompressible/icoFoam/cavity6.4Example1——帶熱源的熱傳導OpenFOAM求解器對應算例構造——帶熱源熱傳導算例6.4Example1——帶熱源的熱傳導OpenFOAM求解器對應算例構造——帶熱源熱傳導算例6.4Example1——帶熱源的熱傳導OpenFOAM求解器對應算例構造算例文件編輯完成后以此執(zhí)行：blockMesh、../app、paraFoam分別進行網(wǎng)格劃分、算例求解及后處理，對應各步驟過程分別如下：6.4Example1——帶熱源的熱傳導OpenFOAM求解器對應算例構造算例文件編輯完成后以此執(zhí)行：blockMesh、../diffusionEquationSolver、paraFoam

分別進行網(wǎng)格劃分、算例求解及后處理，對應各步驟過程分別如下：6.4Example1——帶熱源的熱傳導自編程序簡介——2D有熱源熱傳導通過自編程序方式進行問題求解，程序主要包含模型參數(shù)設置、網(wǎng)格劃分、矩陣系數(shù)計算、線性系統(tǒng)構造、求解及后處理。程序各部分定義參考如下：6.4Example1——帶熱源的熱傳導自編程序簡介——2D有熱源熱傳導通過自編程序方式進行問題求解，程序主要包含模型參數(shù)設置、網(wǎng)格劃分、矩陣系數(shù)計算、線性系統(tǒng)構造、求解及后處理。程序各部分定義參考如下：6.4Example1——帶熱源的熱傳導算例結果在進行后處理時較為方便的繪圖方法是將結果保存在網(wǎng)格節(jié)點上而非網(wǎng)格單元中心，如圖6–10所示，因為邊界上的溫度值是分布在網(wǎng)格節(jié)點上（圖中紅色部分），如果采用單元中心（圖中藍色部分）進行后處理則無法得到邊界上的溫度分布。常見的CFD程序采用單元中心進行計算，因此，內(nèi)部節(jié)點值（圖6–10中以白色顯示）需要通過單元中心點值（藍色顯示）之間的插值計算得到。6.5熱平衡對于2D熱平衡檢驗時，內(nèi)部單元一般有限體積離散格式如下：根據(jù)傅里葉定律通過每個單元面的熱通量為：將上式代入對應離散格式中可以得到：負號是為了確保熱量流動方向與溫度梯度方向相反對于2D四邊形單元而言，熱平衡誤差表示為：6.5熱平衡計算前述算例問題中各單元的熱平衡誤差，結果如下表所示：6.5熱平衡單元熱通量平衡表還可用于計算平板每個面（跨越所有邊界）流出的總熱通量：將每個邊界面（上表中用紅色標出的部分）流出的熱通量相加，結果如表6–8所示：正熱通量表示熱量穿過平板從平板邊界流出，負熱通量則表示熱量通過邊界進入平板內(nèi)。左邊界作為溫度最低的面（????

=100℃），該平面的熱通量最大。相反，上表面作為最熱面（????

=250℃），熱量通過該平面進入平板內(nèi)，從而產(chǎn)生