程林松，時俊杰，曹仁義，楊晨旭，杜旭林

中國石油大學(xué)（北京）石油工程學(xué)院，北京 102249

全球碳中和進(jìn)程的加速，使清潔可再生能源的地位不斷提高，倒逼中國能源體系朝著能源多元化發(fā)展結(jié)構(gòu)持續(xù)轉(zhuǎn)型．干熱巖（hot dry rock，HDR）是一種儲量豐富、分布廣、埋藏深的清潔地?zé)崮茉?，近年來受到廣泛關(guān)注［1-2］．HDR的有效開發(fā)需要依托體積壓裂技術(shù)構(gòu)建的增強(qiáng)型地?zé)嵯到y(tǒng)，這就導(dǎo)致大量溫度低于儲層溫度的流體進(jìn)入到地層中．復(fù)雜的地下儲層環(huán)境導(dǎo)致了復(fù)雜的流動和傳熱過程，精確模擬和描述致密基質(zhì)-復(fù)雜裂縫傳質(zhì)傳熱過程是探究干熱巖開發(fā)過程的關(guān)鍵．

裂縫是地下流體主要的滲流通道，研究這一問題的首要任務(wù)是精確表征復(fù)雜裂縫的分布．雙重介質(zhì)模型是目前商業(yè)數(shù)值模擬中廣泛使用的描述裂縫性儲層的傳統(tǒng)方法［3］，趙延林等［4］建立了雙重介質(zhì)溫度場-滲流場-應(yīng)力場耦合數(shù)學(xué)模型，研究認(rèn)為，變溫條件下溫度場對裂縫滲流場具有顯著的耦合作用．WANG等［5］基于雙重介質(zhì)模型模擬了增強(qiáng)型地?zé)嵯到y(tǒng)中傳質(zhì)傳熱過程，發(fā)現(xiàn)冷水注入使裂縫滲透率提高了7倍．雙重介質(zhì)模型所需的參數(shù)較少，計算效率較高，但其假設(shè)裂縫是規(guī)則分布的，導(dǎo)致無法準(zhǔn)確刻畫裂縫分布形態(tài)．離散裂縫模型（discrete fracture model，DFM）能夠根據(jù)裂縫幾何形狀準(zhǔn)確地捕捉基質(zhì)-裂縫之間的質(zhì)量或熱量傳遞，從而提供更準(zhǔn)確的模擬結(jié)果．SUN等［6］基于離散裂縫模型和有限元方法建立了滿足局部非熱平衡假設(shè)的熱流固（thermal-hydraulic-mechanical，THM）耦合數(shù)學(xué)模型，表明了考慮熱流固耦合在增強(qiáng)型地?zé)嵯到y(tǒng)研究中的重要性．但DFM在計算效率和計算成本方面，仍有一定局限性．LEE等［7］首次提出嵌入式離散裂縫模型（embedded discrete fracture model，EDFM），該方法可以避免非結(jié)構(gòu)網(wǎng)格的剖分，大幅減少網(wǎng)格數(shù)量，從而提高計算效率．EDFM方法已經(jīng)應(yīng)用到各類裂縫性油氣藏模擬的各種復(fù)雜情況中．WANG等［8］利用積分有限差分法，將這種嵌入式裂縫建模方法直接應(yīng)用到地?zé)醿恿芽p系統(tǒng)中，認(rèn)為該方法在真實三維離散裂縫的傳質(zhì)傳熱問題上具有很好的應(yīng)用前景．基于EDFM方法，LI等［9］提出了一種研究裂縫地?zé)醿拥牧鲃雍蛡鳠岬谋菊髡环纸饨惦A模型，對基質(zhì)與裂縫之間的傳質(zhì)傳熱求解進(jìn)行簡化，使計算效率提高了10～15倍．PRADITIA等［10］在EDFM的框架下，提出了一種裂縫性油藏的單相熱流耦合方程的多尺度方法．SHI等［11］針對傳統(tǒng)EDFM方法中對基質(zhì)與裂縫之間的傳質(zhì)傳熱量計算精度不足的問題，提出用局部邊界元方法提高了溫度場和滲流場的計算準(zhǔn)確性．然而，目前基于EDFM方法的傳熱過程研究均是基于局部熱平衡假設(shè)（local thermal equilibrium，LTE），尚無研究考慮局部非熱平衡假設(shè)（local thermal non-equilibrium，LTNE）下的EDFM方法．

關(guān)于局部非熱平衡效應(yīng)對傳熱過程的影響已有較多研究．路朗［12］研究發(fā)現(xiàn)，多孔介質(zhì)入口位置處局部非熱平衡效應(yīng)更顯著．YANG等［13］討論了在考慮局部非熱平衡假設(shè)下，Darcy數(shù)和Biot數(shù)等對多孔介質(zhì)-流體界面熱條件的影響．曲占慶等［14］建立了基于局部非熱平衡的雙重介質(zhì)模型，討論裂縫分布和注采參數(shù)對于干熱巖采熱性能的影響，但未討論局部非熱平衡假設(shè)考慮的必要性．駱雄飛等［15］討論了內(nèi)含熱源的平行平板模型換熱過程中局部非熱平衡模型的必要性，但主要集中在特征參數(shù)的影響上．歐陽小龍［16］用Peclet數(shù)和Stanton數(shù)作為增強(qiáng)型地?zé)嵯到y(tǒng)是否考慮局部非熱平衡假設(shè)的判別依據(jù)，針對模型為一維有限空間而未考慮裂縫存在的情況，目前，尚無針對致密基質(zhì)-復(fù)雜裂縫傳熱過程中局部非熱平衡效應(yīng)影響程度的研究，相應(yīng)的儲層參數(shù)對考慮局部非熱平衡假設(shè)必要性的影響規(guī)律還未明確．

本研究聚焦增強(qiáng)型地?zé)嵯到y(tǒng)的工程背景，基于局部非熱平衡假設(shè)，從微觀孔喉模型和嵌入式離散裂縫模型兩個層面重點探究局部非熱平衡對于致密基質(zhì)-復(fù)雜裂縫系統(tǒng)傳熱過程的影響．明確了壓力梯度、基質(zhì)孔隙度、裂縫孔隙度、基質(zhì)滲透率、裂縫滲透率、對流換熱系數(shù)以及熱擴(kuò)散系數(shù)的大小對于考慮局部非熱平衡假設(shè)必要性的影響，從而為增強(qiáng)型地?zé)嵯到y(tǒng)的數(shù)值模擬工作提供參考．

1 耦合模型建立與求解

1.1 物理模型及控制方程

1.1.1 微觀傳質(zhì)傳熱耦合模型

微觀孔隙尺度的壓差傳質(zhì)可用Navier-Stokes方程進(jìn)行描述，

其中，v為流體速度矢量，單位：m/s；t為流動時間，單位：s；ρf為流體密度，單位：kg/m3；pf為流體壓力，單位：Pa；μf為流體黏度，單位：Pa·s；I為單位張量；F為作用在流體上的體積力矢量，單位：N/m3；上標(biāo)T表示對矩陣進(jìn)行轉(zhuǎn)置．

用對流擴(kuò)散方程描述微觀孔隙尺度下的濃度擴(kuò)散過程為

其中，c為流體的濃度，單位：mol/m3；D為擴(kuò)散系數(shù)，單位：m2/s，R為源匯項，單位：mol/（m3·s）．在巖石區(qū)域只考慮熱傳導(dǎo)作用，能量守恒方程為

在流體區(qū)域考慮熱傳導(dǎo)和熱對流作用的能量守恒方程為

其中，下標(biāo)r和f分別表示巖石和流體；Tr和Tf分別為巖石和流體的熱力學(xué)溫度，單位：K；ρr和ρf分別為巖石和流體的密度，單位：kg/m3；cpr和cpf分別為巖石和流體的比熱，單位：J/（kg·K），λr和λf分別為巖石和流體的導(dǎo)熱系數(shù)，單位：W/（m·K）；q為熱源項，單位：W/m3．

利用COMSOL Multiphysics軟件的PARDISO求解器進(jìn)行傳質(zhì)傳熱過程的順序耦合求解，其優(yōu)點在于占用內(nèi)存小，計算效率高．

通過掃描電鏡圖片提取微觀孔吼結(jié)構(gòu)見圖1，其中，藍(lán)色部分為流體流動通道，灰色部分為巖石骨架，2條紅色虛線之間為主要的滲流通道（裂縫區(qū)域），周圍為近裂縫基質(zhì)區(qū)域．

圖1 孔喉模型輪廓及網(wǎng)格劃分 （a）孔喉模型及主流通道；（b）網(wǎng)格劃分結(jié)果Fig.1 Pore-throat model outline and mesh generation.(a)Pore-throat model and mainstream channel(region between two red lines is the main seepage channel),and(b)results of grid division.

1.1.2 嵌入式離散裂縫模型

在嵌入式離散裂縫模型中，基質(zhì)和裂縫被劃分為兩個相對獨立的系統(tǒng)，系統(tǒng)之間的質(zhì)量和熱量傳遞分別表示為各自方程的源項和匯項，以基質(zhì)系統(tǒng)為例，基于EDFM方法的質(zhì)量守恒方程為

其中，下標(biāo)M、F和W分別指代基質(zhì)、裂縫和井；kM為基質(zhì)滲透率，單位：um2；φM為基質(zhì)孔隙度；V為基質(zhì)網(wǎng)格或裂縫網(wǎng)格內(nèi)部區(qū)域體積；dV為基質(zhì)網(wǎng)格或裂縫網(wǎng)格內(nèi)的體積微元；η為單位轉(zhuǎn)換系數(shù)；Bf為流體體積系數(shù)．式（5）等號左邊第1項為時間累計項；第2項為擴(kuò)散項；等號右側(cè)是源匯項，QMF為單位時間內(nèi)流體通過基質(zhì)網(wǎng)格和裂縫單元間傳遞的流體質(zhì)量，單位：kg/s；QMW為單位時間內(nèi)流體通過基質(zhì)網(wǎng)格和井單元間傳遞的流體質(zhì)量，單位：kg/s．

式（5）基于有限體積法的離散形式為

其中，下標(biāo)i和j分別為相鄰基質(zhì)網(wǎng)格或裂縫單元的編號；n為相鄰基質(zhì)網(wǎng)格邊數(shù)；G為幾何傳質(zhì)系數(shù)，根據(jù)圖2的4種網(wǎng)格連接類型對此物理量進(jìn)行定義．

表1列舉了不同網(wǎng)格連接方式下幾何系數(shù)的計算公式．其中，下標(biāo)g表示相交裂縫處任一裂縫單元的編號；S為網(wǎng)格接觸面的面積，單位：m2；d為裂縫網(wǎng)格中心到交點或交面的距離，單位：m；Δx、Δy和Δz分別為基質(zhì)網(wǎng)格在x、y和z方向的長度，單位：m；γ為儲層厚度，單位：m；L為裂縫單元的長度，單位：m；w為裂縫寬度，單位：m；N為基質(zhì)網(wǎng)格中包含的裂縫網(wǎng)格總數(shù)．

表1 不同網(wǎng)格連接方式下幾何系數(shù)的計算公式Table 1 The calculation formula of geometric coefficient with different grid connection types

在傳統(tǒng)EDFM方法中，對于圖2中所示的第④類連接方式，通常假設(shè)在基質(zhì)網(wǎng)格中壓力沿裂縫的垂直方向線性變化，但在致密基質(zhì)滲透率較低時，壓力傳播速度較慢，基質(zhì)網(wǎng)格中壓力分布沿垂直于裂縫方向呈線性分布的假設(shè)在模擬過程中是不合理的，特別是在生產(chǎn)初期存在一定的計算誤差．本研究采用邊界元法計算基質(zhì)與裂縫之間的傳遞的流體質(zhì)量，可大幅提高計算精度［16］，其矩陣形式為

圖2 EDFM中4種網(wǎng)格連接類型Fig.2 Four types of grid connection in EDFM.

其中，PM為基質(zhì)網(wǎng)格壓力矩陣，PF為裂縫網(wǎng)格壓力矩陣；K為轉(zhuǎn)換矩陣；其余矩陣A、B、C、E、J和H為邊界積分矩陣，其內(nèi)部元素表達(dá)式請掃描論文末頁右下角二維碼，查看補(bǔ)充材料1．井單元在圓柱形儲集層的中心注入或生產(chǎn)流體．邊界條件可以考慮定流量邊界或者定壓邊界．注采過程假定為擬穩(wěn)態(tài)徑向滲流．與井相關(guān)的質(zhì)量源匯項為

其中，rW為井的半徑，單位：m；θ為水平井與裂縫網(wǎng)格平面法向量的夾角，單位：（°）；re為等效供給半徑，單位：m，可采用Peaceman公式［18］進(jìn)行求解．式（10）表示井與基質(zhì)相連，式（11）表示井與裂縫相連．

同樣以基質(zhì)系統(tǒng)為例，給出基于嵌入式離散裂縫模型的基質(zhì)內(nèi)流體能量守恒方程為

巖石能量守恒方程為

其中，a為巖石比表面積，單位：m-1；h為對流換熱系數(shù)，單位：W/（m2·K）；Uf為流體熱傳導(dǎo)產(chǎn)生的源匯項，單位：W；Ur為巖石熱傳導(dǎo)產(chǎn)生的源匯項，單位：W．式（12）和式（13）兩個能量守恒方程的有限體積離散形式分別為

其中，H為幾何傳熱系數(shù)，計算公式見表1．

對于圖2所示的第④類連接的傳熱計算，傳統(tǒng)EDFM沿用與傳質(zhì)過程相似的線性分布假設(shè)，這同樣會導(dǎo)致基質(zhì)與裂縫之間的傳導(dǎo)傳熱量計算出現(xiàn)誤差．對傳統(tǒng)EDFM在處理致密基質(zhì)-裂縫熱流耦合問題上存在的精度問題，誤差一方面來源于第④類連接關(guān)系中基質(zhì)網(wǎng)格的壓力和溫度沿垂直于裂縫面方向呈線性變化的假設(shè)，這會導(dǎo)致基質(zhì)與裂縫之間的傳質(zhì)量和熱傳導(dǎo)傳熱量產(chǎn)生計算誤差．另一方面，基質(zhì)與裂縫之間的傳質(zhì)量又會影響基質(zhì)與裂縫之間熱對流傳熱量的計算，同時，溫度場和滲流場相互耦合又會讓誤差不斷累加．本質(zhì)上需要解決的是EDFM第④類連接關(guān)系中的傳質(zhì)傳熱量的計算．為提高計算精度，本研究采用局部邊界元法對熱傳導(dǎo)傳熱量進(jìn)行求解問題［11］．以流體能量守恒方程中的基質(zhì)與裂縫之間的熱傳導(dǎo)項的計算為例，其矩陣形式為

對于直井和水平井注水井，流體流入或流出基質(zhì)網(wǎng)格，熱對流源匯項為

對于壓裂水平井注水井，流體流入或流出裂縫網(wǎng)格，熱對流源匯項為

熱流耦合求解通常有順序隱式耦合和完全耦合全隱式兩種耦合方法．順序隱式耦合求解方法將整個非線性方程組分成壓力部分和溫度部分并依次求解，主要用于傳質(zhì)和傳熱方程之間耦合不明顯的單相系統(tǒng)［10］．完全耦合全隱式方法是在單個牛頓迭代步內(nèi)同時求解整個非線性方程組．由于流動方程和溫度方程之間的強(qiáng)耦合作用，順序隱式耦合方法通常會由于許多外環(huán)順序迭代而表現(xiàn)不佳．因此，本研究考慮完全隱式耦合，采用牛頓迭代法隱式求解整個方程組．模型的驗證請掃描論文末頁右下角二維碼查看補(bǔ)充材料2．

1.2 熱流耦合機(jī)理

在地下儲層的高壓高溫環(huán)境中，水相的黏度μf、密度ρf、導(dǎo)熱系數(shù)λf和比熱容cpf不再恒定，這將極大地改變流體的質(zhì)量和能量傳遞過程．通過式（19）—式（22）［19-20］計算水相的相關(guān)物理量．

2 LTE和LTNE模擬結(jié)果對比

2.1 微觀傳質(zhì)傳熱模擬結(jié)果

微觀傳質(zhì)傳熱模型仿真模擬參數(shù)見表2．

表2 微觀傳質(zhì)傳熱模型主要物理量及取值Table 2 Values of main physical quantities of the microscopic mass and heat transfer model

圖3為t=0.01 s時考慮局部熱平衡假設(shè)和局部非熱平衡假設(shè)下流體和巖石的溫度分布．其中，從圖3（a）和（b）的局部熱平衡假設(shè)模擬結(jié)果可見，流體與接觸巖石表面溫度相同，即流體流過巖石表面的瞬間達(dá)到了熱平衡狀態(tài)，體現(xiàn)為巖石顆粒在較短的時間內(nèi)被冷卻．圖3（c）和（d）為與局部熱平衡假設(shè)相同模擬條件下，考慮局部非熱平衡假設(shè)時的模擬結(jié)果，可以看到流體與接觸巖石表面溫度并不相同，流體的冷卻速度較局部熱平衡假設(shè)更快，而巖石顆粒冷卻過程更慢，原因是巖石向流體傳遞的能量變少．

圖3 t=0.01 s時的溫度分布（a）LTE流體溫度場；（b）LTE巖石溫度場；（c）LTNE流體溫度場；（d）LTNE巖石溫度場Fig.3 Temperature distribution(t=0.01 s).(a)Fluid temperature field(LTE),(b)rock temperature field(LTE),(c)fluid temperature field(LTNE),and(d)rock temperature field(LTNE).

通過對比上述兩個案例的溫度分布情況，可以看出局部熱平衡假設(shè)低估了在入口處和近裂縫區(qū)域的巖石溫度，而高估了傳熱前緣區(qū)域的流體溫度．

2.2 復(fù)雜裂縫模型模擬結(jié)果

通過一個多裂縫傳質(zhì)傳熱模型（圖4）研究復(fù)雜裂縫條件下的局部非熱平衡假設(shè)的影響．模型大小為100 m×100 m×10 m，流體黏度、密度、導(dǎo)熱系數(shù)和比熱隨溫度變化，原始地層溫度為180℃，在模型左側(cè)模擬水平井定壓30 MPa注入，在模型右側(cè)模擬水平井定壓10 MPa采出．其余主要的模擬參數(shù)見表3．

表3 復(fù)雜裂縫模型主要參數(shù)及取值Table 3 Main parameters and values of the complex fracture model

圖4 復(fù)雜裂縫模型示意圖Fig.4 Schematic diagram of complex fracture model.Red line is the injection end,light blue line is the produced end,and blue lines are the fractures.

模擬結(jié)果如圖5，在局部非熱平衡假設(shè)下，巖石和流體的溫度分布之間長時間尺度上的差別不大，而局部熱平衡模型和局部非熱平衡模型之間存在較大差異，局部熱平衡實際上高估了產(chǎn)出流體溫度，這對于增強(qiáng)型地?zé)嵯到y(tǒng)來說，可能會導(dǎo)致系統(tǒng)的生命周期被高估，從而帶來經(jīng)濟(jì)上的損失．因此，局部非熱平衡模型中巖石和流體溫度之間的微小差異，不一定會導(dǎo)致局部熱平衡模型和局部非熱平衡模型之間的溫度差異小，有必要對局部非熱平衡模型的適用性和敏感性進(jìn)一步研究．

圖5 局部熱平衡和局部非熱平衡條件下的溫度場（t=1000 d）（a）LTE溫度場；（b）LTNE流體溫度場；（c）LTNE巖石溫度場Fig.5 Schematic diagram of multi-fracture model(t=1000 d).(a)LTE temperature field,(b)LTNE fluid temperature field,and(c)LTNE rock temperature field.

3 參數(shù)敏感性分析

3.1 微觀孔喉模型參數(shù)敏感性

圖6為不同壓力梯度下流體和巖石在局部熱平衡和局部非熱平衡假設(shè)條件下的溫度場分布．由圖6可見，局部熱平衡假設(shè)比局部非熱平衡假設(shè)流體冷卻巖石的速度要快．對比來看，局部熱平衡假設(shè)會高估產(chǎn)出流體的溫度，低估巖石的溫度，且壓力梯度越小，局部非熱平衡效應(yīng)越不明顯．這是因為壓力梯度越小，流體流動速度越慢，使得流體和巖石有更多的時間進(jìn)行對流換熱，從而使接觸面上的溫度達(dá)到平衡．

圖6 t=0.01 s時不同壓力梯度下巖石和流體溫度分布（a）和（b）分別為1.0 MPa/m時LTE和LTNE巖石溫度場；（c）和（d）分別為1.0 MPa/m時LTE和LTNE流體溫度場；（e）和（f）分別為0.5 MPa/m時LTE和LTNE巖石溫度場；（g）和（h）分別為0.5 MPa/m時LTE和LTNE流體溫度場；（i）和（j）分別為0.2 MPa/m時LTE和LTNE巖石溫度場；（k）和（l）分別為0.2 MPa/m時LTE和LTNE流體溫度場Fig.6 Temperature distribution of fluid and rock with different pressure gradients(t=0.01 s).(a)1.0 MPa/m,LTE,rock temperature distribution,(b)1.0 MPa/m,LTNE,rock temperature distribution,(c)1.0 MPa/m,LTE,fluid temperature distribution,(d)1.0 MPa/m,LTNE,fluid temperature distribution,(e)0.5 MPa/m,LTE,rock temperature distribution,(f)0.5 MPa/m,LTNE,rock temperature distribution,(g)0.5 MPa/m,LTE,fluid temperature distribution,(h)0.5 MPa/m,LTNE,fluid temperature distribution,(i)0.2 MPa/m,LTE,rock temperature distribution,(j)0.2 MPa/m,LTNE,rock temperature distribution,(k)0.2 MPa/m,LTE,fluid temperature distribution,and(l)0.2 MPa/m,LTNE,fluid temperature distribution.

由于對流換熱系數(shù)h在實際中往往難以準(zhǔn)確測定，本研究對比了不同對流換熱系數(shù)下的局部非熱平衡效應(yīng)．圖7為不同對流換熱系數(shù)下的溫度場分布．對比局部熱平衡結(jié)果可以看出，當(dāng)h較小時，局部熱非平衡效應(yīng)影響較大，隨著h的增大，兩種假設(shè)的溫度場逐漸趨于相似．

圖7 局部熱平衡假設(shè)與局部非熱平衡假設(shè)不同對流換熱系數(shù)下溫度場（t=0.01 s）（a）LTE；（b）h=5×103 W/(m2·K)；（c）h=1×104 W/(m2·K)；（d）h=2×104 W/(m2·K)Fig.7 Comparison of different surface heat transfer coefficient and local thermal equilibrium temperature field(t=0.01 s).(a)LTE,(b)h=5×103 W/(m2·K),(c)h=1×104 W/(m2·K),and(d)h=2×104 W/(m2·K).

3.2 復(fù)雜裂縫模型參數(shù)敏感性

針對復(fù)雜裂縫模型，基于EDFM方法對不同儲層參數(shù)的局部非熱平衡效應(yīng)大小及隨時間變化規(guī)律進(jìn)行分析．模型的基礎(chǔ)參數(shù)如下：基質(zhì)孔隙度為0.1，基質(zhì)滲透率為2×10-4μm2，裂縫孔隙度為0.2，裂縫滲透率為20 μm2，巖石導(dǎo)熱系數(shù)為3 W/（m·K），巖石比熱為850 J/（kg·K），巖石密度為2.6×103kg/m3，ah=500 W/（m3·K），在此基礎(chǔ)上通過控制單一變量進(jìn)行研究．定義一個流體平均溫度誤差值（由式（23）計算）量化整個基質(zhì)裂縫區(qū)域局部非熱平衡效應(yīng)的大小，該值越大表示局部熱平衡假設(shè)對基質(zhì)和裂縫平均流體溫度的高估程度越大．

3.2.1 不同基質(zhì)孔隙度和裂縫孔隙度的影響

圖8和圖9分別為不同基質(zhì)孔隙度和不同裂縫孔隙度相對誤差．從圖8和圖9可見，局部熱平衡假設(shè)對儲層中流體溫度場的高估幅度先上升后下降，原因是低溫區(qū)是逐漸波及又趨于穩(wěn)定的，最終區(qū)域內(nèi)巖石溫度降低至注入流體溫度．同時可以看到，基質(zhì)孔隙度越小，局部非熱平衡假設(shè)的影響越顯著．裂縫孔隙度越高，局部非熱平衡假設(shè)的影響越顯著．可以推斷，當(dāng)基質(zhì)孔隙度和裂縫孔隙度差別越大，越需要考慮局部非熱平衡假設(shè)．總體而言，孔隙度對局部非熱平衡效應(yīng)影響較大，原因是在LTE模型中，孔隙度影響了體積平均參數(shù)的大小，比如巖石流體的平均密度、平均比熱和平均導(dǎo)熱系數(shù)，而由式（12）和式（13）可知，在LTNE模型中，孔隙度對巖石及流體熱量的傳遞均有影響．

圖8 不同基質(zhì)孔隙度相對誤差Fig.8 Relative errors of different matrix porosity.

圖9 不同裂縫孔隙度相對誤差Fig.9 The relative errors of different fracture porosity.

3.2.2 不同基質(zhì)滲透率和裂縫滲透率的影響

圖10和圖11分別為不同基質(zhì)滲透率和不同裂縫滲透率相對誤差．從圖10和圖11可見，兩種模型計算結(jié)果的誤差最大值隨基質(zhì)滲透率和裂縫滲透率的降低而增加．從時間尺度上總體來看，隨著基質(zhì)滲透率和裂縫滲透率的增大，早期局部非熱平衡效應(yīng)越顯著，但在中后期，基質(zhì)滲透率和裂縫滲透率越低局部非熱平衡效應(yīng)越顯著，這是由于隨著滲透率的增大，滲流速度增大導(dǎo)致的熱對流作用增強(qiáng)，使得巖石被冷卻的速度越快，在早期，流速越快巖石和流體換熱越不充分導(dǎo)致兩者差異越大，而隨著時間推移，兩種模型計算結(jié)果的差異性逐漸變小，且局部非熱平衡效應(yīng)對基質(zhì)滲透率的敏感程度大于裂縫滲透率．

圖10 不同基質(zhì)滲透率相對誤差Fig.10 Relative errors of porosity of different matrix.

圖11 不同裂縫滲透率相對誤差Fig.11 The relative errors of different fracture porosity.

3.2.3 不同熱擴(kuò)散系數(shù)和對流換熱系數(shù)的影響

從不同熱擴(kuò)數(shù)系數(shù)相對誤差（圖12）和不同對流換熱系數(shù)相對誤差（圖13）可以看出，巖石熱擴(kuò)散系數(shù)(α=λr/(cpr·ρr))越大，局部非熱平衡效應(yīng)越顯著．對流換熱系數(shù)越小，局部非熱平衡效應(yīng)越顯著，這是由于對流換熱系數(shù)越小，巖石表面與流體的換熱能力越差．因此，當(dāng)巖石的熱擴(kuò)散系數(shù)較大或?qū)α鲹Q熱系數(shù)較小時，需要考慮局部非熱平衡模型．

圖12 不同熱擴(kuò)散系數(shù)相對誤差Fig.12 Relative errors of different thermal diffusivity matrix.

圖13 不同對流換熱系數(shù)相對誤差Fig.13 Relative errors of different convective heat transfer coefficient.

4 結(jié)論

考慮局部非熱平衡假設(shè)，建立了基于有限元方法的微觀傳質(zhì)傳熱的毫米尺度數(shù)值實驗?zāi)Ｐ?，和基于局部邊界元和有限體積方法的宏觀嵌入式離散裂縫傳質(zhì)傳熱模型，形成了致密基質(zhì)-復(fù)雜縫網(wǎng)傳質(zhì)傳熱數(shù)值模擬方法，并分析了局部非熱平衡假設(shè)對傳熱過程的影響．研究可知：

1）所提出的改進(jìn)EDFM方法可以有效模擬考慮局部非熱平衡假設(shè)的致密基質(zhì)-復(fù)雜裂縫的傳質(zhì)傳熱過程，考慮局部非熱平衡假設(shè)能提高基質(zhì)-裂縫傳質(zhì)傳熱過程模擬的準(zhǔn)確性；

2）局部熱平衡假設(shè)會低估入口處和近裂縫區(qū)域的巖石溫度，而高估傳熱前緣區(qū)域的流體溫度；

3）局部非熱平衡對增強(qiáng)型地?zé)嵯到y(tǒng)的影響主要體現(xiàn)在注入早期，而注入強(qiáng)度越高、或基質(zhì)裂縫滲透率越低、或基質(zhì)裂縫孔隙度差異越大、或巖石熱擴(kuò)散系數(shù)越大、或?qū)α鲹Q熱系數(shù)越小，都越需要考慮局部非熱平衡假設(shè).

事實上，影響對流傳熱過程的諸多因素均體現(xiàn)在對流換熱系數(shù)的確定上，模型中將其考慮成常數(shù)存在局限性，未來工作中需進(jìn)一步探究．