地震学报  2008, Vol. 30 Issue (5): 464-473
深部洞室中微小温度年度变化足以造成地应变年度变化
孙玉军1 , 李杰2, 曹建玲1, 卢双苓3, 石耀霖1    
1. 中国北京100049中国科学院研究生院计算地球动力学重点实验室
2. 中国济南250014山东省地震局
3. 中国山东泰安271600山东省泰安基准地震台
摘要:利用轴对称线弹性模型有限元方法, 计算了洞室中温度微小年变化引起的地倾斜和地应变. 结果表明, 在满足地震台站建设规范的情况下(洞室内温度年变化幅度应小于等于0.5℃), 当洞室中温度年变化幅度仅为0.2℃时, 仍会引起量级为10-7 rad的地倾斜和10-7的地应变, 并且在洞室的尽头和拐弯处由于热汇聚和热扩散的作用引起地倾斜和地应变的变化更为突出. 因此今后在观测中, 有必要进一步减小洞室内温度年变化;在台站设计建设和仪器布设中, 对基线式仪器和摆式仪器都可以做一些恰当考虑, 以尽量压制温度年变化带来的影响.
关键词有限单元法    地倾斜    地应变    年变化    
Small variation of annual temperature in deep tunnel can produce annual variation in tilt and strain
Sun Yujun1 , Li Jie2, Cao Jianling1, Lu Shuangling3, Shi Yaolin1    
1. Laboratory of Computational Geodynamics, Graduate University of the Chinese Academy of Sciences, Beijing 100049, China
2. Earthquake Administration of Shangdong Province, Jinan 250014, China
3. The Standard Seismic Station of Taian, Taian 271600, China
Abstract: In this paper, finite element method (FEM) of axisymmetric linear elastic model has been used to calculate the tilt and strain induced by small annual temperature variations in a deep tunnel. The results show that even if the amplitude of the annual variation meets the construction standard of seismic station issued by China Earthquake Administration (the annual temperature variation amplitude in the tunnel is no more than 0.5°C), a small annual temperature variation of amplitude just 0.2°C in the tunnel would produce 10-7 rad changes in tilt and 10-7 changes in strain. Especially, at the end and the corner of the tunnel, changes of tilt and strain can be even larger. Therefore, in the future, it is an important task to reduce the annual temperature variation in the tunnel as far as possible. Within the tunnel, for both baseline instrument and pendulum instrument, the modeling suggests ways of construction of the tunnel and installation of the instrument to decrease the influence of the annual temperature variation.
Key words: finite element method (FEM)    tilt    strain    annual variation    
引言

地壳形变前兆观测中的地应变和地倾斜观测是分析地震前兆异常的主要数据来源之一.洞室内观测条件比较好,温度、 湿度、 气压、 大风、 降雨等气象因素影响较小,地表噪声较小,仪器安置在岩体洞室内,可以进行连续自动观测,精度比较高. 因此,洞室形变观测在地震监测中发挥着重要作用. 但是在洞室中各种观测条件都较好的情况下,究竟哪些干扰因素占主导作用? 各种干扰因素对测量结果定量的影响程度有多少? 这些对于排除干扰因素、 突出地震响应,确定合理的观测条件以及仪器精度的选取等都具有重要的参考价值.

形变测量受到多种因素影响,其中包含了气象变化、 人为作用、 仪器故障及未知的诸多情况,使得干扰排除成为识别异常的一项必不可少的工作. 但在目前的认识水平下,干扰与被视为地震前兆的变化形态有时难以区分,给前兆的判定带来很大困难(王梅等,2004). 这些干扰因素中,气象变化的因素主要有气压、 气温、 降雨、 大风、 雷电等;人为作用主要有荷载、 地表噪声等;仪器故障主要是仪器的零漂. 国内外已开展了不少对各种干扰因素定量的研究工作. 例如,王曰风等(2004)根据怀来台体应变观测资料对气压和零漂的影响因素作了分析研究,并指出了一些剔除和资料处理方法;牟雅元(1999)根据西昌地震台压容应变仪的观测结果对仪器零漂给出了一些数据处理方法;胡卫建等(2002)也对地表荷载对地应力和地应变的影响作了一些试验研究;Zadro等(1987)Moro和Zadro(1998)以及Garavaglia等(2000)根据意大利Friuli地区的观测资料,详细分析了大气压、 降雨量、 温度对地倾斜和地应变的影响. 此外,不同台站对各种影响因素的响应不同,一般对某一台站而言,对某些特定的因素比较敏感. 王梅和季爱东(2000)对山东省内各个台站的异常特征进行了分析总结,如烟台台站受气压和降雨的影响比较明显,而泰安台站则受温度的影响比较明显.

在诸多影响因素中,温度条件的影响无疑是重要因素之一. 前人对温度、 地倾斜和地应变的影响往往仅限于作相关分析,但对这种相关的物理原因则没有进行定量分析. 近年来研究者开始注重这种定量的物理研究. 例如,曹建玲和石耀霖(2005)利用有限单元法,对地形起伏条件下地表温度年变化对地应力和地倾斜的影响进行了定量研究;孙玉军等(2008)则进一步利用泰安台的实际资料,对地表温度年变化对地倾斜和地应变的影响进行了有限单元法计算和分析.

我们注意到,在洞室测量中,依据地震台站建设规范(中国地震局,2004),对洞室的条件和观测设备都有一定的要求. 洞室室温要求年变化幅度小于等于0.5℃,日变化幅度小于等于0.03℃. 但是否满足了规范的这一要求,就可以完全消除了温度的年变化影响呢?本文利用有限元方法,在洞室内不超过规范要求的微小温度年变化条件下,计算了地倾斜和地应变会受到什么影响,并根据这些结果探讨了在观测中应该如何减少温度年变化因素带来的影响.

1 模型的建立和采用的方法 1.1 计算采用的模型

首先,注意到当地表水平、 介质均匀的情况下,如果地表温度变化规律为T=T0Tsin(ωt),所引起的地下温度变化的解析解(Turcotte,Schubert,1986)为

式中,T0为地表平均温度,ΔT为地表温度变化幅度,ω为温度周期性变化角频率,κ为热扩散系数,y为深度. 根据该解析解可知,如果取年周期性变化角频率ω=2×10-7 rad·s-1,热扩散系数κ=1×10-6 m2·s-1,地表温度的变化幅度为10℃时,引起地下50 m深度的温度变化幅度仅为1.37×10-6 ℃. 因此,在足够深的深部岩体中,地表温度年变化可以忽略不计. 这样,我们在下面的模型中不必讨论岩体受地表温度年变化的影响,集中讨论洞室内部温度的微小年变化对地倾斜和地应变的影响问题.

由于实际洞室中洞体结构非常复杂,在考虑洞室中微小温度年变化引起的地倾斜和地应变,进行三维有限元模拟计算时,要想得到需要的精度,需要剖分的网格数目多,计算量非常大,一般的串行程序是无法解决的.因此我们考虑了一个简单的实际情况并进行适 当的假设简化. 我们考虑这样一个简单的洞室,先垂直打一竖井,然后再打一些水平洞室,简化示意图如图 1所示. 竖井的效应在初步研究时不妨予以忽略,这样在下面的模型中仅考虑水平洞室内部温度的影响. 水平洞室的直径为2 m,向左右的进深均为15 m,以竖井的中部为左右对称轴,计算中仅考虑中心向右的15 m洞室,即考虑图 1aAA′线右边为计算对象,左边的情况应与右边对称. 由于实际的三维特征比较复杂,我们进一步把右边的问题简化为一个轴对称的问题,如图 1b所示,即把水平洞室的水平中心轴线(图 1b中的EF)作为模型的对称轴,平面BACDEF绕轴线EF旋转360°得到的空间圆柱体,对应图 1c中的圆柱体部分(为明显起见,我们取了圆柱体的一部分). 当CA足够大时,圆柱体外的部分受洞室温度变化的影响可以忽略不计;当洞室足够深时,地表温度变化的影响也可以忽略不计(或鉴于问题线性可叠加,另行单独考虑). 对于空间轴对称问题,旋转平面BACDEF上的变形情况即代表了整个圆柱体的情况. 这样,我们研究的问题简化为一个二维平面,如图 1d所示,EF为对称轴,模型的水平尺度为50 m,模型的半径也取为50 m. 模型的温度边界条件为: 水平洞室壁(CD, DE)的温度年变化遵从温度边界条件T=0.2sin(ωt),其它边界均为绝热边界. 模型的力学边界条件为: 左右端面边界(AC, BF)的水平位移为0,垂向位移自由;水平洞室壁为位移自由边界;圆柱体模型外缘R向位移自由.

图 1 模型建立示意图 (a) 洞室示意图;(b) 简化的轴对称问题;(c) 网格示意图;(d) 计算采用的模型及尺寸 Fig. 1 Model establishment (a) Schematic diagram of the tunnel; (b) Simplification of the axisymmetric issue; (c) Schematic diagram of the mesh; (d) Calculation model and its sizes
1.2 计算采用的方程

柱坐标r, θ, z下轴对称性的问题(z轴为对称轴),所有量都与坐标θ无关,剪应力τ和τ及环向位移uθ均为零. 因此基本方程可以简化(陆明万, 罗学富,1990). 我们将研究的介质视为线弹性的,采用如下方程进行求解:

热传导方程

平衡方程(无体力)

几何方程

本构方程

水平洞室中温度变化的边界条件设为

式中,T为温度;t为时间;κ为热扩散系数;ur, w分别为r, z方向上的位移; εr, εθ, εz分别为r, θ, z方向上的正应变;γzr为角应变;σr, σθ, σz分别为r, θ, z方向上的正应力;τzr为剪应力;E为杨氏模量;μ为泊松比;α1为热膨胀系数. 计算中这些系数采用的均为国际单位制,κ=1×10-6 m2·s-1,α1=10-5 K-1,E=70 GPa,μ=0.25,ω=2×10-7 rad·s-1.

2 计算结果及讨论

地震台站建设规范(中国地震局,2004)规定, 洞室地倾斜和地应变台站中的仪器主要有两种: 一种是基线式仪器(Baseline instrument). 通过设置于两设定点间的线形装置,以检测两个设定点间相对位置变化为工作原理的测量地倾斜、 地应变的仪器;第二种是摆式仪器(Pendulum instrument). 以检测在一个测点设定的摆的偏转角为工作原理的测量地倾斜的仪器. 根据这样的原理结合我们的模型计算结果,如图 1d所示,由于基线式仪器实际观测时观测点的仪器支墩设立在水平洞室地面下一定距离,因此我们在接近水平洞室内壁0.2 m的水平线上取点MN作为基线式仪器的两个设定点. 摆式仪器则在洞壁表面取3个点(M′,P′,N′),点M′与N′之间的距离为12 m,M′到洞口水平距离为1 m,N′到洞室的尽头水平距离为2 m,P′为M′和N′的中点. 如果我们取点MN作为基线式仪器观测时的两个设定点,那么MN之间的地倾斜和地应变可用以下两式表示:

式中, LM与N之间的距离, uv分别为水平位移和垂直位移.

同样对于摆式仪器的工作原理,我们取一点(例如P′)分析. 如果在该点设定摆式仪器来测定地倾斜,则可以用下式来表示该点的地倾斜:

如果能够在该点局部测量地应变,则应变为

用有限元的方法对以上模型进行计算分析. 采用四结点单元类型,并对接近洞室的部分区域进行网格加密,单元数为149301,结点数为148500. 时间步长取为1天,计算6个周期,即6年,使得结果脱离初始暂态影响而接近稳定的周期性变化,取其中最后的周期进行分析.

设定点M到洞口的水平距离为1 m,改变点N的位置使得MN之间的水平距离分别为6 m、 9 m和12 m. 图 2给出的是以MN为基线,在这3种不同基线长度情况下,用式(7)和式(8)分别计算出的地倾斜和地应变的结果. 可以看出,虽然洞室中的温度年变化幅度很小,只有0.2℃,但是这种变化却可以引起量级为10-7 rad的地倾斜和10-7的地应变,实际的地倾斜和地应变仪器的精度可达10-9 rad和10-9(中国地震局,2004),因此这种变化足以影响到实际仪器的观测结果. 根据本文计算结果,在满足洞室的室温年变化幅度小于等于0.5℃条件下,由于洞室中温度年变化的影响引起的地倾斜和地应变还是很显著的,因此,必须引起足够重视.

图 2MN为基线计算出的地倾斜(a)和地应变(b)随时间步变化图 Fig. 2 Result of tilt (a) and strain (b) variation with time steps based on the baseline MN

另外,从图 2中还可以看出,不同基线长度所计算出的地倾斜和地应变受温度年变化的影响是不同的,因此在实际中我们可以选取适当的长度和设定合适的两个基线点位置来减小温度的年变化影响. 但是如何选取才是合理的呢? 图 3给出了在计算的最后一个周期中选取不同的时间步长,由式(9)和式(10)计算出的地倾斜和地应变随洞室进深(从洞口算起)的变化曲线. 其中选取的时间步分别为1910 d, 1850 d, 1760 d, 1670 d和1580 d. 结合图 2图 3可以看出随着洞室进深的增加,尤其是在接近洞室的尽头处(对地倾斜来说是在洞室进深的11—14 m处,对地应变来说是在洞室进深的14 m处),地倾斜和地应变都出现了突然增大的现象. 地倾斜由10-8—10-7 rad量级增加到10-6 rad量级,地应变也从10-8—10-7量级突增到10-6量级,同样的现象也会出现在洞室的拐弯处. 图 4给出了时间步为1910 d,依据式(9)和式(10)计算出来的地倾斜和地应变分布云图. 由于在远离洞室的部位地倾斜和地应变变化极小,因此我们选取了靠近洞室的一个局部位置. 图 4c为选取的位置和尺寸以及计算时采用的网格剖分图. 从图 4a地倾斜云图和图 4b地应变云图中可以明显看出,随着洞室进深的增加,地倾斜和地应变在变化幅度上都呈增加的趋势. 由于洞室的存在,会导致温度分布的不均匀,使得温度的年变化存在热汇聚和热扩散,因此造成了在同一地层平面上地倾斜和地应变的不同,类似于曹建玲和石耀霖(2005)指出的在山顶和山脚出现热汇聚和热扩散的现象,在洞室的拐弯处和尽头处也会出现热汇聚和热扩散. 而热的汇聚和扩散必然会导致地倾斜和地应变以及地应力分布的不均匀,使得测量的结果变得不稳定. 因此为减小温度年变化带来的影响,在布设仪器的时候也应该尽量避开在洞室的拐弯和尽头处设置仪器支墩. 尤其对于基线式测量,测量室空间应该足够大,布设的两个设定点也应该尽可能设置在洞室的中央,以避开在洞室的拐弯和尽头处由于热汇聚和热扩散作用的影响. 根据地震台站建设规范(中国地震局,2004)的要求: 基线式仪器室,长大于等于10 m,宽2 m,高2.5 m;摆式仪器室,长大于等于3 m,宽2 m,高2.5 m . 另外,考虑到我们建立的模型具有对称性(如图 1b中所示),如果实际中左右两个水平洞室是对称的,那么图 3所得到的结果对于图 1a中左边部分也是适用的. 因此对于基线式仪器测量,可以相对于竖井(图 1a)对称的布置两个设定点.这样对两个设定点来说,由于温度年变化而受热汇聚和热扩散的影响是同等的,结合式(7)和式(8)可以看出,受同样温度年变化引起的地倾斜和地应变在两式相减的时候会消掉. 因此,对于基线式仪器的测量,选择合适的两个设定点,可以有效地压制由于温度年变化而带来的影响.

图 3 不同时间步地倾斜(a)和地应变(b)随洞室进深变化图 Fig. 3 Tilt (a) and strain (b) variation with tunnel length at different time steps

图 4 时间步长为1910 d的地倾斜和地应变分布云图和计算采用的网格剖分图
(a) 地倾斜分布云图; (b) 地应变分布云图;(c) 模型局部网格剖分图
Fig. 4 The contour of tilt and strain at 1910 days and part mesh of the model
(a) The contour of tilt;(b) The contour of strain; (c) Part mesh of the model

图 5给出的是依据摆式仪器原理选取的3个点M′、 P′和N′,由式(9)和式(10)计算出的地倾斜和地应变随时间变化图. 由图 5a可以看出,对于摆式仪器来说,正如以上所讨论的一样,受热汇聚和热扩散的影响更大,放在接近洞口处的点M′地倾斜变化幅度量级为10-8 rad,而在洞室尽头处的点N′量级则可以上升到10-6 rad. 因此摆式仪器的放置地点是一个重要问题. 实际中我们可以选择受温度年变化影响小的地方放置摆式仪器进行测量,以尽量避开热汇聚和热扩散比较明显的洞室拐弯和尽头处. 摆式仪器的灵敏度一般比较高,但是易受地面局部的影响;日本人认为这种仪器只反映测点本身的地面倾斜,因此不稳定(曾融生,1984). 近年来,一些高精度、 智能化、 网络化的仪器(王秀英等,2005肖峻等,2005施志龙,吴书朝,2007)不断产生,观测的准确度和抗干扰度大大提高. 由图 5b也可以看出,如果选择的地点合理,对于局部的地应变受温度年变化的影响相对是比较小的. 实际观测中压容应变仪就可以对这种局部的地应变进行测量.

图 5M′,P′和N′的地倾斜(a)和地应变(b)随时间步变化图 Fig. 5 Tilt (a) and strain (b) of points M′, P′, N′vary with the time steps

对比图 2中基线式计算结果和图 5中摆式计算结果可以发现,对于这两种仪器,年变化幅度可以有一定差别. 对于地倾斜,图 2a中年变化幅度在10-7 rad量级,而图 5a中年变化幅度为10-8—10-6 rad量级;对于地应变,图 2b中变化幅度为10-7量级,而图 5b中变化幅度为10-6—10-7量级. 因此,对于同一类仪器,其放置部位不同,年变化幅度会有所不同. 而对于不同仪器(基线式仪器和摆式仪器),其年变化幅度也会有所不同. 这一点在实际观测中也得到了验证. 图 6为泰安台实测的地倾斜和地应变在一年内的变化曲线. 图中伸缩仪和水管仪都为基线式仪器,其长度均为10.8 m; 石英摆倾斜仪属于摆式仪器. 由于实际泰安台站洞室结构比较复杂,地倾斜和地应变的变化方向与仪器的具体布置位置有密切关系,因此根据本文的研究需要,我们只考虑了地倾斜和地应变的变化幅度,而对其方向性暂不考虑. 对于地倾斜,图 6a图 6b中年变化幅度均为10-6 rad,但具体数值上有一定差别.而对于地应变,图 6c中的年变化幅度在 10-7量级,也有一定差别. 这说明在同一洞室中,同样的影响条件,对于不同仪器其测量结果也有明显区别. 需要指出的是, 图 6中地倾斜和地应变的年变化除洞室中温度年变化的影响外,仪器自身也存在温度系数的年变化. 实际台站中,由于洞室的复杂性及各种条件的影响,尤其是各种不均匀因素的影响作用,如上面提到的热汇聚和热扩散的作用,使得对于不同仪器其测量结果存在一定差别,甚至同一种仪器在不同部位测量结果也会有很大差别. 因此在布设仪器初期,选择合适的布设地点就显得尤为重要.

图 6 泰安台实际观测图
(a) 水管仪实测地倾斜; (b) 石英摆倾斜仪实测地倾斜; (c) 伸缩仪实测地应变
Fig. 6 Observation results at Taian seismic station
(a) Observation results of tilt with water pipe; (b) Observation results of tilt with quartz pendulum; (c) Observation results of strain with extensometer

虽然本文考虑介质为单一均匀的线弹性体,并且在建立模型的时候做了一定的假设,模型相对比较简单,但是从计算结果可以看出,即使洞室内温度年变化幅度很小的情况下(这里年变化幅度仅为0.2℃),仍然可以引起地倾斜和地应变较大的变化,这些变化在实际仪器测量中是能够观测到的. 因此,虽然根据规范要求,洞室中温度年变化幅度小于等于0.5℃是合格的,但在实际观测中应该尽量减小洞室中温度的年变化幅度,尽可能使洞室中保持恒温. 无论是基线式仪器或者是摆式仪器,在洞室的施工和仪器的布置方面都可以做一些详细的工作来减小温度年变化带来的影响.

如果实际情况更复杂,存在不均匀性和地表起伏地形,则温度年变化的干扰的幅度和位相都会更加复杂(曹建玲,石耀霖,2005). 但多种效应的叠加,幅度一般只会比现在计算的更大. 本文结论依然成立.

3 结论

本文利用有限元的方法,计算了洞室内温度年变化引起的地倾斜和地应变. 结果表明,即使当洞室中仅存在微小温度年变化(幅度仅为0.2℃)时,仍可以引起量级为10-7 rad的地倾斜和量级为10-7的地应变;热应变在洞室的不同部位分布是不均匀的,在洞室的拐弯和尽头处存在热汇聚和热扩散作用,因此在这些部位引起地倾斜和地应变的变化更为突出. 基线式仪器和摆式仪器的年变化幅度未必相同,同一类仪器安置部位不同,年变化幅度也会不同.

本研究得到的结果对台站建设和观测有重要的启示. 虽然台站建设规范对洞室内温度年变化提出了一定要求(年变化幅度小于等于0.5℃),但是这不等于说达到了规范规定的指标就可以消除年变化干扰. 台站有必要进一步采取措施,减少洞室内的温度年变化,才能提高观测质量. 在不可避免存在温度微小年变化时,在设计洞室和布设仪器的时候就应该对地热变形问题予以考虑: 水平的两个洞室相对于竖井尽量采取对称的开挖形式,避开在洞室的拐弯和尽头处设置仪器支墩. 对于基线式仪器来说,如果采用将基线式的两个测定点对称地布置在两个水平洞室内,则可以减小温度年变化的影响;对于摆式仪器来说,其放置地点是一个重要问题,同样应该选择受温度年变化影响小的地方,需避开洞室的拐弯和尽头处热汇聚和热扩散明显的地方.

参考文献
[1] 曹建玲, 石耀霖. 2005. 地表温度年变化对地应力和地倾斜的影响[J]. 中国科学院研究生院学报, 22(3): 303-308.(3)
[2] 胡卫建, 张俊山, 谢智, 王志敏, 李安印, 冯建民, 刘丽君. 2002. 荷载对钻孔应变测值影响的实验及力学解析[J]. 地震, 22(3): 95-104.(1)
[3] 陆明万, 罗学富. 1990. 弹性理论基础[M]. 北京: 清华大学出版社: 372-396.(1)
[4] 牟雅元. 1999. 西昌地震台压容应变仪观测数据的处理及效果检验[J]. 西北地震学报, 21(3): 337-339.(1)
[5] 施志龙, 吴书朝. 2007. 垂直摆倾斜仪灵敏度标定方法研究[J]. 测绘科学, 32(4): 63-64.(1)
[6] 孙玉军, 李杰, 曹建玲, 李峰, 石耀霖. 2008. 山东泰安台温度年变化对地应变与地倾斜影响的模拟研究[G]//中国地震局地震预测研究所, 中国地震台网中心, 中国地震学会编. 中国地震预报探索. 北京: 地震出版社: 460-466.(1)
[7] 王曰风, 张秀萍, 李海孝, 王长江, 段立新, 张常慧. 2004. 地应变震兆异常信息提取方法研究[J]. 地震地磁观测与研究, 25(1): 35-40.(1)
[8] 王秀英, 李海亮, 周振安. 2005. 洞室多分量定点高精度地应变测量仪器网络化的实现[J]. 地震地磁观测与研究, 26(5): 121-125.(1)
[9] 王梅, 李峰, 孔向阳, 刘厚明. 2004. 数字化形变观测干扰识别[J]. 大地测量与地球动力学, 24(1): 94-98.(1)
[10] 王梅, 季爱东. 2000. 山东省定点形变倾斜观测异常指标研究[J]. 地震研究, 23(4): 398-404.(1)
[11] 肖峻, 张军, 莫易敏, 解金芳. 2005. 基于空间多体中心构型的高精度垂直摆倾斜仪[J]. 中国机械工程, 16(11): 949-951.(1)
[12] 中国地震局. 2004. 中华人民共和国地震行业标准: 地震台站建设规范地形变台站(DB/T 8.1-2003)[S]. 北京: 地震出版社: 1-11.(4)
[13] 曾融生. 1984. 固体地球物理学导论[M]. 北京: 科学出版社: 351-356.(1)
[14] Garavaglia M, Moro G D, Zadro M. 2000. Radon and tilt measurements in a seismic area: Temperature effects[J]. Phys Chem Earth,, 25(3): 233-237.(1)
[15] Moro G D, Zadro M. 1998. Subsurface deformations induced by rainfall and atmospheric pressure: Tilt/strain measurements in the NE-Italy seismic area[J]. Earth Planet Sci Lett, 164(1): 193-203.(1)
[16] Turcotte D L, Schubert G. 1986. Geodynamics[M]. Second Edition. Cambridge: Cambridge University Press: 125-168.(1)
[17] Zadro M, Plenizio E, Ebblin C, Chiaruttini C. 1987. Influence of groundwater table level variations and of rainfall on tilting in the Friuli area, N-E Italy[J]. Acta Geophysica Polonica,, 35: 323-338.(1)