引入接触单元模拟桩土共同作用

⎤Τ ∂g
⎥ ⎦
Dc ( ∂Pc
)Τ
(9)
Dc
=
⎡ Eτ ⎢⎣ 0
0⎤ E n ⎥⎦
(10)
参变量最小势能原理：在所有满足应变-位移关
系式和几何边界条件式的可能位移场中，真实解使
总势能泛函 ∏[λ(•)] 在接触系统状态方程的控制下
取总体最小值：
∏[λ(•)]
=
∫Ω
1ε 2
Τ DεdΩ
−
⎡ ⎢⎣
∫Ω
λkυk = 0 υk , λk ≥ 0 k =1，2，3
(7)
将 fk 作一阶 Taylor 展开，可以将状态方程写成另一 种表达形式：
f
0 k
+ ωkεc
−
mk λ
+υk
=
0；
λk ⋅υk = 0 υk , λk ≥ 0 k =1，2，3
(8)
其中
f
0 k
为起始时
fk
值，而
mk
=
⎡ ⎢ ⎣
∂f k ∂Pc
128
岩土力学
2005 年
模拟接触的单元，它根据接触条件，把接触点对（在
接触面上坐标相同的节点）的位移和接触力，以单
元的形式进行表示，可以直接向刚度矩阵中组装，
形成的总刚可以进行向接触面上的凝聚，得到在接
触点的经过缩聚的刚度矩阵。
接触面上一点的接触力可分解为法向接触力
Pn 和切向接触力 Pτ ，即 Pc = {Pτ ,Pn }Τ ，接触条件为：
QI Liang-feng1 , Jian Hao2 , TANG Li-yun3
(1. Department of Civil Eng., Shanghai Institute of Technology, Shanghai 200235, China ; 2. Shanghai Geotechnical Engineering Investigation & Design Institute Co.Ltd. ,Shanghai 200002, China ; 3. Department of Civil Eng, Xi’an University of Science and Technology, Xi’an 710054, China)
摘 要：在桩与桩侧土界面引入接触单元模型，并采用参变量变分原理及基于此原理的参变量二次规划法，对接触单元刚度
矩阵进行了推导。模型计算结果与实测比较表明，接触单元能较好地模拟桩与土之间的剪力传递和相对位移。
关 键 词：接触单元；参变量变分原理；桩土共同作用
中图分类号：O 343.3
文献标识码：A
The computing model of adopting contact element to simulate the pile-soil’ s reciprocity
−
u
(2) n
+δ
*
= =
∆uτ ∆u n
⎫
+
δ
*
⎬ ⎭
(2)
式中
u(a) n
，
u(a) τ
分别为桩土域相对接触点局部坐标
系下的切向、法向位移；δ * 为桩土间的初始间隙（如
打入桩存在过盈配合， δ * < 0 ），也可以表示为：
εc
=
ε
e c
+
ε
p c
，ε
e c
为弹性相对位移，即发生接触但未
达到滑动时的相对位移；ε
用分析中，桩采用弹性材料，土采用理想弹塑性材
料，采用 Drucker-Prager 屈服准则，使用相关流动
准则，其屈服面并不随着土的逐渐屈服而改变，因
而没有强化准则。
3 实例分析
陕西信息大厦主楼试桩桩长 83 m，地质资料和 物理力学性质指标见表 1，用本文有限元计算模型 进行计算，其桩侧摩阻力分布如图 3，实测试桩截 面轴力传递[4]如图 4 所示
有限元分析平面点接触单元的刚度矩阵和约束
矩阵为
kce
=
N
e c
Τ
Dce
N
e c
=
⎡Eτ ⎢⎣ 0
0⎤ En ⎥⎦
(15)
Φ
c e
=
N
e c
Τ
ReΤ
=
⎡Eτ ⎢⎣ 0
− Eτ 0
0⎤ En ⎥⎦
(16)
⎡ Eτ µEn ⎤
Cec
=
ωe
N
e c
=
⎢− ⎢⎣
Eτ 0
µEn En
⎥ ⎥⎦
(17)
⎡ Eτ − Eτ µEn ⎤
⎬
(19)
∑8
u = N iui
i=1
⎪ ⎪ ⎪
∑ w =
8 i=1
N i wi
⎪ ⎪⎭
单元刚度矩阵表达式：
[k] e= ∫∫A[B]Τ[D] [B]dA
(20)
2.3 材料特性
在桩土体系中，桩土的弹性模量相差很大，有
2 个量级左右，桩受荷后，一般处于弹性阶段，而
其周边土容易达到塑性状态。因此，在桩土相互作
bΤudΩ
+
∫P Sp
Τ
ud s
⎤ ⎥⎦
+
∫Sc
(
1ε 2
c
Τ
Dcε
c
−
λΤ Rε c )dS
(11)
式中
R
=
(
∂g ∂Pc
)Dc
为常数矩阵；位移
u
为自变量；λ
为不直接参加变分的参变量。因此，接触问题的求
解化为：min. ∏[λ(•)]
s.t.
f (uc , υΤλ =
λ) + 0,υ,
υ λ
=0⎫ ≥ 0⎭⎬
(12)
此问题可以通过参数二次规划法求解，具体解 法参见文献[3]。
将接触物体进行单元划分，在接触边界用接触 单元划分，接触点对的主位移与通常的节点位移一 起编排在总体位移向量 uˆ 中，而接触点对的相对位 移编排在总体位移向量的最前端，离散化后平面弹性 接触问题有限元二次规划求解方程的提法归纳如下：
—
—
—
—
⑥ 粉质粘土 11.80 369.30
24.3
20.1
0.690
0.31
0.17
10.1
⑦ 粉质粘土 6.60
362.70
24.0
20.0
0.688
0.26
0.14
15.0
b
中砂
2.50
360.20
—
—
—
—
—
—
⑧ 粉质粘土 8.70
351.50
23.5
20.3
0.665
0.23
0.14
13.4
(4)
Pτ
=
⎩⎨⎧−
Eτ ετ µPn sign(ε
n
)
(5)
Pτ < −µPn Pτ = −µPn
(6)
其中， En , Eτ → ∞ ，为惩罚因子。 引入原约束松弛变量υk ,令： fk + υk = 0 , k =1，2，3 则接触问题的状态方程为：
f k(∆uτ , ∆un ,λk ) + υk = 0 ;
p c
为滑动相对位移，定义
滑动相对位移
ε
p c
如下：当状态点落在
f1
=
0 时，对
应滑动量为 λ1 ， λ1 > 0 ；当状态点落在 f2 = 0 时，
对应滑动量为 −λ2 ，λ2 > 0 ；当状态点落在 f3 = 0 时，
物体间相互滑脱，对应滑脱量为 λ3 ， λ3 > 0 ；对应
fk ，k =1, 2, 3 定义滑动势函数：
表 1 地基土物理力学性质指标统计表 Table 1 Statistic table of foundation soil’s indexes of physico-mechanical properties
层号
土层 名称
土层厚度 层底标高
/m
/m
ω
γ
%
/kN·m-3
主要物理力学性质指标平均值
α1-2
⑨ 粉质粘土 13.40 338.10
24.2
20.0
0.682
0.32
0.14
12.5Βιβλιοθήκη g1 g2= Pτ + c0 = −Pτ + c0
⎪⎫ ⎬
(3)
g3 = Pn ⎪⎭
c0 为任意常数，则滑动相对位移可表示为：
∑ ε
p c
=
3
λk
k =1
∂g k ∂Pc
，式中 ∂gk ∂Pc
体现相对滑动位移方向，
λk 代表滑动量大小。位移与接触力的关系条件可表 示为：
β (ε n ) =Pnε=n [1E−n βsi(gεnn()ε n )]/ 2⎭⎬⎫
第1期
齐良锋等：引入接触单元模拟桩土共同作用
129
min. ∏[λ(•)] = 1 uˆΤ kuˆ − uˆ(Φλ + Pˆ)
(13)
2
Cuˆ −Uλ − d +υ = 0⎫
υΤ λ = 0,υ ,λ
≥
0
⎬ ⎭
(14)
式中 k, Φ , Pˆ, C, U , d 的矩阵的组装与普通有限单
元法相同。
1 引言
随着计算机的发展，有限元在桩土共同作用分 析中得到了愈来愈多的应用，其一个突出的优点就 是通过对土体区域的划分，解决土体分布的层状性 问题。同济大学沈伟跃、赵锡宏[1]应用有限元法对 桩土进行了分析，得出了一些满意的结果，但这种 分析方法认为，桩侧和桩周土之间没有相对位移。 大量工程实践表明，桩侧和桩周土之间不仅能传递 摩阻力，而且它们之间还有较大的相对位移。文献 [2]在有限元分析的基础上，在桩与土之间引入摩擦 单元（节理单元），本文在桩与土之间引入另一种单 元--接触单元以模拟桩与土之间的剪力传递和相对 位移。