謝一凡 曾祎芃 楊杰 葉逾 吳吉春 魯春輝

摘要：探究水平雙裂隙含水層中的海水入侵機制，對理解裂隙含水層中的海水入侵機制具有重要意義。本文構建了大型室內砂槽模型，并通過離散裂隙基質數(shù)值模型分析了不同分布形式的水平雙裂隙對鹽水楔、混合帶等海水入侵特征的影響。結果表明：① 水平雙裂隙能顯著增強鹽水楔及混合帶的發(fā)展，裂隙間的咸淡水混合尤為顯著;2條主干裂隙會限制彼此的作用，位于近海側、含水層中部的裂隙作用更加顯著。② 位于近海側的裂隙高度越高，則鹽水楔發(fā)展程度和裂隙間混合作用越顯著，而近陸側裂隙高度則僅對鹽水楔尾部有明顯影響。③ 位于含水層中部、底部的裂隙，越靠近內陸則鹽趾長度更大、鹽水體積更小，而遠離鹽水楔時則會使混合帶和鹽趾長度急劇變小。

關鍵詞：海水入侵;地下水動力;水平雙裂隙;室內試驗;數(shù)值模擬

中圖分類號：P641

文獻標志碼：A

文章編號：1001-6791（2023）03-0454-11

收稿日期：2022-11-24;

網絡出版日期：2023-04-14

網絡出版地址：https：∥kns.cnki.net/kcms/detail/32.1309.P.20230413.1832.003.html

基金項目：國家重點研發(fā)計劃資助項目（2021YFC3200500）;國家自然科學基金資助項目（42277190）

作者簡介：謝一凡（1987—），男，江蘇南京人，副教授，博士，主要從事海岸帶水動力方面研究。

E-mail：yfxie@hhu.edu.cn

通信作者：魯春輝，E-mail：clu@hhu.edu.cn

海水入侵是海岸帶一種常見的環(huán)境地質災害，嚴重制約著沿海經濟和社會的可持續(xù)發(fā)展［1-5］。裂隙含水層是全球分布最廣泛的含水層類型之一，也是中國較常見的濱海含水層類型［6］。相比其他類型的含水層，裂隙含水層更易受海水入侵的影響［7-8］，多個國家的原位觀測數(shù)據(jù)和場地案例研究均證實了這一點［5，7，9-10］。目前，裂隙含水層是海水入侵研究領域的前沿方向與研究熱點［8］，但尚在起步階段，主要集中在靜態(tài)邊界下裂隙含水層海水入侵的數(shù)值模擬研究，缺乏基于室內試驗的海水入侵物理機制研究。因此，通過構建室內物理模型深入探索裂隙含水層中的海水入侵過程與規(guī)律，不僅能夠豐富、完善海水入侵機制，也對環(huán)境地質災害防治等實際需求具有重要實踐意義。

雖然裂隙含水層中的海水入侵機制尚不清晰，但其他類型非均質含水層的研究已表明非均質性可影響地下水流動路徑，進而影響海水入侵［1，11］。Lu等［12］研究表明，分層非均質含水層對咸淡水混合帶位置有決定性作用，且流線的折射和分散作用會令高滲透層上方的低滲透層中的混合帶顯著變寬。隨機非均質含水層則會導致混合體積的增加［13］。同時，國內外關于單裂隙、裂隙網絡含水層中水流和溶質運移的研究，揭示了裂隙產生的優(yōu)先流動路徑對于溶質運移的作用，以及裂隙的物理與幾何特征（尺寸、分布、數(shù)量等）對于動力學系統(tǒng)的重要性［14-18］，但大部分研究尚未涉及海岸帶復雜水鹽運移過程。

現(xiàn)有裂隙含水層中海水入侵研究已初步給出研究所用的基本假設與關鍵因子?，F(xiàn)場勘測研究發(fā)現(xiàn)裂隙的位置及其與抽水井、海洋的水力連接決定了海水入侵程度［5，7］。裂隙含水層的現(xiàn)有數(shù)值模型主要有隱式模型和顯式模型2類［19］，分別具有高效率和高精度的特點，2類模型的典型方法分別為等效多孔介質模型（EPM）、雙重孔隙度模型（DPM）和離散裂隙基質模型（DFM）、離散裂隙網絡模型（DFN）［19］。其中，DFM僅保留了部分裂隙，其余裂隙則通過提升尺度技術與基質一同被等效為多孔介質，并同時考慮裂隙與等效多孔介質中的流動［19］。同時，本領域的數(shù)值模擬研究已發(fā)現(xiàn)裂隙的密度、方向等特征對海水入侵有宏觀影響，且顯式模型比隱式模型更適合研究僅含少量主干裂隙的含水層中的海水入侵機制［20-21］?；谶_西定律和DFM，Etsias等［19］于2021年建立了裂隙含水層海水入侵物理模型，并對裂隙特征進行了數(shù)值敏感性分析，發(fā)現(xiàn)裂隙的位置和長度會顯著影響鹽水楔的分布，且對抽水引發(fā)的鹽水上升過程也有重要影響［22］。然而，Etsias等［19］僅考慮了含單條裂隙的含水層，而實際濱海裂隙含水層常發(fā)育多條平行的主干裂隙［23］。

本文建立水平雙裂隙非承壓含水層海水入侵的大型砂槽物理模型，應用SUTRA（saturated-unsaturated transport）構建DFM，以厘清不同位置組合的2條水平主干裂隙對海水入侵的宏觀影響規(guī)律和裂隙間的水鹽運移模式。通過構建4種不同類型的裂隙分布組合（含5組試驗工況，21組數(shù)值模擬工況），詳細研究不同水平雙裂隙分布對11個典型海水入侵特征值的影響。本研究有助于理解裂隙含水層中海水入侵的基本機制，為實際濱海裂隙含水層的管理、利用和保護提供重要科學依據(jù)。

1?模型設置

1.1?室內物理試驗設置

試驗在長4 m、高0.8 m、寬0.02 m的大型砂槽中進行，砂槽系統(tǒng)如圖1（a）所示，砂槽左右兩側設置了鹽水和淡水水槽，鹽水水槽內設的柵欄式溢流槽以及與淡水水槽相接的蠕動泵可分別精確控制海、陸邊界所需的水位（0.54 m）與流量（0.62 mL/s）。參考Etsias等［19］的試驗設置，主干裂隙由長1 m、高0.01 m、寬0.02 m且表面分布0.8 mm孔徑細孔的中空不銹鋼沖孔管模擬（圖1（b））;等效多孔介質由粒徑為1 mm的石英砂模擬。經過測定，試驗砂的飽和滲透率為1.33×10-9 m2，孔隙度為0.45，與相關研究的參數(shù)值相近［19，22］。試驗砂的縱向彌散度（αL）為0.005 m，橫向彌散度（αT）設為αL的1/10［12］，即0.000 5 m。裂隙的孔隙度設為1，其滲透率則通過數(shù)值敏感性分析確定［19］，為多孔介質滲透率的300倍，滿足Matthi等［24］工作中的裂隙-基質滲透率比范圍（102～104）。

配置密度為1 025 kg/m3的鹽水模擬海水并加入誘惑紅作為染色劑［19］。試驗期間每30 min使用Atago PAL-SALT鹽度計監(jiān)測鹽度，并緩慢加入濃縮鹽水以確保試驗期間海水邊界的鹽度和水位恒定。當30 min內的鹽水鹽度和鹽趾位置變化不超過1%時，認為試驗已經達到穩(wěn)態(tài)。砂槽裝填高度為0.7 m，并在海水邊界設置長1.2 m，坡度為1∶4的海灘斜坡。本文共進行5組試驗工況，工況1為均質情況;以砂槽左下角為原點，工況2—5中2條裂隙的底部中點坐標分別為（1.0 m，0.3 m）和（2.0 m，0.1 m）、（1.0 m，0.2 m）和（2.0 m，0.2 m）、（1.7 m，0.3 m）和（2.0 m，0.1 m）、（1.0 m，0.3 m）和（1.5 m，0.1 m）。

1.2?數(shù)值模型設置

采用由美國地質調查局開發(fā)的SUTRA進行二維數(shù)值模擬。SUTRA是基于有限差分法和有限元法模擬變密度、變飽和度水流運動及溶質運移的程序，已廣泛應用于海水入侵模擬研究［19，22，25］，所用控制方程和van Genuchten經驗公式參見文獻［25］。本文假設流體密度（ρ）隨質量濃度（c）線性變化，變化率ρ/c=714.3 kg/m3［25］，并控制海、陸邊界的水位和流量以保證水流為層流［19］。模擬的21個工況中，最大流速為0.005 7 m/s，雷諾數(shù)為75.49，遠小于Etsias等［19］工作中層流的雷諾數(shù)值（256.3），驗證了達西定律在本文工作中的適用性。

數(shù)值模型示意如圖1（c），模型不考慮降水和蒸散發(fā)，所有設置均與室內試驗一致（表1）。研究區(qū)被剖分為500×120個網格，共60 621個節(jié)點，其中裂隙被顯式剖分為二維網格。離散網格的最大Péclet數(shù)Pe=1.2，滿足Pe≈ΔL/αL≤4的數(shù)值穩(wěn)定性判據(jù)，其中ΔL為沿地下水流動方向單元兩側的傳輸距離［25］。邊界ABC上，對流入含水層的節(jié)點設為海水質量濃度（cs=35 kg/m3），而對流出含水層的節(jié)點設為零濃度梯度［25］。模擬時間步長為1 s、總時間為1 d，所有模擬工況均可達到穩(wěn)態(tài)。

2?結果與分析

2.1?室內試驗與數(shù)值模擬結果

圖2展示了5組試驗工況達到穩(wěn)態(tài)的結果和對應的數(shù)值模擬結果，顯示裂隙改變了其附近流體流動方向、加速了其內的流體流動，并在其上、下方形成低流速區(qū)。工況1為均質情況。工況2中鹽水楔占領了裂隙下方的低流速區(qū)，且裂隙中的高速水流限制了鹽水楔的高度，故鹽水楔呈現(xiàn)階梯狀；同時，咸淡水在低速區(qū)流速接近，近海側裂隙的上、下方均形成較大的混合區(qū)。工況3中2條裂隙形成1條長裂隙，促進了鹽趾長度。對比工況2、工況3的結果，可以發(fā)現(xiàn)裂隙高度是鹽水楔高度的重要影響因素。工況4中含水層中部裂隙對鹽水楔形態(tài)的影響較小，但其下方低流速區(qū)使鹽水突破了含水層底部裂隙的限制，在裂隙間形成溶質濃度較小的混合區(qū)。工況5中含水層底部裂隙完全處于鹽水楔內部，故鹽水楔特征與Etsias等［19］工作中的單裂隙情況類似。同時，圖2也顯示數(shù)值結果與試驗結果吻合較好，僅有的細微差別可能來源于填砂均勻程度、裂隙位置等試驗中難以準確控制的部分。

2.2?敏感性分析

根據(jù)均質情況的海水入侵水平和垂直范圍，構建了4組數(shù)值模擬敏感性分析，共21個工況（表2），以研究水平雙裂隙含水層中裂隙的高度、水平位置對海水入侵的影響。所分析的海水入侵特征包括鹽趾趾長（海水50%濃度等值線與含水層底部的交點位置［19］）、鹽水楔質心位置、鹽水體積、混合體積（以海水10%和90%濃度等值線為界限［20］）、混合帶質心位置、混合帶垂直厚度和水平寬度。

2.2.1?近海側裂隙的高度變化對海水入侵特征的影響

圖3展示了水平雙裂隙含水層中近海側裂隙位于不同高度時的海水入侵，相關數(shù)值總結于表3。

如圖3（a）所示，近海側裂隙對淡水的匯聚作用隨其高度上升變強，裂隙匯聚了大量淡水，減少了淡水對鹽水楔的擠壓，促進了鹽水楔的擴張。同時，當近海側裂隙升高時，近陸側裂隙對淡水的匯聚隨著其與淡水接觸面積的減小而遞減，對底層鹽水的抑制作用逐步減弱。相較均質情況，A1—A5的鹽趾長度分別增加了29.11%、30.61%、32.57%、36.01%和40.52%，且增量逐漸遞增（圖3（b）、表3）;鹽水楔質心位置逐漸向右上方遷移。其中，A1—A3中質心垂直位置略低于均質情況，展現(xiàn)了裂隙對鹽水楔高度的限制作用。此外，5組工況的鹽水體積均大于均質情況。其中，A2、A3中鹽水體積較小，此時因近海側裂隙匯聚淡水而減少的鹽水體積（裂隙上方）大于因其高度上升而增加的鹽水體積（裂隙下方）。隨后，鹽水體積隨近海側裂隙高度上升遞增。

結合圖3和表3可以看到，A1—A5中混合體積以及混合帶的厚度、寬度均大于均質情況，裂隙間產生了大面積的混合帶，其混合程度隨近海側裂隙高度上升而增強。A1—A5中混合體積的變化趨勢為先增大后減小;峰值于A3發(fā)生（約為均質情況的5.14倍），此時近海側裂隙位于淡水排泄出口和近陸側裂隙中間，在其上、下方均有較大面積的低流速區(qū)，為混合帶的產生和發(fā)展提供了條件。隨著近海側裂隙的位置更高，主要混合帶逐漸從其上部轉移到下部，表現(xiàn)為混合帶質心水平位置的右移。A1—A3中最大混合帶垂直厚度出現(xiàn)在近海側裂隙上方，隨著上方混合帶的變小而遞減（A2—A3）;A4、A5中則轉移至裂隙之間，并因高濃度鹽水體積擴張而遞減。A1、A2中最大混合帶水平寬度出現(xiàn)在近陸側裂隙下方，A3—A5中則出現(xiàn)在近海側裂隙下方?；旌蠋骄穸茸兎淮?，但其平均寬度隨近海側裂隙高度上升遞增。

2.2.2?近陸側裂隙的高度變化對海水入侵特征的影響

圖4展示了水平雙裂隙含水層中近陸側裂隙位于不同高度時的海水入侵，相關數(shù)值總結于表4。

近陸側裂隙的主要影響范圍是鹽水楔尾部（圖4（a）），故其對鹽水楔形態(tài)及溶質混合程度的影響小于近海側裂隙。海水入侵程度隨著近陸側裂隙高度上升呈先增大后減小的趨勢，具體體現(xiàn)在鹽水楔尾部形狀、鹽趾長度和鹽水體積的變化上。B2中2條裂隙首尾相連，此時鹽趾長度、鹽水體積達到峰值，分別為均質情況的1.48倍和1.57倍（表4）。B1、B2中近陸側裂隙的高度分別小于、等于近海側裂隙，故鹽水楔尾部的擴張程度隨其位置更高而加劇。B3—B5中近陸側裂隙不再與鹽水楔接觸，此時近海側裂隙限制了鹽水楔尾部的擴張高度。隨著近陸側裂隙的位置升高，含水層下部淡水逐漸從被近陸側裂隙匯聚轉變?yōu)楸桓叨容^低的近海側裂隙匯聚，淡水對鹽趾長度的抑制作用增加（圖4（b））。此外，鹽水楔質心高度對近陸側裂隙高度變化的敏感性較小，這一方面是由于近海側裂隙高度對鹽水楔尾部垂直擴張程度的限制，另一方面是由于近海側裂隙的位置決定了鹽水楔中大部分鹽水的分布。

受近海側裂隙的限制，近陸側裂隙對混合帶質心的影響主要針對其水平位置。B1中最大混合帶垂直厚度（26.08 cm）遠高于其他工況（表4）。這是由于此時近陸側裂隙靠下，其匯聚的水大部分又被近海側裂隙匯聚，為近海側裂隙上部咸淡水混合帶的發(fā)展提供了條件。B2中2條裂隙連接，出現(xiàn)最大混合帶水平寬度的峰值為1.631 5 m（表4），大于A1中的1.293 4 m（表3），這是由于B2中裂隙與淡水的接觸面積更大;B3中近陸側裂隙仍能吸引下方鹽水，形成較寬的混合帶;但隨著近陸側裂隙與鹽水楔的距離更遠（B4—B5），其對鹽水楔的影響陡然減小，表現(xiàn)為最大混合帶水平寬度的突降。5組工況中的混合體積、混合帶平均厚度和平均寬度均高于均質情況。

2.2.3?含水層中部裂隙的水平位置變化對海水入侵特征的影響

圖5展示了水平雙裂隙含水層中含水層中部裂隙位于不同水平位置時的海水入侵，相關數(shù)值總結于表5。

圖5和表5顯示，隨著含水層中部裂隙更靠近內陸，鹽趾長度總體呈現(xiàn)增加的趨勢、而鹽水體積總體呈現(xiàn)減小的趨勢。這是由于隨著含水層中部裂隙更靠近內陸，其對淡水的匯聚逐步增強，排泄區(qū)淡水流速增大而對鹽水楔上部的鹽水形成更強的擠壓和稀釋作用，鹽水體積從而減小。同時，含水層底部裂隙對淡水的匯聚因含水層中部裂隙匯聚更多淡水而減弱，其對鹽趾長度的抑制作用減小。圖5（b）中鹽水楔質心的水平位置先增大后減小，對應了含水層中部裂隙對淡水的匯聚增強使鹽水楔更伸向內陸（C1—C3）、又在遠離鹽水楔時排放淡水而對鹽水楔上部形成擠壓的過程（C4—C6）。鹽水楔質心高度與均質情況接近，變化不大。

如圖5和表5所示，混合體積、混合帶厚度和寬度總體均呈現(xiàn)先增大后減小的變化趨勢。這是由于隨著含水層中部裂隙逐漸靠近鹽水楔中部（C1—C4），吸引了底部裂隙中的混溶水，主混合區(qū)逐漸從含水層中部裂隙的上方和含水層底部裂隙附近轉移至裂隙間，并隨著含水層中部裂隙右移而擴張;隨著含水層中部裂隙更靠近內陸（C5—C6），主混合區(qū)轉移至含水層底部裂隙附近（圖5（a））。因此，混合帶質心位置先向海洋和含水層頂部變化，再向內陸和含水層底部變化（圖5（b））。

2.2.4?含水層底部裂隙的水平位置變化對海水入侵特征的影響

圖6展示了水平雙裂隙含水層中含水層底部裂隙位于不同水平位置時的海水入侵，相關數(shù)值總結于表6。

圖6（b）顯示鹽趾長度的變化趨勢為先增大后減小，且總大于均質情況。這是由于含水層底部裂隙在逐步脫離鹽水楔的過程中（D1—D4），其對鹽水的匯聚和吸引的增強令鹽趾長度變長;隨著遠離鹽水楔（D5），其左端排放的淡水使鹽水楔尾部被推向海洋，令鹽趾長度陡降0.706 2 m（表6），顯示兩者的相對位置決定了裂隙對鹽水楔的作用（吸引或抑制）。鹽水楔質心高度變化不大，水平位置則先逐步減小然后突變增大再逐步減小，對應了底部裂隙逐漸脫離鹽水楔再逐漸遠離鹽水楔的過程。同時，鹽水體積也因含水層底部裂隙匯聚更多淡水而擠壓變?。▓D6（b）、表6）。

結合圖6和表6來看，D2中混合體積達到峰值同時鹽水體積陡降，此時含水層底部裂隙右側開始匯聚內陸淡水，又被中部裂隙下方的低流速區(qū)所影響，2條裂隙之間產生了較大的混合帶;隨著含水層底部裂隙遠離鹽水楔，其左側開始排放大量淡水擠壓鹽水楔尾部，也令混合帶減小。圖6（b）顯示混合體積、混合帶厚度和寬度的變化趨勢均為先增大后減小，對應于含水層底部裂隙水平位置向內陸變化時，2條裂隙對咸淡水混合的促進作用先增大后減小?；旌蠋恢玫淖兓厔菀脖砻?條裂隙之間和含水層底部裂隙附近的溶質混合先加強后減弱。

3?結?論

本文通過室內試驗及數(shù)值模擬研究了不同水平雙裂隙分布含水層中的海水入侵特征與機制，主要結論如下：

（1） 位于海水入侵范圍的水平雙裂隙能夠增強鹽水楔以及混合帶的發(fā)展，但水平雙裂隙對海水入侵的影響會彼此限制，位于近海側和含水層中部的裂隙作用更明顯。裂隙的匯聚淡水作用會在其上、下方產生低流速區(qū)，故裂隙在鹽水楔附近時會增強咸淡水混合，且裂隙間的混合帶尤為顯著;裂隙遠離鹽水楔時則會促進內陸淡水向海洋排放，抑制海水入侵。

（2） 水平雙裂隙含水層中近海側裂隙高度越高，鹽水楔更向內陸和含水層頂部發(fā)展，混合帶逐漸從其上部轉移到裂隙間且混合作用逐漸加劇，并始終高于均質情況。近陸側裂隙的主要影響范圍是鹽水楔尾部，雖受到近海側裂隙高度限制而對海水入侵影響相對較弱，但其高度變化對鹽趾長度、鹽水體積、混合帶范圍等仍有重要作用，且可決定鹽水楔質心和混合帶質心的水平位置。

（3） 水平雙裂隙含水層中含水層中部裂隙水平位置越靠近內陸，鹽趾長度更大、鹽水體積減小;當其遠離鹽水楔時，含水層中部的混合帶急劇變小。含水層底部裂隙水平位置越靠近內陸，鹽趾長度更大、鹽水體積減小，而當其遠離鹽水楔時，鹽趾長度陡降;其對混合帶參數(shù)的影響相對較小。

參考文獻：

［1］WERNER A D，BAKKER M，POST V E A，et al.Seawater intrusion processes，investigation and management：recent advances and future challenges［J］.Advances in Water Resources，2013，51：3-26.

［2］楊蘊，吳劍鋒，林錦，等.控制海水入侵的地下水多目標模擬優(yōu)化管理模型［J］.水科學進展，2015，26（4）：579-588.（YANG Y，WU J F，LIN J，et al.A multi-objective simulation-optimization model constrained by the potential seawater intrusion［J］.Advances in Water Science，2015，26（4）：579-588.（in Chinese））

［3］陳飛，徐翔宇，羊艷，等.中國地下水資源演變趨勢及影響因素分析［J］.水科學進展，2020，31（6）：811-819.（CHEN F，XU X Y，YANG Y，et al.Investigation on the evolution trends and influencing factors of groundwater resources in China［J］.Advances in Water Science，2020，31（6）：811-819.（in Chinese））

［4］王焰新，甘義群，鄧婭敏，等.海岸帶海陸交互作用過程及其生態(tài)環(huán)境效應研究進展［J］.地質科技通報，2020，39（1）：1-10.（WANG Y X，GAN Y Q，DENG Y M，et al.Land-ocean interactions and their eco-environmental effects in the coastal zone：current progress and future perspectives［J］.Bulletin of Geological Science and Technology，2020，39（1）：1-10.（in Chinese））

［5］李闖，倪金，宋慶春，等.大連典型海水入侵區(qū)地下水供水適宜性及合理利用對策［J］.節(jié)水灌溉，2021（7）：65-70.（LI C，NI J，SONG Q C，et al.Suitability and rational utilization strategy of groundwater supply in a typical seawater intrusion area of Dalian［J］.Water Saving Irrigation，2021（7）：65-70.（in Chinese））

［6］董健，曾獻奎，吳吉春.不同類型海岸帶海水入侵數(shù)值模擬研究進展［J］.高校地質學報，2018，24（3）：442-449.（DONG J，ZENG X K，WU J C.Advances in numerical simulation of seawater intrusion in different coastal zones［J］.Geological Journal of China Universities，2018，24（3）：442-449.（in Chinese））

［7］PARK H Y，JANG K，JU J W，et al.Hydrogeological characterization of seawater intrusion in tidally-forced coastal fractured bedrock aquifer［J］.Journal of Hydrology，2012，446/447：77-89.

［8］GIESE M，BARTHEL R.Review：saltwater intrusion in fractured crystalline bedrock［J］.Hydrogeology Journal，2021，29（7）：2313-2328.

［9］MACALLISTER D J，JACKSON M D，BUTLER A P，et al.Remote detection of saline intrusion in a coastal aquifer using borehole measurements of self-potential［J］.Water Resources Research，2018，54（3）：1669-1687.

［10］SUN N B，ZHUANG H B，WANG S Y，et al.Study on the causes of seawater intrusion in coastal tuff hill area［J］.Water Resources and Hydropower Engineering，2021，52（S1）：298-305.

［11］譚博，劉曙光，代朝猛，等.濱海地下水交互帶中的膠體運移行為研究綜述［J］.水科學進展，2017，28（5）：788-800.（TAN B，LIU S G，DAI C M，et al.A review of colloid transport behaviors based on coastal groundwater interaction［J］.Advances in Water Science，2017，28（5）：788-800.（in Chinese））

［12］LU C H，CHEN Y M，ZHANG C，et al.Steady-state freshwater-seawater mixing zone in stratified coastal aquifers［J］.Journal of Hydrology，2013，505：24-34.

［13］MAHMOODZADEH D，KARAMOUZ M.Seawater intrusion in heterogeneous coastal aquifers under flooding events［J］.Journal of Hydrology，2018，568：1118-1130.

［14］MEDICI G，SMERAGLIA L，TORABI A，et al.Review of modeling approaches to groundwater flow in deformed carbonate aquifers［J］.Ground Water，2021，59（3）：334-351.

［15］ZHU W W，KHIREVICH S，PATZEK T W.Impact of fracture geometry and topology on the connectivity and flow properties of stochastic fracture networks［J］.Water Resources Research，2021，57（7）：e2020WR028652.

［16］鄧亞平，張燁，施小清，等.非均質裂隙介質中重非水相流體運移［J］.水科學進展，2015，26（5）：722-730.（DENG Y P，ZHANG Y，SHI X Q，et al.Study on the migration of dense non-aqueous phase liquids in heterogeneous fractured media［J］.Advances in Water Science，2015，26（5）：722-730.（in Chinese））

［17］羅明明，季懷松.巖溶管道與裂隙介質間溶質暫態(tài)存儲機制［J］.水科學進展，2022，33（1）：145-152.（LUO M M，JI H S.Mechanism of solute transient storage between Karst conduit and fissures［J］.Advances in Water Science，2022，33（1）：145-152.（in Chinese））

［18］成建梅，羅一鳴.巖溶多重介質地下水模擬技術及應用進展［J］.地質科技通報，2022，41（5）：220-229.（CHENG J M，LUO Y M.Simulation technology and application progress of Karst multi-media groundwater［J］.Bulletin of Geological Science and Technology，2022，41（5）：220-229.（in Chinese））

［19］ETSIAS G，HAMILL G A，CAMPBELL D，et al.Laboratory and numerical investigation of saline intrusion in fractured coastal aquifers［J］.Advances in Water Resources，2021，149：103866.

［20］SEBBEN M L，WERNER A D，GRAF T.Seawater intrusion in fractured coastal aquifers：a preliminary numerical investigation using a fractured Henry problem［J］.Advances in Water Resources，2015，85：93-108.

［21］KOOHBOR B，F(xiàn)AHS M，ATAIE-ASHTIANI B，et al.Uncertainty analysis for seawater intrusion in fractured coastal aquifers：effects of fracture location，aperture，density and hydrodynamic parameters［J］.Journal of Hydrology，2019，571：159-177.

［22］ETSIAS G，HAMILL G A，THOMSON C，et al.Laboratory and numerical study of saltwater upconing in fractured coastal aquifers［J］.Water，2021，13（23）：3331.

［23］宋明春，林少一，楊立強，等.膠東金礦成礦模式［J］.礦床地質，2020，39（2）：215-236.（SONG M C，LIN S Y，YANG L Q，et al.Metallogenic model of Jiaodong Peninsula gold deposits［J］.Mineral Deposits，2020，39（2）：215-236.（in Chinese））

［24］MATTHI S K，BELAYNEH M.Fluid flow partitioning between fractures and a permeable rock matrix［J］.Geophysical Research Letters，2004，31（7）：L07602.

［25］SHEN C J，F(xiàn)AN Y，WANG X B，et al.Effects of land reclamation on a subterranean estuary［J］.Water Resources Research，2022，58（6）：e2022WR032164.

Abstract：Exploring the seawater intrusion mechanism in coastal aquifer with horizontal double fractures is critical for understanding the seawater intrusion mechanism in coastal region.Sand-tank experiments and?numerical simulation based on discrete fracture matrix model were performed.The effects of horizontal double fractures with different spatial distribution on the seawater intrusion were analyzed in terms of the saltwater wedge and mixing zone.Results showed that：① The development of saltwater wedge and mixing zone can be enhanced significantly by the horizontal double fractures.The seawater-freshwater mixing between the fractures is more pronounced.The individual effects of the two fractures are limited to each other.The effects are more significant for the fracture that is horizontally close to the offshore and vertically more close to the middle of the aquifer.② A higher offshore-side fracture can exacerbate the development of the saltwater wedge and the mixing between fractures.However，the height of land-side fracture has a significant influence only on the tail of the saltwater wedge.③ When the upper or lower fracture is horizontally closer to the inland，the saltwater wedge becomes wider with the total volume being shrank.However，once they are distributed so far to the inland that being out of the range of the saltwater wedge，the mixing zone area and the width of the saltwater wedge can drop sharply.

Key words：seawater intrusion;groundwater dynamics;horizontal double fractures;laboratory experiment;numerical simulation