2021年Fagradalsfjall火山喷发前的变形和地震活动性下降gydF4y2Ba

世界各地的火山观测站旨在向平民、航空当局和其他利益攸关方及时提供火山爆发预警，以防止生命损失和基础设施破坏。要做到这一点，正确理解火山喷发前兆的模式是至关重要的gydF4y2Ba^{4gydF4y2Ba，gydF4y2Ba5gydF4y2Ba，gydF4y2Ba6gydF4y2Ba}在火山爆发前，火山喷发的速度通常会逐渐加快gydF4y2Ba^{2gydF4y2Ba，gydF4y2Ba3.gydF4y2Ba}．对这种行为进行量化的开创性工作是由博伊特完成的gydF4y2Ba^1gydF4y2Ba他提出了一种解释火山中观测到的前兆活动(如应变或地震活动)的上升速率的方法，以找到失效和喷发开始的时间。这种材料失效预测方法，或其改进版本，已被广泛用于事后预测火山爆发的时间，在某些情况下接近实时gydF4y2Ba^{7gydF4y2Ba，gydF4y2Ba8gydF4y2Ba，gydF4y2Ba9gydF4y2Ba}．然而，有些火山爆发表现出不同的行为，甚至在喷发前的前兆活动减少。2010年3月，冰岛Eyjafjallajökull火山的侧翼喷发发生在山顶爆发之前，经过几个月的活跃活动后，地震活动和变形率有所下降gydF4y2Ba^10gydF4y2Ba地震开始前没有立即发出短期警告gydF4y2Ba^11gydF4y2Ba．在一些层状火山潜水爆炸或喷发开始之前，地震活动也有所减少，例如，在1989-1990年和2009年Redoubt火山活动之前gydF4y2Ba^{12gydF4y2Ba，gydF4y2Ba13gydF4y2Ba}在此之前，1999年在特利卡(Telica)。gydF4y2Ba^14gydF4y2Ba)．在某些情况下，这归因于气体运移通道的密封，从而导致系统的增压gydF4y2Ba^15gydF4y2Ba导致地表变形增加(隆起)gydF4y2Ba^16gydF4y2Ba．然而，我们在Fagradalsfjall观察到的前兆活动包括地震活动的下降和变形。gydF4y2Ba

雷克雅内斯半岛斜裂带是北美-欧亚伸展板块边界的一部分，它出现在冰岛的海平面以上(图2)。gydF4y2Ba1gydF4y2Ba和扩展数据图。gydF4y2Ba1gydF4y2Ba)．板材铺展19mm yrgydF4y2Ba^{−1gydF4y2Ba}在N104°E左右方向上，雷克雅内斯半岛的中轴线高度倾斜，雷克雅内斯半岛的中轴线位于N77°E左右，与开放区相比，剪切作用更大gydF4y2Ba^{17gydF4y2Ba，gydF4y2Ba18gydF4y2Ba}．沿着雷克雅内斯半岛东北至西南排列的梯队排列的裂缝群形成了将该地区划分为火山系统的基础(相当于在海洋中脊上扩展的中心)，每个火山系统都包括一个裂缝群和一个高温地热区gydF4y2Ba^19gydF4y2Ba．在非爆发期，沿板块边界释放剪应力的南北向走滑断层阵列主要分布在裂隙群之间gydF4y2Ba^20.gydF4y2Ba．在2021年3月19日在Geldingadalir的Fagradalsfjall火山爆发之前，雷克雅内斯半岛已经有大约800年没有爆发过了。过去3000年的火山活动以几百年的喷发期为特征，间隔800 - 1000年的无喷发期gydF4y2Ba^{19gydF4y2Ba，gydF4y2Ba21gydF4y2Ba}．20世纪90年代以来的大地测量结果表明，沿板块边界的应变和应力积累与板块运动模型一致，其锁定深度约为5-8 km，在此深度以下，韧性变形占主导地位，构造板块受地震和地热变形扰动而自由滑动gydF4y2Ba^{22gydF4y2Ba，gydF4y2Ba23gydF4y2Ba，gydF4y2Ba24gydF4y2Ba}(扩展数据图。gydF4y2Ba1gydF4y2Ba)．2019年12月中旬，雷克雅内斯半岛开始发生动荡，Fagradalsfjall山南缘发生了持续一周的强烈地震群，随后雷克雅内斯半岛大部分地区的地震活动加剧。与反复的岩浆侵入、地震触发和构造活动有关的高强度群期发生在不同的地区，沿板块边界的三个地区发现了几次膨胀gydF4y2Ba^{25gydF4y2Ba，gydF4y2Ba26gydF4y2Ba}．在这里，我们描述了从2021年2月24日起的活动，当时发生了动荡事件中最大的地震(震级)gydF4y2Ba米gydF4y2Ba_WgydF4y2Ba5.64)，直到2021年3月19日爆发。gydF4y2Ba

的gydF4y2Ba米gydF4y2Ba_WgydF4y2Ba2021年2月24日5.64地震之前，在即将到来的堤坝位置发生了3小时的强烈集中微地震活动，每分钟在近7公里深的狭窄深度范围内发生几次，这表明岩浆侵入开始。在沿着未来堤坝的大致方向向东北传播了几百米之后，地震群达到了顶点gydF4y2Ba米gydF4y2Ba4走滑事件，在蜂群区域内触发。在17秒内gydF4y2Ba米gydF4y2Ba_WgydF4y2Ba5.64级地震发生在东南方1公里处，靠近板块边界的中轴线(图2)。gydF4y2Ba1gydF4y2Ba)．在接下来的4小时里，所有的活动都局限在板块边界的中轴线上，地震两侧10公里长的板块被激活了60多次，其中大部分是走滑，gydF4y2Ba米gydF4y2Ba3次地震。在前两个半小时内，八个gydF4y2Ba米gydF4y2Ba_WgydF4y2Ba> 4.0事件沿这些段发生，表明板块边界约20公里处发生滑动。gydF4y2Ba1gydF4y2Ba和扩展数据图。gydF4y2Ba2gydF4y2Ba)．该地区的变形模式在2月24日之后发生了变化，表明岩浆持续流入脆性地壳中的垂直岩脉，并在板块边界的中轴线上发生了一些滑移。通过全球导航和卫星系统(GNSS)的大地测量和合成孔径雷达图像(InSAR)的干涉分析，可以很好地记录堤防就位的地面变形。在大堤东南方的KRIV GNSS站显示了一个同震跳gydF4y2Ba米gydF4y2Ba_WgydF4y2Ba在接下来的几天里，地震在N148°E附近以大约每天10毫米的速度发生位移，并随着时间的推移逐渐减小(图2)。gydF4y2Ba1 bgydF4y2Ba)．在3月3日之前，大堤西北方的LISK GNSS站在N305°E方向上以每天约4毫米的速度移动，在类似的方向上加速到每天约14毫米，但此后随着时间的推移逐渐减少。由Sentinel-1数据形成的为期6天的干涉图，揭示了2021年2月和3月从卫星到地面的视线(LOS)位移如何随时间变化。从InSAR中还可以明显看出，变形速率随着时间的推移而下降，最慢的变形速率发生在喷发前的最后几天(图2)。gydF4y2Ba2gydF4y2Ba和扩展数据图。gydF4y2Ba2gydF4y2Ba)．因此，GNSS和InSAR的数据都揭示了地震后的高变形率gydF4y2Ba米gydF4y2Ba_WgydF4y2Ba5.64级地震，然后在喷发开始时逐渐减少到几乎为零。gydF4y2Ba

后,gydF4y2Ba米gydF4y2Ba_WgydF4y2Ba2021年2月24日5.64级地震，9次以上地震gydF4y2Ba米gydF4y2Ba_WgydF4y2Ba4.8级地震发生在雷克雅内斯半岛，其中最大的一次是a级gydF4y2Ba米gydF4y2Ba_WgydF4y2Ba2021年3月14日5.33事件(gydF4y2Ba方法gydF4y2Ba及补充表gydF4y2Ba1gydF4y2Ba)．地震震级和地震矩释放率在3月19日晚爆发前一直在下降(图19日晚)。gydF4y2Ba1 bgydF4y2Ba)，上面只有一项活动gydF4y2Ba米gydF4y2Ba_wgydF4y2Ba4.0 (gydF4y2Ba米gydF4y2Ba_wgydF4y2Ba3月15日4.20)。喷发前两天，地震活动主要集中在两个小地震群，其中一个在即将爆发的地点。从2月24日至3月19日爆发，自动探测到53,000多次地震(gydF4y2Ba方法gydF4y2Ba)．这种地震活动的一个特点是地震的数量异常多(64)gydF4y2Ba米gydF4y2Ba_WgydF4y2Ba≥4.0级地震，比该地区过去20年的记录多了近20次。史无前例的地震活动，其特征是重叠的主震-余震序列，在冰岛首都地区连续三周昼夜都能感觉到，因为沿半岛引发了南北方向的断层，库仑破坏应力增加了高达4兆帕，以响应沿堤坝的开放(扩展数据图)。gydF4y2Ba1gydF4y2Ba)．通过波形互相关获得的堤防地震活动的相对重定位(gydF4y2Ba方法gydF4y2Ba)，显示出北段向N45°E，南段向N25°E的垂直两段岩脉(扩展数据图)。gydF4y2Ba3.gydF4y2Ba)．2月24日至3月3日，大堤北段发生地震活动。此后，地震活动向堤坝南段和向西移动，在接下来的4天里激活了大约15公里的板块边界段。3月7日，这一趋势停止，活动从堤坝中心沿南堤段向西南方向迁移。地震活动速率受岩脉和板块边界活动总体演化的影响，在火山爆发前几天整体呈下降趋势。与堤区相比，板块边界区开始下降的时间提前了3天。gydF4y2Ba1 cgydF4y2Ba)．在火山爆发前15小时，地震信号在较浅震源的特征中高频衰减，在约20个浅层(<1 km深度)记录到，小于gydF4y2Ba米gydF4y2Ba即将爆发的地点将发生2次地震。gydF4y2Ba

地表变形的来源是使用改良版的GBIS大地贝叶斯反演软件进行建模推断的gydF4y2Ba^27gydF4y2Ba．同震变形被模拟为剪切运动在矩形平面上的结果，以及在平面上有开口和剪切的堤坝，嵌入均匀弹性半空间(gydF4y2Ba方法gydF4y2Ba)．观测到的地表变形也类似于地震活动，揭示了在研究期间沿板块边界中轴线的活动。类似的，但更小的变形，也出现在动荡之前的几年(扩展数据图。gydF4y2Ba4gydF4y2Ba)，可以用沿板块边界趋势的亚垂直平面上的剪切模型来解释，我们的模型中也包括了这个模型。沿这一平面的剪切作用和岩脉平面的剪切作用结合在一起，再现了该地区分布剪切作用的影响，这种影响可能发生在该地区映射的一系列紧密分布的南北向走滑断层上gydF4y2Ba^{20.gydF4y2Ba，gydF4y2Ba21gydF4y2Ba}．由于整个事件过程中压力点源逐渐收缩，这里还模拟了靠近堤坝中心的局部下沉区。模拟结果表明，大部分的变形与垂直堤坝的侵位有关，该地区发生了大部分较小规模的地震。初始模型考虑具有均匀开剪的单段堤坝(扩展数据图。gydF4y2Ba5gydF4y2Ba)．在此基础上，考虑了符合重定位地震活动性、开口和滑移分布(gydF4y2Ba方法gydF4y2Ba和扩展数据图。gydF4y2Ba6gydF4y2Ba)．北段走向N45°E，南段走向N23.5°E。gydF4y2Ba3gydF4y2Ba)．日体积变化(岩浆流入速率)求解为(gydF4y2Ba方法gydF4y2Ba)，计算出岩浆侵位的平均深度(图;gydF4y2Ba3 b, cgydF4y2Ba)．这些数据总体上受到GNSS和InSAR联合数据的约束，其中GNSS数据约束每日变化。推断岩浆流入岩脉的速率(体积变化)是波动的。最高值出现在堤坝形成初期的2月24日至3月3日(30-35米)gydF4y2Ba^3.gydF4y2Ba年代gydF4y2Ba^{−1gydF4y2Ba})，取值范围为10-20米gydF4y2Ba^3.gydF4y2Ba年代gydF4y2Ba^{−1gydF4y2Ba}在3月3日至15日期间，最小值(<10 mgydF4y2Ba^3.gydF4y2Ba年代gydF4y2Ba^{−1gydF4y2Ba})发生在火山爆发前几天(3月15日至19日)。岩浆侵位的日平均深度也有波动，但在3月11日以后逐渐变浅。岩浆流速与平均侵位深度呈线性关系，岩浆侵位越浅，岩浆流速越低(图2)。gydF4y2Ba3 cgydF4y2Ba)．这与由沿堤防下方的管道向上流动而无需沿堤防侵入的显著压力所控制的流量是一致的(gydF4y2Ba方法gydF4y2Ba)．相关预测初步喷发速率约为7米gydF4y2Ba^3.gydF4y2Ba年代gydF4y2Ba^{−1gydF4y2Ba}，与观察结果一致gydF4y2Ba^28gydF4y2Ba．根据地球化学观测结果，它还可以估算出供给大堤的水源深度为19公里gydF4y2Ba^29gydF4y2Ba．据推测，连接岩浆源和堤坝底部的管道的横截面积约为几平方米。gydF4y2Ba

模拟变形源的综合效应(图;gydF4y2Ba4gydF4y2Ba)是释放板块运动积累的应力。在该模型中，估计3.5 km深度处的扩张应力变化为数十兆帕(扩展数据图。gydF4y2Ba7gydF4y2Ba)，而以前由于板块运动引起的每年的扩张应力变化要小3个数量级或每年几十千帕(扩展数据图。gydF4y2Ba1gydF4y2Ba)．因此，这些事件中的应力变化可能与约800年前雷克雅内斯半岛上的火山爆发以来长期积累的大量构造应力相对应。岩脉侵入整个时期北向东位移的导数，可以直接与InSAR观测结果进行比较(图2)。gydF4y2Ba4 b, cgydF4y2Ba)，很好地说明了板块边界中轴线对剪切应力的释放也起着重要作用。在这种斜向扩张的背景下，这种剪切作用与堤坝开口共同释放了构造应力(参见扩展数据图)。gydF4y2Ba7gydF4y2Ba)．由于累积地震滑动造成的地表位移的空间格局，不包括两个最大的地震事件，与由于其他震源的综合影响造成的地表位移的空间格局类似(扩展数据图)。gydF4y2Ba9gydF4y2Ba)，表明较小的地震也有助于压力的释放。gydF4y2Ba

我们的观测表明，构造应力的释放伴随着变形和地震活动率的下降可能是某类火山爆发预期的典型前兆活动。我们的观测结果与以下几个演化阶段是一致的:(1)岩浆进入岩脉的过程伴随着特征变形和地震活动，以及周边地区的变形和地震活动的触发。在活动前积累的构造应力在高应变和应力释放期得到缓解。(2)构造应力释放的主要阶段是变形率和地震活动速率的下降，两者都与岩浆流入速率的下降有关。随着应力在脆性地壳底部逐渐释放，岩浆在侵入岩脉之前被迫向更高的地方流动。这导致驱动压力下降更大，从而导致岩浆流动速率下降。(3)一旦积累的构造应力沿岩脉释放，岩浆以相对平静的方式冲出地表，在没有明显地震能量释放或变形的情况下启动喷发。即使堤坝需要突破最上面的地壳来引发喷发，火山裂谷中最上面的1公里地壳是断裂的，而且很弱，所以信号不容易被探测到。没有明显的前兆活动就可能发生喷发。这是一种非常不同的前兆活动模式，与物质失效预测方法所预测的喷发前活动升级的情况截然不同。 Although this method predicts escalating rates of precursory activity, it is known that the inflation rate and the rate of magma inflow into a shallow magma body may slow down before eruptions as was the case for Krafla volcano, North Iceland^{30.gydF4y2Ba，gydF4y2Ba31gydF4y2Ba}．由于岩浆流速率与通道内的压力梯度成正比，因此，在连接深部恒压源和压力上升的浅层岩浆体的岩浆管道中，岩浆流速率将呈指数级下降gydF4y2Ba^32gydF4y2Ba．如果岩浆密度小于周围的寄主岩石，那么岩浆浮力在推动岩浆在这样的通道中向上起着重要作用gydF4y2Ba^32gydF4y2Ba．gydF4y2Ba

我们的发现与被动堤坝侵入和裂谷的观察有相似之处，例如，在夏威夷基拉韦厄火山，由于基拉韦厄南翼滑动，在先前的高应力和应变积累的地区gydF4y2Ba^{33gydF4y2Ba，gydF4y2Ba34gydF4y2Ba}．在1999年9月侵入之前，没有观测到先兆的膨胀，地下岩浆流速也没有增加，表明岩浆管道系统内没有压力增加的先兆。在东非，是被动裂谷作用的前兆gydF4y2Ba^{35gydF4y2Ba，gydF4y2Ba36gydF4y2Ba，gydF4y2Ba37gydF4y2Ba}受高张应力的影响。我们观察到的前兆活动也类似于堪察加Tolbachick火山1975-1976年裂缝爆发前的观测，在那里，地震群在新裂缝开放前1-2天急剧减少gydF4y2Ba^38gydF4y2Ba．我们观察到的前兆行为可能与一般火山裂谷中的喷发有关，这是在从深层岩浆体进入预应力地壳的不断演变的堤坝侵入形成之后。在这样的构造环境下，地壳的非应力在喷发前是可以预料到的。例如，这可能会对夏威夷莫纳罗亚火山的下一次喷发产生重大影响，那里一直存在火山动荡gydF4y2Ba^39gydF4y2Ba．如果前兆活动将与Fagradalsfjall相似，那么在下一次喷发之前可能不会发生不断升级的变形和地震活动。监测和了解火山爆发的前兆对于预测潜在的相关危险和减少对人类和基础设施的潜在影响至关重要gydF4y2Ba^40gydF4y2Ba．因此，在通常由火山观测站获取的地球物理数据集中确定不太常见的喷发前趋势是支持民事当局及时决策和采取风险缓解措施的一个重要步骤。gydF4y2Ba

常规地震位置和地震活动性gydF4y2Ba

该地区的地震活动由冰岛气象局(IMO)监测，地震记录在SIL国家地震台网上。地震的P波和S波到达时间由SIL分析系统自动检测，该系统定位事件并分配初步震级。然后根据需要手动检查、修改和更新事件。这些常规分析使用一个通用的一维速度模型(SIL模型)来分析冰岛的所有事件gydF4y2Ba^{42gydF4y2Ba，gydF4y2Ba43gydF4y2Ba}．在2021年2月24日至3月14日期间，SIL系统自动检测到的地震总数约为45,000次，其中15%为人工检测。在随后的3月15日至19日期间，大约7500次自动定位地震中的19%已被审查。为了改善对雷克雅内斯半岛动荡的地震监测，由捷克科学院(CAS)拥有、冰岛地质调查局(ÍSOR)运营的半岛上另外8个监测站与国际海事组织的监测系统相连，并用于自动和审查的标准位置。其中三个站(lsf、iss和lag)于2020年加入系统，另外四个站(lat、ash、moh、faf)于2021年3月初加入系统(前三个站于3月5日，第四个站于3月11日)。3月5日，剑桥大学的一个站(odf)也被包括在常规的自动和手动分析中。雷克雅内斯半岛最近增加的台站密度提高了完整性的量级(gydF4y2Ba米gydF4y2BaC)该地区的地震活动性，从2002年至2013年之前0.5或更好的估计值(参考文献)。gydF4y2Ba^44gydF4y2Ba)．然而，在强震群和大地震期间，自动定位系统饱和，背景噪声大大增加，导致灵敏度显著降低。在这里分析的激烈时期，gydF4y2Ba米gydF4y2Ba雷克雅内斯半岛的c估计在gydF4y2Ba米gydF4y2Ba1、因此地震活动性分析仅限于gydF4y2Ba米gydF4y2Ba1次地震。地震活动性(图;gydF4y2Ba1 cgydF4y2Ba)分别估算了堤坝和板块边界区域的地震强度，方法是将堤坝周围2公里宽的条带和板块边界轴周围2.5公里宽的条带的所有地震强度相加(重叠区域的地震被划分在堤坝和板块边界区域)。数字gydF4y2Ba1 cgydF4y2Ba还报告了其他地区的地震数量。由于人工检查了地震活动性的一小部分，地震活动性分析是基于自动定位的。gydF4y2Ba

地震震源机制及震级确定gydF4y2Ba

在人工检查地震位置后，计算出震源机制和震级，并将事件添加到SIL目录中。计算是基于对所有可能的走向、倾角和前倾角进行网格搜索，其中第一运动P波极性匹配，P波和s波振幅在观测到的可接受范围内gydF4y2Ba^45gydF4y2Ba．SIL分析系统中的震级计算往往低估事件的震级gydF4y2Ba米gydF4y2Ba_WgydF4y2Ba> 3.5(参考。gydF4y2Ba^46gydF4y2Ba)．通过单独的实时分析自动计算出这些事件更好的震级估计，该分析连续评估每个地震台站在不同频段观测到的峰值地面速度(PGV)gydF4y2Ba^47gydF4y2Ba．使用这些值和从SIL系统自动定位事件gydF4y2Ba米gydF4y2Ba> 3.5矩的大小，gydF4y2Ba米gydF4y2Ba_WgydF4y2Ba，是利用地面运动速度预测方程为每个站计算的gydF4y2Ba^48gydF4y2Ba．最后一个gydF4y2Ba米gydF4y2Ba_WgydF4y2Ba每个事件的估计值是中位数附近计算值的三分之一的稳健平均值，从裁剪范围外(约15公里)到距离事件230公里的站点gydF4y2Ba^49gydF4y2Ba．这个值是稳健的，并与国际上确定的矩幅值相比较。补充表gydF4y2Ba1gydF4y2Ba列出所有确定的震级gydF4y2Ba米gydF4y2Ba_WgydF4y2Ba2021年2月至4月雷克雅内斯半岛发生≥4.0次事件。这些震级用于计算动荡时期最大事件释放的累积力矩。gydF4y2Ba

地震相对移位gydF4y2Ba

2021年2月24日至3月20日期间的地震活动位于岩脉侵入体和由侵入体引发的许多大约南北向走滑断层上及其周围。数字gydF4y2Ba1gydF4y2Ba及补充表gydF4y2Ba1gydF4y2Ba的15个走滑事件gydF4y2Ba米gydF4y2Ba≥4.5级地震在堤岸附近触发，导致其断平面上多次余震，震级范围在附近gydF4y2Ba米gydF4y2Ba1.5.因此，为了绘制出岩脉侵入的整体特征和尺度，并尽量减少大约南北向走滑断层余震的污染，选择了较小震级事件的数据集进行高精度相对重定位。该数据集包括迄今为止常规定位和审查的所有1321次地震，震级范围在0.5级以下。gydF4y2Ba米gydF4y2Ba< 1.4。gydF4y2Ba

重新定位方法使用了“选择”的绝对到达时间和改进的相对到达时间测量的双重差异，这些测量是通过不同地震的P波和S波的相互关联获得的gydF4y2Ba^50gydF4y2Ba．该方法迭代反演P波和S波绝对到达时间差的加权平方和，以及(1)P波和S波的绝对到达时间、(2)P波和S波的相对到达时间和(3)S - P波的相对到达时间的双重差异。由此产生的事件位置分布具有较高的内部/相对定位精度，因此非常适合绘制活动地下断层和岩脉侵入的详细信息gydF4y2Ba^{51gydF4y2Ba，gydF4y2Ba52gydF4y2Ba}．雷克雅内斯半岛的速度结构与常规分析中使用的标准sil -速度模型有很大不同，这尤其影响了震源深度，常规分析中震源深度往往大约深了1公里。因此，本文利用体波沿雷克雅内半岛传播的走时和相速度分析，建立了拟合较好的一维速度模型gydF4y2Ba^53gydF4y2Ba雷克扬斯-冰岛地震实验地震折射剖面gydF4y2Ba^54gydF4y2Ba．搬迁使用了距离地震50公里范围内的15个地震台站，包括所有IMO台站和2020年连接到系统的3个CAS/ÍSOR台站。排除最靠近堤坝的站点是有意的，以满足重新定位方法正常工作的要求，因为它假设来自距离较近的事件的波通过相同的速度结构传播。由此得到的地震震源深度参考了最近台站的平均海拔高度，大约在海平面以上80米。gydF4y2Ba

由此产生的重定位事件分布(图。gydF4y2Ba3.gydF4y2Ba和扩展数据图。gydF4y2Ba3.gydF4y2Ba)显示了长约9公里的堤防侵入的两个主要段，其走向和深度范围不同。北段至少4.5公里长，更偏东，倾角约45°，深达8公里。南段4.5公里，倾角25°，深达6公里。在两段的交点处深度变化相当剧烈(扩展数据图B)。gydF4y2Ba3.gydF4y2Ba)．扩展数据图。gydF4y2Ba3.gydF4y2Ba还显示了爆发区0.5公里深度附近的事件集群(图中星形和白线)。扩展数据图中事件的颜色编码。gydF4y2Ba3.gydF4y2Ba揭示了活动在此期间的整体时间演变。地震活动开始于靠近堤坝中心的深度(扩展数据图B)。gydF4y2Ba3.gydF4y2Ba)并在2月24日至3月3日期间激活北段，然后在3月4日至3月19日左右向南传播激活南段。交叉口周围区域(扩展数据图B)。gydF4y2Ba3.gydF4y2Ba)在整个时间段内保持活跃。gydF4y2Ba

大地造型gydF4y2Ba

我们使用贝叶斯方法和改进版的GBIS软件对大地测量数据进行反演gydF4y2Ba^27gydF4y2Ba得到能解释观测到的形变的源参数的范围。利用覆盖2021年2月23日至3月19日整个喷发前时期的18个站的GNSS观测(运动和连续测量的结合)和两张Sentinel-1干涉图(从下降轨道T155覆盖2021年2月23日至3月19日和上升轨道T16覆盖2021年2月19日至3月21日)进行了初步反演。这个贝叶斯反演软件利用马可夫链蒙特卡罗方法和大都会黑斯廷斯算法gydF4y2Ba^{55gydF4y2Ba，gydF4y2Ba56gydF4y2Ba}．通过运行200万次迭代，得到了各种源参数的联合概率分布函数。在我们最初的反演中，我们用矩形位错来模拟岩脉gydF4y2Ba^57gydF4y2Ba，允许开放和剪切，以及在板块边界中轴线上观察到的带有第二次位错的运动，也允许开放和剪切。我们将两次最大的构造地震导致的地表位移建模为单独的矩形位错，只考虑剪切运动。靠近喷发地点的下沉信号被建模为放气点源gydF4y2Ba^58gydF4y2Ba简单地说，尽管这种收缩可能是与岩脉侵入有关的正常断层活动的结果。我们假设泊松比为0.27。所有的源参数都允许变化，并预先定义了宽的边界。gydF4y2Ba

该反演的输出参数证实了岩脉沿地震活动性确定的垂直平面侵入。在第二次反演中，我们用单个岩脉位错代替了两个相连的位错，这些位错的走向固定在地震活动转移所指示的走向上(北段为N45°E，南段为N23.5°E)。我们允许位置发生变化，认识到地震活动转移的绝对位置并不准确，并发现首选位置在地震活动转移以西400米左右，并穿过喷发裂缝。gydF4y2Ba

对于我们的最终反演，我们将堤坝的位置固定在第二次反演的最佳位置，其最大长度(每段4.5公里)和深度(北部7.5公里，南部6公里)来自地震活动。我们还从第二次反演中确定了板块边界位错的位置、倾角(垂直)和走向(N65.4°E)。在第二次反演中，这一段的开口很小(9厘米)，因此我们只允许剪切运动。考虑到堤坝和板块边界段的垂直性质，我们只考虑走滑意义上的走滑剪切，不考虑倾角运动。我们将堤和板块边界段划分为750 × 750 m的斑块gydF4y2Ba^2gydF4y2Ba在3公里以上和1.5 × 1.5公里gydF4y2Ba^2gydF4y2Ba在3公里以下深度，并在每个补丁上独立求解开口和/或滑移。我们允许放气点源的位置变化，我们的模型将其放置在1480 - 1670米的深度，体积变化范围为−2.6至−3.4 × 10gydF4y2Ba^6gydF4y2Ba米gydF4y2Ba^3.gydF4y2Ba．gydF4y2Ba

为了约束沿堤坝的日增量开口，我们使用了约束线性最小二乘方法。我们从最终反演的整个喷发前期的后验分布中选取了100个样本，以限制每个堤坝斑块的总开口，并仅使用GNSS求解增量开口。为了减少GNSS时间序列中的噪声，在计算增量位移之前，我们首先使用2天三角形移动窗口对其进行滤波。我们估计每个补丁的每日开口，限制只允许正面开口。为了考虑影响整个GNSS网络的系统误差，我们同时求解了东、北、向上的每日参考系调整。然后我们将每个patch的日体积变化率相加，得到日体积变化率，并计算出这个日体积变化率的平均深度(重心)。gydF4y2Ba

通过模拟试验研究了岩浆侵位深度与体积日变化之间的平衡关系。对模拟数据进行了增量式堤开度反演，以检验所观察到的日解深度和体积变化之间的权衡(图2)。gydF4y2Ba3.gydF4y2Ba)影响深度与体积变化的整体关系。我们将我们的分布式开放反演的最大后知概率解的总开放划分为28天的等体积变化，平均深度变浅。我们使用这些每日开堤模型来模拟GNSS站点的位移，并添加了具有代表性的噪声。然后，我们以与之前相同的方式对模拟的GNSS数据进行反转。结果如扩展数据图所示。gydF4y2Ba8gydF4y2Ba．尽管在个别日子的深度和体积之间进行了权衡，但模拟的恒定体积变化率被检索出来，而不依赖于深度。gydF4y2Ba

有记录的64次地震造成的地面变形gydF4y2Ba米gydF4y2Ba_wgydF4y2Ba≥4，除2月24日和3月14日两个最大事件外，采用以下程序进行评估。根据主导断层机制，每一次地震都假定发生在一条南北走向的垂直走滑矩形断层上。每个断层的中心位置和深度都是从国际海事组织震源目录中获得的。对于所有小于gydF4y2Ba米gydF4y2Ba5、长2 × 2 km的矩形断层gydF4y2Ba^2gydF4y2Ba被认为;为gydF4y2Ba米gydF4y2Ba≥5的长度为4公里，高度为2公里的地震。对于非常浅的震心(深度小于1.1公里)的测绘，在中心深度上增加了1公里，以避免人工断层从地面伸出来的人工制品。矩震级转换为地震矩gydF4y2Ba米gydF4y2Ba_0gydF4y2Ba使用关系gydF4y2Ba米gydF4y2Ba_0gydF4y2Ba= 10gydF4y2Ba^{1．5gydF4y2Ba米gydF4y2Ba+ 9.05gydF4y2Ba}(ref。gydF4y2Ba^59gydF4y2Ba)，剪切模量为gydF4y2BaμgydF4y2Ba= 30 GPa以获得平均成绩gydF4y2Ba年代gydF4y2Ba从地震时刻开始，gydF4y2Ba米gydF4y2Ba_0gydF4y2Ba=gydF4y2Baμ为gydF4y2Ba,在那里gydF4y2Ba一个gydF4y2Ba是断层区。所有62个断层(也就是说，不包括两次最大地震)的累积力矩相当于一次gydF4y2Ba米gydF4y2Ba_WgydF4y2Ba5.9级地震。库仑软件gydF4y2Ba^60gydF4y2Ba在密(0.1 × 0.1 km ?gydF4y2Ba^2gydF4y2Ba)网格。得到的位移场(扩展数据图。gydF4y2Ba9gydF4y2Ba)与利用沿堤、板块边界段和两条主要断裂的分布滑移反演的大地模型具有相同的特征，但量级要小得多。考虑到上述反演过程中较小断层的贡献并没有提高与大地测量数据的拟合，因此这里不考虑这些。较小断层的影响可以预期在主反转中已经得到解释。gydF4y2Ba

管道流动模型gydF4y2Ba

我们模拟了通过圆柱形导管流向堤坝底部的水流。在现实中，导管可能是不同的形状，在形成时更像裂缝。然而，它是在地壳锁定深度以下的韧性下地壳中形成的，在那里，板块构造应力不像在最上面的弹性地壳中那样以相同的方式积累。此外，在有限尺寸的裂缝中，岩浆流动可能会迅速集中到裂缝最宽的部分，有效流动路径可能会变得更圆柱形gydF4y2Ba^{32gydF4y2Ba，gydF4y2Ba61gydF4y2Ba，gydF4y2Ba62gydF4y2Ba}．在任何情况下，非圆形截面都不会显著影响我们的结论(见下文)。假定不可压缩流体在管道中层流，由岩浆源区和岩浆头的压力差驱动，并将堤坝视为对流动阻力可以忽略不计的储层，gydF4y2Ba

在哪里gydF4y2BaPgydF4y2Ba是岩浆源的压力，gydF4y2BaρgydF4y2Ba_米gydF4y2Ba是岩浆的密度，gydF4y2Ba\ (g \)gydF4y2Ba是重力加速度，gydF4y2BalgydF4y2Ba而且gydF4y2BargydF4y2Ba分别为导管的长度和半径，gydF4y2BahgydF4y2Ba是岩浆侵入的堤坝底部以上的高度，gydF4y2BaνgydF4y2Ba是动力粘度和gydF4y2BaugydF4y2Ba是导管内的平均岩浆流速。左边一项为驱动压力，右边第一项为层流(Hagen-Poiseuille)流动造成的粘性损失，右边第二项为动压力损失。利用体积流量与速度的关系，gydF4y2Ba\ (Q = {\ rm{\π}}{r} ^ {2} u \)gydF4y2Ba，gydF4y2Ba

在这种情况下gydF4y2Ba\(L\gg {\rho}_{\rm{m}}}{\rm{Q}}/16\pi \nu \)gydF4y2Ba．注意，在这种情况下gydF4y2BahgydF4y2Ba＝gydF4y2Ba0gydF4y2Ba而且gydF4y2Ba\ (P ={\ρ}_ {{\ rm {c}}} gh \)gydF4y2Ba,在那里gydF4y2BaρgydF4y2Ba_cgydF4y2Ba地壳以上的平均密度是不是岩浆源区，这就进一步减小到比较熟悉的关系了gydF4y2Ba\ (Q = {{\ rm{\π}}r} ^{4} \离开({\ρ}_ {{\ rm {c}}} -{\ρ}_ {{\ rm {m}}} \右)g / 8 \ν\)gydF4y2Ba(ref。gydF4y2Ba^63gydF4y2Ba)．对(2)求微分gydF4y2BahgydF4y2Ba给了gydF4y2Ba

在我们的概念模型中，我们假设使堤坝保持开放的应力是由上地壳中的偏应力提供的，并且堤坝实际上是一个由进入其底部的管道填充的容器。分gydF4y2Ba\ (P \)gydF4y2Ba进入堤坝底部上下地壳的重量所产生的压力，并考虑进入堤坝底部的流速，gydF4y2Ba问gydF4y2Ba_DgydF4y2Ba，由(2)gydF4y2Ba

在哪里gydF4y2BaρgydF4y2Ba_ugydF4y2Ba是上地壳(堤坝底部以上)的平均密度，gydF4y2BaDgydF4y2Ba到堤坝底部的深度是Δ吗gydF4y2BaρgydF4y2Ba是下地壳和岩浆密度的平均差值。gydF4y2Ba

结合(3)和(4)得到gydF4y2Ba

因此，如果我们能测量出堤坝底部的流入速度，以及流入速度在堤坝内随高度的变化，就可以确定管道的长度:gydF4y2Ba

假设的值gydF4y2BaρgydF4y2Ba_米gydF4y2Ba= 2,700 kg mgydF4y2Ba^{−3gydF4y2Ba}，gydF4y2BaρgydF4y2Ba_ugydF4y2Ba= 2,700 kg mgydF4y2Ba^{−3gydF4y2Ba},ΔgydF4y2BaρgydF4y2Ba= 300 kg mgydF4y2Ba^{−3gydF4y2Ba}而且gydF4y2BaDgydF4y2Ba= 6 km，值由我们的数据估计gydF4y2Ba问gydF4y2Ba_DgydF4y2Ba= 32米gydF4y2Ba^3gydF4y2Ba年代gydF4y2Ba^{−1gydF4y2Ba}而且gydF4y2Ba\ \(压裂{{rm \ d {}} Q} {{rm \ d {}} h} \)gydF4y2Ba=−0.0043 mgydF4y2Ba^2gydF4y2Ba年代gydF4y2Ba^{−1gydF4y2Ba}给了gydF4y2BalgydF4y2Ba= 13公里，也就是说震源深度估计为19公里。注意，这个结果与粘性压力损失项无关，因此不依赖于管道是圆柱形的(两者都是圆柱形的)gydF4y2Ba问gydF4y2Ba和dgydF4y2Ba问gydF4y2Ba/ dgydF4y2BahgydF4y2Ba取决于粘性项，当取比值时，粘性项会被抵消)。还要注意，模型不需要这样做gydF4y2BaρgydF4y2Ba_米gydF4y2Ba等于gydF4y2BaρgydF4y2Ba_ugydF4y2Ba．假设gydF4y2BaggydF4y2Ba= 9.8米gydF4y2Ba^{−2gydF4y2Ba}而且gydF4y2BaνgydF4y2Ba= 100 Pa s然后给出一个估计值gydF4y2Ba\ (r \)gydF4y2Ba0.9米。因此，我们推断导管的横截面积约为几平方米。如果有效的流动路径不能与圆柱形导管相比，那么它的横截面积可以预期是相似的大小，尽管粘性阻力取决于导管的形状gydF4y2Ba^58gydF4y2Ba．gydF4y2Ba