基于malab的土压力问题的分析与求解

基于malab的土压力问题的分析与求解

0土压力及其作用点位置coulumb土壤压力理论是法国科学家coulumb在1773年的力学文章中提出的“建筑力学的最大和小规律应用”中提出的。这是200多年前提出的，在重力墙的结构设计中仍然很常见。Coulomb土压力理论虽然解决了主动土压力和被动土压力的大小问题，但对力的作用点位置问题依然采用了土压力强度沿墙高线性分布的假定。然而，大量的现场观测和室内测试资料表明：土压力沿墙高并非线性分布，而是某种曲线形式；土压力的作用点位置与Coulomb理论的结果不一致，并非总是作用在1/3墙高处，而是与墙体变位模式的类型等多种因素有关，是个变化的值。另外，Coulomb土压力理论是以土体平面滑动破坏假定为基础的，此假定的物理背景仍需在理论上进行分析。1土壤闭合负荷的基本方程1.1物权粘聚力分析建模计算分析的基本假定：(1)所研究的问题为平面应变问题；(2)墙后土体为Coulomb材料，可用其物性参数粘聚力c和内摩擦角ϕ表征；(3)墙后土体产生主动或被动土压力时，土体形成滑动楔体，其滑裂面通过墙踵；(4)挡土结构刚性且其运动不受限制，运动位移与墙高相比可忽略不计；(5)墙土之间的摩擦可用摩擦角δ表示。1.2土体平衡方程图1和图2分别表示墙后土体处于主动和被动极限状态的情况，AC表示挡土结构，AB为土层的水平面，BC为土层的破坏面，W为ABC部分土体的重量，Pa、Pp分别表示作用在墙背上的主动和被动土压力，δ为墙土间的摩擦角，并且假定其小于土体的内摩擦角ϕ，d为力的作用点距墙底的垂直距离，h为墙高，σ、τ分别为作用在滑动面上的法向、切向应力。取ABC部分的土体为隔离体，可以建立其处于平衡状态时的2个力的平衡方程。由ΣX=0可得式中Pn为土体的极限荷载。当n=1时，Pn为主动土压力Pa，当n=-1时，Pn为被动土压力Pp;y˙=dy/dx。由ΣY=0可得式中γs为土体的重度。设在破裂面上法向应力σ与切向应力τ服从Mohr-Coulomb破坏准则，即式中c、ϕ分别为土体的粘聚力和内摩擦角，代入式（3），式（1）、式（2）变为将式（4）代入式（5）可得至此，推导出了滑动土体ABC处于极限平衡状态的2个力的平衡方程，式（4）、式（6）。1.3lagrange函数式（4）、式（6）表示土压力Pn为含有2个自变量函数y(x）、σ（x）的泛函，且为边界待定的条件变分极值问题。滑裂面曲线的起始点为（墙踵处）C点，坐标为xc=0,yc=0，终点为坡面上的B点，坐标为x=xB,xB待定，yB=h。依据泛函在约束条件下的变分法，用Lagrange乘子法构造如下的泛函J*，从而使上述条件极值问题转化为无约束的极值问题：式中λ1为Lagrange乘子；其余符号的意义同上。滑动面的方程y(x）及沿滑动面的y(x）分布的法向应力σ（x）须满足以下条件：(1)辅助函数F的Euler微分方程(2)积分约束方程：式（6）。(3)两类边界条件1)可动边界点处的变分边界条件—横截条件：式中δ为变分算子；其余符号的意义同上。2)固定边界点条件：1.4土体平动破坏模型求解Euler方程，不难证明沿平面滑动是墙后土体可能的破坏机制。由式（11）及式（8）、式（9）、式（10）、式（14）可得由式（12）、式（13）及式（8）、式（9）、式（10）可得至此，已经求出了土体平动破坏时滑动平面的方程—式（15）和沿滑动面的法向应力分布—式（16）。1.5土体极限荷载的推导当土体的内摩擦角ϕ为常数时，式（15）是一条直线方程。从而可知，当墙后土体达到临界状态时，其破坏形式是沿平面滑动。下面推导临界状态下土体极限荷载的表达式。首先将滑动面BC的方程改写为式中α为滑裂面与水平面的夹角，如图3、图4所示。由式（15）和式（17）可得另外，由图3、图4可得下面的关系式：式中xB为滑动体右上角点B的X坐标。将式（17）、式（18）、式（20）及式（16）代入式（1）和式（2），化简整理可得至此，推导出了土体极限荷载Pn的表达式—式（21）和式（22）。1.6土压力pne土体极限荷载的作用点位置在确定了土体极限荷载的大小后，根据滑动体ABC静力平衡的力矩式求出。此时，作用在滑动体ABC上的外力有主动土压力和被动土压力Pn的极值Pne、滑动面BC上作用的法向力σ和切向力τ、滑动体ABC受到的重力W对C点的力矩为MP、Mσ、Mτ、MW。因τ的作用线过C点，所以Mτ=0。力系对C点的力矩平衡可写为参照图3、图4容易写出由式（16）可知，土压力Pn取极值Pne时滑动面上的从而同时，由图3、图4可知将式（24）、式（25）、式（27）、式（28）、式（29）代入式（23），可得式中Pne为土体极限荷载Pn的极值；αe、λ1e、xBe分别为土体极限荷载Pn取极值Pne时相应滑动面与水平方向的夹角、Lagrange乘子、滑动体上B点x坐标。确定Pne作用点的计算步骤如下：首先利用式（21）、式（22）确定αe、Pne，接下来利用式（19）、式（20）确定λ1e、xBe，最后把求得的参数代入式（30）就可求出Pne的作用点到墙踵的距离d与墙高h的比值。2土压力函数的求解在利用式（21）和式（22）确定土体极限荷载的大小时，涉及到求解带有约束的一元函数的极值问题。需借助于某种优化方法确定其极值。经过分析，本文选取了Matlab6.1优化工具箱提供的fmincon函数。fmincon函数的功能为求多变量有约束非线性函数的最小值。因而，在求解主动土压力问题时，需转化为求取-Pa的最小值。对于本文求解问题，相应的约束为：(1)角度α限定的区间为[0.π/2]，α为滑动面与水平方向的夹角；(2)滑动破裂面上的切向应力τ大于0；(3)式（21）×sinδ—式（22）×cosδ=0。利用fmincon函数，编制了各种情况下的土体极限荷载计算程序的m文件。3土体极限荷载计算前文将土体极限荷载的泛函极值问题转化为一元函数的极值问题。可利用fmincon函数确定土体极限荷载的大小，同时也确定了临界滑动土体的几何尺寸和位置，最后利用滑动体静力平衡的力矩式求出土体极限荷载作用点在墙体的相对位置d/h。3.1最佳土压力的确定为了验证本文计算方法的合理性及正确性，下面先将c=0情况下的部分计算结果与Coulomb(1776）的解答进行了对比。此时，无量纲化的土压力P为式中K为土压力系数。表1、表2分别为主动土压力、被动土压力计算结果的对比，可以看出，土压力大小的计算结果完全一致的，但力作用点位置并非总是作用在墙高的1/3处。对于土体极限荷载作用点不是总在墙高的1/3处，相关测试已经证明。3.2土体粘聚力的极限荷载当粘聚力c≠0时，本文计算了其变化对土体极限荷载大小的影响，见图5。可见，无量纲化的土体极限荷载大小Pne/（γsh2）与无量纲化后的土体粘聚力c/（γsh）成线性关系。考虑了土体粘聚力c的影响，土体的极限荷载的大小Pne在形式上可以表达为式中Kϕ为与土体的内摩擦角ϕ相关的土压力系数，可通过式（31）计算；Kc为与土体的粘聚力c相关的土压力系数。(2）土体粘聚力c的影响利用滑动体ABC满足的力矩平衡方程—式（30）确定土体极限荷载的作用点在墙高的相对位置d/h。土体的粘聚力c=0时，可以推导出土体极限荷载作用点在墙高的相对位置d/h的解析表达式。如图6，此时沿滑动面BC作用的法向应力σ在墙高的方向成线性分布，墙顶处为0，因而，滑动体受到的反力R作用在滑动面BC上距C点h/（3sinαe）长度处。根据平面汇交力系三力（重力W、反力R和土体极限荷载Pne）共点的特性，可以推导出其中n=1时，d为主动土压力的作用点距墙底的高度；n=-1时，d为被动土压力的作用点；其余符号的意义同上。土体的粘聚力c=0时，土体极限荷载作用点在墙高的相对位置d/h的计算结果见表1、表2。可以发现，.在主动极限状态下，随着土体内摩擦角ϕ和墙土之间摩擦角δ的增大，土压力作用点相对位置d/h在减少；在被动极限状态下，随着土体内摩擦角ϕ和墙土之间摩擦角δ的增大，土压力作用点相对位置d/h在增大。需要说明的是，墙后土体平动破坏机制的滑动体并不是在所有情况下都满足力矩平衡的条件的，这可从式（33）看出。显然，d/h的取值要大于0小于1，否则就没有实际意义。对于一定的内摩擦角ϕ和墙土之间的摩擦角δ，式33就会得出小于0或大于1的结果。合理的解释只能是，土体发生平动破坏是有条件的，对于一定的范围的内摩擦角ϕ和墙土之间的摩擦角δ，墙后土体在临界状态时不是沿面滑动破坏的，或者对于给定的内摩擦角ϕ，墙土之间取过大的摩擦角δ是不合理的，也或者是，利用平面假定计算土体极限荷载有时是不满足力矩平衡的条件的。因而，本文舍去了d/h不合理的计算结果。当土体粘聚力c≠0时，可以计算各种情况下粘聚力c变化对d/h的影响。可以得出，随着土体的粘聚力c的增大，主动土压力作用点相对位置d/h的计算结果在减小；被动土压力作用点相对位置d/h的计算结果在增大；对于墙后土体的两种极限状态，土压力作用点相对位置d/h与土体的内摩擦角ϕ、墙土之间的摩擦角δ及土体的粘聚力c都有关系。土压力的作用点并非总是作用在墙高的1/3处。只有当墙背光滑，土体的粘聚力c=0且平动破坏的情况下，力的作用点才作用在墙高的1/3处。4粘聚力c=0时，土压力系数及kc、c、h(1)以滑动体静力平衡的力的极限平衡方程出发是可以推导出其沿平面滑动破坏的结论的，证明了Coulomb土压力理论平面滑动破坏假定的合理性。(2)当粘聚力c=0时，土压力大小的结果与Coulo