石窟寺表面泛鹽病害熱學響應的離散元模擬

關鍵詞：石窟寺；離散元；紅外成像；表面泛鹽；溫度分布；升溫指數

石窟寺是中華優秀傳統文化的重要載體，以其分布廣泛、規模宏大、體系完整的特點，集中體現了建筑、雕塑、壁畫、書法等多種藝術的完美融合，深刻見證了中華民族的歷史背景與文化內涵[1?2]。石窟寺大多依山而建，在使用和保存過程中長期遭受自然與人為因素的破壞，出現了各種不同類型的病害[3]。其中，鹽害是一種破壞性較大的病害，不僅可對石質文物造成直接損害，還可能進一步誘發其他病害[4]。鹽害主要是由于鹽溶液在巖石孔隙中的遷移與沉淀，隨著水分蒸發而析出晶體，從而對石窟寺巖體造成機械及化學性損傷[5]。鹽害的演化一般分為4個階段，表面泛鹽[6]是鹽害發育的第一階段（見圖1），若不及時治理，當其發育到一定程度后，可使巖體進一步誘發粉化[7]、空鼓[8]、剝蝕[9]等病害，直至最終破壞[10]。

目前，無損檢測技術（如紅外熱成像檢測技術[11]、近景攝影[12?13]、探地雷達[14]、三維激光掃描[15?16]等）已被廣泛應用于文物保護工作中[17]。其中，紅外熱成像檢測技術獲取被測物體內、外部溫度信息，具有檢測速度快、范圍廣、安全性高等特點[18]，可較好地適應石窟寺特殊的內部環境。RuizValero等[19]利用紅外熱成像檢測技術檢測了教堂壁畫的潮濕現象，并根據所獲熱成像結果確定了壁畫表面的潮濕程度。Bisegna等[20]使用被動紅外熱成像檢測技術檢測了教堂內壁，論證了紅外檢測技術在文物保護工作中的有效性。楊文峰[21]使用主動式紅外熱成像檢測技術對慶陽某石窟寺開展了全面檢測，針對自然冷卻過程提出了降溫指數，并成功將該指數應用于石窟寺剝落病害的定量評估中。上述成功案例已充分證明了紅外熱成像檢測技術在石窟寺保護工作中的可行性。針對鹽害，國內外學者已開展了較為豐富的試驗研究與理論分析。Martínez-Martínez等[22]發現鹽類結晶是圣母瑪利亞大教堂中石質文物發生風化的重要誘因，從而提出了以預防鹽類結晶為首要任務的保護策略。Hu等[23]通過對莫高窟中收集的顆粒樣品進行微觀分析，揭示了鹽分隨毛細現象而向上遷移的機理。Zhao等[24]使用低場核磁共振、X射線斷層掃描等技術明確了氯化鈉溶液在多孔二氧化硅和多孔硅酸鹽材料中的時空演化規律，進而提出鹽溶液對硅酸鹽材料的破壞機理。Michette等[25]對已受鹽害影響的萊蓋特石樣品進行了不同環境條件下的加速老化試驗，并指出不同鹽結晶可能會導致不同的鹽害特征。張虎元等[26]進行了不同濕度條件下的硫酸鈉溶液毛細上升試驗，建立了巖石表面蒸發速率與鹽害特征的內在聯系。在石窟寺保護工作中，表面泛鹽階段的量化評估是至關重要的，不僅是制定修復方案的重要基礎，更是控制病害進一步發育的重要節點。上述研究成果盡管厘清了鹽害的發生機理與破壞模式，卻未能充分結合現有檢測技術提出指導工程實踐的定量化評估指標。

針對上述瓶頸，以慶陽北石窟寺為例開展離散元模擬研究，在充分考慮巖體與析出晶體顆粒特性的前提下構建表面泛鹽的熱學模型，并進一步結合主動式紅外熱成像檢測技術獲取析出晶體與巖體在不同病害程度下的熱學響應，以此提出表面泛鹽的量化評估指標。

1 主動紅外熱成像檢測技術

紅外熱成像檢測技術主要包括主動式和被動式紅外熱成像檢測方法。被動式方法依賴接收被測物體輻射出的紅外線進行檢測；而主動式方法則利用熱激勵裝置施加外部熱源（見圖2），強迫被測物體發生溫度變化從而完成檢測。對于內部環境較為封閉的石窟寺而言，其窟內自然環境的溫度變化較小，被動式方法往往難以使鹽害區域與周圍巖體出現差異明顯的熱學響應特征。因此，采用主動式方法，通過施加外部熱源打破析出晶體與巖體的熱平衡狀態，從而獲取兩者的升溫曲線，完成表面泛鹽厚度的定量化評估。

2 離散元熱學模型

2.1 熱學計算原理

研究對象主要由粗、細顆粒組成，故采用離散元方法獲取顆粒尺度上泛鹽病害與周圍巖體的不同熱學響應特征。顆粒流程序（particleflowcode，PFC）是一款支持并行運算的離散元模擬軟件[27]，可準確模擬顆粒尺度上巖土體的非線性力學行為，已在巖土工程領域得到了廣泛應用[28?29]，其熱學計算原理如下。

對于連續體而言，假設應變變化對溫度的影響忽略不計，其熱傳導方程為

PFC進行熱學計算時，需要設置4種參數，即密度ρ、熱阻η、比熱C、線膨脹系數α。實際上，材料的熱阻與其形狀、尺寸等因素有關，在PFC2D中，采用如下的熱阻計算公式。

式中：k為材料的導熱系數；Vb為顆粒b的體積；lp為接觸p的長度；n為孔隙率。

2.2 研究方案

采用PFC2D軟件建立簡化的二維表面泛鹽熱學模型。在使用主動熱激勵裝置以前，石窟寺內巖壁及表面泛鹽病害處于熱平衡狀態，因此，模型整體的初始溫度按現場實測溫度設置（26.7℃）。考慮到主動激勵輻照面積較小，故模型邊界溫度均設置為26.7℃，用于模擬巖體無窮遠處的溫度。根據泛鹽程度選定5個特定區域記錄其溫度變化曲線，具體位置及泛鹽厚度如圖3所示。為實現主動熱激勵功能，在監測區域左側均勻布置具有一定熱功率的顆粒。在實際檢測中，熱激勵裝置與巖壁間有一定間隔，實際傳輸到巖壁的熱流有一定損失，因此，需根據實測的砂巖升溫曲線反演熱源總功率。

2.3 模型參數的選取

對于具有一定連續性的砂巖巖體，采用平行膠結模型[30]，模型細觀力學參數參考了葉永芃[31]的研究成果，其余參數則通過現場試驗直接測定[21]。根據曹張喆等[4]的研究，北石窟寺中鹽害的主要成為硫酸鹽，硫酸鹽晶體的模型參數參考羅崇亮等[32]的研究成果。在PFC2D中建立表面泛鹽病害模型主要包括2個步驟：

1）砂巖模型。試樣寬為320mm，試樣高為640mm，砂巖顆粒最大半徑和最小半徑分別為1.2mm和0.8mm，滿足均勻分布，最大最小粒徑比為3：2。試樣按照孔隙率0.01%生成密實樣，包含約6.39×104個顆粒。

2）泛鹽病害。使用硫酸鹽顆粒在砂巖表面形成1個弓形區域，模擬實際的表面泛鹽病害（見圖4）。為減少邊界效應影響，弓形區域高度為砂巖試樣高度的1/3，其最大厚度設置為7mm。硫酸鹽顆粒半徑按照砂巖顆粒半徑的1/6設置，孔隙率按照0.001%設置，包含約1.26×104個顆粒。

因主要探究表面泛鹽厚度與溫度響應間的關系，只需進行熱學計算，不考慮由熱應變導致的力學響應，2種顆粒的線膨脹系數均設為0K-1。具體模型參數如表1～2所示。

圖5展示了實測[21]與模擬的升溫曲線，可以看出兩者的變化規律基本吻合且數值差距較小。當升溫結束時，實測砂巖表面的溫度為31.3℃，數值模擬的溫度為31.6℃，兩者的溫度僅相差0.3℃，同時R2為0.996。以上說明采用的計算原理與模型參數能基本再現實測結果，此時熱源的總功率為115W。

模擬結果和實測結果出現誤差的原因可能有2點：1）巖體內部的溫度可能低于26.7℃，而數值模擬中所有顆粒的初始溫度均統一設為26.7℃；2）北石窟寺的巖體是非均質的，主要由石英、長石等礦物構成[33]，而數值模擬中的簡化模型將顆粒處理為均質砂巖。

3 結果與分析

3.1 數值模擬結果分析

圖6是不同加熱時間下的溫度云圖，在初始狀態時，所有顆粒的溫度是一致的，均為石窟寺內部實測環境溫度。從圖6（a）～（g）中的云圖變化可以看出，在升溫900s的過程中，熱量逐漸向沒有直接接觸熱源的上下兩側及內部傳導。在同一水平位置，泛鹽病害后側的砂巖溫度普遍較低，說明在泛鹽病害影響下的砂巖升溫速度較慢。圖6（h）為升溫900s后表面泛鹽位置的局部溫度云圖，不難發現表面泛鹽病害的表面溫度遠高于砂巖表面，這說明泛鹽使熱量在表面形成堆積并阻礙熱量進一步向內傳播。

圖7為砂巖表面（測點a）及不同泛鹽厚度（測點b、c、d、e）處的溫度變化曲線。升溫900s后，具有泛鹽病害的4個位置的最終溫度依次為31.8、32.3、33.0、33.7℃，比砂巖表面溫度分別高了0.2、0.7、1.4、2.1℃。這表明在相同時間內，泛鹽厚度與最終溫度呈現明顯正相關性，同時泛鹽越厚溫度上升越快。

當砂巖表面的鹽害物質與其熱學特性不一致時，熱量從鹽害處進入砂巖內部，與直接從砂巖表面進入其內部的過程相比，會表現出不同的熱學響應特征。在本文中，鹽顆粒的導熱系數為0.140W/（m·℃），砂巖的導熱系數為0.995W/（m·℃），鹽顆粒較小的導熱系數使其產生比砂巖更好的隔熱性能，在相同熱激勵條件下，熱量難以通過表面泛鹽向內傳遞，從而體現出較高的熱學響應。對于泛鹽中心位置和邊緣位置，中心位置的泛鹽厚度比兩側邊緣更大；同時，邊緣位置的局部熱阻受到泛鹽和砂巖2種材料的影響會比泛鹽中間位置的局部熱阻低。因此，在相同時間內，泛鹽病害邊緣傳遞的熱量比泛鹽病害中間傳遞的熱量更多，導致病害邊緣表面的溫度會比病害中心處的溫度低。由此可以推斷，泛鹽越厚的位置，表面溫度也越高。

3.2 升溫指數

為消除時間和溫度的不同量綱影響，使時間和溫度處于同一數量級，將時間和溫度按照式（8）進行歸一化處理，即進行最大最小標準化處理。

式中：tn與Tn分別是時間和溫度采用最大最小標準化進行歸一化處理后的數據。ti為初始時刻，通常情況下為0；Ti為初始時刻目標表面的溫度；tf為結束時刻；Tf為結束時刻目標表面的溫度；t為時間；T為t時刻下目標表面的溫度。

以時間歸一化后的值tn為橫坐標，溫度歸一化后的值Tn為縱坐標，建立Tn與tn的關系，如圖8所示。

從圖8可以看出，在恒定熱功率的情況下，Tn和tn呈現出明顯的非線性關系，其形狀與冪函數相似。故采用Allometricl函數（y=AxB）進行擬合，使用Levenberg-Marquardt優化算法進行迭代。由于對數據進行了歸一化處理，則函數必過（0，0）和（1，1）這2點，則控制擬合函數的系數A值為1.0，擬合結果如表3所示。

每條曲線的擬合的R2均為0.989以上，說明擬合效果良好。從擬合結果可以看出，泛鹽厚度與指數B存在著明顯的負相關關系。以B為橫坐標，泛鹽厚度H為縱坐標進行擬合，得到圖9的結果。可見，指數B與泛鹽厚度H具有良好的線性關系。

根據升溫指數的定義可知，升溫指數Bw即為表3中擬合曲線的指數B。通過圖9可以發現在數值模擬中，隨著泛鹽厚度的增加，升溫指數逐漸變小，兩者呈明顯的線性負相關關系。由此可見，通過記錄主動式紅外檢測的升溫過程，使用升溫指數即可實現表面泛鹽病害厚度的定量化評估。

4 結論

基于離散元數值模擬方法，構建了表面泛鹽病害模型，并采用主動式紅外熱成像檢測技術原理，獲取了表面泛鹽不同病害位置的熱學響應。研究結論如下：

1）因鹽顆粒的低熱導系數使其成為天然隔熱層，故在升溫過程中泛鹽病害處表現出更高的熱學響應。

2）對時間和溫度進行歸一化處理后，兩者呈現冪函數關系。擬合結果顯示，泛鹽厚度與冪函數指數呈負相關線性關系，因此將冪函數指數作為升溫指數，用以定量化評估泛鹽病害發育程度。