直角坐标系黏弹性层状地基荷载作用下的解析刚度矩阵解
寇磊，白云
(同济大学 地下建筑与工程系，上海，200092)
摘要：基于直角坐标系下黏弹性力学的基本控制方程，通过Fourier-Laplace变换及矩阵理论，推导出三维空间问题和平面应变问题在积分变换域的解析解，进而得到相应问题的精确单元刚度矩阵；然后根据对号入座原则组装得到总体刚度矩阵；通过求解总体刚度矩阵形成的代数方程，得到层状地基相应问题在积分变换域内的解答；应用Fourier-Laplace逆变换技术，得到其物理域内的解。求解黏弹性问题退化的弹性问题并与已有解答进行比较，验证本文计算方法的正确性，并分析黏弹性地基黏滞系数对沉降的影响。研究结果表明：黏滞系数越大，土体的蠕变越明显，地基达到最终沉降的时间将越长。
关键词：直角坐标系；黏弹性层状地基；刚度矩阵法；Fourier-Laplace变换；矩阵理论
中图分类号：TU433 文献标志码：A 文章编号：1672-7207(2014)07-2346-07
Analytical stiffness matrix of viscoelastic layered foundation under loading in Cartesian coordinate system
KOU Lei, BAI Yun
(Department of Geotechnical Engineering, Tongji University, Shanghai 200092, China)
Abstract: Based on the basic viscoelastic equations of Cartesian coordinate system, the analytical solutions of three-dimensional and plane strain in the integral transform domain were obtained by using the Fourier-Laplace transform and the matrix theory, and then the corresponding element stiffness matrix was derived. The global stiffness matrixes were assembled by matrix matching method, and the solutions for the corresponding problem of multilayered foundation in the transform domain were obtained by solving the algebraic equations of the global stiffness matrix. The solutions in the physical domain were acquired by inverting the Fourier-Laplace transform. The validity of the proposed method was examined by comparing the results of viscoelastic problem reducing to elastic problem with existing solutions, and the influences of viscosity parameter of viscoelastic foundation on the settlement were analyzed. The result shows that the soil creep increases with the increase of the viscosity parameter, and the duration of the final foundation settlement is longer.
Key words: Cartesian coordinate system; viscoelastic layered foundation; stiffness matrix method; Fourier-Laplace transform; matrix theory
土体承受荷载是结构、地下、水利和道路路基等实际土木工程中所关注的一类问题，天然地基在沉积形成中呈现出层状特征；在一般情况下，各层内土体性质均匀，而各层间差异较大；同时，荷载作用下的天然地基具有明显随时间变化的流变现象：因此，采用层状黏弹性地基模型描述地基更加符合实际状况。
目前，针对层状地基的研究大多集中在求解极坐标系下相关问题和弹性相关问题[1-10]。层状地基的求解方法主要包括解析法[1-7]、半解析半数值法[8-11]和数值方 法[12]；解析法主要有传递矩阵法[1-3]和刚度矩阵法[4-7]；半解析半数值法主要为有限层法[8-10]及样条半解析 法[11]等。鉴于解析法不仅具有精确求解的优点，而且其解答可以作为其他求解方法的基础，因此，研究直角坐标系黏弹性层状地基的解析求解十分有必要；同时，三维空间问题简化为平面应变问题同样具有实际应用意义。在现有解析求解层状地基的研究中，传递矩阵法在计算指数函数乘积时存在溢出问题[3]；钟阳[4-5]推导得到的弹性三维空间问题的单元刚度矩阵为6行6列，且矩阵元素表达式复杂不便实际应用。本文作者从直角坐标系黏弹性力学的基本控制方程入手，推导得到黏弹性层状地基三维空间问题和平面应变问题均为4行4列的解析单元刚度矩阵，并集成层状地基总体刚度矩阵；根据边界条件和层间连续条件求解总体刚度矩阵形成的代数方程得到相应问题的解答，为实际工程中的土体受荷分析提供理论基础。
1 单元刚度矩阵的推导
1.1 三维空间单元刚度矩阵
黏弹性力学中几何方程、积分型本构方程和不计体力的静力平衡方程为
(1)
(2)
(3)
式中：
为应变分量；ui为位移分量；
为应力分量；G(t)为剪切松弛函数；
；K(t)为体积松弛函数，
；
；Y(1)(t)为剪切松弛模量；Y(2)(t)为体积松弛模量。
将几何方程式(1)代入积分型本构方程式(2)，然后代入平衡方程式(3)，得到用位移表示的黏弹性三维空间问题的静力平衡方程，即
(4)
(5)
(6)
式中：
为体积应变；
为Laplace算子。
将积分型本构方程式(2)展开得到黏弹性三维空间问题的物理方程为
(7)
(8)
(9)
对式(4)，(5)和(6)进行x和y方向的双重Fourier变换和对时间t进行Laplace变换，得到
(10)
(11)
(12)
式中：
；
；
；
为便于推导，作以下定义：
；
；
则得到
；
(13)
将式(13)代入式(10)，(11)和(12)得到
(14)
(15)
(16)
将式(14)和(15)写成齐次常微分方程组的矩阵形式，即
(17)
令
且为可逆矩阵，
，
，
，
，
，
，则方程(17)变化为
(18)
式中：
；
，为系数矩阵。求得系数矩阵S的特征值为
，
得到矩阵S的Jordan标准形为
设可逆矩阵P=(P1, P2, P3, P4)使
，求得
根据微分方程组理论[13]，得到方程(18)的通解为
则得到
(19)
对式(7)，(8)和(9)进行x和y方向的双重Fourier变换和对时间t进行Laplace变换，得到
(20)
(21)
(22)
为便于推导，作以下定义：
，
，
由式(13)，(19)和(20)得到：
(23)
(24)
(25)
由式(19)得到位移分量与待定常数Ci (i=1, 2, 3, 4)之间的矩阵关系式，即
(26)
由式(23)和式(24)得到应力分量与待定常数Ci (i=1, 2, 3, 4)之间的矩阵关系式，即
(27)
因U和V均可以由待定常数Ci (i=1, 2, 3, 4)表达，可进一步得到之间的关系为
。其中：
，为三维空间问题在积分变换域内的单元刚度矩阵。
1.2 平面应变单元刚度矩阵
根据平面应变的条件，即
，得到用位移表示的黏弹性平面应变问题的静力平衡方程为
(28)
(29)
式中：
；
为Laplace算子。
平面应变问题的物理方程为
(30)
(31)
对式(28)，(29)，(30)和(31)进行x方向的Fourier变换和对时间t进行Laplace变换，得到
(32)
(33)
(34)
(35)
将式(32)和(33)写成齐次常微分方程组的矩阵形式，即
(36)
不难发现，方程(36)和方程(17)的表达形式是一样的，根据1.1节的求解过程得到
(37)
由式(34)，(35)和(37)得到
(38)
(39)
由式(37)得到位移分量与待定常数Ci (i=1, 2, 3, 4)之间的矩阵关系式，即
(40)
由式(38)和式(39)得到应力分量与待定常数Ci (i=1, 2, 3, 4)之间的矩阵关系式，即
(41)
因U1和V1均可以由待定常数Ci (i=1, 2, 3, 4)表达，可进一步得到之间的关系为
。其中：
，为平面应变问题的在变换空间内的单元刚度矩阵，且
，即三维空间问题和平面应变问题的单元刚度矩阵均为4行4列。
2 层状地基的刚度矩阵求解
2.1 三维空间问题求解
黏弹性层状地基三维空间问题的计算模型如图1所示，根据天然土层及计算点将土层划分为n个计算层，第i个计算层底部距地表的距离为Hi(i=1, 2, …, n)，第i个计算层的厚度为
，在第i计算层处作用荷载
。
将1.1节推导的三维空间问题的单元刚度矩阵应用于各计算层，并结合层间连续条件，按照对号入座原则集成整个层状地基的总体刚度矩阵，得到
(42)
式中：K(i)为第i层刚度矩阵；
和
为作用在Hi深度处的荷载经过Fourier-Laplace变换后的形式。
图1 黏弹性层状地基三维空间问题
Fig. 1 Three-dimensional of viscoelastic layered foundation
结合已知边界条件求解方程(42)，可得到层状地基三维问题在积分变换域内的解，对其按照文献[14]中的方法进行Fourier-Laplace逆变换即可得到真实物理域内的解。
2.2 平面应变问题求解
黏弹性层状地基平面问题的计算模型如图2所示，由于平面应变问题作为三维空间问题的1个特例，层状地基平面问题的求解同样采用2.1节的过程。
3 计算结果
由于研究直角坐标系下黏弹性问题的文献不多，为验证本文方法的正确性，将黏弹性问题退化为弹性问题。在弹性问题中，
，
图2 黏弹性层状地基平面应变问题
Fig. 2 Plane strain of viscoelastic layered foundation
式中：E为弹性模量；
为泊松比。对积分变换域内的解进行Fourier逆变换即可得到真实物理域内的解。
3.1 平面应变问题验证
采用本文方法与文献[15]中分析的有限深土层表面作用竖向均布荷载时沿深度方向位移的变化进行对比。本文方法计算时Fourier逆变换采用分段32 点Guass 积分实现，步长取为π；当段数取到20时计算结果趋于稳定。当=0.3和=0.4，H=10a时，
的计算结果如图3所示。
图3 均布荷载作用下有限深土层的表面变形
Fig. 3 Surface displacement of finite soil layer subjected to uniform loading
从图3可知：本文计算结果与文献[15]的解析解相当吻合。
3.2 三维空间问题验证
计算弹性半无限体表面作用竖向集中荷载Fz=1 kN时土体沉降随深度的变化。计算中，E=2 000 kPa，=0.3，计算点分别为A点x=1.0 m，y=0；B点x=2.0 m，y=0；C点x=3.0 m，y=0。由表1可知：本文计算结果和Boussinesq解的计算结果相当吻合。
3.3 黏弹性层状地基分析
探讨均布矩形荷载作用下黏弹性黏滞系数变化对多层黏弹性地基变形的影响。计算时选取4层黏弹性土作为地基模型进行研究，如图1所示，各层土参数见表2。均布矩形荷载作用在地基表面，为20 kPa，荷载作用面积为a×b=2 m×2 m=4 m2，土体体积变形符合弹性性质，剪切变形采用Kelvin模型，即
；
当第1层土的黏弹性黏滞系数从4 MPa·s变化至20 MPa·s时，地基内土体位移随时间变化情况如图4所示。
表1 A点、B点和C点随深度变化变形值
Table 1 Displacement along depth of A, B and C mm
表2 黏弹性层状地基的土体参数
Table 2 Soil parameters of viscoelastic layered foundation
图4 黏弹性系数变化时土层变形随时间的变化
Fig. 4 Soil time-displacement under different viscoelastic coefficients
由图4可知：当第1层土体的黏弹性系数变化时，随着黏弹性黏滞系数的增大，土体的蠕变性越明显，地基达到最终沉降的时间将变长，但沉降的最终值保持不变。
4 结论
(1) 基于直角坐标系黏弹性力学问题的控制方程，通过Fourier-Laplace变换及矩阵理论推导得到三维空间问题和平面应变问题均为4行4列的单元刚度矩阵，根据边界条件和层间连续条件求解层状黏弹性地基在荷载作用下的总体刚度矩阵代数方程，通过积分逆变换得到真实物理域内的解。
(2) 由于刚度矩阵的对称性及矩阵元素只存在负指数，从而避免传递矩阵法中正指数存在导致的计算溢出问题，提高计算精度和稳定性。
(3) 对于黏弹性层状地基，随着黏滞系数的增大，土体的蠕变性越明显，地基达到最终沉降的时间将变长，但变形的最终值保持不变。
参考文献：
[1] 张子明. 用传递矩阵法求解多层地基的轴对称问题[J]. 工程力学, 1986, 3(2): 96-103.
ZHANG Ziming. Solving axis-symmetry problems of layers by recurrence matrix method[J]. Engineering Mechanics, 1986, 3(2): 96-103.
[2] 钟阳, 王哲人, 郭大智. 求解多层弹性半空间轴对称问题的传递矩阵法[J]. 土木工程学报, 1992, 25(6): 38-43.
ZHONG Yang, WANG Zheren, GUO Dazhi. The transfer matrix method for solving axisymmetric problems in muti-layered elastic half space[J]. China Civil Engineering Journal, 1992, 25(6): 38-43.
[3] 艾智勇, 吴超. 三维直角坐标系下分层地基的传递矩阵解[J]. 重庆建筑大学学报, 2008, 30(2): 43-46.
AI Zhiyong, WU Chao. Transfer matrix solutions for multi-layered soils in rectangular coordinate system[J]. Journal of Chongqing Jianzhu University, 2008, 30(2): 43-46.
[4] 钟阳, 耿立涛. 多层弹性平面问题解的精确刚度矩阵法[J]. 岩土力学, 2008, 29(10): 2829-2832.
ZHONG Yang, GEN Litao. Explicit solution of multiplayer elastic plane by exact stiffness matrix method[J]. Rock and Soil Mechanics, 2008, 29(10): 2829-2832.
[5] 钟阳. 多层弹性半空间问题解的精确刚度矩阵法[J]. 应用力学学报, 2008, 25(2): 316-319.
ZHONG Yang. Explicit solution of multilayered elastic half space by exact stiffness matrix method[J]. Chinese Journal of Applied Mechanics, 2008, 25(2): 316-319.
[6] 艾智勇, 曾凯, 曾文泽. 层状地基三维问题的解析层元解[J].岩土工程学报,2012,34(6):1154-1158.
AI Zhiyong, ZENG Kai, ZENG Wenze. Analytical layer- element solution for three-dimensional problem of multilayered foundation[J]. Chinese Journal of Geotechnical Engineering, 2012, 34(6): 1154-1158.
[7] 艾智勇, 苏辉, 成怡冲. 求解层状地基平面应变问题的解析层元法[J]. 岩土工程学报, 2011, 33(11): 1797-1799.
AI Zhiyong, SU Hui, CHENG Yichong. Analytical layer element method for solving plane strain problem of multi- layered soils[J]. Chinese Journal of Geotechnical Engineering, 2011, 33(11): 1797-1799.
[8] Small J C, Booker J R. Finite layer analysis of layered elastic materials using a flexibility approach: Part 1. Strip loadings[J]. International Journal for Numerical Methods in Engineering, 1984, 20(6): 1025-1037.
[9] Small J C, Booker J R. Finite layer analysis of layered elastic materials using a flexibility approach: Part 2. Circular and rectangular loadings[J]. International Journal for Numerical Methods in Engineering, 1986, 23(5): 959-978.
[10] 张问清, 赵锡宏, 宰金珉. 任意力系作用下的层状弹性半空间的有限层分析方法[J]. 岩土工程学报, 1982, 3(2): 27-42.
ZHANG Wenqing, ZHAO Xihong, ZAI Jinmin. Finite layer analysis for layered and elastic half-space under arbitrary force system[J]. Chinese Journal of Geotechnical Engineering, 1982, 3(2): 27-42.
[11] 姬亦工, 王复明, 栾茂田. 层状黏弹性地基动力响应的样条半解析法[J]. 岩土工程学报, 2002, 24(1): 35-37.
JI Yigong, WANG Fuming, LUAN Maotian. A simplified method for vibration analysis of layered visco-elastic foundation[J]. Chinese Journal of Geotechnical Engineering, 2002, 24(1): 35-37.
[12] Runesson K, Booker J R. Finite element analysis of elastic- plastic layered soil using discrete fourier series expansion[J]. International Journal for Numerical Methods in Engineering, 1983, 19(4): 473-478.
[13] 严国政. 常微分方程[M]. 北京: 科学出版社, 2012: 32-36.
YAN Guozheng. Ordinary differential equations[M]. Beijing: Science Press, 2012: 32-36.
[14] Small J C, Booker J R. A Fourier-Laplace transform finite element method for the analysis of contaminant transport in porous media[J]. International Journal for Numerical and Analytical Methods in Geomechanics, 1996, 23(14): 1763-1796.
[15] 赵光恒, 张子明. 有限深弹性层上基础梁的计算[J]. 华东水利学院学报, 1984(2): 32-44.
ZHAO Guangheng, ZHANG Ziming. Computation of foundation beams on elastic layers of finite depth[J]. Journal of East China Institute of Water Conservancy, 1984(2): 32-44.
(编辑 杨幼平)
收稿日期：2013-07-08；修回日期：2013-10-16
基金项目：教育部长江学者和创新团队发展计划项目(IRT1029)
通信作者：寇磊(1983-)，男，河南许昌人，博士研究生，从事岩土力学及在隧道与地下工程中应用的研究；电话：021-65980401；E-mail: klyhe@163.com
