2. 中国科学院研究生院, 北京 100049;
3. Department of Geoscience, University of Calgary, Calgary T2N 1N4, Canada
2. Graduate University of Chinese Academy of Sciences, Beijing 100049, China;
3. Department of Geoscience, University of Calgary, Calgary T2N 1N4, Canada
末次冰期始于约110kaBP,在约18kaBP 的盛冰期,北美劳伦地区、北欧芬诺斯坎底亚和极地发育巨厚的冰盖,在大约6—8kaBP完成了消融过程[1].相对海平面(RSL)变化是海平面相对地壳的上升或下降,自盛冰期以来,RSL 被海岸和海洋各类沉积物记录下来,成为研究这段历史RSL 变化最重要的数据集,为研究海岸、海洋环境变化和地壳运动等提供重要证据.研究表明,盛冰期以来RSL 变化的机制主要有:(1)冰融水流进了海洋,由于重力场的区域差异,其分布也有区域变化;(2)冰盖消融和海水增加改变了地表负荷,固体地球要进行均衡调整,其直接的响应是地壳的垂直运动,在近场(冰盖及其周围)较强,在远场(例如东亚地区)较弱;(3)沉积(或剥蚀)作用[2];(4)沉积压实作用[3];(5)局部和区域构造运动[3].本文关注前两个因素,即冰盖消融和地壳均衡调整,可广义地称为冰川均衡调整(GIA)[1],我们将从全球GIA 模型出发,模拟和给出东亚历史相对海平面变化的结果,并结合RSL 观测进行对比分析研究.长期以来,地球物理学家在研究GIA 对RSL的贡献时,主要精力放在近场,对处在远场的东亚地区讨论较少[4].地质学家在研究东亚历史RSL 变化时,需要考虑GIA 的影响,采用了冰融等效海平面(IESL),即冰盖消融成的海水体积除以全球海洋面积[2]或者整个地区统一的模型[4],没有考虑冰融水分布或者地壳均衡垂直运动的区域性差异.事实上,地球物理学家不断完善冰盖消融模型和地幔黏滞度模型(例如Peltier的模型系列[5-6]),可以根据严密的负荷理论和海平面方程[7-8],模拟出GIA 对东亚历史RSL 变化的影响,但是对于东亚地区一直没有详细的结果和分析研究.
我们最近利用近十年空间大地测量观测数据的约束,包括北欧和北美GPS网络监测的地壳运动速率[9-10]、在北美陆地用卫星测高[11-12]、在五大湖地区用验潮站结合卫星测高[13-14]解算的地壳隆升速率和GRACE 卫星重力计划发布的时变重力观测结果,建立了基于横向非均匀的岩石圈和地幔流变的GIA 模型RF3L20(β=0.4)+ICE-4G[15-16].该模型已用于模拟对东亚现今重力场变化、空间大地测量监测地表现今(陆地水和海水)质量变化的影响[17].本文将利用该模型,进一步模拟东亚盛冰期以来的RSL 变化,评估模拟的不确定性,将模拟结果与观测[18-19]进行比较分析.
2 模型与算法在我们发展的GIA 模型RF3L20(β=0.4)+ICE-4G 中,横向非均匀地球模型RF3L20(β=0.4)包含厚度有横向变化的岩石圈(图 1)和黏滞度有横向变化的四层地幔(图 2)[7-8],盛冰期以来冰盖消融模型采用Peltier的ICE-4G 模型[5].
岩石圈厚度是根据地震波速转换估计的[20].在东亚(图 1)沿海地区、近海和较深海区,岩石圈厚度为50km;向东亚大陆西部岩石圈增厚,大部分地区岩石圈厚度为65km,在青藏高原,岩石圈厚度增加到90km;向日本以东太平洋海域岩石圈厚度为65km.
四层地幔包括上地幔(UM)、过渡带(TZ)、下地幔浅部(LM1)、下地幔深部(LM2),参考模型RF3给出黏滞度为0.6×1021 Pa·s(UM+TZ)、3.0×1021 Pa·s (LM1)、6.0×1021 Pa·s (LM2),在此(对数黏滞度)基础上叠加对数黏滞度的横向扰动,就得到三维地幔对数黏滞度.对数黏滞度的横向扰动是根据地震剪切波速异常转换的(图 2),所用线性比例因子β=0.4[15-16].由图 2 可知,在东亚地区较大的地幔黏滞度横向非均匀出现在上地幔,最大超过一个数量级,其次在过渡带,最大也接近一个数量级,而横向非均匀在下地幔的两个分层中则不明显.地球密度和弹性参数仅随深度变化,详见文献[15-16].
使用了耦合拉普拉斯方程的有限元算法进行GIA 模拟[7-8, 16].有限元模型分辨率在浅层为2°×2°,向核幔边界逐渐减少到8°×8°.在5次迭代计算过程中,利用有限元位移结果,通过解海平面方程计算新的RSL,从而不断修正海水负荷.
为了评估RSL 模拟结果的不确定性,对当前参考模型黏滞度和冰盖消融模型认识差异的影响进行计算,取两方面结果的平方根作为RSL 估计的不确定性.
对参考模型黏滞度,考虑两个代表性模型结果的差异.一是Peltier(1998)的VM2模型[5],即本文简化的RF3模型,其黏滞度为0.6×1021 Pa·s(UM+TZ)、3.0×1021 Pa·s(LM1)、6.0×1021 Pa·s(LM2);另一是Mitrovica和Forte(1997)[21]根据地幔对流和GIA联合反演的结果,这里用RF2模型简化表示,给出黏滞度为0.7×1021 Pa·s (UM+TZ)、1.0×1022 Pa·s(LM1+LM2),可见下地幔较上地幔黏滞度有大幅度的升高.
对冰盖消融模型,考虑有代表性的ICE-4G[5]及其升级版ICE-5G[6]模型的差异.当前这两个模型使用最为广泛,不过其差异也是明显的:冰盖的厚度有显著差异,特别在北美哈德逊湾以西地区,ICE-5G模型推测有巨厚的冰盖,而ICE-4G 则没有;冰盖负荷总质量也有明显差异.由于IESL 变化与冰盖负荷总质量成反比,因此图 3 所示两模型的IESL 随时间上升,反映了随时间冰盖负荷总质量的减少.ICE-5G 在26 kaBP 的盛冰期冰盖开始融化,而ICE-4G 则推迟至18kaBP.在12kaBP 以前ICE-4G 的IESL 较高,说明冰盖负荷总质量相对较小.12kaBP以后两模型的IESL和冰盖负荷总质量相近.
本节给出了东亚RSL 模拟结果,评估结果的不确定性,最后结合图 1测点实测数据进行对比分析.
3.1 东亚RSL的模拟结果根据RF3L20(β=0.4)+ICE-4G 模型[15-16],模拟出东亚盛冰期(18kaBP)以来每隔2ka共8幕的RSL 变化图像,结果如图 4 所示.前6 幕RSL 均是负值(图 4a—图 4f),表明整个东亚海平面相对固体地壳降低;后2 幕有负也有正(图 4g—图 4h),说明海平面相对固体地壳在一些区域降低、另一些区域升高.由图可知,8 幕的RSL 变化结果分别是:-141~-105m、-117~ -81 m、-80~ -52 m、-56~-33m、-32~ -15 m、-13~ -1 m、-5~3 m、-1.8~1.9m,可见盛冰期以来海平面随着冰盖消融逐渐上升.这反映冰融海水是RSL 上升最重要的控制因素,直接影响了RSL 变化的主体,但对于冰盖消融结束后(图 4f—图 4h),RSL 变化主要与地壳黏性均衡调整引起的垂直运动有关.由图还可以看出,8幕RSL 具有显著的区域性差异,其最大最小幅度差为28、23、17、12、8、3.7 m.一般而言在水浅区域RSL 高于IESL,而在较深水区域RSL 则在IESL之下,这反映了不同区域地壳对海水负荷的均衡作用差异,部分反映不同区域全球冰盖消融的远场均衡效应的差异,地壳均衡下沉则RSL 就增高,地壳均衡隆升则RSL 就下降.
RSL 模拟的不确定性来自参考模型黏滞度差异和冰模型差异的影响.
在计算参考模型黏滞度差异对模拟RSL 的影响时,考虑当前对参考模型黏滞度认识的差异,针对模型RF3 和RF2 分别进行有限元模拟,冰模型都采用ICE-4G[5],根据有限元输出的位移分别计算RSL 变化,然后计算二者之间的差值.图 5a—5d给出了18kaBP以来每隔4ka共4幕RSL 变化差异的图像,可以看出,每幕参考模型黏滞度差异对RSL 的影响都不大,差异基本上在2m 以内;此外,从RSL 差异平面分布看,发现均具有以冲绳以东太平洋负值异常和沿东亚大陆近海正异常梯度带为特征,与图 2所示4层地幔黏滞度横向非均匀图像相似,推测由参考黏滞度不同和黏滞度横向非均匀影响耦合引起的.
关于冰盖消融模型差异对RSL 的影响,在已有RF3L20(β=0.4)+ICE-4G 模型[15-16]结果的基础上,进一步采用ICE-5G[6]冰盖消融模型,用有限元法进行模拟和RSL 变化的计算,地球模型仍采用RF3L20(β=0.4)[15-16],然后计算ICE-4G 与ICE-5G 预测RSL 变化的差值.图 5e—5h 给出了18kaBP 以来每隔4ka 共4 幕的结果.由图可知,18kaBP、14kaBP、10kaBP、6kaBP的结果分别为2~13m、8~17m、-6~-1m、0.1~1.3m,正值或负值表示ICE-4G给出的RSL较高或较低.ICE-4G在12kaBP前的RSL 较高,由图 4 可知,与这一时期冰盖负荷总质量相对较小或IESL 相对较高有关;12kaBP以后,两模型的冰盖负荷总质量或IESL 几乎相同,其RSL 差异则与12kaBP 前冰盖负荷作用的差异有关.每幕结果都呈现区域性变化,其变化规律较复杂,即在18kaBP、14kaBP 向西南海域和海岸冰模型差异的影响减少,向东北海域则影响增强,而10kaBP、6kaBP的情形则正好相反,这与冰模型负荷时空演化差异以及与岩石圈厚度横向变化、黏滞度横向非均匀影响耦合等因素有关.
计算参考模型黏滞度差异和冰模型差异对RSL 影响的平方根,以此评估RSL 模拟结果(图 4)的不确定性,图 6给出了18kaBP以来每隔2ka共8幕的结果.计算表明,在18—12kaBP 冰模型差异的影响较大(图 5e—5g),因此成为图 6a—6d 所示不确定性的控制因素,而在10—4kaBP,参考模型黏滞度差异与冰模型差异的影响量级逐渐相当,因此图 6e—6h所示不确定性为两因素共同作用的结果.8 幕的RSL 变化不确定性范围分别是:1.0~7.0m、8.5~14.5 m、4.0~8.0 m、0.4~3.4 m、0.6~2.8m、1.0~1.9m、0.1~0.7m、0.0~0.3m,其幅度基本上随时间流逝而减少;RSL 变化不确定性分布也呈现显著的区域性变化,例如在18kaBP、16kaBP、14kaBP 最大最小幅度差分别为6.0 m、4.0m、3.0m.一般而言,每幕不确定性从东亚沿海向东北向近海、太平洋逐渐增加,但10kaBP 和8kaBP则情形正好相反;从图 7 可知,RSL 变化不确定性一般在观测的误差范围内.
本节将模拟计算的RSL 变化结果与东亚现有观测进行对比和研究.如图 1 所示,RSL 观测在中国沿海和台湾包括测点(编号):扬子三角洲(1)、福建—台湾海峡(2)、韩江三角洲(3)、广东东部(4)、珠江三角洲(5)、广东西部(6)[18]、天津(635)和台湾北端(649)[22];在国外包括俄罗斯纳霍德卡(632)、日本北海道(637)、日本东京湾(639)、日本鹿儿岛(644)[22].中国测点主要采用Zong(2004)收集和重新审定的结果[18],C14年龄已根据与树年轮、珊瑚年龄的关系进行了校正[23];天津(635)、台湾北端(649)和国外的4 个测点的RSL 数据采用Tushingham和Peltier(1991)[22]收集整理的结果,但我们利用Fairbanks等(2005)[24]给出的C14与原生珊瑚年龄关系,对C14年龄进行校正.
图 7给出了东亚测点根据RF3L20(β=0.4)+ICE-4G 模型[15-16]的RSL 模拟结果和实测结果.从图可以看出,除日本3个测点和中国台湾北端(649),东亚海洋统一的IESL 曲线显然不能反映大多数测点模拟的RSL 区域变化,因此利用统一的IESL[2]或统一的GIA RSL[4]研究东亚RSL 变化机制是不合适的.还可看出,在扬子三角洲(1)、天津(635)、台湾北端(649)、纳霍德卡(632)、北海道(637)、东京湾(639)、鹿儿岛(644),RSL的GIA理论预测结果与观测较符合,说明这些地区的RSL 变化基本受GIA控制,其他因素影响非常小.图 8 给出了RSL 实测值经过GIA 改正的结果,在华南沿海地区的测点,RSL 的GIA 理论预测结果与观测有较大差异,这是由非GIA因素引起的,下面进行解释.
图 8中剩余RSL 负值表示古地壳高于现今位置的高度.因此,第一阶段在8kaBP以前,地壳经历下沉和隆升(测点2、4、6),第二阶段在8kaBP 以后一段时间,基本以较稳定的速率下沉,最后是第三阶段则以较小速率下沉或隆升.每个测点三个阶段的区间不同,第二第三阶段地壳变化速率不同,表现出显著区域性差异.对测点2、3、4、5,第二阶段区间和地壳下沉速率分别为6—8kaBP 和-4.2 mm/a、2.8—7kaBP和-2.1mm/a、3—8kaBP和-1.5mm/a、2.8—8kaBP和-1.0mm/a;对测点2、3、5,随后的第三阶段下沉速率分别为-0.1mm/a、-0.3mm/a、-0.3mm/a,而测点4地壳则隆升,速率为+0.2mm/a.广东西部(6)情况较复杂,近7ka以来,地壳先隆升后有下沉,速率分别为+0.5 mm/a、-0.6 mm/a.因此,在华南测点RSL 观测经过GIA 模型改正后,揭示了较复杂的地壳运动,推测与东南部的菲律宾板块俯冲有关,具体机制有待进一步研究.在东亚的其他地区测点,RSL 的GIA 理论预测与观测结果相当,说明构造运动对RSL 变化影响较小.
我们将图 8华南地壳垂直运动长期历史趋势与大地测量观测结果对比.由图可以看出,用最近1.4ka的2个数据,计算出粤西海岸和粤东海岸地壳垂直运动速率为+0.6mm/a,用最近1ka的2个数据计算出珠江三角洲结果为+2.0mm/a.卢汝圻(1997)根据30多年的精密水准复测数据,揭示粤西海岸段平均以+0.6 mm/a速率上升、粤东海岸段平均以-0.7mm/a速率下降、珠江三角洲以-1.8 mm/a速率下沉[19].由此可见,在粤西海岸,现今地壳垂直运动与本文揭示的长期地壳垂直运动基本吻合,但在粤东海岸和珠江三角洲则二者存在运动方向相反的根本差异,可能与人类活动导致的沉降有关.
4 结 论利用新的冰川均衡调整(GIA)模型,模拟了末次盛冰期以来东亚RSL 变化,评估了其不确定性,将模拟结果与观测数据进行了比较分析.主要结论如下:
(1) 给出了东亚RSL 时空分布模拟结果,其不确定性一般小于历史RSL 的观测误差,该模型结果对研究东亚历史RSL 具有重要意义.
(2) 模拟的RSL 在8kaBP 前随时间的上升与全球冰盖消融有关,而后期的变化则主要受地壳黏性均衡调整控制.每个时间的RSL 均具有显著的区域性差异,与局部地壳对海水负荷的均衡作用及远场均衡效应差异有关.RSL 的不确定性主要受当前冰盖消融模型差异的影响,其次是参考黏滞度的影响.
(3) 东亚GIA 相关的RSL 有显著的区域性变化,因此用IESL 或代表性观测站的GIA 结果研究RSL 变化机制是不合适的.
(4) 东亚历史海平面变化主要受全球GIA 控制.在扬子三角洲、天津、台湾北端、纳霍德卡、北海道、东京湾、鹿儿岛测点,RSL 的GIA 理论预测较好地解释了观测结果.
(5) 在华南海岸测点,RSL在排除GIA影响后,清楚地显示全新世以来的地壳垂直运动,发现早期地壳经历下沉和隆升,中期地壳以较稳定的速率(1~4mm/a)下沉,后期地壳以较小速率下沉或隆升,推测可能与东南部菲律宾板块的俯冲有关.
(6) 在粤东海岸和珠江三角洲,现今地壳垂直运动是下沉的,与本文揭示的近千年地壳运动隆升趋势不同,推测与人类活动导致的沉降有关.
[1] | 姚衍桃, HarffJ, MeyerM, 等. 南海西北部末次盛冰期以来的古海岸线重建. 中国科学D辑:地球科学 , 2009, 39(6): 753–762. Yao Y T, Harff J, Meyer M, et al. Reconstruction of paleocoastlines for the northwestern South China Sea since the Last Glacial Maximum. Sci. China Ser. D-Earth Sci. (in Chinese) , 2009, 39(6): 753-762. |
[2] | 辛立国, 李广雪, 李西双, 等. 中国东海2万年来海平面变化分析. 中国海洋大学学报(自然科学版) , 2006, 36(5): 699–704. Xin L G, Li G X, Li X S, et al. Analysis of sea-level variation in the East China Sea since the Last Glacial Maximum. Periodical of Ocean University of China (in Chinese) , 2006, 36(5): 699-704. |
[3] | 汪汉胜, WuP, 许厚泽. 冰川均衡调整(GIA)的研究. 地球物理学进展 , 2009, 24(6): 1958–1967. Wang H S, Wu P, Xu H Z. A review of research in glacial isostatic adjustment. Progress in Geophys. (in Chinese) , 2009, 24(6): 1958-1967. |
[4] | 陈伟光, 张虎男, 李子权. 闽粤沿海全新世垂直构造运动的速率. 地震地质 , 1986, 8(2): 33–42. Chen W G, Zhang H N, Li Z Q. The velocity of Holocene vertical tectonic motion in the coastal area of Fujian and Guangdong. Seismology and Geology (in Chinese) , 1986, 8(2): 33-42. |
[5] | Peltier W R. Postglacial variations in the level of the sea: implications for climate dynamics and solid-earth geophysics. Rev. Geophys. , 1998, 36(4): 603-689. DOI:10.1029/98RG02638 |
[6] | Peltier W R. Global glacial isostasy and the surface of the ice-age earth: the ICE-5G (VM2) Model and GRACE. Annual Review of Earth and Planetary Science , 2004, 32(1): 111-149. DOI:10.1146/annurev.earth.32.082503.144359 |
[7] | Wu P. Using commercial finite element packages for the study of earth deformations, sea levels and the state of stress. Geophys. J. Int. , 2004, 158(2): 401-408. DOI:10.1111/gji.2004.158.issue-2 |
[8] | Wang H S, Wu P. Effects of lateral variations in lithospheric thickness and mantle viscosity on glacially induced relative sea levels and long wavelength gravity field in a spherical, self-gravitating Maxwell Earth. Earth Planet. Sci. Lett. , 2006, 249(3-4): 368-383. DOI:10.1016/j.epsl.2006.07.011 |
[9] | Lidberg M, Johansson J M, Scherneck H, et al. An improved and extended GPS-derived 3D velocity field of the glacial isostatic adjustment (GIA) in Fennoscandia. J. Geod. , 2007, 81(3): 213-230. DOI:10.1007/s00190-006-0102-4 |
[10] | Sella G F, Stein S, Dixon T H, et al. Observation of glacial isostatic adjustment in "stable" North America with GPS. Geophys. Res. Lett. , 2007, 34: L02306. |
[11] | Lee H, Shum C K, Kuo C Y, et al. Application of TOPEX altimetry for solid earth deformation studies. Terr. Atmos. Ocean. Sci. , 2008, 19(1-2): 37-46. DOI:10.3319/TAO.2008.19.1-2.37(SA) |
[12] | Lee H, Shum C K, Yi Y C, et al. Laurentia crustal motion observed using TOPEX/POSEIDON radar altimetry over land. J. Geodyn. , 2008, 46(3-5): 182-193. DOI:10.1016/j.jog.2008.05.001 |
[13] | Kuo C Y, Shum C K, Braun A, et al. Vertical crustal motion determined by satellite altimetry and tide gauge data in Fennoscandia. Geophys. Res. Lett. , 2004, 31: L01608. |
[14] | Kuo C Y, Shum C K, Braun A, et al. Vertical motion determined using satellite altimetry and tide gauges. Terr. Atmos. Ocean. Sci. , 2008, 19(1-2): 21-35. DOI:10.3319/TAO.2008.19.1-2.21(SA) |
[15] | Wang H S, Wu P, van der Wal W. Using postglacial sea level, crustal velocities and gravity-rate-of-change to constrain the influence of thermal effects on mantle lateral heterogeneities. J. Geodyn. , 2008, 46(3-5): 104-117. DOI:10.1016/j.jog.2008.03.003 |
[16] | 汪汉胜, WuP, vander Wal W, 等. 大地测量观测和相对海平面联合约束的冰川均衡调整模型. 地球物理学报 , 2009, 52(10): 2450–2460. Wang H S, Wu P, van der Wal W, et al. Glacial isostatic adjustment model constrained by geodetic measurements and relative sea level. Chinese J. Geophys. (in Chinese) , 2009, 52(10): 2450-2460. |
[17] | 汪汉胜, 贾路路, WuP, 等. 冰川均衡调整对东亚重力和海平面变化的影响. 地球物理学报 , 2010, 53(11): 2590–2602. Wang H S, Jia L L, Wu P, et al. Effects of global glacial isostatic adjustment on the secular changes of gravity, and sea level in East Asia. Chinese J. Geophys. (in Chinese) , 2010, 53(11): 2590-2602. |
[18] | Zong Y Q. Mid-Holocene sea-level highstand along the Southeast Coast of China. Quaternary International , 2004, 117(1): 55-67. DOI:10.1016/S1040-6182(03)00116-2 |
[19] | 卢汝圻. 广东沿海地区现代地壳垂直运动研究. 华南地震 , 1997, 17(1): 25–33. Lu R Q. Study on the modern crustal vertical movement in Guangdong coast. South China Journal of Seismology (in Chinese) , 1997, 17(1): 25-33. |
[20] | Wang H S, Wu P. Effects of lateral variations in lithospheric thickness and mantle viscosity on glacially induced surface motion on a spherical, self-gravitating Maxwell Earth. Earth Planet. Sci. Lett. , 2006, 244(3-4): 576-589. DOI:10.1016/j.epsl.2006.02.026 |
[21] | Mitrovica J X, Forte A M. Radial profile of mantle viscosity: results from the joint inversion of convection and postglacial rebound observables. J. Geophys. Res. , 1997, 102(B2): 2751-2769. DOI:10.1029/96JB03175 |
[22] | Tushingham A M, Peltier W R. Ice-3G: a new global model of late Pleistocene deglaciation based upon geophysical predictions of post-glacial relative sea level change. J. Geophys. Res. , 1991, 96(B3): 4497-4523. DOI:10.1029/90JB01583 |
[23] | Stuiver M, Bard E, Beck J W, et al. INTCAL98 radiocarbon age calibration, 24, 000-0 cal BP. Radiocarbon , 1998, 40(3): 1041-1083. DOI:10.1017/S0033822200019123 |
[24] | Fairbanks R G, Mortlocka R A, Chiu T C, et al. Radiocarbon calibration curve spanning 0 to 50, 000 years BP based on paired 230Th/234U/238U and 14C dates on pristine corals. Quaternary Science Reviews , 2005, 24(16-17): 1781-1796. DOI:10.1016/j.quascirev.2005.04.007 |