張帥領(lǐng),張春艷,江承陽,楊亞東

(1.中國電建集團(tuán)河南省電力勘測(cè)設(shè)計(jì)院有限公司,河南 鄭州 450007;2.華北水利水電大學(xué) 地球科學(xué)與工程學(xué)院,河南 鄭州 450046)

0 引言

裂隙作為巖溶含水系統(tǒng)的重要組成部分,對(duì)巖溶水流運(yùn)動(dòng)起到重要作用[1]。裂隙水流運(yùn)動(dòng)規(guī)律的研究對(duì)巖溶水流運(yùn)動(dòng)規(guī)律、溶質(zhì)運(yùn)移規(guī)律以及裂隙巖體滲流理論的發(fā)展起到重要作用[2-4]。裂隙水流運(yùn)動(dòng)一般可以用納維-斯托克斯方程(簡(jiǎn)稱NS方程)來描述[5]。但是由于NS方程解的復(fù)雜性,限制了方程的實(shí)際應(yīng)用。眾多研究者按照水流運(yùn)動(dòng)實(shí)際情況,對(duì)NS方程進(jìn)行了不同程度的簡(jiǎn)化,第一層將NS方程化為斯托克斯(Stokes)方程;第二層將Stokes方程簡(jiǎn)化為雷諾方程;第三層將雷諾方程簡(jiǎn)化為立方定律(CL),由于立方定律的簡(jiǎn)潔性,得到了廣泛的引用。文章以具有四個(gè)裂隙面的矩形狹縫裂隙為研究對(duì)象(如圖1),探討立方定律在矩形狹縫裂隙中的適用性以及矩形狹縫裂隙水流運(yùn)動(dòng)特性。

圖1 矩形狹縫平行裂隙示意圖

NS方程化為Stokes方程的簡(jiǎn)化條件為:慣性力可以忽略,即流速較小,粘性力占主導(dǎo)地位,雷諾數(shù)Re較小。研究者主要從雷諾數(shù)出發(fā),研究了Stokes方程的適用條件。Sharp[6],Iwai[7]以及Schrauf[8]利用野外裂隙樣本通過室內(nèi)試驗(yàn)分析得到,當(dāng)雷諾數(shù)大于1到10的某個(gè)數(shù)值時(shí),慣性力不可忽略。Zimmerman[9]提出當(dāng)Re(b)/∧b(波長(zhǎng))小于8時(shí),慣性力對(duì)滲流量修正系數(shù)的影響最高達(dá)10%。與NS方程相比,Stokes方程有了很大程度的簡(jiǎn)化,但是三維Stokes方程的解仍非常復(fù)雜。

對(duì)LCL進(jìn)行簡(jiǎn)化,簡(jiǎn)化條件為兩裂隙面光滑且平行分布(即uy=uz=0),則三維的流體流動(dòng)可概化為一維流動(dòng),設(shè)流體沿x軸流動(dòng),則沿x軸方向的裂隙滲流量為:

Q=Cb3J

(1)

式中:Q為裂隙間的流量,m3·s-1;b為裂隙開度,m;J為水力坡度,無量綱;C為常數(shù),1 s-1。

方程(1)由于流量Q與裂隙開度b的三次方呈正比,故稱為立方定律(Cubic Law,CL)(Boussinesq[16])。對(duì)于直線流,常數(shù)C可表達(dá)為:

(2)

聯(lián)合方程(1)和(2),可得裂隙滲流量方程為:

(3)

式中:W為裂隙寬度,m;L為裂隙長(zhǎng)度,m;ρ為流體密度,kg·m-3;g為重力加速度,m·s-2;μ為動(dòng)力黏滯系數(shù),N·s·m-2。

將流量方程Q=Wbv以及水力坡度與水頭損失之間的關(guān)系J=hf/L代入式(3),得到立方定律的另一種表達(dá)形式:

(4)

式中:hf為水流通過裂隙時(shí)產(chǎn)生的水頭損失,m;v為水流流速,m·s-1。

式(4)表明,恒溫條件下,水流流經(jīng)光滑裂隙時(shí)的水頭損失與水流流速之間呈線性關(guān)系。

由于立方定律的簡(jiǎn)潔性,立方定律得到了廣泛應(yīng)用[17],但由于立方定律的局限性,國內(nèi)外許多學(xué)者主要從裂隙開度[18-19]、裂隙壁面粗糙性[20]、節(jié)理粗糙度系數(shù)修正[21]、面積接觸率[22]等方面進(jìn)行了研究,并基于各自的研究對(duì)象提出了立方定律的修正形式[23]。綜合國內(nèi)外研究成果,絕大多數(shù)研究者注重于開放裂隙中立方定律的應(yīng)用,在實(shí)際野外地區(qū),完全開放的裂隙較少存在,矩形狹縫裂隙廣泛存在,針對(duì)這一問題,以過水?dāng)嗝鏋榫匦蔚莫M縫裂隙(如圖1所示)為研究對(duì)象,通過數(shù)值模擬,對(duì)比NS的解,對(duì)立方定律的有效性進(jìn)行研究,提出以極限流速及極限雷諾數(shù)來衡量立方定律是否適用,提出適用于矩形狹縫裂隙的修正立方定律,并對(duì)裂隙流的水力特性進(jìn)行了研究,該研究不僅對(duì)以往的理論研究進(jìn)行一定的補(bǔ)充,而且為巖溶含水系統(tǒng)水流運(yùn)動(dòng)特征的研究提供理論支撐。

1 矩形狹縫裂隙數(shù)值模型的建立

1.1 控制方程

非穩(wěn)定流條件下,對(duì)于不可壓縮流體,NS方程[24]為:

(5)

式中:f為體積力,N;p為壓力,Pa;t為時(shí)間,s。

穩(wěn)定流情況下,對(duì)于不可壓縮流體,NS方程為:

ρ(v)v=μ2v-p

(6)

對(duì)于兩光滑平板間的流體,其運(yùn)動(dòng)符合泊肅葉定律,設(shè)裂隙兩側(cè)壁分布在z=±b/2處,那么流速分布表達(dá)式[25]為:

(7)

1.2 邊界條件

對(duì)不同尺寸的裂隙進(jìn)行二維和三維模擬,不同模擬情景下,邊界條件及裂隙尺寸的設(shè)置見表1。

表1 不同模擬情景下邊界條件、裂隙尺寸設(shè)置

1.3 模型驗(yàn)證(S0模擬結(jié)果)

對(duì)尺寸為W=2 cm,b=0.840 0 mm,L=20 cm的裂隙中的水流運(yùn)動(dòng)進(jìn)行數(shù)值模擬,數(shù)值模型設(shè)置溫度T=13℃,水的密度ρ=999.4 kg/m3,動(dòng)力粘度μ=0.001 2 Pa·s,重力加速度g=9.8 m/s2,邊界條件的設(shè)置見表1中模擬情景S0。

水力坡度試驗(yàn)測(cè)量值與數(shù)值模擬值對(duì)比如圖2所示。經(jīng)計(jì)算,試驗(yàn)測(cè)量值與數(shù)值模擬值相對(duì)誤差((數(shù)值模擬值-試驗(yàn)測(cè)量值)/試驗(yàn)測(cè)量值)在10%以內(nèi)。

圖2 試驗(yàn)測(cè)量值與數(shù)值模擬值對(duì)比圖

2 模型結(jié)果與分析

2.1 立方定律適用性驗(yàn)證(S1模擬結(jié)果)

數(shù)值模擬得到7個(gè)不同裂隙寬度(W=0.5,1,2,3,4,5,6 cm),6個(gè)不同裂隙開度(b=0.05,0.06,0.07,0.08,0.09,0.1 cm),共42個(gè)不同尺寸的裂隙水流模擬結(jié)果。以數(shù)值模擬所得水力坡度與立方定律計(jì)算所得水力坡度的比值(J-NS/J-CL)為研究對(duì)象(J-NS與J-CL對(duì)比圖見附件1),來確定立方定律的適用范圍?，F(xiàn)將不同裂隙寬度(W),不同裂隙開度(b)條件下,J-NS/J-CL隨流速變化作對(duì)比分析,如圖3所示(以W=1,2 cm為例,其余結(jié)果見附件2)。結(jié)果表明,流速較小時(shí),J-NS與J-CL重合度較高,表明流速較小時(shí),可用立方定律來描述裂隙水流,然而隨著流速的增大,J-CL與J-NS之間的差值越來越大,且J-CL均小于J-NS,表明流速較大時(shí),立方定律計(jì)算值與實(shí)際值偏離越大,立方定律不可以再用來描述裂隙水流運(yùn)動(dòng)。

圖3 不同裂隙寬度及開度條件下,水力坡度NS值與CL值之比J-NS/J-CL隨流速變化

結(jié)果表明,在本次研究中(裂隙寬度分布于0.5～6 cm之間,裂隙開度分布于0.05～0.1 cm之間),J-NS/J-CL分布于0.9～2之間,J-CL與J-NS的差值隨流速的增大而增大,隨裂隙開度的增大而增大,隨裂隙寬度的增大而減小,因?yàn)榱魉僭酱?慣性力的作用越不可忽略。為確定立方定律的適用范圍,定義相對(duì)偏差為((CL計(jì)算值-NS方程模擬值)/NS方程模擬值),相對(duì)偏差達(dá)到10%時(shí)對(duì)應(yīng)的流速稱為極限流速(vlim),對(duì)應(yīng)的雷諾數(shù)稱為極限雷諾數(shù)(Relim)。將不同裂隙寬度條件下,極限流速以及極限雷諾數(shù)隨裂隙開度的變化作對(duì)比分析,分別如圖4(a)(b)所示?？傮w來講,極限流速或極限雷諾數(shù)隨裂隙寬度的增大而增大,隨裂隙開度的增大而減小,說明裂隙寬度越大,開度越小,立方定律越適用。本次研究中,極限流速的最大值為30.08 cm/s,極限雷諾數(shù)的最大值為75.2。

(a)極限流速 (b)極限雷諾數(shù)

2.2 修正立方定律的提出

(8)

式中的n值可由S1模擬結(jié)果擬合確定(NS方程模擬值與修正立方定律計(jì)算值對(duì)比結(jié)果見附件3)。不同裂隙寬度條件下,n值隨裂隙開度的變化如圖5所示。結(jié)果表明,裂隙開度越大,n值越大;總體來講,裂隙寬度越大,n值越小,本次研究中,n值分布于1.4～1.9。

圖5 不同裂隙寬度條件下,n值隨裂隙開度的變化

2.3 矩形狹縫平行裂隙滲透系數(shù)分析

根據(jù)修正立方定律和達(dá)西定律,矩形狹縫平行裂隙滲透系數(shù)可以用下式計(jì)算:

(9)

分析矩形狹縫裂隙滲透系數(shù)隨裂隙寬度、裂隙開度以及裂隙寬度與裂隙開度比值的變化(分別如圖6(a)(b)(c)所示),同一裂隙開度條件下,滲透系數(shù)隨裂隙寬度的增大而增大,但增大幅度相對(duì)較小。隨裂隙開度的增大而增大,增大幅度相對(duì)較大。滲透系數(shù)隨裂隙寬度與裂隙開度比值的增大而增大,通過曲線擬合得到二者之間呈指數(shù)關(guān)系。本次研究中,矩形狹縫裂隙滲透系數(shù)分布于10cm·s-1-70 cm·s-1。

圖6 滲透系數(shù)隨裂隙寬度、裂隙開度以及裂隙寬度與裂隙開度比值的變化圖

2.4 流速分布剖面分析

2.4.1 裂隙開度對(duì)流速分布剖面影響分析(S2模擬結(jié)果)

不同裂隙開度條件下,X=10 cm處流速分布如圖7所示(以b=0.05,0.1 cm為例,其余結(jié)果見附件4)。從圖中可以看出,二維裂隙面流速剖面形狀變化與理想泊肅葉流體流速變化形狀基本一致,數(shù)值上有所差異。流速的最大值位于裂隙中心線處。流速最小值位于裂隙壁處且為零。分析不同開度條件下,流速模擬值與泊肅葉值之差(絕對(duì)偏差)及相對(duì)偏差沿質(zhì)點(diǎn)位置(Z軸)的變化(如圖8所示),絕對(duì)偏差隨裂隙開度的增大而增大,且離裂隙中心線(X軸)越近,偏差越大。相對(duì)偏差沿質(zhì)點(diǎn)位置(Z軸)呈波動(dòng)變化,變幅很小;相對(duì)偏差隨裂隙開度的增大而增大,表明裂隙開度越大,流速分布越不符合標(biāo)準(zhǔn)拋物線分布。

圖7 不同裂隙開度條件下,流速剖面分布

圖8 不同開度條件下,流速模擬值與泊肅葉值的絕對(duì)偏差及相對(duì)偏差沿質(zhì)點(diǎn)位置(Z軸)的變化圖

2.4.2 流速大小對(duì)流速分布影響分析(S3模擬結(jié)果)

不同流速條件下,流速分布剖面如圖9所示。從圖中可以看出,無論流速大小,模擬所得流速分布大致符合拋物線分布,流速最大值位于裂隙中心線處(即X軸上),離裂隙壁越近,流速越小,在裂隙壁處流速為零。分析不同流速條件下,流速模擬值與泊肅葉值之差(絕對(duì)偏差)及相對(duì)偏差變化(如圖10所示),流速越大,絕對(duì)偏差越大,在裂隙中心線處最大,離裂隙壁越近,偏差越小。而相對(duì)偏差隨流速的變化幅度相對(duì)較小,表現(xiàn)為集中型。

圖9 不同流速條件下,流速剖面分布對(duì)比(NS值與泊肅葉值對(duì)比圖)

圖10 不同流速條件下,流速模擬值與泊肅葉值的絕對(duì)偏差及相對(duì)偏差沿質(zhì)點(diǎn)位置(Z軸)的變化圖

2.4.3 裂隙寬度對(duì)流速分布影響分析(S4模擬結(jié)果)

取XY面中心線處Z=b/2,沿X=L/2質(zhì)點(diǎn)流速為研究對(duì)象,不同裂隙寬度條件下,流速分布如圖11所示,從圖中可以看出流速值沿Y軸波動(dòng)較小,在裂隙壁處驟減為0,隨著裂隙寬度的增大,流速波動(dòng)幅度變小,裂隙寬度大于等于3 cm(對(duì)應(yīng)的裂隙寬度與開度之比W/b≥40)的四條曲線重合度較高,表明隨著裂隙寬度及裂隙寬度與開度之比的增大,裂隙寬度的大小對(duì)流速分布曲線的影響減弱。

圖11 不同裂隙寬度條件下,XY面中心線處Z=b/2,沿X=L/2質(zhì)點(diǎn)流速分布圖

3 結(jié)語

以矩形狹縫平行單裂隙為研究對(duì)象,對(duì)不同尺寸的矩形狹縫平行裂隙進(jìn)行NS方程模擬計(jì)算,并將模擬結(jié)果與CL計(jì)算結(jié)果進(jìn)行對(duì)比,首先對(duì)立方定律在矩形狹縫平行裂隙中的適用性進(jìn)行驗(yàn)證,發(fā)現(xiàn)流速較小時(shí),CL值與NS值吻合較好,流速越大,立方定律越不適用,進(jìn)而提出極限流速和極限雷諾數(shù)的概念來確定立方定律在矩形狹縫裂隙中的適用范圍,在本次研究中,極限流速的最大值為30.08 cm/s,極限雷諾數(shù)的最大值為75.2;基于修正立方定律和達(dá)西定律,提出矩形狹縫裂隙中立方定律的修正形式,式中的參數(shù)n值分布與1.4～1.9之間,矩形狹縫裂隙滲透系數(shù)隨裂隙寬度的變化幅度較小,隨裂隙開度的增大而增大,隨裂隙寬度與裂隙開度比值的增大而增大。分析流速剖面分布隨裂隙開度、裂隙寬度以及流速的變化規(guī)律,結(jié)果表明裂隙開度越大,流速越大,流速分布越不符合標(biāo)準(zhǔn)拋物線分布;裂隙寬度越大,對(duì)流速分布曲線的影響越弱。