章文博，張南，李昊遠(yuǎn)

北京大學(xué)地球與空間科學(xué)學(xué)院，北京 100871

0 引言

月幔被認(rèn)為由大碰撞后廣泛熔融的月球巖漿洋(Lunar Magma Ocean,LMO)結(jié)晶形成.巖漿洋逐漸冷卻凝固，結(jié)晶出不同的礦物，使得月幔自下而上形成了橄欖石、單斜輝石、斜方輝石的化學(xué)分層結(jié)構(gòu).最后階段，月幔會(huì)形成致密的富含鈦鐵礦等含鐵礦物的堆晶層 (以下簡(jiǎn)稱鈦鐵礦層)，同時(shí)具有不相容性的KREEP也會(huì)逐漸在結(jié)晶的殘余液體中富集(Snyder et al.,1992;Elardo et al.,2011;Elkins-Tanton et al.,2011;Lin et al.,2017,2020;Rapp and Draper,2018).由于鈦鐵礦層密度遠(yuǎn)大于下覆的月幔橄欖巖密度，因此堆晶層會(huì)發(fā)生重力失穩(wěn)下墜，此過程被稱為月幔倒轉(zhuǎn)(Lunar Cumulate Mantle Overturn)，許多天體演化的過程中都會(huì)經(jīng)歷類似的階段，比如地球和火星(Hess and Parmentier,1995;Parmentier et al.,2002).

巖漿洋結(jié)晶形成固態(tài)橄欖巖月幔和斜長(zhǎng)巖月殼，之后月球繼續(xù)經(jīng)歷地質(zhì)演化.一方面，月表早期經(jīng)歷了頻繁的撞擊事件導(dǎo)致月表部分熔融，進(jìn)而引發(fā)了許多巖漿事件(Norman,2009;Ding and Zhang,2018).另一方面，月幔部分熔融所產(chǎn)生的巖漿在表面噴發(fā)形成了月海玄武巖 (主要填充在地勢(shì)較低的區(qū)域，這些區(qū)域因?yàn)槌尸F(xiàn)暗色故被稱為月海)、月球玻璃等巖石(Delano,1986;Shearer et al.,2006).此外，盡管還有爭(zhēng)議(Mighani et al.,2020;Tarduno et al.,2021)，月球古強(qiáng)磁場(chǎng)的觀測(cè)建議月球內(nèi)部活動(dòng)至少持續(xù)到1.5 Ga(Mighani et al.,2020)，在這之后月球內(nèi)部活動(dòng)仍然持續(xù)，但強(qiáng)度已經(jīng)迅速衰減，其活動(dòng)證據(jù)為月球深震所記錄(Khan et al.,2014;Weiss and Tikoo,2014).

月幔演化主要的地質(zhì)和動(dòng)力學(xué)約束來自于月海玄武巖.月海玄武巖具有高FeO，低Al2O3，部分富含TiO2的巖石化學(xué)特征，根據(jù)TiO2的含量，被分為極低鈦(Very-low-Ti/VLT basalts,<1 wt.% TiO2)、低鈦(Low-Ti basalts,1～5 wt.% TiO2)以及高鈦玄武巖 (High-Ti basalts,>5 wt.% TiO2).空間上月海玄武巖的面積大約占月表面積的17%，主要分布在近地側(cè) (Head,1976).噴發(fā)時(shí)間主要集中于3.9～3.1 Ga(Nyquist et al.,2001;Hiesinger et al.,2011)，最年輕的月海玄武巖樣品年齡為1.96 Ga(圖1，Che et al.,2021).月海玄武巖的形成機(jī)制主要與月幔倒轉(zhuǎn)后富鈦物質(zhì)的上涌有關(guān)(Hess and Parmentier,1995;Zhong et al.,2000;Zhang et al.,2013a;Li et al.,2019).致密的鈦鐵礦層可以下沉進(jìn)入月幔，甚至到達(dá)月球的核幔邊界(Core Mantle Boundary,CMB).隨著鈦鐵礦層中攜帶的產(chǎn)熱元素的衰變，發(fā)生倒轉(zhuǎn)的鈦鐵礦層會(huì)隨時(shí)間逐漸積累熱量并升溫.一方面，鈦鐵礦層加熱過程會(huì)影響月球核幔邊界的熱流，從而激發(fā)月核的能量變化，進(jìn)而參與月球發(fā)電機(jī)的演化(Christensen et al.,2009;Weiss and Tikoo,2014)；另一方面，高溫的鈦鐵礦層會(huì)產(chǎn)生熱不穩(wěn)定性，在某些條件下會(huì)形成上涌的熱柱.在熱的鈦鐵礦物質(zhì)上升過程中，月幔也會(huì)被加熱熔融，鈦鐵礦物質(zhì)與月幔發(fā)生熔融與同化混染形成巖漿噴發(fā)到月表，可以解釋月海玄武巖和月球玻璃的形成機(jī)制(Elkins-Tanton et al.,2003;Grove and Krawczynski,2009;Barr and Grove,2013;Brown and Grove,2015).

圖1 Apollo玄武巖樣品,高地巖石樣品,玄武質(zhì)和斜長(zhǎng)質(zhì)月球隕石樣品的年齡的總結(jié).定年方法在圖例中給出月海玄武巖的年齡主要來自ApolloBasalt DB_V2(Cone et al.,2020)以及Snape等(2019).玄武質(zhì)隕石包括Dhofar 287,NWA 032,NWA 773,KAL 009,MIL 13317,Y 793169,A 881757,MET 01210,MIL 05035,EET 96008,LAP 02205,NEA 003A和NWA 4734.年齡來自于Merle等(2020)和文中的參考文獻(xiàn).月球玻璃包括Apollo 15綠色玻璃(Shih et al.,2001),Apollo 17橙色玻璃(Shih et al.,2001)和Apollo 17 VLT(Heiken et al.,1991).高地巖石年齡來自于Elkins-Tanton等(2011),Carlson等(2014),Borg等(2015)以及文中的參考文獻(xiàn).高地隕石包括 Y 86032(Nyquist et al.,2006)和Dhofar 025(Leont′eva et al.,2005).CSFD (Crater Size-Frequency Distribution)方法的結(jié)果來自于Hiesinger等(2011).1.96 Ga為嫦娥五號(hào)著陸區(qū)玄武巖樣品U-Pb年齡(Che et al.,2021).Fig.1 Compilation of ages from Apollo mare basalts,highland rocks,basaltic and anorthositic meteorites.Analytical methods are shown Ages of mare basalts are mainly from ApolloBasalt DB_V2 (Cone et al.,2020)with the addition of Snape et al.(2019).Basaltic meteorites include Dhofar 287,NWA 032,NWA 773,KAL 009,MIL 13317,Y 793169,A 881757,MET 01210,MIL 05035,EET 96008,LAP 02205,NEA 003A and NWA 4734.Their ages are from Merle et al.(2020)and reference therein.Lunar glasses include Apollo 15 green glass (Shih et al.,2001),Apollo 17 orange glass (Shih et al.,2001)and Apollo 17 VLT (Heiken et al.,1991).Summary of highlands rock ages are from Elkins-Tanton et al.(2011),Carlson et al.(2014),Borg et al.(2015)and reference therein.Highland meteorites consist of Y 86032 (Nyquist et al.,2006)and Dhofar 025 (Leont′eva et al.,2005).The result of CSFD (Crater Size-Frequency Distribution)is from Hiesinger et al.(2011).1.96 Ga is the U-Pb age of mare basalts from Chang-E′5 landing site (Che et al.,2021).

除了月海玄武巖的分布以外，月球表面KREEP的空間分布也具有不對(duì)稱性.遙感觀測(cè)表明，月表K、Th、U等元素主要在風(fēng)暴洋克里普巖地體(Procellarum KREEP Terrane,PKT)富集，而該區(qū)域與月海玄武巖的分布大致重合，因此兩者的空間分布很可能存在相關(guān)的作用機(jī)制(Lucey et al.,1998;Prettyman et al.,2006).

前人對(duì)月海玄武巖形成的內(nèi)動(dòng)力過程研究主要集中于對(duì)鈦鐵礦層下沉的月幔倒轉(zhuǎn)以及下沉鈦鐵礦層受熱發(fā)生上涌?jī)蓚€(gè)主要過程分別的數(shù)值模擬和參數(shù)探索，而將上述兩個(gè)過程統(tǒng)一起的數(shù)值模擬尚未被探索研究.

鈦鐵礦層的倒轉(zhuǎn)過程容易受到該層厚度和密度的影響(Parmentier et al.,2002;Li et al.,2019).這兩種性質(zhì)都依賴于巖漿洋的凝固過程.巖石學(xué)實(shí)驗(yàn)和計(jì)算模擬實(shí)驗(yàn)均給出了鈦鐵礦層的體積和結(jié)晶礦物組合的約束，前者對(duì)應(yīng)于鈦鐵礦層厚度，后者對(duì)應(yīng)于鈦鐵礦層密度.然而，由于這些實(shí)驗(yàn)對(duì)巖漿洋的初始化學(xué)組成和深度有不同的假設(shè)，它們得出的鈦鐵礦層的厚度和密度會(huì)有一定程度的變化(Snyder et al.,1992;Elardo et al.,2011;Elkins-Tanton et al.,2011;Lin et al.,2017,2020;Rapp and Draper,2018).

內(nèi)部生熱在模型演化中也起著至關(guān)重要的作用.235U,238U,232Th和40K的放射性衰變是熱化學(xué)演化的主要熱源.同時(shí)，這些元素都具有不相容性，所以會(huì)在月球巖漿洋結(jié)晶的殘余液體中富集，也就是KREEP.盡管現(xiàn)在關(guān)于全硅酸鹽月球(Bulk Silicate Moon,BSM)的產(chǎn)熱元素豐度方面仍然存在分歧，但月球上U和Th元素的豐度被認(rèn)為更可能接近于地球的豐度，而月球相對(duì)于地球更可能虧損鉀元素(Taylor,1982;Taylor and Wieczorek,2014).

前人的研究揭示了鈦鐵礦層的流變學(xué)性質(zhì)也對(duì)理解月球的熱化學(xué)演化有必不可少的幫助(Zhang et al.,2013a,b,2017;Dygert et al.,2016;Li et al.,2019).Zhang等(2013a)的數(shù)值模擬結(jié)果表明，如果鈦鐵礦層具有更高的活化能，那么被產(chǎn)熱元素加熱后，其黏度會(huì)更小，從而增強(qiáng)倒轉(zhuǎn)后的鈦鐵礦層內(nèi)部對(duì)流，反過來抑制熱量積累.因此，倒轉(zhuǎn)的鈦鐵礦層會(huì)更難產(chǎn)生熱不穩(wěn)定性，也就是說更難形成上涌熱柱.Zhang等(2017)考慮了鈦鐵礦的流變性質(zhì)，發(fā)現(xiàn)鈦鐵礦的軟化效應(yīng)會(huì)產(chǎn)生類似于高活化能的效果，這也不利于上涌熱柱的形成.Li等(2019)將化學(xué)組分和溫度相關(guān)的流變性質(zhì)應(yīng)用于月幔倒轉(zhuǎn)的數(shù)值模擬過程中，其結(jié)果表明鈦鐵礦對(duì)鈦鐵礦層黏度的減少會(huì)有利于月幔倒轉(zhuǎn)過程的進(jìn)行.

本文的工作主要是同時(shí)結(jié)合了月幔倒轉(zhuǎn)和鈦鐵礦層上涌的過程，并探索鈦鐵礦層厚度、內(nèi)部生熱、瑞利數(shù)、成分和溫度相關(guān)的流變學(xué)的影響，進(jìn)一步討論月海玄武巖的形成機(jī)制和相應(yīng)的生熱元素等的參數(shù)空間.

1 研究方法

1.1 控制方程與有限元數(shù)值模擬

月球熱化學(xué)演化的實(shí)現(xiàn)采用了三維球面幾何數(shù)值模型，模型假設(shè)無限大的普郎特?cái)?shù)(Prandtl Number)，使用了布辛涅斯克近似(Boussinesq Approximation)，化學(xué)成分以及溫度相關(guān)的流變學(xué)和參數(shù)化的月核演化方程(Zhong et al.,2008;Zhang et al.,2013a,b,2017;Li et al.,2019).

無量綱形式的模型控制方程為：

(1)

(2)

(3)

(4)

(5)

(6)

其中α,ρ0,g,ΔT,R0,Rc,η0,κ和ρIBC分別是熱膨脹系數(shù)、月幔密度、月表重力加速度 (重力加速度隨深度增加而減小)、上下邊界的初始溫度差、月球半徑、月核半徑、參考黏度、熱擴(kuò)散系數(shù)以及鈦鐵礦層密度.Ra和Rb可以用來定義浮力數(shù)B：

(7)

利用阿倫尼烏斯定律(Arrhenius Law)，我們計(jì)算鈦鐵礦層和月幔的混合密度：

(8)

其中ηIBC,E,R分別是鈦鐵礦在參考溫度下的黏度、有效活化能以及氣體常數(shù).

我們使用有限元計(jì)算軟件CitcomS來數(shù)值求解三維球殼中的不可壓縮流方程(1)—(4)，用來模擬月幔固結(jié)后內(nèi)部物質(zhì)的重力倒轉(zhuǎn)和加熱產(chǎn)生熔體上涌噴發(fā)的過程.球殼從下到上分別為核幔邊界、月幔層、鈦鐵礦層、月殼層以及月表，三維球殼被均分成12個(gè)區(qū)域，每個(gè)區(qū)域有65×65×65結(jié)點(diǎn).為了提高分辨率，鈦鐵礦層的結(jié)點(diǎn)在徑向上被進(jìn)一步細(xì)分，根據(jù)其厚度，結(jié)點(diǎn)數(shù)量從5～15不等.使用自由滑動(dòng)條件作為上下邊界的力學(xué)邊界條件.依照Z(yǔ)hang等(2013a)的處理方式，月表溫度的無量綱值恒為1，CMB的初始無量綱溫度為1，之后隨著月核與月幔之間的能量守恒關(guān)系演化.

1.2 參數(shù)空間

在我們的模型中，有四個(gè)關(guān)鍵參數(shù)需要被探索，分別是：鈦鐵礦層厚度、內(nèi)部生熱、參考黏度 (瑞利數(shù))以及鈦鐵礦層黏度.其他參數(shù)在表1中給出.

表1 模型參數(shù)Table 1 Model parameters

我們的模型研究了50～150 km的鈦鐵礦層厚度.對(duì)于1000 km深度的LMO深度，通過質(zhì)量守恒和完全分離結(jié)晶的相平衡關(guān)系得到的鈦鐵礦層厚度范圍是18～50 km，對(duì)應(yīng)的鈦鐵礦層中鈦鐵礦體積占比為16 vol.%到3.6 vol.%(Snyder et al.,1992;Elardo et al.,2011;Elkins-Tanton et al.,2011;Lin et al.,2017,2020;Rapp and Draper,2018;Li et al.,2019).然而，考慮到在LMO完全固結(jié)之前月幔倒轉(zhuǎn)就已經(jīng)啟動(dòng)了，因此我們引入固定鈦鐵礦總量但是對(duì)體積進(jìn)行一定程度稀釋的鈦鐵礦層.基于上述文獻(xiàn)提到的鈦鐵礦層厚度以及鈦鐵礦含量的上下界，我們考慮厚度和鈦鐵礦含量的算數(shù)平均值，即20 km的鈦鐵礦層厚度和11 vol.%的鈦鐵礦含量，并對(duì)其進(jìn)行稀釋 (圖2a).

圖2 (a)鈦鐵礦含量與鈦鐵礦層厚度的關(guān)系,修改自 Li等(2019).淺藍(lán)色實(shí)線代表50 km 和3.6 vol.%的下界.淺藍(lán)色虛線代表18 km和16 vol.%的上界.深藍(lán)色實(shí)線為算數(shù)平均值20 km和11 vol.%.(b)月球內(nèi)部生熱的四套方案.兩條實(shí)線分別為Turcotte和Schubert(2002)、McDonough和Sun(1995)的結(jié)果,兩條虛線為對(duì)應(yīng)的忽略鉀元素的結(jié)果.Taylor′s Moon的產(chǎn)熱量用黑色點(diǎn)劃線表示，與Turcotte和Schubert (2002)無鉀方案相當(dāng)，即紅色虛線Fig.2 (a)The derived relation between the ilmenite volume fraction and IBC thickness,from Li et al.(2019).Solid light blue line is the lower bound with 50 km and 3.6 vol.%.Dashed light blue line is the upper bound with 18 km and 16 vol.%.Solid dark blue line is the average value with 20 km and 11 vol.%.(b)Four schemes of lunar internal heating.Solid lines are original estimates from Turcotte and Schubert (2002)and McDonough and Sun (1995),respectively.Dash lines are corresponding K-free results.Dash-dotted black line is heat-producing rate of Taylor′s Moon,which is equivalent to the K-free scheme of Turcotte and Schubert (2002),that is,the dashed red line

正如之前提到的，鈾、釷元素作為難熔元素，其在BSM的含量與全硅酸鹽地球(Bulk Silicate Earth,BSE)的含量接近，但是作為揮發(fā)性元素的鉀元素在月球中的含量則有很大的不確定度(Taylor and Wieczorek,2014).因此，對(duì)鈾、釷元素，我們采用地球自然豐度.但是目前對(duì)全硅酸鹽地球的放射性同位素含量估計(jì)也存在一些分歧，于是我們考慮兩套方案：鈾、釷元素豐度的下界來自于McDonough和Sun(1995)的地幔巖模型，而上界則來自于Turcotte和Schubert(2002)的地球物理模型(表2).相對(duì)于地球，月球很可能高度虧損鉀元素，因此我們?nèi)サ羯鲜鰞蓚€(gè)模型的鉀元素豐度，從而導(dǎo)出兩個(gè)無鉀的方案.為了便于討論內(nèi)部生熱的效應(yīng)，我們將以上四套產(chǎn)熱方案統(tǒng)一到鈾、釷元素豐度分別為20.3 ng·g-1和79.5 ng·g-1的情況，相應(yīng)的K/U比值分別為0，8500，11800和23500(圖2b).前人的數(shù)值和化學(xué)模型所使用的豐度含量大多落在K/U值為0～8500的范圍內(nèi)(Zhong et al.,2000;Zhang et al.,2013a,b,2017;Laneuville et al.,2018;Li et al.,2019).此外，我們也考慮了月球相對(duì)地球更富集鈾、釷元素的模型，即 Taylor′s Moon(Taylor et al.,2006)，發(fā)現(xiàn)其產(chǎn)熱量與K/U=8500的產(chǎn)熱方案相當(dāng)(圖2b).被鈦鐵礦攜帶的KREEP含量目前仍然缺少實(shí)驗(yàn)約束，因此我們采用了Zhong等(2000)和Zhang等(2013a)的假設(shè)，即模型開始時(shí)所有的KREEP與鈦鐵礦層完全混合.產(chǎn)熱元素的衰變常數(shù)依據(jù)Turcotte和Schubert(2002).

表2 U,Th和K的預(yù)計(jì)濃度Table 2 Estimated concentrations of U,Th and K

我們將參考溫度Tref設(shè)為1300 ℃，對(duì)應(yīng)的橄欖巖黏度η0范圍為5×1018到1×1021Pa·s(Hirth and Kohlstedt,2004).考慮1×1020Pa·s,5×1020Pa·s和1×1021Pa·s的參考黏度，對(duì)應(yīng)的熱瑞利數(shù)為7.6×106,1.5×106和7.6×105.鈦鐵礦層的黏度由鈦鐵礦含量、含水量以及物理狀態(tài)比如熔融和溫度決定 (Dygert et al.,2016).取鈦鐵礦層相對(duì)于橄欖巖的黏度為1×10-1到1×10-2.此外，考慮橄欖巖和鈦鐵礦層具有相同的活化能100 kJ·mol-1，雖然鈦鐵礦和橄欖石的流變學(xué)實(shí)驗(yàn)得到的活化能更大(Karato and Wu,1993;Hirth and Kohlstedt,2004;Dygert et al.,2016)，但那是在非線性流變的情況下得到的值，Christensen(1984)認(rèn)為在假設(shè)應(yīng)變速率不變的情況下，將實(shí)驗(yàn)得到的活化能乘以0.3，可以將非線性流變的效果轉(zhuǎn)化為線性流變的效果，約為100 kJ·mol-1，并由此定義出有效黏度.

1.3 初始條件

模型開始的初始年齡為 4.42 Ga，盡管最近的研究表明LMO完全結(jié)晶的年齡可能更晚，但是我們的模型仍然將月球斜長(zhǎng)巖樣品的平均年齡作為L(zhǎng)MO結(jié)晶末期年齡，即4.42 Ga(Elkins-Tanton et al.,2011;Borg et al.,2015,2019;Maurice et al.,2020).我們的模型采用了和Zhang等(2013a)一樣的橄欖巖固相線溫度作為初始溫度剖面(圖3c)，模型初始密度分布在圖3d中給出.算例的命名方式會(huì)給出鈦鐵礦層的厚度、相對(duì)黏度、參考黏度以及K/U比值.比如對(duì)于算例“D50v1e-1V1e20K0”，即代表初始的鈦鐵礦層厚度為50 km，相對(duì)黏度為1×10-1，參考黏度為1×1020Pa·s，K/U比值為0.

圖3 (a)鈦鐵礦層厚度為50 km和75 km的算例初始黏度剖面.(b)鈦鐵礦層厚度為100 km和150 km的算例初始黏度剖面.(c)初始溫度剖面以及四個(gè)算例即將產(chǎn)生上升熱柱時(shí)的溫度剖面.(d)四種鈦鐵礦層厚度的初始密度剖面Fig.3 (a)Initial viscosity profiles for cases of 50 km and 75 km.(b)Initial viscosity profiles for cases of 100 km and 150 km.(c)Initial temperature profile and temperature profiles in the beginning of upwelling for four cases.(d)Initial density profiles for four sets of IBC thickness

1.4 玄武巖上涌時(shí)間定義

為了研究發(fā)生月幔倒轉(zhuǎn)后的內(nèi)部演化，我們計(jì)算了演化過程中，在深度為450 km處的球面上，TiO2含量大于0.46 wt.%的結(jié)點(diǎn)占比，我們認(rèn)為結(jié)點(diǎn)占比超過15%是成功的上涌熱柱的標(biāo)志，基于此，成功上涌的時(shí)間被定義為：

tupwelling=t[ratio=15%]-4.42 Ga,

(9)

我們使用該簡(jiǎn)化準(zhǔn)則來判斷是否形成了上涌熱柱，理由如下：(1)月海玄武巖來源于更為分異的巖漿，而與玄武巖共存的月球玻璃則被認(rèn)為代表了更為原始的巖漿成分(Shearer et al.,2006；Xiao et al.,2021).本模型不考慮復(fù)雜的化學(xué)演化，因此將鎂鐵質(zhì)月球玻璃的組分作為玄武巖的替代指標(biāo).月球玻璃的多相飽和壓強(qiáng)(Multiple Saturation Pressures)可以通過高溫高壓實(shí)驗(yàn)確定，前人的實(shí)驗(yàn)表明，月球玻璃的最大深度約為450 km，其對(duì)應(yīng)的熔體的TiO2含量為0.46 wt.%(圖4，Delano,1986;Elkins-Tanton et al.,2003;Grove and Krawczynski,2009).因此，我們采用0.46 wt.%作為模型計(jì)算的玄武巖TiO2含量.(2)對(duì)月表玄武巖的觀測(cè)表明，月海玄武巖的面積占月表面積的16%～17%，并且其噴發(fā)年齡主要分布于3.8～3.3 Ga之間(Head,1976;Whitten and Head,2015).為了模擬月海玄武巖母巖漿的生成，我們將月海玄武巖所占面積的比例大致對(duì)應(yīng)于在源區(qū)450 km深度時(shí)的面積比例，即在450 km深度達(dá)到(1)中的TiO2含量的結(jié)點(diǎn)比例.

圖4 綠色玻璃的密度、多相飽和壓強(qiáng)以及相應(yīng)TiO2含量，圖修改自 Grove和Krawczynski (2009)，TiO2含量來自Delano(1986)Fig.4 Green glasses′ density,multiple saturation pressure from Grove and Krawczynski (2009)and contents of TiO2 from Delano (1986)

2 結(jié)果

2.1 鈦鐵礦層厚度的影響

首先，我們研究了鈦鐵礦層厚度對(duì)上涌時(shí)間的影響.當(dāng)參考黏度為1×1020Pa·s,K/U比為0,鈦鐵礦層相對(duì)黏度為1×10-1時(shí)，分別指定鈦鐵礦層厚度為50 km、75 km、100 km和150 km(圖5a).

圖5 不同鈦鐵礦層厚度，K/U 比值，參考黏度下的上涌時(shí)間圓圈、叉號(hào)以及三角形符號(hào)分別代表成功上涌的算例、無法形成上涌熱柱的算例以及根據(jù)規(guī)律推測(cè)的無法形成上涌的參數(shù)空間.四條豎線分別代表圖2b中的四種生熱方案.(a),(b),(c)分別為不同參考黏度下的結(jié)果.(b)和(c)中的粗虛線是能否成功形成上升熱柱的分界線.Fig.5 Upwelling time with respect to the IBC thickness,K/U ratio,reference viscosity Circles,crosses and triangles are cases of successful upwelling,failed upwelling and deduced failed upwelling according to upwelling patterns,respectively.Four vertical lines correspond to four schemes in Fig.2b.(a),(b),(c)are results of different reference viscosity.Thick dashed lines in (b)and (c)are dividing lines between successful upwelling and that of failed upwelling.

標(biāo)準(zhǔn)算例為D50v1e-1V1e20K0.熱化學(xué)演化包括三個(gè)階段.開始階段，致密的冷的鈦鐵礦層下沉(圖6a,6b).鈦鐵礦層將一直下沉直到核幔邊界深度，最終經(jīng)過約200 Ma后，48.5 vol.%的鈦鐵礦層會(huì)進(jìn)入下月幔.第二個(gè)階段，在底部聚集的鈦鐵礦層會(huì)停留約100 Ma，期間鈦鐵礦層所攜帶的產(chǎn)熱元素會(huì)不斷衰變并加熱鈦鐵礦層(圖6c，6d).在約4.12 Ga時(shí)，被加熱的鈦鐵礦層產(chǎn)生了足夠的浮力，開始形成上涌熱柱(圖6e，6f).在接下來的100 Ma間，底部的鈦鐵礦層會(huì)不斷地被熱能驅(qū)動(dòng)，形成很強(qiáng)的上涌熱柱，熱柱的三維形態(tài)經(jīng)過三維球諧函數(shù)分解，主要的球諧階數(shù)為3階(圖6g，6h).

圖6 標(biāo)準(zhǔn)算例D50v1e-1V1e20K0的三維熱化學(xué)結(jié)構(gòu)以及相應(yīng)的二維溫度剖面(a)和(b)是倒轉(zhuǎn)階段.(c)和(d)是加熱階段.(e)和(f)是上涌的啟動(dòng)階段.(g)和(h)是上涌階段.三維結(jié)構(gòu)中，透明的綠色等值面代表化學(xué)組分為0.14的面.紅色和藍(lán)色等值面分別代表溫度比球面平均高250 ℃，低250 ℃的面.白色的球面代表月核.在(h)中,箭頭代表流場(chǎng)速度.三維和二維圖中的藍(lán)色圓代表450 km深度.Fig.6 3-D thermochemical structures and corresponding 2-D cross sections for the standard case D50v1e-1V1e20K0(a)and (b)are overturned phase.(c)and (d)are heating phase.(e)and (f)are onset of upwelling.(g)and (h)are upwelling phase.In the 3-D view,transparent green iso-surfaces represent the chemical composition at a contour level of 0.14.Red and blue iso-surfaces represent the residue temperature with the contour levels at 250 ℃,-250 ℃.White sphere stands for the lunar core.In (h),the arrows show the flow velocity.Blue circles in both 3-D and 2-D graphs stand for the 450-km depth.

與D50v1e-1V1e20K0相比，更厚的鈦鐵礦層，即算例D75v1e-1V1e20K0,D100v1e-1V1e20K0,D150v1e-1V1e20K0有著相似的演化模式.但是，也存在一些值得關(guān)注的區(qū)別.由上涌時(shí)間的變化可以看出，其明顯與鈦鐵礦層厚度負(fù)相關(guān).隨著鈦鐵礦層厚度從50增大到150 km，上涌時(shí)間從387 Ma單調(diào)縮短到354 Ma.對(duì)上涌時(shí)間的影響來源于月幔倒轉(zhuǎn)、加熱和上涌過程各自的貢獻(xiàn).一方面，更厚的鈦鐵礦層有利于月幔倒轉(zhuǎn)過程的進(jìn)行，這不僅會(huì)增加進(jìn)入下月幔鈦鐵礦層的比例，形成更厚的底部鈦鐵礦層，還會(huì)縮短發(fā)生月幔倒轉(zhuǎn)的時(shí)間尺度(圖7a).另一方面，更厚的鈦鐵礦層會(huì)導(dǎo)致更小的化學(xué)密度，而這有利于底部鈦鐵礦層更快地產(chǎn)生上涌熱柱.

2.2 內(nèi)部生熱的影響

其次，我們探索了內(nèi)部生熱的影響.與2.1節(jié)類似，固定鈦鐵礦層厚度為50 km，鈦鐵礦層相對(duì)黏度為1×10-1，參考黏度為1×1020Pa·s，分別指定K/U比值為0、8500、11800和23500(圖5a).

K/U比值的增加等效于內(nèi)部生熱的增加，因此有利于形成上涌熱柱.上涌時(shí)間從387 Ma縮短到184 Ma.除了演化時(shí)間以外，其他指標(biāo)如發(fā)生倒轉(zhuǎn)的鈦鐵礦層比例則幾乎不受內(nèi)部生熱影響.因此，增加內(nèi)部生熱僅僅會(huì)加速加熱階段，而對(duì)其他階段的影響很小(圖7b，圖8).

圖7 下月幔中鈦鐵礦層比例隨時(shí)間的演化垂直虛線代表上涌時(shí)間點(diǎn)，下同.(a)代表不同鈦鐵礦層厚度的四個(gè)算例.(b)代表不同內(nèi)部生熱的四個(gè)算例.Fig.7 Evolutions of overturned IBC fraction in the lower lunar mantleVertical dashed lines stand for the upwelling time and are the same as following figures.(a)is for four cases of different IBC thickness.(b)is for four cases of different internal heating.

圖8 內(nèi)部生熱對(duì)平均溫度的影響的演化圖Fig.8 Evolutions of the azimuthally averaged temperature with respect to internal heating

2.3 月幔參考黏度的影響

參考黏度/瑞利數(shù)對(duì)熱化學(xué)演化也有顯著的影響.考慮兩組算例，第一組算例為D150v1e-1V1e20K0,D150v1e-1V5e20K0和D150v1e-1V1e21K0；第二組算例為D50v1e-1V1e20K23500,D50v1e-1V5e20K23500和D50v1e-1V1e21K23500.

對(duì)第一組算例，可以看出上涌時(shí)間與參考黏度呈正相關(guān)關(guān)系.從1×1020Pa·s到1×1021Pa·s，參考黏度增加了一個(gè)數(shù)量級(jí)，相應(yīng)的上涌時(shí)間增加了一倍，即從354 Ma增加到了793 Ma(表3,圖9a、b、c).因此，增加參考黏度會(huì)極大推遲上涌熱柱的形成.

更重要的是，第二組算例證明了增加參考黏度甚至?xí)柚篃嶂纬?由于K23500是我們考慮的最大產(chǎn)熱方案，但僅有算例D50v1e-1V1e20K23500在184 Ma后形成了上升熱柱(表3,圖9d、e、f).實(shí)際上，增加參考黏度等效于減小瑞利數(shù)，而瑞利數(shù)是對(duì)流系統(tǒng)活躍程度的表征.除了上涌時(shí)間以外，倒轉(zhuǎn)的鈦鐵礦層比例也在很大程度上受到參考黏度的影響.對(duì)于以上考慮的兩組算例，黏度從1×1020Pa·s,5×1020Pa·s到1×1021Pa·s變化時(shí)，倒轉(zhuǎn)的鈦鐵礦層比例也減小為原來的50%(表3,圖9).最大的倒轉(zhuǎn)比例發(fā)生在1×1020Pa·s的算例中，分別是69.3%和48.7%，而最小的比例則發(fā)生在1×1021Pa·s的算例中，對(duì)應(yīng)比例分別減小到35.5%和22.0%(表3，圖9).因此，當(dāng)瑞利數(shù)減小時(shí)上涌熱柱更快形成的結(jié)果是符合預(yù)期的.

圖9 參考黏度對(duì)平均密度的影響的演化圖Fig.9 Evolutions of azimuthally averaged density with respect to reference viscosity

2.4 鈦鐵礦層黏度的影響

最后，我們分析了鈦鐵礦層黏度的影響效果，考慮了四個(gè)算例：D75v1e-1V5e20K11800,D75v1e-2V5e20K11800,D100v1e-1V5e20K8500和D100v1e-2V5e20K8500.

與鈦鐵礦層厚度、內(nèi)部生熱以及參考黏度不同，鈦鐵礦層黏度對(duì)于上涌時(shí)間的影響更為復(fù)雜.對(duì)于D75的兩個(gè)算例，鈦鐵礦層黏度從1×10-1Pa·s到1×10-2Pa·s的減小會(huì)導(dǎo)致上涌失敗(表3).相反，在D100的兩個(gè)算例，當(dāng)鈦鐵礦層黏度同樣減小一個(gè)數(shù)量級(jí)時(shí)，上涌時(shí)間從449 Ma減小到325 Ma.

表3 關(guān)鍵輸入輸出參數(shù)Table 3 Key input and output parameters

解釋上述效應(yīng)需要考慮兩個(gè)因素，分別是底部鈦鐵礦層的厚度以及化學(xué)組分的混合效應(yīng).首先，更小的鈦鐵礦層黏度有利于倒轉(zhuǎn)過程的進(jìn)行，因?yàn)殁佽F礦層黏度越小，在受重力影響進(jìn)入下月幔的過程中會(huì)受到更小的阻力，使得倒轉(zhuǎn)過程更早啟動(dòng)，并且有更多的鈦鐵礦層會(huì)發(fā)生倒轉(zhuǎn)(圖10).對(duì)于更厚的底部鈦鐵礦層，熱量會(huì)更容易聚集起來，從而有利于上涌熱柱的形成(圖3c).其次，化學(xué)組分的混合效應(yīng)也對(duì)上涌機(jī)制有著重要的影響.但鈦鐵礦層的相對(duì)黏度為1×10-1時(shí)，我們注意到鈦鐵礦層和月幔發(fā)生了更高程度的混合，而相對(duì)黏度為1×10-2時(shí)，鈦鐵礦層與月?；旌铣潭葎t更低(圖11b，d，f，h).發(fā)生混合的鈦鐵礦層比純的鈦鐵礦層密度更小，因此在加熱過程中更容易獲得足夠的浮力.所以，在某種程度上，鈦鐵礦層混合越充分，越容易形成上涌熱柱.

圖10 鈦鐵礦層黏度對(duì)進(jìn)入下月幔的鈦鐵礦層比例的影響的演化圖Fig.10 Evolutions of overturned IBC fraction in the lower lunar mantle with respect to IBC viscosity

圖11 鈦鐵礦層黏度的影響中的4個(gè)算例在上涌之處的三維形態(tài)結(jié)構(gòu)圖和對(duì)應(yīng)二維化學(xué)組分圖三維圖(a),(c)和(f)的綠色透明等值面代表化學(xué)組分為0.21的面，在(e)中為0.27的等值面.其他均與圖6相同.Fig.11 3-D thermochemical structures and corresponding 2-D chemical cross sections for 4 cases in the beginning of upwelling to test the effect of IBC viscosityIn the 3-D view,transparent green iso-surfaces in (a),(c)and (f)represent the chemical composition at a contour level of 0.21 and that in (e)is at a contour level of 0.27.Other color schemes are the same as those in Fig.6.

總而言之，鈦鐵礦層黏度更小，會(huì)導(dǎo)致底部鈦鐵礦層厚度更大，但是化學(xué)混合程度更弱.這兩個(gè)效應(yīng)對(duì)于上涌時(shí)間是相互對(duì)抗的.對(duì)于D100的兩個(gè)算例，鈦鐵礦層黏度對(duì)于底部厚度的影響效果超過了化學(xué)混合的效應(yīng)(圖11e,f,g和h)，而對(duì)于D75的兩個(gè)算例，兩者效果的相對(duì)權(quán)重則發(fā)生了顛倒(圖11a,b,c和d).

3 討論

月海玄武巖主要集中在近地側(cè)區(qū)域，這種不對(duì)稱性需要長(zhǎng)波的上涌，即主導(dǎo)球諧階數(shù)最小(Zhong et al.,2000;Zhang et al.,2013a,2017).基于當(dāng)前已知的鈦鐵礦的黏度，很遺憾我們暫時(shí)無法產(chǎn)生一個(gè)包含倒轉(zhuǎn)和上涌過程的球諧一階結(jié)構(gòu)，最接近長(zhǎng)波上涌的算例為D150v1e-1V1e21K0，其上升熱柱的形態(tài)經(jīng)過球諧函數(shù)分解主要為二階(表3,圖12).

圖12 最長(zhǎng)波長(zhǎng)的算例的三維形態(tài)結(jié)構(gòu)圖以及二維溫度剖面三維圖中綠色透明等值面代表化學(xué)組分為0.4的面.其余均與上同.Fig.12 The 3-D thermochemical structures and 2-D cross sections for the case of the largest upwelling wavelengthIn the 3-D view,transparent green iso-surfaces represent the chemical composition at a contour level of 0.4.Other schemes are the same as those above.

重力失穩(wěn)的解析解也可以通過瑞利-泰勒不穩(wěn)定分析得到(圖13a)，其結(jié)果表明一階的上涌要么需要更小的底部鈦鐵礦層厚度(至多200 km)和鈦鐵礦層黏度(至多10-3的相對(duì)黏度)，要么需要黏度更大的鈦鐵礦層(至少101的相對(duì)黏度).因此，D150v1e-1V1e21K0落在了一階上涌的參數(shù)空間外(圖13a,b).如果要將鈦鐵礦層黏度減小到10-3，那么根據(jù)之前的結(jié)果，很難形成上涌熱柱.而如果要提高鈦鐵礦層黏度，我們需要將其他因素納入流變學(xué)性質(zhì)的考慮中，最新的巖石力學(xué)實(shí)驗(yàn)發(fā)現(xiàn)鈦鐵礦有較大的活化體積(Tokle et al.,2021)，這意味著隨著深度(壓強(qiáng))增加，鈦鐵礦的黏度會(huì)大大增加，希望這能夠幫助我們產(chǎn)生更加長(zhǎng)波的月幔熱柱上涌.

圖13 (a)通過瑞利-泰勒不穩(wěn)定分析得到的在不同鈦鐵礦層黏度和底部鈦鐵礦層厚度下發(fā)生上涌時(shí)的主導(dǎo)球諧階數(shù).(b)兩個(gè)最長(zhǎng)波長(zhǎng)算例在即將形成上升流時(shí)的黏度剖面Fig.13 (a)The predicted dominant spherical harmonic degree of overturn as a function of the IBC viscosity contrast and the thickness of bottom IBC.(b)Viscosity profiles in the beginning of upwelling for two cases of the largest upwelling wavelength

熔融對(duì)模型演化的影響主要體現(xiàn)在兩個(gè)方面，在力學(xué)效應(yīng)方面，鈦鐵礦層受熱發(fā)生部分熔融，產(chǎn)生熔體增加巖石格架的孔隙率，并會(huì)極大的降低鈦鐵礦層的黏度；在熱效應(yīng)方面，鈦鐵礦層熔融時(shí)，熔體以潛熱形式吸收熱量，從而影響模型的熱演化過程.巖石學(xué)實(shí)驗(yàn)表明對(duì)于部分熔融的橄欖巖，其應(yīng)變速率與eα φ成正比，其中φ是部分熔融程度，α是部分熔融系數(shù)，根據(jù)形變機(jī)制不同，取值在21～45之間(Hirth and Kohlstedt,2004;Mei et al.,2002;Scott and Kohlstedt,2006)，據(jù)此約10%的部分熔融可以導(dǎo)致黏度減小一個(gè)數(shù)量級(jí).但是以上估計(jì)是基于熔體不發(fā)生運(yùn)移的情況，實(shí)際情況下，固體熔融產(chǎn)生熔體后，固體格架會(huì)發(fā)生改變，同時(shí)熔體與固體存在密度差，因此產(chǎn)生的液相熔體會(huì)立刻發(fā)生運(yùn)移和分隔，這使得，想要定量研究熔融對(duì)黏度的影響，需要利用考慮液相熔體的雙相流模型(例如Zhang et al.,2020).前人主要考慮熔融對(duì)熱演化的影響(Laneuville et al.,2018;Li et al.,2019)，我們也嘗試了一個(gè)只考慮熔融對(duì)熱演化影響的算例，為D150v1e-1V1e20K0M2(表3)，參照Li等(2019)，只改變控制方程中的能量方程，其余不變：

(10)

(11)

圖14 熔融(D150v1e-1V1e20K0M2)對(duì)進(jìn)入下月幔的鈦鐵礦層比例的影響的演化圖Fig.14 Evolutions of overturned IBC fraction in the lower lunar mantle with respect to melting (D150v1e-1V1e20K0M2)

目前現(xiàn)存的月球巖石樣品、隕石樣品以及隕石坑數(shù)量定年方法均認(rèn)為月球火山作用是一個(gè)長(zhǎng)期過程，歷時(shí)超過10億年(Nyquist et al.,2001;Hiesinger et al.,2011;Joy and Arai,2013).阿波羅任務(wù)獲得的月球玄武巖樣品年齡主要分布于3.1～3.9 Ga，而玄武質(zhì)的月球隕石樣品年齡則擴(kuò)展至4.37 Ga(KAL 009)到2.7 Ga(NWA 032,773)(圖1)(Nyquist et al.,2001;Fagan et al.,2002;Terada et al.,2007;Joy and Arai,2013;Snape et al.,2019)，中國(guó)嫦娥5號(hào)任務(wù)的著陸區(qū)的玄武巖樣品年齡為1.96 Ga(Che et al.,2021).然而，通過隕石坑數(shù)量定年得到的玄武巖年齡會(huì)更為年輕，可達(dá)1.2 Ga(Hiesinger et al.,2011,Qian et al.,2021).月球火山作用的時(shí)間分布特征需要一個(gè)能夠長(zhǎng)期驅(qū)動(dòng)月幔熔融的有效機(jī)制.我們的模型可以得到上涌熱柱的幕式行為(圖15).從3.9 Ga到3.3 Ga，450 km深度到初始鈦鐵礦層底部之間的鈦鐵礦層比例三次達(dá)到了極值，表示在這期間有三次主要的上涌，從而可能形成三次月幔熔融和玄武巖噴發(fā)事件(Qin and Zhong,2014).但我們的模型產(chǎn)生的最后一次上涌為3.3 Ga，這個(gè)時(shí)間比觀測(cè)到最新的月海玄武巖(1.96 Ga)還有很大的時(shí)間差別，所以可能需要更多的熱源來驅(qū)動(dòng)月球內(nèi)部的熱演化.此前的熱演化研究中，月核與月幔差異運(yùn)動(dòng)導(dǎo)致的耗散功率常常被忽略(Williams et al.,2001;Dwyer et al.,2011)，Dwyer等(2011)估計(jì)了月球早期核幔差異運(yùn)動(dòng)的耗散功率，得到在4.0 Ga左右其耗散功率的量級(jí)為1×1012W.但是他們?cè)谟懻撝姓J(rèn)為可能高估了耗散功率，因此又定性地討論了核幔耦合作用下，月球演化早期差異運(yùn)動(dòng)摩擦產(chǎn)生的功率為3×1011W.通過對(duì)本文中的四套生熱方案的單位生熱率(圖2b)可以估計(jì)得到四套方案的總生熱功率分別為18×1011W，14×1011W，12×1011W，7×1011W.所以核幔摩擦這一熱源對(duì)月球內(nèi)部更長(zhǎng)期的熱演化有一定作用，我們認(rèn)為他們的估計(jì)還不夠精確，我們將在以后的工作中精確評(píng)估核幔差異運(yùn)動(dòng)的效果，來對(duì)比和生熱元素貢獻(xiàn)的大小.

圖15 2000 Ma間在450 km深和初始鈦鐵礦層層底部之間的鈦鐵礦層比例隨時(shí)間的演化圖橙色陰影區(qū)標(biāo)記了主要的上涌階段.Fig.15 Evolutions of upwelling IBC fraction between the 450-km depth and the bottom of initial IBC layer for 2000 MaOrange shadow regions denotes major upwelling periods with peaks of IBC fraction.

通過改變生熱元素豐度和參考黏度，模型給出的上涌時(shí)間可以用來限定這兩者的參數(shù)空間.我們的模型用到的參考黏度范圍為1×1020Pa·s到1×1021Pa·s，對(duì)地球上地幔的實(shí)驗(yàn)結(jié)果給出黏度的范圍為1×1018Pa·s到1×1021Pa·s(Hirth and Kohlstedt,2004)，而前人模型中考慮的參考黏度范圍為5×1019Pa·s到1×1022Pa·s(Zhang et al.,2013a;Laneuville et al.,2018;Qin et al.,2018;Li et al.,2019;Zhao et al.,2019).我們的結(jié)果表明當(dāng)參考黏度為1×1021Pa·s時(shí)，只有位于粗虛線上方的參數(shù)空間才能發(fā)生上涌，進(jìn)一步來說，如果要月球主要玄武巖來自鈦鐵礦層倒轉(zhuǎn)又加熱上升的機(jī)制，而且集中在3.9 Ga大規(guī)模噴發(fā)，需要極厚的初始鈦鐵礦厚度(>100 km)和較大的生熱元素豐度，至少K/U為8500，逼近地球生熱元素豐度(圖5c).而如果參考黏度達(dá)到1×1022Pa·s(Qin et al.,2018)，對(duì)我們模型改變參考黏度得到上涌時(shí)間進(jìn)行外推，可知該參考黏度下無法在3.9 Ga產(chǎn)生大規(guī)模玄武巖巖漿，Qin等(2018)也認(rèn)為他們基于月球赤道突起而約束得到的月幔黏度1×1022Pa·s實(shí)際上對(duì)黏度范圍1×1021～1×1022Pa·s沒有分辨率，因此我們認(rèn)為月幔參考黏度最大應(yīng)為1×1021Pa·s.總而言之，我們的模型表明為了鈦鐵礦和KREEP能重新噴發(fā)，月幔的參考黏度上限應(yīng)在1×1021Pa·s這個(gè)量級(jí)，此時(shí)月球生熱元素的豐度至少為Taylor′s Moon的程度，接近McDonough和Sun(1995)的地球豐度；當(dāng)月幔參考黏度更小時(shí)，允許月球的生熱元素豐度小于McDonough和Sun(1995)的地球豐度.

4 結(jié)論

總而言之，我們使用了三維的有限元模型來研究月幔倒轉(zhuǎn)和倒轉(zhuǎn)之后的熱化學(xué)演化得到了如下結(jié)論.我們發(fā)現(xiàn)內(nèi)部生熱對(duì)玄武巖產(chǎn)生的影響較大，對(duì)于相同的鈦鐵礦層厚度、黏度以及月幔參考黏度，內(nèi)部生熱越強(qiáng)，熱柱上涌產(chǎn)生玄武巖巖漿所用的時(shí)間越短，反之若內(nèi)部生熱越弱，則所用時(shí)間越長(zhǎng).

我們通過月球在3.9 Ga集中噴發(fā)月海玄武巖的觀測(cè)結(jié)論約束得到了月幔參考黏度的上限為1×1021Pa·s，此時(shí)對(duì)于較厚的初始鈦鐵礦層，月球生熱元素的豐度至少應(yīng)該達(dá)到Taylor′s Moon的程度，如果月球初始鈦鐵礦層較薄 (<100 km)，則需要更多的生熱元素豐度，即McDonough和Sun(1995)的地球豐度或更高.

我們考慮了熔融對(duì)月球內(nèi)部熱演化的影響，結(jié)果表明，如果鈦鐵礦層發(fā)生熔融，那么將要求月球內(nèi)部有更高的生熱元素豐度或者其他貢獻(xiàn)較大的熱源，才能滿足3.9 Ga大規(guī)模玄武巖噴發(fā)的觀測(cè)約束.

本文還得到了月海玄武巖周期性的噴發(fā)模式以及玄武巖的不對(duì)稱分布，但是與實(shí)際月海玄武巖的時(shí)間和空間分布還有一定差距，未來的工作將繼續(xù)探索參數(shù)空間來接近實(shí)際情況.

致謝本研究的資金來源于CNSA D020205.感謝2位匿名評(píng)審人的詳細(xì)意見.其中建議的文獻(xiàn)提高了本文的視角.本文的大型并行計(jì)算在北京大學(xué)高性能計(jì)算平臺(tái)上完成.本文通訊作者感謝本科時(shí)期王仁老師的言傳身教，他在科研工作上的耐心和仔細(xì)給學(xué)生時(shí)期的作者留下深刻的印象.