俯沖板塊脫水相變對(duì)板塊俯沖過(guò)程的影響

王瑞澤, 王建超, 金宗瑋, 劉駿標(biāo), 王雪純, 靳國(guó)棟, 邢會(huì)林, 2*

俯沖板塊脫水相變對(duì)板塊俯沖過(guò)程的影響

王瑞澤1, 王建超1, 金宗瑋1, 劉駿標(biāo)1, 王雪純1, 靳國(guó)棟1, 邢會(huì)林1, 2*

(1. 深海圈層與地球系統(tǒng)教育部前沿科學(xué)中心, 海底科學(xué)與探測(cè)技術(shù)教育部重點(diǎn)實(shí)驗(yàn)室, 中國(guó)海洋大學(xué) 海洋地球科學(xué)學(xué)院, 山東 青島 266100; 2. 嶗山實(shí)驗(yàn)室, 山東 青島 266237)

板塊俯沖是地球上最重要的動(dòng)力學(xué)過(guò)程之一, 關(guān)于俯沖帶相變對(duì)地幔楔和上覆板塊影響前人已開展大量研究, 但關(guān)于相變反應(yīng)對(duì)俯沖板塊的動(dòng)力學(xué)影響研究較少。本研究在拉格朗日積分點(diǎn)有限元代碼的基礎(chǔ)上, 實(shí)現(xiàn)了俯沖板塊的相變, 建立了二維自由俯沖模型, 研究了俯沖過(guò)程中俯沖板塊上礦物脫水相變對(duì)俯沖速度、俯沖板塊形態(tài)、板塊上的應(yīng)力變化的影響。同時(shí)我們研究了薄弱帶的黏度變化對(duì)俯沖起始與俯沖過(guò)程的影響。數(shù)值模擬實(shí)驗(yàn)結(jié)果顯示: ①俯沖板塊脫水相變現(xiàn)象的發(fā)生降低了板塊俯沖的速度, 改變了俯沖板塊的形態(tài), 在俯沖板塊的相變界面處發(fā)生應(yīng)力累積; ②相較于礦物相變導(dǎo)致的黏度增大, 俯沖起始角度仍是影響俯沖板塊俯沖速度的主要因素; ③薄弱帶黏度越高越不利于俯沖板塊與上覆板塊的解耦, 抑制了俯沖過(guò)程的進(jìn)行。該研究為探究板塊回卷的控制因素以及俯沖板塊應(yīng)力狀態(tài)變化的影響因素提供了新的視角。

俯沖板塊; 俯沖過(guò)程; 礦物相變; 有限元數(shù)值模擬

0 引 言

俯沖帶是地幔對(duì)流單元的下降分支, 是地球內(nèi)部占主導(dǎo)地位的物理和化學(xué)系統(tǒng), 也是地球上最大的“回收系統(tǒng)”(Stern, 2002)。俯沖板塊將遠(yuǎn)洋陸源沉積物、蝕變的和新鮮的玄武質(zhì)洋殼物質(zhì)以及地幔巖石圈等“原材料”運(yùn)送到俯沖帶, 使海洋巖石圈、沉積物和海水與周圍的地幔物質(zhì)發(fā)生熔融, 重新平衡(Tatsumi, 2005)。

隨著俯沖深度的增加, 俯沖帶中溫度和壓力不斷增大。俯沖帶熱力學(xué)模擬是研究俯沖帶熱結(jié)構(gòu)最好的方法, 通過(guò)計(jì)算發(fā)現(xiàn)板塊在俯沖過(guò)程中被高溫的地幔加熱(Peacock and Wang, 1999; Gerya et al., 2002; Gorman et al., 2006; Syracuse et al., 2010)。同時(shí), 研究證明, 在高溫高壓條件下, 俯沖板塊中礦物主要發(fā)生兩類相變反應(yīng): 一類為俯沖板塊上深度較淺的玄武巖?榴輝巖系列相變(Schmidt and Poli, 1998; Fryer et al., 1999; Forneris and Holloway, 2003; Forneris and Holloway, 2004); 另一類則是發(fā)生在俯沖地幔以及俯沖板塊上覆地幔楔中的橄欖巖?蛇紋巖?橄欖巖系列相變。蛇紋石化是一個(gè)重要的地球動(dòng)力學(xué)過(guò)程, 它極大地改變了地幔和洋殼的物理和力學(xué)性質(zhì)。蛇紋巖是由超基性巖和基性巖在俯沖帶淺部、低溫(100 ℃)?中溫(700 ℃)條件下水合而成(Bose and Ganguly, 1995;Ulmer and Trommsdorff, 1995; Kawamoto and Holloway, 1997; Grove et al., 2006; Till et al., 2012; Green et al., 2014)。隨著俯沖深度的增加, 在高溫高壓條件下, 蛇紋巖失穩(wěn)發(fā)生脫水相變反應(yīng), 形成橄欖石與其他礦物組合(Ulmer and Trommsdorff, 1995; Wunder and Schreyer, 1997; Evans, 2004; Perrillat et al., 2005; Schwartz et al., 2013; Guillot et al., 2015)。

由于俯沖帶巖石樣品獲取難度較大, 而實(shí)驗(yàn)又難以解決俯沖過(guò)程中復(fù)雜的地質(zhì)問(wèn)題。隨著計(jì)算機(jī)技術(shù)的進(jìn)步, 數(shù)值模擬可以有效地解決地質(zhì)問(wèn)題中存在的時(shí)空尺度局限性問(wèn)題, 因此許多學(xué)者通過(guò)數(shù)值模擬手段對(duì)俯沖帶的相變反應(yīng)及其影響進(jìn)行了研究。Gerya et al. (2002)結(jié)合前人巖石學(xué)成果, 基于有限差分方法, 針對(duì)俯沖通道的發(fā)展和地幔楔水化作用這一問(wèn)題進(jìn)行了模擬研究; Arcay et al. (2005)通過(guò)有限元數(shù)值模擬的方法, 研究俯沖板片脫水對(duì)地幔楔的影響; Wada et al. (2008)通過(guò)有限元數(shù)值模擬, 研究上覆板塊與俯沖板塊之間薄弱帶以及板塊脫水對(duì)地幔楔的影響; Gerya and Meilick (2011)通過(guò)有限差分程序, 研究俯沖帶熔體和流體對(duì)流變性質(zhì)的影響, 以及其對(duì)俯沖動(dòng)力學(xué)的影響; Faccenda et al. (2012)基于有限差分方法, 研究板塊相變脫水對(duì)中?深地震活動(dòng)、板塊弱化和深部水循環(huán)的影響; Liao et al. (2017)利用同樣的方法, 研究水化作用引起的上地幔弱化對(duì)克拉通巖石圈動(dòng)力學(xué)的影響。

盡管前人對(duì)俯沖板塊的相變過(guò)程以及相變脫水對(duì)地幔物質(zhì)的影響進(jìn)行了大量研究, 但主要側(cè)重于相變反應(yīng)對(duì)地幔楔以及上覆板塊的影響, 而很少考慮俯沖過(guò)程中礦物脫水相變對(duì)俯沖板塊動(dòng)力學(xué)行為的影響(Schellart et al., 2010; Sandiford et al., 2019)。為了充分研究板塊俯沖至平衡深度時(shí)發(fā)生蛇紋巖巖石圈地幔相變對(duì)俯沖速度、俯沖板片形態(tài)的影響, 以及對(duì)俯沖板塊應(yīng)力應(yīng)變狀態(tài)的影響, 本次研究基于有限元數(shù)值模擬方法, 通過(guò)構(gòu)建俯沖板塊自由俯沖模型, 對(duì)俯沖板塊中各種材料參數(shù)進(jìn)行設(shè)置, 對(duì)板塊俯沖過(guò)程進(jìn)行模擬, 并結(jié)合前人研究成果, 進(jìn)一步討論礦物相變對(duì)板塊俯沖的影響。

1 控制方程

1.1 有限元部分

通過(guò)求解穩(wěn)態(tài)不可壓縮蠕動(dòng)流的斯托克斯方程進(jìn)行數(shù)值模擬。

其中質(zhì)量守恒方程為:

動(dòng)量守恒方程為:

對(duì)于不可壓縮流, 可通過(guò)公式(3)計(jì)算偏應(yīng)力張量:

本次研究采用拉格朗日積分點(diǎn)有限元法進(jìn)行有限元數(shù)值模擬計(jì)算。拉格朗日積分點(diǎn)有限元法主要優(yōu)點(diǎn)在于: ①能在不顯著改變精度的情況下繼續(xù)進(jìn)行極大變形; ②能夠跟蹤材料歷史和界面, 其精度高; ③在整個(gè)運(yùn)行過(guò)程中保留規(guī)則網(wǎng)格, 可以使用快速數(shù)值求解器。求解域由歐拉網(wǎng)格和嵌入的拉格朗日追蹤點(diǎn)或粒子表示, 未知變量在網(wǎng)格節(jié)點(diǎn)上計(jì)算, 拉格朗日粒子在變形過(guò)程中攜帶歷史變量。該方法非常適合于模擬地質(zhì)環(huán)境中經(jīng)常遇到的連續(xù)介質(zhì)固體的類流體行為, 能夠?qū)Πǖ蒯?duì)流、板塊俯沖和大陸碰撞在內(nèi)的構(gòu)造過(guò)程進(jìn)行二維/三維熱?變形?化學(xué)反應(yīng)建模(Moresi et al., 2003, 2007)。

1.2 相變部分

葉蛇紋石是俯沖帶高溫高壓蛇紋石的變種, 是深度大于100 km俯沖帶中最主要的含水礦物, 俯沖帶中主要脫水相變過(guò)程(圖1)(Guillot et al., 2015)為

圖1 根據(jù)實(shí)驗(yàn)數(shù)據(jù)和天然蛇紋石巖石學(xué)觀測(cè)得出的蛇紋石穩(wěn)定場(chǎng)(Perrillat et al., 2005; Guillot et al., 2015)

Fig.1 Serpentine stability fields based on experimental data and natural serpentine lithological observations

據(jù)Carminati et al. (1999)和Schwartz et al. (2001)的研究, 俯沖帶蛇紋巖脫水相變會(huì)導(dǎo)致巖石流變強(qiáng)度變大。本次研究基于有限元理論、固定的計(jì)算網(wǎng)格和拉格朗日積分點(diǎn)法的程序Underworld, 通過(guò)判斷板塊上的追蹤粒子是否發(fā)生相變進(jìn)而改變材料黏度, 在程序中實(shí)現(xiàn)了俯沖板片上蛇紋巖脫水相變這一過(guò)程(圖2)。

圖2 相變部分算法計(jì)算流程圖

Fig.2 The algorithm flow chart of phase transition

圖3 不同俯沖起始角度模型幾何示意圖

2 模型設(shè)置

2.1 模型的幾何和流變學(xué)

在海洋板塊彎曲前端與大陸板塊之間設(shè)置寬度約10 km的薄弱帶, 初始模型中黏度設(shè)置為上地幔黏度的十分之一, 上覆板塊厚度設(shè)置為90 km, 與俯沖板塊俯沖起始深度相同。本次模擬實(shí)驗(yàn)只考慮主要的相變過(guò)程(即蛇紋巖中葉蛇紋石失穩(wěn)脫水相變, 形成橄欖石為主以及少量蛇紋石的巖石組分)對(duì)俯沖板塊的流變學(xué)性質(zhì)的影響, 模型中各種材料參數(shù)如表1所示。

表1 模型使用材料參數(shù)表

2.2 邊界條件

本研究假設(shè)成熟洋板塊和上地幔之間的重力不穩(wěn)定性導(dǎo)致了俯沖板塊的下沉, 俯沖是由板塊前端的彎曲尖端觸發(fā)。在模型的邊界施加自由滑動(dòng)速度邊界。因此, 只有負(fù)浮力(板塊拉張)才能驅(qū)使海洋板塊俯沖, 誘發(fā)地幔流動(dòng), 沒(méi)有速度或應(yīng)力施加到俯沖板塊或上覆板塊上。

2.3 模型設(shè)計(jì)

本模擬設(shè)置了俯沖起始角度為30°的對(duì)照組模型1～2, 不同俯沖起始角度的實(shí)驗(yàn)組模型3～6(表2)。對(duì)照組為在重力作用下板塊發(fā)生俯沖, 俯沖過(guò)程中俯沖板塊、上地幔等不同材料的黏度保持不變; 實(shí)驗(yàn)組為在重力作用下板塊發(fā)生俯沖, 在俯沖至固定深度時(shí), 俯沖板塊會(huì)發(fā)生相應(yīng)的相變, 從而導(dǎo)致黏度變化。本次研究主要探究板塊俯沖過(guò)程中相變對(duì)俯沖板塊內(nèi)應(yīng)力場(chǎng)的變化、俯沖板塊的動(dòng)力學(xué)行為的影響。同時(shí)根據(jù)模擬結(jié)果, 增設(shè)兩組模型以探討不同的薄弱帶黏度對(duì)俯沖板塊俯沖的影響(模型7、8), 以及兩組模型以探討不同薄弱帶寬度對(duì)俯沖板塊俯沖過(guò)程的影響(模型9、10)。

表2 模型參數(shù)設(shè)置表

3 結(jié)果分析

3.1 俯沖過(guò)程分析(模型1、2)

板塊俯沖運(yùn)動(dòng)是由板塊自身重力和克服板塊底部黏滯阻力的能力決定的, 其中下降板的負(fù)浮力是板塊運(yùn)動(dòng)的主要?jiǎng)恿?。板塊的自由俯沖模型結(jié)果顯示俯沖具有兩個(gè)階段。

階段一: 由于地幔和大洋板塊之間的重力不穩(wěn)定, 大洋板塊從預(yù)先設(shè)置的俯沖彎曲前端擾動(dòng)處開始俯沖(圖4a、e)。隨著板塊俯沖, 上覆板塊持續(xù)向海溝方向后退。上地幔在俯沖作用下發(fā)生擾動(dòng), 以平行于俯沖板塊下方海溝的方向運(yùn)動(dòng)(圖4b、d、f、g)。俯沖板塊后撤引起的水平向地幔角流既驅(qū)動(dòng)上覆板塊向海溝側(cè)運(yùn)動(dòng), 同時(shí)也拖拽俯沖板塊持續(xù)向下俯沖(Zhang et al., 2017)。階段二: 在俯沖進(jìn)行一段時(shí)間后, 自由俯沖模型中的俯沖板塊在重力作用下以接近90°向地幔深處俯沖(圖4d、h)。

(a)～(d) 為板塊未發(fā)生相變俯沖至50 Ma結(jié)果圖; (e)～(h) 為板塊發(fā)生相變俯沖至50 Ma結(jié)果圖。白色箭頭標(biāo)記為均方根速度(Vrms)。

當(dāng)俯沖起始角度為30°時(shí), 巖石圈地幔未發(fā)生相變時(shí), 整個(gè)俯沖板塊黏度均勻, 因此整個(gè)板塊俯沖趨勢(shì)保持一致。俯沖彎曲前端在自身重力的作用下, 俯沖角度由起始時(shí)的30°逐漸變?yōu)?0°垂直向下俯沖。若在模型中添加俯沖過(guò)程中蛇紋石化巖石圈地幔的相變反應(yīng)后, 俯沖過(guò)程則變得有所不同。由于大洋巖石圈在俯沖過(guò)程中發(fā)生了相變, 導(dǎo)致黏度變大, 進(jìn)而使得俯沖更難以進(jìn)行, 因此俯沖角度從30°至90°所需時(shí)間明顯少于未發(fā)生相變的模型。

3.2 蛇紋石化巖石圈地幔相變對(duì)俯沖速度的影響

為了研究不同俯沖起始角度對(duì)俯沖速度的影響, 以及蛇紋石化巖石圈地幔發(fā)生相變對(duì)俯沖速度帶來(lái)的影響, 本次研究設(shè)置了模型3～6。首先, 設(shè)定俯沖起始角度為30°, 俯沖至同一時(shí)刻發(fā)生相變與未發(fā)生相變的模型(模型1、2), 并進(jìn)行了均方根速度成圖(圖5)。結(jié)果顯示, 板塊俯沖至相同時(shí)刻下, 未發(fā)生相變的模型中俯沖板塊處的均方根速度明顯大于發(fā)生相變模型。由于板塊俯沖是在板塊的負(fù)浮力以及板塊與地幔之間的黏滯阻力共同作用下完成的, 而發(fā)生相變俯沖板塊底部需要克服的黏滯阻力變大, 導(dǎo)致俯沖板塊俯沖速度變小。

圖5 俯沖起始角度為30°模型俯沖至20.3 Ma時(shí)均方根速度圖

模擬中通過(guò)追蹤不同俯沖起始角度以及是否發(fā)生相變影響下, 各組模型俯沖板塊最底部粒子所在深度, 得到模型1～6俯沖深度與時(shí)間的關(guān)系。通過(guò)對(duì)未發(fā)生相變、俯沖起始角度不同的模型1、3、5進(jìn)行分析(圖6), 結(jié)果顯示,在相同俯沖時(shí)間內(nèi), 俯沖起始角度為30°的模型俯沖速度最快, 隨著俯沖起始角度的增加, 俯沖速度相應(yīng)變慢。在相同的俯沖時(shí)間內(nèi), 俯沖起始角度越小, 俯沖板塊俯沖所能達(dá)到的深度越深。

圖6 不同俯沖起始角度未發(fā)生相變模型俯沖速度圖

發(fā)生相變、俯沖起始角度不同的模型2、4、6俯沖深度隨俯沖時(shí)間變化的結(jié)果(圖7)顯示, 板塊相變?cè)谝欢ǔ潭壬嫌绊懥税鍓K俯沖的速度, 導(dǎo)致板塊俯沖速度變慢。但是隨著俯沖起始角度的增大, 板塊相變對(duì)俯沖速度的影響逐漸減小, 因此, 俯沖起始角度是控制俯沖板塊俯沖速度的主要因素。

圖7 不同俯沖起始角度發(fā)生相變模型俯沖速度圖

3.3 俯沖過(guò)程應(yīng)力變化分析

圖8 不同俯沖起始角度發(fā)生相變與未發(fā)生相變模型板塊上最大主應(yīng)力的相對(duì)變化

3.4 薄弱帶作用結(jié)果

為了研究薄弱帶的黏度和寬度在俯沖過(guò)程中的作用, 本次研究設(shè)置了模型7～10。當(dāng)模型中薄弱帶固定寬度為10 km, 薄弱帶黏度為低于俯沖板塊與上覆板塊的黏度1×1016Pa·s時(shí), 俯沖板塊彎曲前端可以較快速的與上覆板塊發(fā)生解耦, 在自身重力的作用下, 耗時(shí)約23 Ma進(jìn)入俯沖角度近乎90°的垂直俯沖階段, 并在俯沖后期發(fā)生板塊回卷現(xiàn)象(圖4h)。而增設(shè)兩組更大黏度的模型7和模型8則表明, 隨著薄弱帶流變強(qiáng)度增大, 俯沖角度從起始的30°變化至90°所需的時(shí)間變長(zhǎng)。其中薄弱帶黏度為1×1019Pa·s的模型7(圖9a、b); 即薄弱帶不存在, 俯沖洋殼與流變強(qiáng)度相同的上覆地殼直接接觸時(shí), 耗時(shí)約28 Ma變?yōu)楦_角度接近90°的垂直向下俯沖; 而薄弱帶黏度為1×1022Pa·s的模型8(圖9c、d), 則在俯沖角度俯沖至約45°時(shí), 俯沖趨勢(shì)就被上部高黏度的薄弱帶所抑制, 黏度比俯沖板塊和上覆板塊高的薄弱帶不利于俯沖板塊與上覆板塊之間的解耦, 抑制了俯沖過(guò)程的進(jìn)行。

圖9 不同薄弱帶黏度設(shè)置下板塊俯沖模型黏度場(chǎng)與速度場(chǎng)

為探究不同的薄弱帶寬度對(duì)俯沖過(guò)程的影響, 本次研究增設(shè)了薄弱帶寬度為5 km的模型9和薄弱帶寬度為20 km的模型10。結(jié)果顯示, 在薄弱帶黏度固定為1×1016Pa·s時(shí), 薄弱帶寬度增大時(shí), 板塊俯沖速度明顯增大(圖10)。通過(guò)對(duì)比不同薄弱帶寬度模型中俯沖起始角度由30°變?yōu)榇怪备_階段所用的時(shí)間, 顯示薄弱帶寬度為20 km時(shí), 所需的時(shí)間約為15 Ma; 薄弱帶寬度為10 km時(shí), 所需時(shí)間約為23 Ma; 薄弱帶寬度為5 km時(shí), 所需時(shí)間約為27 Ma。因此, 薄弱帶黏度固定且為較低黏度時(shí), 薄弱帶寬度越大, 俯沖板塊與上覆板塊之間解耦越快。

圖10 不同薄弱帶寬度下發(fā)生相變模型俯沖速度圖

4 討 論

本次研究在拉格朗日積分點(diǎn)有限元方法的基礎(chǔ)上, 添加了與俯沖板塊蛇紋巖脫水相變相關(guān)的代碼, 研究了礦物相變對(duì)俯沖過(guò)程產(chǎn)生的一系列影響。但是由于本次研究的模型中只考慮材料黏度變化對(duì)俯沖板塊的影響, 因此只考慮了俯沖板塊中蛇紋巖脫水相變導(dǎo)致的黏度差異對(duì)俯沖過(guò)程的影響, 而認(rèn)為俯沖板塊中大洋中脊玄武巖系列的相變反應(yīng)對(duì)俯沖板塊的影響較小。

4.1 板塊俯沖形態(tài)

在前人研究的基礎(chǔ)上(Zhang et al., 2017), 選擇的基礎(chǔ)模型為自由俯沖模型, 沒(méi)有額外恒定的邊界條件, 所以板塊不能保持傾角恒定的持續(xù)俯沖, 表現(xiàn)出俯沖角度從30°逐漸向90°增加。而在模型2中, 由于俯沖板塊120 km以下板塊發(fā)生相變, 俯沖角度變化更大, 出現(xiàn)類似于Perchuk et al. (2020)模擬結(jié)果中的板塊回卷現(xiàn)象(圖4h)。前人研究認(rèn)為榴輝巖化導(dǎo)致的重力異常是引起板塊回卷的主要原因之一(孫圣思和嵇少丞, 2011; 皇甫鵬鵬等, 2016; Dai et al., 2020)。而在本次模擬中由于設(shè)定密度固定不變, 只有板塊俯沖至一定深度時(shí)黏度發(fā)生變化, 俯沖板塊才會(huì)發(fā)生了類似板塊回卷現(xiàn)象, 所以俯沖板塊相變導(dǎo)致板塊黏度變大可能也是導(dǎo)致板塊回卷現(xiàn)象的原因之一。

4.2 板塊應(yīng)力變化

模擬結(jié)果顯示, 板塊俯沖至不同深度時(shí), 在板塊發(fā)生相變的界面處出現(xiàn)應(yīng)力集中導(dǎo)致應(yīng)力陡變。應(yīng)力集中出現(xiàn)的原因是在俯沖過(guò)程中, 俯沖板塊相變導(dǎo)致黏滯力增大, 因此板塊自相變界面以下部分俯沖速度降低, 發(fā)生了短暫的停滯。Kita et al. (2006)發(fā)現(xiàn)日本東北部下方70～100 km深度處雙地震帶上平面存在帶狀地震活動(dòng)集中區(qū), 并將其解釋為與Hacker et al. (2003)模擬的相邊界附近的脫水反應(yīng)有關(guān); Tsuji et al. (2008)利用雙差層析成像技術(shù), 估算出日本東部太平洋板塊下一個(gè)厚度約為10 km的含水低速帶在70～90 km處逐漸發(fā)生脫水消失。本次模擬結(jié)果也證實(shí)了相變界面處發(fā)生了應(yīng)力集中。因此, 俯沖板塊黏度變化導(dǎo)致的應(yīng)力集中可能是導(dǎo)致雙地震帶上形成帶狀地震活動(dòng)集中區(qū)的原因之一。

4.3 薄弱帶的影響

俯沖帶的數(shù)值模擬實(shí)驗(yàn)中, 俯沖洋殼和上覆板塊之間的薄弱帶黏度的大小會(huì)影響板塊速度(Androvi?ová et al., 2013; Holt et al., 2015), 較高黏度的薄弱帶會(huì)抑制板塊回卷現(xiàn)象的發(fā)生(Arredondo and Billen, 2017)。因此, 俯沖洋殼和上覆板塊的解耦對(duì)俯沖板塊深部變形具有很強(qiáng)的敏感性。如圖4h、9b、9d所示, 較低的薄弱帶黏度導(dǎo)致板塊回卷較容易發(fā)生, 而較強(qiáng)的薄弱帶黏度導(dǎo)致板塊回卷現(xiàn)象難以發(fā)生。同時(shí), 薄弱帶的寬度也會(huì)對(duì)俯沖板塊與上覆板塊的解耦產(chǎn)生影響, 在薄弱帶黏度較低時(shí), 越寬的薄弱帶越有利于俯沖板塊與上覆板塊的解耦。

5 結(jié) 論

本次研究通過(guò)在拉格朗日積分點(diǎn)有限元代碼中添加俯沖板塊礦物相變部分, 建立俯沖板塊自由俯沖相變模型, 模擬了相變對(duì)俯沖板塊的動(dòng)力學(xué)行為與應(yīng)力狀態(tài)的影響, 以及不同的俯沖角度與薄弱帶對(duì)俯沖過(guò)程的影響, 得到以下幾點(diǎn)認(rèn)識(shí):

(1) 礦物相變導(dǎo)致了俯沖板塊上發(fā)生相變的部分黏度變大, 進(jìn)而使板塊俯沖需要克服的板塊底部黏滯阻力變大, 從而降低了板塊俯沖的速度, 改變了俯沖板塊的俯沖形態(tài), 同時(shí)在俯沖板塊的相變界面處發(fā)生應(yīng)力累積。

(2) 雖然俯沖板塊上相變降低了其俯沖速度, 但俯沖起始角度才是影響俯沖板塊俯沖速度的主要因素。

(3) 當(dāng)薄弱帶寬度一定時(shí), 黏度較低的薄弱帶有利于上覆板塊與俯沖板塊之間的解耦, 薄弱帶黏度越高越不利于俯沖板塊與上覆板塊的解耦, 在一定程度上抑制俯沖過(guò)程的進(jìn)行; 當(dāng)薄弱帶黏度(較低黏度)一定時(shí), 薄弱帶的寬度越大越有利于俯沖板塊與上覆板塊的解耦, 反之則會(huì)對(duì)俯沖過(guò)程起到抑制作用。

致謝:中國(guó)科學(xué)院廣州地球化學(xué)研究所郭鋒研究員及另一名匿名審稿人對(duì)本文提出了建設(shè)性修改意見, 在此致以特別感謝。

皇甫鵬鵬, 王岳軍, 范蔚茗, 李忠海, 王喻鳴, 周永智. 2016. 大洋平板俯沖的數(shù)值模擬再現(xiàn): 洋?陸匯聚速率影響. 大地構(gòu)造與成礦學(xué), 40(3): 429–445.

李三忠, 王光增, 索艷慧, 李璽瑤, 戴黎明, 劉一鳴, 周潔, 郭玲莉, 劉永江, 張國(guó)偉. 2019. 板塊驅(qū)動(dòng)力: 問(wèn)題本源與本質(zhì). 大地構(gòu)造與成礦學(xué), 43(4): 605–643.

孫圣思, 嵇少丞. 2011. 大洋板塊俯沖帶地震波各向異性及剪切波分裂的成因機(jī)制. 大地構(gòu)造與成礦學(xué), 35(4): 628–647.

Androvi?ová A, ?í?ková H, van den Berg A. 2013. The effects of rheological decoupling on slab deformation in the Earth’s upper mantle., 57(3): 460–481.

Arcay D, Tric E, Doin M. 2005. Numerical simulations of subduction zones: Effect of slab dehydration on the mantle wedge dynamics., 149(1–2): 133–153.

Arredondo K M, Billen M I. 2017. Coupled effects of phase transitions and rheology in 2-D dynamical models of subduction.:, 122(7): 5813–5830.

Bose K, Ganguly J. 1995. Experimental and theoretical studies of the stabilities of talc, antigorite and phase A at high pressures with applications to subduction processes., 136(3–4): 109–121.

Carminati E, Giunchi C, Argnani A, Sabadini R, Fernandez M. 1999. Plio-Quaternary vertical motion of the Northern Apennines: Insights from dynamic modeling., 18(4): 703–718.

Cloos M. 1993. Lithospheric buoyancy and collisional orogenesis: Subduction of oceanic plateaus, continental margins, island arcs, spreading ridges, and seamounts., 105(6): 715–737.

Dai L M, Wang L L, Lou D, Li Z H, Dong H, Ma F F, Li F K, Li S Z, Yu S Y. 2020. Slab rollback versus delamination: Contrasting fates of flat-slab subduction and implicationsfor South China evolution in the Mesozoic.:, 125(4), e2019JB019164.

Duretz T, Gerya T V, May D A. 2011. Numerical modelling of spontaneous slab breakoff and subsequent topographic response., 502(1): 244–256.

Evans B W. 2004. The serpentinite multisystem revisited: Chrysotile is metastable., 46(6): 479–506.

Faccenda M, Gerya T V, Mancktelow N S, Moresi L. 2012. Fluid flow during slab unbending and dehydration: Implications for intermediate-depth seismicity, slab weakening and deep water recycling.,,, 13(1), Q01010.

Forneris J F, Holloway J R. 2003. Phase equilibria in subductingbasaltic crust: Implications for H2O release from the slab., 214(1–2): 187–201.

Forneris J F, Holloway J R. 2004. Evolution of mineral compositions during eclogitization of subducting basaltic crust., 89(10): 1516–1524.

Fryer P, Wheat C G, Mottl M J. 1999. Mariana blueschist mud volcanism: Implications for conditions within the subduction zone., 27(2): 103–106.

Gerya T V, Meilick F I. 2011. Geodynamic regimes of subduction under an active margin: Effects of rheologicalweakening by fluids and melts., 29(1): 7–31.

Gerya T V, St?ckhert B, Perchuk A L. 2002. Exhumation of high-pressure metamorphic rocks in a subduction channel: A numerical simulation., 21(6), 1056.

Gerya T V, Yuen D A. 2003. Rayleigh-Taylor instabilities from hydration and melting propel ‘cold plumes’ at subduction zones., 212(1): 47–62.

Gorman P J, Kerrick D M, Connolly J. 2006. Modeling open system metamorphic decarbonation of subducting slabs.,,, 7(4), Q04007.

Green D H, Hibberson W O, Rosenthal A, Kovács I, Yaxley G M, Falloon T J, Brink F. 2014. Experimental study of the influence of water on melting and phase assemblagesin the upper mantle., 55(10): 2067– 2096.

Grove T L, Chatterjee N, Parman S W, Médard E. 2006. The influence of H2O on mantle wedge melting., 249(1–2): 74–89.

Guillot S, Schwartz S, Reynard B, Agard P, Prigent C. 2015. Tectonic significance of serpentinites., 646: 1–19.

Hacker B R, Peacock S M, Abers G A, Holloway S D. 2003. Subduction factory 2. Are intermediate-depth earthquakes in subducting slabs linked to metamorphic dehydration reactions?:, 108(B1), ESE11.

Holt A F, Becker T W, Buffett B A. 2015. Trench migration and overriding plate stress in dynamic subduction models., 201(1): 172–192.

Kawamoto T, Holloway J R. 1997. Melting temperature and partial melt chemistry of H2O-saturated mantle peridotite to 11 gigapascals., 276(5310): 240–243.

Kita S, Okada T, Nakajima J, Matsuzawa T, Hasegawa A. 2006. Existence of a seismic belt in the upper plane of the double seismic zone extending in the along-arc direction at depths of 70–100 km beneath NE Japan., 33(24), L24310.1- L24310.5.

Korenaga T, Korenaga J. 2016. Evolution of young oceanic lithosphere and the meaning of seafloor subsidence rate.:, 121(9): 6315–6332.

Liao J, Wang Q, Gerya T, Ballmer M D. 2017. Modeling craton destruction by hydration-induced weakening of the upper mantle.:, 122(9): 7449–7466.

Moresi L, Dufour F, Mühlhaus H. 2003. A Lagrangian integrationpoint finite element method for large deformation modeling of viscoelastic geomaterials., 184(2): 476–497.

Moresi L, Quenette S, Lemiale V, Meriaux C, Appelbe B, Mühlhaus H B. 2007. Computational approaches to studying non-linear dynamics of the crust and mantle., 163(1–4): 69–82.

Peacock S M, Wang K, 1999. Seismic consequences of warm versus cool subduction metamorphism: Examples from southwest and northeast Japan., 286(5441): 937–939.

Perchuk A L, Gerya T V, Zakharov V S, Griffin W L. 2020. Building cratonic keels in Precambrian plate tectonics., 586(7829): 395–401.

Perrillat J P, Daniel I, Koga K T, Reynard B, Cardon H, Crichton W A. 2005. Kinetics of antigorite dehydration: A real-time X-ray diffraction study., 236(3–4): 899–913.

Sandiford D, Moresi L, Sandiford M, Yang T. 2019. Geometric controls on flat slab seismicity., 527: 115787.

Schellart W P, Stegman D R, Farrington R J, Freeman J, Moresi L. 2010. Cenozoic tectonics of western North America controlled by evolving width of Farallon slab., 329(5989): 316–319.

Schmidt M W, Poli S. 1998. Experimentally based water budgets for dehydrating slabs and consequences for arc magma generation., 163(1–4): 361–379.

Schwartz S, Allemand P, Guillot S. 2001. Numerical model of the effect of serpentinites on the exhumation of eclogitic rocks: Insights from the Monviso ophiolitic massif (Western Alps)., 342(1–2): 193–206.

Schwartz S, Guillot S, Reynard B, Lafay R, Debret B, Nicollet C, Lanari P, Auzende A L. 2013. Pressure- temperature estimates of the lizardite/antigorite transition in high pressure serpentinites., 178: 197–210.

Stern R J. 2002. Subduction zones., 40(4): 3–1.

Syracuse E M, van Keken P E, Abers G A. 2010. The global range of subduction zone thermal models., 183(1–2): 73–90.

Tatsumi Y. 2005. The subduction factory: How it operates in the evolving Earth., 15(7): 4.

Till C B, Grove T L, Withers A C. 2012. The beginnings of hydrous mantle wedge melting., 163(4): 669–688.

Tsuji Y, Nakajima J, Hasegawa A. 2008. Tomographic evidence for hydrated oceanic crust of the Pacific slab beneath northeastern Japan: Implications for water transportation in subduction zones., 35(14), L14308.

Ulmer P, Trommsdorff V. 1995. Serpentine stability to mantle depths and subduction-related magmatism., 268(5212): 858–861.

Wada I, Wang K, He J, Hyndman R D. 2008. Weakening of the subduction interface and its effects on surface heat flow, slab dehydration, and mantle wedge serpentinization.:, 113(B4), B04402.

Wunder B, Schreyer W. 1997. Antigorite: High-pressure stability in the system MgO-SiO2-H2O (MSH)., 41(1–3): 213–227.

Zhang Q W, Guo F, Zhao L, Wu Y M. 2017. Geodynamics of divergent double subduction: 3-D numerical modeling of a Cenozoic example in the Molucca Sea region, Indonesia.:, 122(5): 3977–3998.

Influence of Dehydration Phase Transition of Subducting Plate on the Process of Plate Subduction

WANG Ruize1, WANG Jianchao1, JIN Zongwei1, LIU Junbiao1, WANG Xuechun1, JIN Guodong1, XING Huilin1, 2*

(1. Frontiers Science Center for Deep Ocean Multispheres and Earth System, MOE Key Lab of Submarine Geosciences and Prospecting Techniques, College of Marine Geosciences, Ocean University of China, Qingdao 266100, Shandong, China; 2. Laoshan Laboratory, Qingdao 266237, Shandong, China)

Plate subduction is one of the most important dynamic processes on the earth. The influences of phase change in subduction zones on the mantle wedge and overriding plates have been studied extensively, but researches on the influences of phase change on the dynamics of subducting plates are rare. Based on the Lagrangian integration point finite element code, we implemented the phase change of the plate into the code and performed 2D free subduction modelling to simulate the process of subduction. We investigated the influence of dehydration-related phase change on the subduction speed, shape of the subducting plate, and change of the stress field of the plate during the subducting process. Meanwhile we studied the influence of the viscosity change of the weak zone on the initiation and process of subduction. Based on these experiments we concluded that: (1) The phase change caused by mineral dehydration of the subducting plate reduces the speed of plate subducting and changes the shape of the subducting plate, and leads to stress accumulation at the phase change interface of the subducting plate. (2) Compared with the increase in viscosity caused by mineral phase transitions, initial subducting angle is still the main factor affecting the subduction velocity of subducting plates. (3) The higher the viscosity of the weak zone, the more difficult the decoupling of the subducting plate and the overlying plate, and thus inhibiting the progress of the subduction process. The results of simulation provide a new perspective for us to explore the controlling factors of rollback and the influencing factors of changes in the stress state of the subducting plate.

subducting plate; subducting process; phase change; finite element simulation

2021-12-05;

2022-03-28

國(guó)家自然科學(xué)基金重大計(jì)劃重點(diǎn)支持項(xiàng)目(92058211)、國(guó)家自然科學(xué)基金面上項(xiàng)目(52074251)、中央高?；究蒲袠I(yè)務(wù)經(jīng)費(fèi)(842012003)和高等學(xué)校學(xué)科創(chuàng)新引智計(jì)劃項(xiàng)目(B20048)聯(lián)合資助。

王瑞澤(1997–), 男, 碩士研究生, 海洋地質(zhì)專業(yè)。E-mail: ruize1997@163.com

邢會(huì)林(1965–), 男, 教授, 博士生導(dǎo)師, 主要從事超級(jí)計(jì)算地球科學(xué)理論、軟件研發(fā)及其應(yīng)用等研究。E-mail: h.xing@ouc.edu.cn

P541

A

1001-1552(2023)06-1232-010

10.16539/j.ddgzyckx.2023.06.003