劉鳳臣，于學(xué)峰，賈超，杜圣賢，宋香鎖，楊斌，張尚坤

(1.山東省地質(zhì)科學(xué)研究院，國土資源部金礦成礦過程與資源利用重點(diǎn)實(shí)驗(yàn)室，山東省金屬礦產(chǎn)成礦地質(zhì)過程與資源利用重點(diǎn)實(shí)驗(yàn)室，山東 濟(jì)南 250013；2.山東大學(xué)土建與水利學(xué)院，山東 濟(jì)南 250013)

0 引言

古生物化石保護(hù)得到人們?cè)絹碓蕉嗟年P(guān)注，恐龍化石在古生物化石中最為重要[1]，然而它同樣經(jīng)歷著大自然的風(fēng)化破壞，研究發(fā)現(xiàn)，溫度和水是其中的兩大主要因素[2-4]。

巖體溫度場(chǎng)與滲流場(chǎng)耦合分析是巖體水力學(xué)科中重要的研究內(nèi)容，在大量的人類工程及開采活動(dòng)中，如地下水開采、石油開采、地?zé)衢_發(fā)、地質(zhì)災(zāi)害治理、大型水利工程等領(lǐng)域得到廣泛的應(yīng)用[5-18]?？铸埢皣鷰r受日照影響，在其表面和內(nèi)部形成溫度場(chǎng)，同時(shí)化石及圍巖內(nèi)部本身存在一定孔隙水，形成滲流場(chǎng)。該文將巖體水力學(xué)中溫度場(chǎng)與滲流場(chǎng)的耦合關(guān)系運(yùn)用到恐龍化石及圍巖中，研究恐龍化石及圍巖中溫度場(chǎng)與滲流場(chǎng)如何相互影響、相互作用對(duì)化石及圍巖造成風(fēng)化破壞。

1 溫度場(chǎng)與滲流場(chǎng)耦合基本方程

該文滲流場(chǎng)采用多孔介質(zhì)布里克曼方程[19-24]，與溫度場(chǎng)相結(jié)合。

(1)多孔介質(zhì)模型傳熱的數(shù)學(xué)模型為：

(1)

(ρCp)eq=θpρpCρ·p+(1-θp)ρCρ

(2)

keq=θp+(1-θp)k

(3)

式中：ρ—液體密度(kg/m3)；Cp—液體等壓熱容(J/K)；k—熱傳導(dǎo)系數(shù)(W/(m2K))；u—速度場(chǎng)；Y—比熱比，等壓熱容與等容熱容之比；(ρCp)eq—等效體積等壓熱容；θp—材料的體積分?jǐn)?shù)；kp—固體材料的熱傳導(dǎo)系數(shù)；ρp—固體材料密度(kg/m3)；Cρ·p—固體材料的比熱容(J/(kg·K))。

(2)多孔介質(zhì)流動(dòng)的數(shù)學(xué)模型是由連續(xù)性方程和動(dòng)量平衡方程組合而成，它們一起構(gòu)成了Brinkman方程：

(4)

(5)

式中：μ—流體的動(dòng)態(tài)粘滯系數(shù)(Pa·s)；u—流體的速度；ρ—流體的密度(kg/m3)；p—壓力(Pa)；εp—孔隙率；k—滲透系數(shù)(m/d)；Qbr—質(zhì)量源；F—體積力。

(3)溫度場(chǎng)與滲流場(chǎng)的耦合：在溫度場(chǎng)的計(jì)算中，速度采用布里克曼方程計(jì)算得出的速度；而在滲流場(chǎng)計(jì)算中，體積力計(jì)算公式為：

F=ρgβ(T-Tc)

(6)

式中：T—溫度場(chǎng)計(jì)算得出的溫度；Tc—初始溫度。

2 溫度場(chǎng)與滲流場(chǎng)耦合數(shù)值模型建立與條件設(shè)定

根據(jù)恐龍化石的簡化模型和圍巖的簡化模型，建立熱流耦合數(shù)值模型如圖1所示[2]。其邊界條件為：內(nèi)部初始溫度采用諸城地區(qū)多年平均氣溫12℃，左右兩邊界對(duì)稱，前后兩邊界對(duì)稱，上表面裸露于空氣中。

圖1 熱流耦合數(shù)值模型 (據(jù)杜圣賢等，2015)

隨著時(shí)間變化，外界溫度變化過程如圖2所示。其中包括溫度上升和下降兩個(gè)階段。0～30000s為溫度上升階段，30000～60000s為溫度下降階段，之后溫度又上升。

圖2 時(shí)間-溫度變化曲線

熱流耦合模型參數(shù)如表1所示。

3 溫度場(chǎng)與滲流場(chǎng)耦合數(shù)值模型結(jié)果分析

3.1 隨外界溫度變化化石及圍巖表面滲流場(chǎng)變化

隨著外界溫度升高，分別模擬了0.5h(1800s)，1h(3600s)，2h(7200s)，3h(10800s)，6h(21600s)，8h(28800s)的變化情況(圖3)。圖中等值線為滲流場(chǎng)流速梯度的變化情況，矢量圖為速度矢量的變化情況。可以看出，在外界溫度變化的前8個(gè)小時(shí)(0～28800s)，隨著外界溫度不斷升高，表面流場(chǎng)速度梯度迅速變化，速度矢量也不斷變化。在0.5個(gè)小時(shí)(1800s)時(shí)，化石表面流場(chǎng)速度梯度最大，流速也是化石兩端最大，而此時(shí)圍巖表面流速梯度普遍變化較小。這是因?yàn)楫?dāng)外界溫度變化時(shí)，由于化石裸露面積較大，化石內(nèi)部孔隙水運(yùn)動(dòng)劇烈導(dǎo)致的；第1個(gè)小時(shí)(3600s)時(shí)，明顯可以看出，在恐龍化石和圍巖接觸處，流場(chǎng)流速梯度和流速均明顯大于其他部分，這是由于化石和圍巖的各力學(xué)性質(zhì)和水力學(xué)性質(zhì)不同，導(dǎo)致它們之間存在相互作用，故在外界條件變化時(shí)，化石和圍巖接觸處受影響最大，此處風(fēng)化速度會(huì)比其他地方快。

隨著外界溫度溫度開始降低時(shí)，分別模擬了13h(46800s)，15h(54000s)，17h(61200s)，24h(86400s)的變化情況(圖4)?？梢钥闯?，外界溫度開始降低后，最先受影響的也是化石裸露區(qū)，然后擴(kuò)展到圍巖部分。15(54000s)～17h(61200s)，通過流速梯度和流速大小變化圖形可以看出，外界溫度降低的情況下，化石兩端和內(nèi)部速度梯度和速度大小非常不均勻，這種情況更會(huì)加速化石的風(fēng)化。根據(jù)總體變化趨勢(shì)可以發(fā)現(xiàn)，圍巖部分內(nèi)部空隙水隨外界溫度變化而變化持續(xù)的時(shí)間更長，故當(dāng)外界條件，特別是水和溫度持續(xù)變化時(shí)，圍巖的風(fēng)化速率可能會(huì)大于化石?？傊?，隨外界溫度的不斷變化，化石和圍巖的表面流場(chǎng)劇烈變化，孔隙水的運(yùn)動(dòng)會(huì)破壞化石和圍巖的內(nèi)部結(jié)構(gòu)，導(dǎo)致其強(qiáng)度降低；水的作用會(huì)帶動(dòng)內(nèi)部化學(xué)離子的運(yùn)動(dòng)，加速化學(xué)風(fēng)化的進(jìn)行。

表1 熱流耦合模型參數(shù)

圖3 溫度上升時(shí)恐龍化石及圍巖表面滲流場(chǎng)變化圖

圖4 溫度降低時(shí)恐龍化石及圍巖表面滲流場(chǎng)變化圖

3.2 隨外界溫度變化化石及圍巖內(nèi)部滲流場(chǎng)變化

隨著外界溫度升高時(shí)，分別模擬了1h(3600s)，3h(10800s)的變化情況(圖5)。可以看出，外界溫度變化1h(3600s)時(shí)，內(nèi)部流場(chǎng)在恐龍化石和圍巖接觸處變化最大，在圍巖內(nèi)部也有幾處變化較大區(qū)域，但比例較?。浑S著外界溫度的不斷升高，到第3個(gè)小時(shí)(10800s)時(shí)，內(nèi)部流場(chǎng)在化石和圍巖中流速均較大，即內(nèi)部空隙水運(yùn)動(dòng)劇烈。

隨著外界溫度開始降低直至再次升高時(shí)，分別模擬了12h(43200s)，18h(64800s)的變化情況(圖6)?？梢钥闯觯S著外界溫度的降低，第12個(gè)小時(shí)(43200s)時(shí)，恐龍化石及圍巖內(nèi)部流速都有所降低，但部分區(qū)域仍有較大的流速梯度；同時(shí)，在由溫度升高到溫度降低的過程中，其內(nèi)部孔隙水的運(yùn)動(dòng)方向發(fā)生了變化。之后，隨外界溫度又回升，內(nèi)外溫差較大，孔隙水運(yùn)動(dòng)更加劇烈，到第18個(gè)小時(shí)(64800s)時(shí)，整個(gè)化石內(nèi)部空隙水流動(dòng)劇烈，圍巖內(nèi)孔隙水運(yùn)動(dòng)也開始深入內(nèi)部。

3.3 隨外界溫度變化化石整體內(nèi)部滲流場(chǎng)變化

隨外界溫度變化化石內(nèi)部滲流場(chǎng)變化，分別模擬了9h(32400s)，24h(86400s)的變化情況(圖7)?？梢钥闯?，隨著外界溫度變化，化石內(nèi)部流速梯度較大，在化石裸露區(qū)兩側(cè)，空隙水向外界流出；而在化石與圍巖接觸區(qū)，孔隙水有穿過圍巖區(qū)向上表面運(yùn)動(dòng)的趨勢(shì)，而這種運(yùn)動(dòng)會(huì)加速化石與圍巖接觸區(qū)化石和圍巖的風(fēng)化；在化石中間部分也出現(xiàn)了較大的流速區(qū)，表明化石內(nèi)部隨著外界溫度的變化，內(nèi)部空隙水流動(dòng)會(huì)出現(xiàn)集中區(qū)，循環(huán)往復(fù)，這部分化石的強(qiáng)度會(huì)降低，出現(xiàn)微小的裂縫，隨后裂縫部分水的運(yùn)動(dòng)會(huì)更加集中，風(fēng)化速度加快。

圖5 溫度升高時(shí)恐龍化石及圍巖內(nèi)部滲流場(chǎng)變化圖

圖6 溫度降低至升高時(shí)恐龍化石及圍巖內(nèi)部滲流場(chǎng)變化圖

圖7 化石內(nèi)部滲流場(chǎng)分布圖

圖8 化石及圍巖熱流耦合內(nèi)部溫度及流速變化圖

3.4 化石及圍巖溫度場(chǎng)和滲流場(chǎng)耦合內(nèi)部溫度及流速變化

圖8為恐龍化石及圍巖熱流耦合的二維切片，分別模擬了1h(3600s)，20h(72000s)的變化情況?？梢钥闯?，初始階段溫度場(chǎng)和滲流場(chǎng)相互作用，使得化石和圍巖表面產(chǎn)生孔隙水運(yùn)動(dòng)；隨著外界溫度的持續(xù)變化，化石和圍巖內(nèi)部孔隙水運(yùn)動(dòng)范圍加大，但可以看出，在化石和圍巖接觸處孔隙水運(yùn)動(dòng)速度遠(yuǎn)大于其他部分；隨著外界溫度的升降變化，化石和圍巖內(nèi)部空隙水運(yùn)動(dòng)的方向會(huì)發(fā)生變化；由于孔隙水的運(yùn)動(dòng)，影響了溫度分布，20h(72000s)時(shí)的溫度云圖非常明顯，這會(huì)更加加劇化石和圍巖內(nèi)部的溫度差，使得內(nèi)部產(chǎn)生溫度應(yīng)力，加速風(fēng)化進(jìn)程；同時(shí)溫度的分布也影響了內(nèi)部空隙水的運(yùn)動(dòng)速度。

4 結(jié)論

(1)從恐龍化石及圍巖的內(nèi)部結(jié)構(gòu)出發(fā)，以恐龍化石及圍巖作為巖體，介紹了恐龍化石風(fēng)化過程中溫度場(chǎng)與滲流場(chǎng)耦合作用的基本理論，并且建立了化石及圍巖熱流耦合作用的多孔介質(zhì)模型，模擬了溫度場(chǎng)與滲流場(chǎng)在恐龍化石及圍巖內(nèi)部及表面的相互作用，描述了恐龍化石在溫度和水的相互作用下的破壞機(jī)制。

(2)實(shí)驗(yàn)結(jié)果顯示，化石及圍巖受日照影響產(chǎn)生的溫度會(huì)影響其內(nèi)部孔隙水的變化，影響滲流場(chǎng)的變化；同時(shí)，化石及圍巖的孔隙水也會(huì)參與到化石及圍巖內(nèi)部的熱量交換中，影響化石及圍巖的溫度場(chǎng)的分布。兩者相互作用，共同對(duì)恐龍化石的破壞產(chǎn)生影響。

(3)通過對(duì)恐龍化石及圍巖進(jìn)行熱流耦合作用的模擬研究，對(duì)恐龍化石的保護(hù)提供了一定的科學(xué)依據(jù)。今后可借助巖石力學(xué)、流體力學(xué)、傳熱學(xué)等相關(guān)學(xué)科，開展化石及圍巖溫度場(chǎng)—滲流場(chǎng)—應(yīng)力場(chǎng)三場(chǎng)之間耦合作用研究，建立三場(chǎng)耦合作用的數(shù)學(xué)模型，為化石保護(hù)建立完整的物理場(chǎng)防護(hù)體系。