滲流模型參數(shù)標(biāo)定與實(shí)操指南)
簡介這是一份面向計算流體動力學(xué)與流體仿真使用者的FLOW-3D多孔介質(zhì)滲流模型講解PPT適合從事油氣藏、水文地質(zhì)、環(huán)境污染模擬的工程師與科研人員。內(nèi)容系統(tǒng)覆蓋達(dá)西定律、拖曳力模型、飽和與非飽和多孔介質(zhì)模型、拖曳力系數(shù)與滲透率關(guān)系等核心理論并結(jié)合軟件界面與公式逐步說明模型參數(shù)設(shè)置思路。資源僅1個pptx文件整體約2.43MB便于直接打開學(xué)習(xí)。已有197人學(xué)習(xí)下載適合快速掌握FLOW-3D多孔介質(zhì)建模的基本框架理解滲透率、孔隙率、阻力系數(shù)等關(guān)鍵概念為開展多孔材料內(nèi)部流動仿真打下基礎(chǔ)。講解兼顧公式推導(dǎo)與工程應(yīng)用尤其適合需要系統(tǒng)入門滲流模擬的初學(xué)者及希望梳理理論脈絡(luò)的中級用戶。 第一次用FLOW3D做滲流項(xiàng)目是在一個堤防加固校核里。甲方要的是浸潤線位置和滲流量我當(dāng)時滿腦子只記得“FLOW3D有多孔介質(zhì)模型”于是急吼吼地把堤身設(shè)成多孔、孔隙率填了個0.3然后開始算。結(jié)果連續(xù)三輪浸潤線都比實(shí)測高不少流量也對不上開會時被問得臉都發(fā)燙。后來回頭查模型定義才發(fā)現(xiàn)我搞混了孔隙率和阻力系數(shù)這兩個完全不同的概念。這篇東西的緣起就是手頭那份《FLOW3D多孔介質(zhì)模型滲流模型》的技術(shù)整理我把多孔介質(zhì)參數(shù)、流動狀態(tài)判斷、飽和非飽和處理、實(shí)操流程和踩過的坑重新梳理了一遍希望能幫正在用或準(zhǔn)備用FLOW3D做滲流模擬的同行少走點(diǎn)彎路。不管你是做巖土、水利、地下水還是砂濾系統(tǒng)這套邏輯基本通用。1. 多孔介質(zhì)模型算的到底是什么——達(dá)西定律背面的流動狀態(tài)判斷1.1 達(dá)西定律的適用范圍很多教程一上來就擺出達(dá)西定律滲流量Q與水力梯度成正比即Q K·A·(ΔH/L)。這個關(guān)系在大學(xué)課本里是從細(xì)砂、土柱里總結(jié)出來的它成立的前提是流動足夠慢慣性力小到可以忽略阻力幾乎全部來自流體與孔隙壁面的黏性剪切。這個狀態(tài)用孔隙雷諾數(shù)判斷更可靠一般Re_p 1~10都算達(dá)西流公式是Re_p ρvd / μ其中v是孔隙內(nèi)的實(shí)際流速d是特征顆粒直徑。注意這里用的是顆粒粒徑不是宏觀幾何尺寸。舉個例子水流在粉細(xì)砂里的孔隙流速通常只有毫米到厘米每秒量級算出來Re_p遠(yuǎn)小于1用達(dá)西完全沒問題??梢坏┙橘|(zhì)變成礫石、拋石、堆石體流速上到每秒分米甚至米級Re_p過百慣性損失就開始占據(jù)主導(dǎo)地位達(dá)西定律預(yù)測的流量會明顯偏小。1.2 FLOW3D如何處理超出達(dá)西范圍的流動FLOW3D在多孔區(qū)域并不強(qiáng)行假設(shè)達(dá)西流而是在動量方程里直接加入一個拖曳力源項(xiàng)形式是“線性項(xiàng) 二次項(xiàng)”。線性項(xiàng)對應(yīng)黏性損失二次項(xiàng)對應(yīng)慣性/形狀阻力。也就是說同一個模型里如果要同時模擬細(xì)砂滲流和粗骨料通道流只要拖曳系數(shù)給得對流動狀態(tài)是模型自己算出來的不需要人肉切換公式。我通常用一個生活里的場景去理解這事在鵝卵石路上騎車特別慢的時候覺得顛簸和輪胎黏滯感強(qiáng)這是線性阻力騎快了風(fēng)阻和碰撞感撲面而來這就是二次阻力。FLOW3D的多孔介質(zhì)拖曳項(xiàng)就是同時裝了這兩種“阻尼”。這里要特別強(qiáng)調(diào)一個認(rèn)知多孔介質(zhì)模型不是“給個低孔隙率就能產(chǎn)生阻力”孔隙率決定的是介質(zhì)里能裝多少水真正控制水怎么流過去的是拖曳系數(shù)??紫堵?.8的泡沫金屬如果拖曳系數(shù)做得很大照樣能憋出很大的壓降孔隙率0.2的碎石堆如果參數(shù)設(shè)不對也可能算得比實(shí)際通暢得多。這點(diǎn)后面展開講。2. FLOW3D參數(shù)界面背后那些物理量怎么填才不瞎填2.1 以組件為中心的多孔介質(zhì)屬性設(shè)置FLOW3D設(shè)置多孔介質(zhì)核心操作是在幾何組件上“掛屬性”比在網(wǎng)格單元上指定物理區(qū)要直觀得多。你先用CAD導(dǎo)入或軟件自帶的primitive命令把堤身、砂濾層、堆石體這些需要當(dāng)多孔的區(qū)域建成獨(dú)立組件然后右鍵激活Porous Media選項(xiàng)接下來需要填的字段大致如下字段物理含義取值依據(jù)Porosity孔隙率流體可占據(jù)的體積分?jǐn)?shù)試驗(yàn)實(shí)測、材料手冊典型砂土0.3~0.45礫石0.25~0.4Linear Drag Coefficient線性拖曳系數(shù)對應(yīng)達(dá)西阻力滲透試驗(yàn)反算或Kozeny-Carman公式估算Quadratic Drag Coefficient二次拖曳系數(shù)對應(yīng)慣性阻力粗顆粒材料壓降試驗(yàn)或Ergun公式Initial Saturation初始飽和度現(xiàn)場含水率/地下水位Residual Saturation殘余飽和度土水特征曲線決定不可動水比例Direction參數(shù)各向同性或各向異性地層分層方向主滲流方向很多版本里線性drag系數(shù)是μ/k的形式單位是Pa·s/m2二次drag系數(shù)是βρ的形式單位是kg/m?。版本不同菜單名稱可能有差異但物理量就這兩個你在用戶手冊里搜Darcian和non-Darcian通常能找到對應(yīng)位置。2.2 拖曳系數(shù)到底怎么來這是新手最容易卡殼的地方。我推薦先用手邊數(shù)據(jù)做估算再用簡單算例標(biāo)定。比如砂礫石取特征粒徑d1mm孔隙率φ0.4用Kozeny-Carman公式估算滲透率k (d2/180) × φ3/(1-φ)2代入數(shù)字d21×10?? m2φ3/(1-φ)2 0.064/0.36 ≈ 0.178于是k ≈ (1×10??/180)×0.178 ≈ 9.9×10?1? m2。水的動力黏滯系數(shù)μ≈1×10?3 Pa·s線性拖曳系數(shù)就是μ/k ≈ 1.0×10? Pa·s/m2。如果材料更粗比如豆礫石d5mm線性系數(shù)會降到10?量級但二次系數(shù)會明顯抬升這時候就不能忽略Forchheimer項(xiàng)了。二次項(xiàng)可用Ergun公式估算β 1.75(1-φ)/(φ3·d)同一批砂樣φ0.4、d1mm時β≈1.75×0.6/(0.064×0.001) ≈ 1.64×10? m?1再乘水的密度1000 kg/m3二次系數(shù)約1.64×10? kg/m?。就這個數(shù)量級算下來當(dāng)流速到0.01 m/s時二次損失已經(jīng)跟線性損失差不多大足以見得高速區(qū)二次項(xiàng)不能省。如果你有現(xiàn)成的室內(nèi)滲透試驗(yàn)數(shù)據(jù)更直接的辦法是反算用不同水力梯度下的流速-壓降曲線擬合出線性系數(shù)和二次系數(shù)比任何估算公式都準(zhǔn)。我沒試驗(yàn)數(shù)據(jù)時至少也會做一組單箱流水算例跟理論解對比后再上完整模型。3. 飽和到非飽和殘余濕度和毛細(xì)壓力在FLOW3D里怎么體現(xiàn)3.1 默認(rèn)情況下的“飽和流”陷阱FLOW3D多孔介質(zhì)模型如果不額外設(shè)置默認(rèn)是飽和流假定孔隙全部被水充滿壓力水頭直接傳遞。這個假定對堤防穩(wěn)定滲流階段沒問題但一旦涉及降雨入滲、庫水位驟降、包氣帶水分運(yùn)移飽和流假定就會出現(xiàn)一個尷尬局面——干燥區(qū)域要么完全不進(jìn)水要么水一下子灌滿整個介質(zhì)跟實(shí)際情況差出幾條街。問題出在非飽和區(qū)的負(fù)壓和殘余水量上。天然土壤里水不是靠重力流就能走遍所有孔隙的部分水被毛細(xì)力吸在細(xì)孔隙里、吸附在顆粒表面這部分水在重力條件下不參與流動就是殘余飽和度。如果沒有定義殘余飽和度模型會把最后那點(diǎn)水也算成可動水導(dǎo)致排水過程偏慢或偏快。3.2 FLOW3D怎么設(shè)置非飽和流動操作上要做兩件事第一在多孔介質(zhì)屬性里指定初始飽和度和殘余飽和度第二啟用毛細(xì)壓力模型常用的是Van Genuchten模型需要填α、n、m這幾個參數(shù)。沃倫Van Genuchten的參數(shù)本質(zhì)是描述“飽和度-負(fù)壓”曲線形態(tài)α控制進(jìn)氣值大小n控制曲線陡緩m通常取1-1/n。砂土的α常在0.1~1/m量級黏土的α?xí)『芏嘁驗(yàn)檫M(jìn)氣值大。如果你手頭有土水特征曲線測試報告直接用實(shí)測點(diǎn)導(dǎo)入比盲填參數(shù)穩(wěn)妥得多。實(shí)際工程里非飽和滲流出現(xiàn)最多的情況是堤防初蓄水階段、邊坡降雨入滲和尾礦庫干灘區(qū)域。我做過一個降雨工況的邊坡滲流分析剛開始沒啟用毛細(xì)模型雨滴落在坡面直接形成地表徑流坡體含水量幾乎不動后來啟用了Van Genuchten模型、把殘余飽和度設(shè)成0.15左右水分才按預(yù)期的浸潤鋒速度向下遷移。另外有個必須注意的坑非飽和模型下流速顯示的是達(dá)西通量宏觀面流速不直接等于孔隙內(nèi)的真實(shí)流速。兩者的關(guān)系是q φvφ是孔隙率。如果你要算粒子示蹤或者污染物運(yùn)移停留時間必須用物理速度v去算否則會差一個孔隙率的倍數(shù)。這點(diǎn)在報告里不寫清楚專家評審一定會挑。4. 用FLOW3D做滲流的完整操作路線——一個二維堤壩算例的服役過程4.1 能直接用到的模型設(shè)定步驟我不整虛的直接以二維土堤滲流算例走一遍。目標(biāo)是算穩(wěn)態(tài)下的浸潤線位置和單寬滲流量介質(zhì)取中砂粒徑0.5mm孔隙率0.35。實(shí)測滲透系數(shù)K≈5×10?? m/s水深差2m堤身長度10m。第一步建幾何。在FLOW3D里用Geometry模塊拉出一個10m×4m的矩形堤身區(qū)域把它設(shè)成一個獨(dú)立組件。上游和下游各留出空槽作為水體區(qū)域堤身組件底部落在邊界上。第二步定義流體和物理模型。激活水Newtonian fluid勾選Porous Media多孔介質(zhì)選項(xiàng)。這一步不做后面組件屬性里再多孔選項(xiàng)也沒用。第三步給組件掛多孔屬性。選擇堤身組件激活Porous Media填孔隙率0.35。拖曳系數(shù)我用滲透率反算K5×10?? m/sμ1×10?3 Pa·s所以k Kμ/ρg (5×10?? × 1×10?3)/(1000×9.81) ≈ 5.1×10?12 m2線性拖曳系數(shù)≈1.96×10? Pa·s/m2這個數(shù)量級跟細(xì)砂匹配。二次系數(shù)顆粒太細(xì)可以不填或填一個很小值。第四步邊界和初始條件。上游邊界選壓力邊界靜水壓力對應(yīng)2m水深下游邊界設(shè)壓力出口對應(yīng)0m水位。初始水位從0開始讓水自己往堤身里推進(jìn)如果要快速收斂也可以直接給定初始水位場。第五步劃分網(wǎng)格。網(wǎng)格尺寸取0.1m堤身厚度10m方向上有100個網(wǎng)格對二維穩(wěn)態(tài)問題足夠。FLOW3D用FAVOR方法處理邊界網(wǎng)格夠細(xì)才能準(zhǔn)確捕捉多孔介質(zhì)區(qū)域的進(jìn)入界面不然孔隙率和邊界體積會被低估。第六步運(yùn)行和后處理。跑穩(wěn)態(tài)看壓力云圖、飽和度等值線取0.5飽和度作為浸潤線近似再從后處理里統(tǒng)計上游進(jìn)流量。我用一維達(dá)西理論解驗(yàn)證Q K×A×ΔH/L 5×10?? × 1 × 2/10 1×10?? m3/s單寬模擬值通常在0.9×10??~1.1×10??之間浮動算對上。如果你算出來差好幾倍先別調(diào)軟件回頭看拖曳系數(shù)單位對不對或者網(wǎng)格是不是太粗。4.2 為什么這個流程不建議跳過任何一步這套流程里最容易讓人“覺得自己會了但其實(shí)沒會”的是跳過第二步直接填組件屬性或者不驗(yàn)證流量直接上復(fù)雜三維模型。多孔介質(zhì)模擬的成敗90%在參數(shù)標(biāo)定剩下10%在邊界和初始條件匹配。所以任何人問我要模板我都建議他先用一個最簡單的二維/一維模型把模擬結(jié)果跟理論解或?qū)崪y數(shù)據(jù)對一遍確認(rèn)參數(shù)和流動狀態(tài)匹配對了再往實(shí)際工程規(guī)模上搬。5. 第一次用多孔介質(zhì)模型最容易翻車的五個細(xì)節(jié)5.1 FAVOR邊界與孔隙率疊加造成的體積丟失FLOW3D用FAVOR方法把幾何體映射到網(wǎng)格多孔介質(zhì)區(qū)域的孔隙率會與網(wǎng)格體積分?jǐn)?shù)疊加如果網(wǎng)格太粗邊界上單元的流動面積會被“雙殺”出現(xiàn)進(jìn)口壓力偏大或流量偏小。這個問題的排查方法是單獨(dú)建一個只有多孔介質(zhì)區(qū)域、沒有其余障礙物的模型對比體積積分是不是接近幾何理論值如果差超過5%就應(yīng)該加密網(wǎng)格。5.2 初始飽和度沒設(shè)對水位推進(jìn)過程嚴(yán)重失真這個問題在非飽和模型里最突出。不設(shè)初始飽和度等于假設(shè)整個多孔區(qū)域一開始全干水淹進(jìn)堤身時會把空氣往外擠但FLOW3D如果沒處理好兩相壓力耦合就會出現(xiàn)一段短暫的高壓異常有時候還表現(xiàn)為浸潤線抬升過頭。實(shí)際工程里上游水位變化前堤身通常有初始含水率在模型里至少給一個與地下水位對應(yīng)的初始飽和度分布結(jié)果會穩(wěn)得多。5.3 二次拖曳系數(shù)缺失高速區(qū)流量嚴(yán)重偏大粗顆粒介質(zhì)或高水頭差場景里如果不填二次系數(shù)等于只算了黏性損失慣性損失被忽略結(jié)果流量會偏大有時候能大出30%~50%。判斷標(biāo)準(zhǔn)很簡單看孔隙雷諾數(shù)Re_p超過10就應(yīng)該考慮二次項(xiàng)。我曾經(jīng)算一個堆石排水體只填線性系數(shù)時流量比實(shí)測大40%加上Ergun估算的二次項(xiàng)后中心通道的流速立刻降下來跟壓力傳感器實(shí)測基本吻合。5.4 各向異性方向定義錯滲流場整個反轉(zhuǎn)層狀地層的水平向和垂直向滲透系數(shù)可以差一個數(shù)量級以上FLOW3D里各向異性參數(shù)是綁定坐標(biāo)系的。如果地層走向跟全局坐標(biāo)系有夾角必須在組件里正確設(shè)定各向異性方向或者分多個組件分別定義。我見過最離譜的案例是把水平向拖曳系數(shù)填到了垂直方向結(jié)果本應(yīng)沿層理方向的水流被憋住垂直方向卻異常通暢滲流場整個歪掉。5.5 后處理速度單位與物理速度混淆FLOW3D后處理里顯示的流速多孔介質(zhì)區(qū)域往往是達(dá)西通量不是孔隙物理流速。做粒子追蹤或溶質(zhì)運(yùn)移時如果不做q/φ修正停留時間會算短。要驗(yàn)證這個坑可以看下游流量與多孔介質(zhì)邊界面積的商再跟后處理里的面平均速度對比兩者應(yīng)該一致如果一致了還覺得粒子跑得太快那就是該做物理速度修正了。最后再分享一個我個人習(xí)慣每次建模前先問自己三個問題——流動是達(dá)西流還是非達(dá)西流介質(zhì)是飽和還是非飽和各向異性方向是不是跟真實(shí)地層對齊。三個問題回答完參數(shù)和模型結(jié)構(gòu)基本就定了。多孔介質(zhì)模擬的難點(diǎn)從來不在軟件操作而在對流動機(jī)制的敬畏心。把孔隙率當(dāng)成阻力參數(shù)把二次拖曳系數(shù)當(dāng)擺設(shè)模型就會在關(guān)鍵時刻給你顏色看。本文還有配套的精品資源點(diǎn)擊獲取