[發明專利]一種基于界帶有限元和拉格朗日坐標的流體仿真方法

有效

• 申請號: 201410483986.2
• 申請日: 2014-09-19
• 公開（公告）號: CN104317985B
• 公開（公告）日: 2017-05-24
• 發明（設計）人: 吳鋒;徐小明;陳飆松;鐘萬勰
• 申請（專利權）人: 大連理工大學
• 主分類號: G06F17/50
• 分類號: G06F17/50
• 代理公司: 大連理工大學專利中心21200
• 代理人: 趙連明,梅洪玉

技術領域

本發明涉及流體仿真分析技術，構建了一種基于界帶有限元和拉格朗日坐標的流體仿真方法。

背景技術

流體的仿真分析廣泛應用于各個不同工程領域，如水利工程、航空飛行、高速列車等。流體的仿真分析，從理論上看主要有歐拉坐標和拉格朗日坐標兩種，從數值分析手段來說，主要有限差分法和有限單元法兩種。其中有限差分方法使用規則的矩形網格，因此對于在不規則體中的流體仿真分析需要加密網格，這導致計算量大幅增加。有限單元法在固體結構仿真分析中得到廣泛運用，該方法的特點是適用于不規則問題，但是需要變分原理。

目前，流體的仿真分析主流是基于歐拉坐標，在歐拉坐標系可以導出不可壓縮流體的運動微分方程，即著名的納維-斯托克斯方程，對于該方程的計算最常用的方法是有限差分法。但是基于歐拉坐標的流體仿真分析在處理自由面流體問題時，分析困難，計算格式復雜。由于自由面隨時間不斷變化，其計算區域也不再規則，因此用有限差分分析困難，而且基于歐拉坐標的流體方程，很難建立變分原理，從而導致有限元分析也很困難。目前僅能對某些特殊的流體運動問題，建立變分，使用有限元分析。在分析不可壓縮流體時，采用普通有限元，還存在體積閉鎖等問題，影響計算精度。

如果運用拉格朗日坐標，則可以很容易得到流體的動能表達式和勢能表達式，并可以利用哈密爾頓變分原理導出拉格朗日坐標系的流體動力微分方程。根據動能表達式和勢能表達式，還可以利用有限元進行分析。但是在拉格朗日坐標下，以流體中質點的位移為基本未知量，以位移為基本未知量建立的傳統有限單元不能完全滿足不可壓縮條件，計算精度不好。實際上在固體不可壓縮材料的仿真分析中，有學者提出界帶有限元方法，通過引入流函數來代替位移作為未知量，并利用界帶單元進行分析，所構造的單元完全滿足不可壓縮條件。界帶單元可以在傳統有限元網格上進行分析，因此可以利用現有有限元網格剖分技術，而且具有計算精度高、所需的自由度小等優點。目前還沒有關于流體仿真分析中采用界帶有限元方法的應用。

發明內容

本發明針對現有技術不足，提出一種分析二維不可壓縮流體的仿真分析方法，該方法將拉格朗日坐標方法與界帶有限元方法結合來解決不可壓縮流體的運動仿真問題，目的是利用界帶有限元精度高和拉格朗日坐標下邊界處理方便，通用性好的優勢，提高分析的計算效率和精度。

為此，本發明提供了一種基于界帶有限元和拉格朗日坐標的流體仿真方法，將二維不可壓縮流體的計算域Ω按照傳統有限元網格剖分成N_e個單元，每個單元為Ω_i；構造單元Ω_i的位移插值場，根據位移插值場構造流體的動力微分方程，求解該動力微分方程得到流體的各種物理參數，從而進行流體的運動分析；本發明用界帶有限單元法構造位移插值場；并基于拉格朗日坐標描述法得到流體的動力微分方程；具體方法如下：

(a)將二維不可壓縮流體的計算域Ω采用傳統有限元網格剖分成N_e個單元，每個單元為Ω_i；

(b)在單元Ω_i上建立流函數ψ(x,y)的插值場時，以單元Ω_i為本體，將Ω_i周邊的單元視為Ω_i的界帶，將這些單元的節點合作進行插值，構造單元Ω_i上的插值函數ψ(x,y)；

(c)將步驟(b)中所述的流函數表達式求偏導，得到流體中在坐標(x,y)處質點的位移表達式；

(d)根據上一步的位移場表達式，得到每個單元上流體的質量矩陣M_i，并把所有單元的質量矩陣計算出來后，通過累加得到總體質量矩陣M；

(e)根據步驟(c)的位移場表達式，得到流體的剛度矩陣：

K(ψ)＝K_l+K_n(ψ)

其中K是剛度矩陣，ψ是所有節點的流函數值所組成的向量，K_l是線性剛度矩陣，K_n是非線性剛度矩陣；

(f)流體的邊界包括兩種，一種為自由面Γ_f，一種為不可穿過邊界Γ_n，根據不可穿過邊界Γ_n上所有有限元單元節點的編號，把剛度矩陣和質量矩陣中相應編號的行和列劃去，得到描述流體運動的非線性微分方程：

(g)利用非線性微分方程求解軟件求解上述非線性微分方程，得到不同時間點上的流函數向量ψ；