自Vening Meinesz采用走航式三摆重力仪进行海洋重力测量和均衡学研究(Vening-Meinesz, 1948), 在20世纪50年代即已认识到岩石圈的弹性薄板模型是足够精确的(Gunn, 1949), 依据弹性力学薄板理论(Timoshenko, 1958)和岩石圈均衡响应理论(Jeffreys, 1945; Gunn, 1949), 伴随20世纪70、80年代至今的海洋和大陆岩石圈重力均衡和挠曲研究(Walcott, 1970; Watts et al, 1974; McNutt et al, 1978; Watts, 2001), 岩石圈在火山岛屿、俯冲带、挠曲盆地等处的挠曲变形可以用弹性薄板模型来近似研究的认识已经获得广泛的承认(Watts, 2001; Turcotte et al, 2014)。

目前岩石圈挠曲研究中采用的是克希霍夫(Kirchhoff)弹性薄板小挠度弯曲平衡方程(钱伟长等, 1956;柳春图等, 2001):

其中, 是薄板中面上任一点的横向(垂直于板面方向, 向下z轴方向)位移, 即挠度; 式中是双调和/双谐算子; q是薄板受到的横向总负荷; D是板的挠曲刚度。

其中, t是板的厚度, E是板材的杨氏模量, μ是泊松比(因下文中板的挠曲位移采用v, 本文用μ记泊松比)。

克希霍夫方程由于其数学形式上的简洁而被广泛应用, 但在推导过程中采用的各项假设中, 有一个关于薄板内垂向应力₂^σ的假设是

它是基于克希霍夫薄板理论的一个主要假设——薄板的厚度比薄板的其它尺寸小, 可以把垂向应力忽略不计(钱伟长等, 1956), 它的物理意义是平行于板面的各个面在挠曲中不相互挤压。此项假设, 式(3), 不仅在物理上非常牵强, 而且在推导式(1)的过程中还与其它假设不相协调。因此假设式(3)是克希霍夫方程推导过程中的一个问题, 也是一个缺陷。

在岩石圈挠曲的研究中仍然采用这一假设, 并因为薄板表面垂向应力为0, 而且薄板相对很薄, 可以取垂向应力在薄板中都为0(特科特等, 1986; Turcotte et al, 2014)。但很难设想岩石圈在施加负载时其内部应力没有变化, 这一假设是许多学者不赞同的(金煜等, 2002)。可是鉴于克希霍夫方程数学上的简洁和国际惯例, 目前仍普遍接受采用这一假设得到的克希霍夫方程。

本文针对克希霍夫方程中的这一问题, 提供了一个不需要假设式(3)的推导过程, 建立了一个新的弹性薄板挠曲方程:

其中D_FF是新的挠曲刚度,

但前人在建立式(1)中所用的其他5个假设, 即下文中式(6)、(9)、(12)、(26)及(35), 在本文建立式(4)的过程中仍保留使用, 但并不增加任何新的假设。

1 推导过程

仍旧把直角坐标系O-xyz之原点取在板之中面内, 且中面选为xOy平面, z轴则垂直板面向下, 如图 1所示。

板内任意一点(x, y, z)在挠曲时的位移记为(u, v, w)。一般而言, u, v及w均是(x, y, z)之函数。

采用z轴方向应变为0的假设

则描述弹性体形变的几何方程(Cauchy方程)组

其中, 第三式给出w与z无关, 仅是x与y的函数, 即有:

采用假设

可由式(7)之第五、第六式给出:

将式(10)两边对z积分, 可得:

采用假设(中面内各点无面内位移):

得及, 于是, 式(11)变为:

将式(13)代入式(7), 就有

及

依弹性体变形的物理方程(胡克定律), 有

其中, σ_x、σ_y、σ_z分别表示各个正应力, τ_xy、τ_yz、τ_xz分别表示各个剪应力, G是板材的切变弹性模量, 即

现在, 本文摒弃假设式(3), 保留σ_z, 进行如下恒等变换。

首先由式(16)中第三式解出σ_z, 得到

将式(18)代入式(16)中第一、二式, 有

视此二式为σ_x及σ_y的二元一次代数方程, 求解有

现在, 引用假设式(6), 并注意到式(14), 此二式变为

又将式(15)及式(17)代入式(16)之第四式得到式(23), 待用。

依照弹性体的静力平衡方程组(Navier方程), 有

及

其中, (X, Y, Z)是弹性体在(x, y, z)处的单位体积所受到的外力(如重力)。

仍旧采用假设

则式(24)化为

注意到式(25), τ_xy= τ_yx, 则式(27)可写成

将式(21)及(23)中的σ_x及τ_xy代入式(30)得到

同样, 将式(22)及(23)代入式(28)可得到

将式(31)及(32)的两边对z积分, 并注意到式(8), 有

及

这里,

而F₁(x, y)及F₂(x, y)是任意函数, 待定。

依旧采用克希霍夫方程推导中采用的假设:

即在板的上、下表面上无剪切应力, 则可由式(36)、(33)与(34)解出F₁(x, y)及F₂(x, y), 从而得出

根据式(25)、(37)、(38), 式(29)变成为

如果弹性薄板受到的体力只有重力, 那么Z=ρg, 这里ρ是薄板的材料密度, 而g是重力加速度。但我们这里不做此限定, 而是一般地, 容许有可变体力

为了由方程(39)推导出薄板挠曲的微分方程, 我们引入一个辅助量P如下:

于是,

仍用q₁及q₂记薄板顶面及底面上受到的横向载荷, 见图 2。

则有如下的边界条件成立(依惯例, 张应力为正):

仍旧用q表示薄板受到的横向总载荷, 即

则有如下边界公式成立:

将(42)式代入(39)式, 有

将此式两边对z积分, 得到

由(46)及(49)可求出F₃(x, y),

再将F₃(x, y)代回(49), 有

将条件(47)代入式(50), 有

这就是弹性薄板小挠度弯曲中面微分方程。它是在没有采用假设(3)的情况下由弹性体几何方程、物理方程和静力平衡方程推导成立的。

为了与经典的克希霍夫方程(1)及(2)式相对照, 引入符号D_FF代表新的挠曲刚度,

则挠曲方程(51)可写成

这个结果显示新建立的方程(52)与克希霍夫方程(1)具有同等的数学简洁性。至于就数值影响而言, 试取泊松比μ=0.25, 则有

对此, 新挠曲方程(52)中的挠曲刚度(因子D_FF)要比经典克希霍夫挠曲方程的D值大12.5%。

2 讨论 2.1 本文推导的弹性薄板小挠度弯曲方程无需垂向应力为0的假设

在弹性薄板小挠度弯曲方程的推导中假设垂直薄板方向的应力为0(钱伟长等, 1956), 鉴于薄板厚度相对很小, 而且获得的克希霍夫挠曲方程具有简介的数学形式, 至今仍是广泛采用的经典公式(柳春图等, 2001;刘人通, 2002)。但在本文的上述推导表明, 这个假设可以省去, 而且本文公式(43)及(44)还表明假设条件式(3)也并不成立。在岩石圈挠曲的研究中, 虽然很早就确认了岩石圈是在冰川、火山岛屿等负载的作用下发生挠曲的, 其挠曲可以近似地采用弹性薄板来拟合(Gunn, 1949; Walcott, 1970), 但在其挠曲方程的推导中仍继承了克希霍夫挠曲方程, 这就默认了垂向应力为0的假设, 并认为薄板的表面无应力、而且薄板厚度很薄, 因此可以假设薄板上下垂向应力均为0(特科特等, 1986; Turcotte et al., 2014)。也曾有人对岩石圈挠曲中采用垂向应力为0的假设提出疑问(金昱等, 2002), 但未解决该问题。而本文推导的岩石圈挠曲方程则不需要这一假设, 并且解决了这一问题。

2.2 在岩石圈挠曲研究中的应用

在俯冲带海沟、大洋岛屿等相对简单的负载下, 可以直接采用岩石圈挠曲方程和海底地形来计算挠曲参数, ρ_m和ρ_w分别为地幔和海水的密度(Walcott, 1970; Turcotte et al, 2014), 进而获得岩石圈的挠曲刚度D和弹性厚度T_e。在更广泛的地质背景下, 则是采用重力异常与地形的导纳来计算岩石圈有效弹性厚度(McKenzie et al, 1976; McNutt et al, 1978; Watts, 1978; Forsyth, 1985; Jin et al, 1994; Watts, 2001)。为了研究区域内岩石圈有效弹性厚度的变化, Kalnins等(2009)和Pérez-Gussinyé等(2004)提出了滑动窗口导纳技术(MWAT), 获得了广泛的应用(杨安等, 2016; Yang et al, 2018)。导纳法也是由获得的挠曲刚度计算岩石圈有效弹性厚度值。

若采用本文推导的岩石圈挠曲新方程, 则由其挠曲刚度D_FF之值得到的岩石圈弹性厚度值T_FF与克希霍夫挠曲方程计算得到的弹性厚度值T_e的关系式为。

上述表明, 虽然本文推导的岩石圈挠曲新方程的挠曲刚度值D_FF与“经典的”D值有明显的不同, 但从D_FF及D值进一步计算的岩石圈有效弹性厚度值之间相差很小, 仅有不到4%的差别。但鉴于新挠曲方程无需垂向应力为0的假设, 这一假设在岩石圈动力学研究中是无法得到合理解释的, 采用新挠曲方程计算的T_FF显然更合理, 而且只需将已有的岩石圈T_e值做约4%的校正即可, 并不要求推翻已有的岩石圈T_e研究中有关的方法。

实际上岩石圈有效弹性厚度是为了研究和比较岩石圈在不同构造背景(即不同应力状态和负载)下挠曲响应的参数, 进而分析其动力学性质, 在物质上是不存在的, 因此本文提出的岩石圈T_FF值虽然比T_e值更合理, 但不会在实际应用中给出显著不同的结论。

2.3 新挠曲方程有待于三维数值模拟或光弹性法实验的验证

本文推导的弹性薄板小挠度弯曲的新方程在弹性力学中也可得到应用, 例如对于同样弹性模量的材料、同样负载、同样挠曲(或弯曲)的情况下, 本方程得到的材料厚度比相对经典克希霍夫方程得到的材料的厚度可以降低约4%。

与岩石圈挠曲中的应用不同, 本文推导的新挠曲方程可以通过三维数值模拟或光弹性法模拟(曾佐勋等, 1992)来验证, 分析垂向应力的影响, 这有待于业内专家的研究工作。

3 结论

与岩石圈挠曲——弹性薄板小挠度弯曲中垂向应力为0的假设不同, 本文在不需要假设岩石圈垂向应力为0的假设, 推导了新的岩石圈－弹性薄板挠曲方程, 在数学形式上与经典的克希霍夫挠曲方程一样, 但新的挠曲刚在μ=0.025时比传统的挠曲刚度大12.5%。因此对于同样数值的挠曲刚度, 新方程计算的T_FF略低于经典方程计算的岩石圈有效弹性厚度值T_e, T_FF=0.9615T_e, 两者仅有不到4%的差值。本方程也可以在弹性力学中得到应用。

致谢 本文是国家自然科学青年基金项目(40206009)的后续研究成果。谨以此文祝贺中国科学院海洋研究所成立七十周年。