【doc】求解几何非线性桩-土耦合系统的微分求积单元法

【doc】求解几何非线性桩-土耦合系统的微分求积法 求解几何非线性桩-土耦合系统的微分求积 单元法 第29卷第2期 2008年6月 固体力学 CHINESEJOURNAL0FSOLIDMECHANICS Vo1.29NO.2 June2OO8 求解几何非线性桩一土耦合系统的微分求积单元法 胡育佳朱媛媛.程昌钧 (上海大学,上海市应用数学和力学研究所,力学系,上海,200072)(.上海师范大学,计算机科学与技术系,上海,200234) 摘要将桩一土系统看成在土层中嵌入了一根等圆截面桩的空间轴对称弹性体,在几何非线性的条件下建 立了具有间断性条件的桩一土系统的非线性控制方程,并运用微分求积方法(DQEM)来求解了该问题.提出了利用 DQEM求解非线性空间轴对称问题中处理单元之间连接条件(包括间断性条件)及边界条件的离散化方法,最终得 到了一组离散化的非线性DQEM代数方程,运用Newton-Raphson迭代方法求解非线性代数方程组可以得到每个 节点处的位移,进一步可以得到系统的应力和应变.给出了两个数值算例,并与有限元解进行了比较,它们是非常 吻合的.将看到,由于在采用DQEM求解时只布置了较少的节点,因此,该文方法具有较小的计算工作量,较高的 精度,良好的收敛性以及应用广泛等优点.该文提出的处理连接条件的方法是一个一般的方法,由于它在数学上遵 循了求解边值问题的思路,因此,数学上也是严谨的. 关键词桩一土耦合系统,几何非线性,间断性条件,微分求积方法 0引言 Bellman和Casti?1于1971—1972年提出了 微分求积方法(differentialquadraturemethod即 DQM),并已成功地应用于许多领域.许多研究者对 DQM的应用和发展做出了重要的贡献?4].对规则 域问题和分布式参数系统,DQM的工作量远少于 差分方法,有限元方法和边界元方法,因此,Bert等 人[3认为DQM可望发展成与差分方法,有限元方 法一样的求解微分方程的强有力的数值方法.但对 于不规则的区域和具有材料和几何不连续或分界界 面的问题,传统的DQM在应用上受到比较大的限 制.仲政课题组采用小波和微分求积混合法研究了 在集中载荷作用下结构的响应?;Sun和Zhu嘲在 对驱动方框粘弹性流体的流动问题和入口流动问题 的计算时,对微分求积法进行了一系列的改造,提出 了一种具有迎风机制的混合型微分求积法和能处理 非区域问题的局部微分求积法.Liew,Liu, Chen,Striz等...通过区域分解技术将DQM推广 为微分求积单元法(differentialquadratureelement method即DQEM),并在一定程度上处理了具有问 断性条件的结构力学问题,拓宽了DQM的运用.但 是,在现有DQM和DQEM应用中仍然有许多问题 难以解决,例如本文中涉及的空间轴对称桩一土系统 问题,它不仅具有奇异性(在对称轴上),而且具有各 种间断性,例如材料间断性,载荷间断性和几何间断 性,还具有非线性等因素. 桩基是建筑物的一种重要的基础形式,但是,由 于桩一土之间的相互作用,载荷的传递方式等非常复 杂,本质上是一个非线性系统,建立比较合理的桩, 土一桩系统的数学模型是系统力学行为分析的一个 关键问题,对于这样的系统,得到解析解一般也是非 常困难的,通常只能采用数值方法来求数值解.虽 然,有限元方法仍然是求解这种复杂系统的主要数 值方法,但是对于桩一土一桩这类复杂系统采用有限 元方法来进行数值计算的工作量非常大,特别是非 线性动力学分析.因此,寻求其他的高精度的数值方 法也显得非常必要. 本文的目的是:(1)在几何非线性的条件下,利 用弹性力学的理论建立桩一土系统力学行为分析的 数学模型,其中,具有半径r的圆截面桩基看成是 被嵌入一个无限大的圆柱形弹性介质中,该弹性介 质的半径为r..,在桩一土系统的连接面上位移和应 力满足适当的连接条件.因此,系统是一个具有界面 连接条件和其他间断性条件的空间轴对称问题.(2) *上海市浦江人才(07pj14073),市重点学科建设项目(Y0103)和上海高校选拔 培养优秀青年教师科研专项基金资助 2007—04—26收到第1稿,2008—05—12收到修改稿. **通讯作者.E,mail:chjcheng@mail.shu.edu.cn. ? 142?固体力学 提出了求解具有间断性条件的空间轴对称桩一土系 统的控制方程的DQEM;特别是提出了在使用 DQEM求解空间轴对称问题时处理单元之间的连 接条件,对称轴处奇性条件和边界条件的离散化方 法,并得到了非线性控制方程的DQEM离散化形 式,并用Newton—Raphson迭代方法求解了离散化 非线性代数方程组.(3)给出了两个数值算例,其一 是作为模型的退化,考察了各向同性介质中的空间 轴对称问题的DQEM解,并与有限元方法的结果进 行比较.另一个算例是研究了几何非线性条件下,层 状介质中桩土系统问题的DQEM解,并由此推广 DQEM在岩土力学这种具有大范围计算区域问题 中的应用. 1问题的数学描述 1.1物理模型和基本假设 采用空间柱坐标系(r,,),将桩一土系统看成 在土层中嵌入了桩基的空间轴对称弹性体(图1), 其中桩基看成是由线性弹性材料组成的圆柱体,其 半径为r,桩长为h,剪切模量为G,泊松比为, 桩基顶面受到分布力q(r)的作用.同时将土看成是 弹性空间轴对称圆筒体,其内外半径分别为r,r.., 弹性模量为G,泊松比为土层厚度为h.还假设 桩一土之间没有剥离. 图1桩一土系统物理模型 Fig.1Physicalmodelofpile-soilsystem 1.2桩一土系统的控制方程 1.2.1基本方程 对于空间轴对称问题,在空间柱坐标系(r,,) 中,设,7.22是径向和轴向位移,,,,,,, '和',',,r三'分别为相应的非零应 变分量和应力分量,其中上标()一(p)表示与桩有 关的量,()一(s)表示与土有关的量. 引入无量纲变量和参数如下: IR一,z一,H一一,H一,R一 lR..一一一警 一 ,w一一 一 ,一,Q一号 其中,G为各项同性弹性土体的剪切模量. 根据连续介质力学,可得到如下的微分方程,连 接条件和边界条件. 空间轴对称问题的平衡微分方程: f3兰R+:: R 二::+生 3 翌Z一.l.. }+警+一. 几何关系:在考虑几何非线性的条件下,有下面 的几.何关系. 一 +().+丢(). 一 +丢(). 一 +丢()+丢(). lv()一++一1+—3W—(o)—3W—(o).3Z.3R.23R3Z.23R3Z f一2G'dc'口'[(1一'),+'(,+,)] {鬟:22GGCc三二:三;:;c4I一'[(1一),+(,+,)] lr一Grl: 2c(+去一)++ rc4-f+:::+ .3RaZ.3R\3Z..R. 去+cc)一+ 丽~2W(oO+ 第2期胡育佳等:求解几何非线性桩土耦合系统的微分求积单元法?143? (1TC(一.(5a) (+去)+(+去)+ z(O2W(a>+(去+ )+(2C?c1--1J(~>++ 20a'1j(a)U',lOW(,1OW'}a.W}aW',【 ROZ\ROR.aZ2.OR/. C?—OW— (o)塑一0(5b) OROROZ 1.2.2桩一土系统的边界条件 假设土层的底部固定,桩头受到竖直向下的静 力Q(R)的作用,土表不受力,并R..在处被固定.由 于对称性,在桩基的轴线上,径向位移U—o,剪力 r一o.于是有如下桩和土体的对称性条件和边界 条件. 对称轴上的对称性条件: R一0,H一H.?Z?H:U一0,r一0 (6a) R一0,0?Z<H一H.:U"一0,r===0 (6b) 底端固定的边界条件: 0?R?R..,Z一0:U''一0,W"===0(6c) 土体侧表面处固定的边界条件: R—R..,0?Z?H:U"一0,W"'一0(6d) 顶端边界条件: 0?R?R.,Z—H:一一Q(R),r'一0 (6e) R.<R?R..,Z—H:一0,r窭'一0(6f) 1.2.3桩一土系统的连接条件 在桩一土系统的连接面上位移和应力应该满足 连续性条件,即在R—R.,H一H.?Z?H上,有 连续性条件: U'一L,",W一w",一,一r (7a) 同时在0?R?R.,Z—H一H.上,有连续性条 件: Up'一U"',W'声一W",一,r一r譬 (7b) 因此,问题归结为求满足平衡微分方程方程 (5),边界条件(6)和连续性条件(7)的解.可见,得到 问题的解析解是非常困难的.本文中将建立相应的 DQEM,并用它来进行求解.将看到,这种方法既保 留了DQM计算的高精度,高效率,又成功地解决了 传统DQM不能解决的处理具有间断性条件的桩一 土系统问题. 2控制方程的DQEM离散化形式与连接条 件的处理 对于具有不连续性条件(包括材料的不连续性 等)的非线性结构力学问题,采用传统的DQM求解 会遇到困难.本文中将采用DQEM来求解,其实质是在DQM的基础上结合单元分割近似原理而形成 的一种数值计算方法,它首先将求解区域分割成若 干个子区域或单元,这些单元在节点处进行连接,并 建立各个子区域的微分方程,利用DQM将这些微 分方程离散为代数方程.然后根据连接条件,将所有 单元的离散化的代数方程集合起来,形成一组整体 代数方程,从而求解得各节点的未知量.这种方法既 保持了DQM精度高,收敛性好,计算效率高的特 性,又大大拓宽了DQM的运用领域,特别是对于非 规则区域和具有各种不连续性条件的情况.因而, DQEM的求解关键在于处理单元的连接条件和边 界条件.本文中将提出一种数学严谨,操作方便并且 精度高的方法来处理单元的连接条件和边界条件. 2.1单元的划分与节点布置 设将系统所占有的空间轴对称区域划分为K 个子单元,为了便于说明问题,现取K一4(图2).对 整个结构进行统一布点和编号,在每个单元内部的 和Z一方向采用非等间距的切比雪夫零点的配点 图2单元划分 Fig.2Elementdivision ? 144?固体力学 方法,并在方向和Z方向分别布置N,和N个 节点,其中,N一1-+-3—1,N一2-+-4—1,1,3 为R_方向的单元节点数,.,为Z_方向的单元节 点数(图2).可以看到,这样的布点方法保证了接点 处位移的连续性条件. 2.2控制方程的DQEM离散化形式 从图2中取出任意一个单元i来考虑,在此单 元内,R方向和Z_方向分布的节点数分别为和 .个,对于如图2所示的单元1,,—,.一.,对 于单元2,—723,一2,对于单元3,,—721,.一 ,对于单元4,一.,一.虽然,由于单元i满 足空间轴对称的控制方程(5ab),可以直接运用 DQM离散方程(5a_h),但是由于连接条件(7)和边 界条件(6),这种处理是十分复杂的,本文采用下面 的方法来进行离散. 根据DQM,在单元i内部对(2)进行离散, 可以得到如下的DQM离散化形式: n r n z f?'z(切+?鼢(r)+l一1=1 I二::一n R} IA()+B;()}+芸一o (8a) 其中,上标i表示与单元有关的量.一2,3,…,,一 1,卵一2,3,…,一1分别为在方向和Z_方向分 布的节点数.A"',Bn"'分别为单元i在R和Z_方 向的一阶权系数矩阵. 由式(3)(4)可以得到单元i内部的应力为: (勤一2G"C"[(1一")(,抑-+- "((,;.)勋-+-(,勋)] (.)勃一2G"C"[(1一")(,)新-+- "((,)白+(,勋)](8b) ()勋一2G"C"[(1一")(,)勋+ "((,)新+(,白)] (r')白一G"() 此处,一1,2,…,,叩一1,2,…,,G",C"'一1/(1 — 2v"')分别为单元i的无量刚剪切模量和系数.注 意到(3),有下面的几何非线性方程的DQM离散化 形式: (,)勖一?A(u?)+ 『(?A(u?))+(?A.(w?)).] (,如一?鼢(w?)+ 『(?B(w?)).-+-(?B()).] (,5勖一+1(). (r)如一?B(u?)-+-?A.(w?)+ 告『?A'(u?)×?B(u")+一 1一1 r ?A(w?)×?B(?)] (8c) 需要指出,当i单元中的点的位于对称轴R一0上 时,由于当R—O时是奇异的,这种奇异性是由坐 标变换引起的,可以采用如下的近似的表达式: 一 U 十 1 IU)?d3Uq_1\d3U),R---~O一十IJ?d\dJ, (8d) 并可得到(8d)的离散化公式: (,5)一?A'(?)-+- (?A(u?)).,当R一0(8e) 其中,卵一1,2,…,.这样我们可以得到单元i内部 节点的离散化方程(8a—e). 2.3边界条件的处理 根据边界条件(6)可得到它们的离散形式. (1)对称轴(R一0,O?Z?H)上对称性条件的 离散方程为: U' r U' f2 7—0,一1,2, 一0,卯一2,3, 7—0,一1,2, 一O,?一2,3, 其中,图2中单元1与单元3的连接a处满足的关 系,可以根据力的平衡条件(图3)得到: (r)11-+-(r)1一0(9b) (2)底端(Z一0,O?R?R..)处的边界条件的离 ) a 9 ll 一, ??44 ,,,, 第2期胡育佳等:求解几何非线性桩一土耦合系统的微分求积单元法 散方程为: 1厂] (r)11IlI单元1l i]IU单元3 图3处受力 Fig.3Forcesatdl U.')71—0, W.')71—0, U)71—0, W')71—0, (3)侧表面处(R—R.. 件离散方程为: U'')3 一O, W') . —O, U'). 一O, W)3一O, 7/=2,3,…,nl 叩一1,2,…, (9c) "一2,3,…,n3 田2,3,…,n. ,O?Z?H)处的边界条 (4)顶端(O?R?R..,z—H)处边界条件离散 方程为: f()}扎,一2G''C'?[(1一u)(s")2+ I'"((s;")+(,:"))]+Q—o I(r)一Gn(y)一o,一2,3,…,nl一1 {()弘.=2G?C趵[(1一')(,)+ {((s;')f.+(,)f)]一o I(r)一G(y)一0,一2,3,…,一1 (9e) 由力的平衡条件,单元1与2的连接处z的边界条 件为(如图4): f?_厂r一(r)+(r).+(").一 {(:).一0(gD j?一()+(')+(r').一 I(r)十Q(R)一0 将(8b—c)代人(9a—f)即可得到边界条件的离散化形 式. 2.4连接条件的处理 (1)单元1和单元2连接处(R—R,H一H <z<H)连接条件的处理:由两单元界面处的 方向的正应力和剪应力相等可以得到: (r') 单元1 图4 Fig.4 (r)In2 单元2 7/=2,3,…,n2--1(1Oa) 同理,可以得到单元3和单元4的连接处(R—R,0 <Z<H一H)连接条件为: 田一2,3,…,n4—1(1Ob) (2)单元1和单元3连接处(O<R<Rp,Z—H 一 H)连接条件的处理:由两单元界面处的方向 的正应力和剪应力相等可以得到: 一 2,3,…,l,1(10c) 同理,可以得到单元2和单元4的连接处(R<R< R..,z—H一H)连接条件为: 'f一' ' 一2,3,…,n.一1(1od)L一厶,J,.,,'31u l(r)fl:(Z-(4一), (3)单元1,2,3和4的连接OL3处的连接条件为 (图5): f?_厂r一()+(;.')一(), }()一(r)一(r)+(r)+ I(r)In—a0 1?一(.)+()一(")一 l()+(r)一(f)+(r)一 l(r2')一0 (1Oe) 将(Sbc)代入(1Oa-e)即可得到连接条件的离散化形式. 这样,微分方程,边界条件和连接条件的离散化 形式(8)一(10)构成了求解在几何非线性条件下空问 轴对称桩一土系统的DQEM离散化的非线性代数方 程组. 2.5DQEM离散化方程组的求解 运用Newton—Raphson迭代方法直接求解离散 化的非线性方程组(8)一(10),可以得到结构在每个 节点处的无量纲位移u勤,w勤.计算中,取迭代初值 ll \/\/ 盯r /,一一 ,/,/ ?" }r盯r /, ,?l ,, \/,/ 44 //一一 n \/,/ 33 // ,??,,? l 44 \,\ ?心 盯r //一一 l1 }f \,\ n船盯r // ,??l ,, d 9 ,'\ ll一一 2244 ,?,?,, ???? ???? ,,,, 2233 ,,,, l122 一一一一 7 ? 146?固体力学 单元1 () 图5a处的受力 Fig.5Forcesata3 单元2 ')11 ( 单元4 R 图7z—H处的R一和方向的位移 Fig.7DisplacementsintheR—andZ-directionsatZ—H 为U勤一W勖一0,其中一1,2,…,N,叩一1,2,…, N;并取误差精度,一1×10.表1各位置处方向位移 Table1DisplacementsintheZ-directionatthedifferentpostion 3数值算例与结论 3.1各向同性介质空间轴对称问题的DQEM解 为了说明方法的正确性,考虑小变形问题,并且 令G"===G===1(即桩和土的介质相同),并取H=== 1,H===0.5,R..一1,R===0.5,Q===0.1,"一'==: 0.3.分别采用DQEM和有限元方法(FEM)对该问 题进行求解,其中DQEM采用4个单元,每个单元 按照切比雪夫规则非均匀分布9×9个节点,FEM 采用(4O×4O)个单元,如图6所示.土体表面处的 R方向和Z_方向的位移如图7所示.可以看到,它 们是完全一致的.表1给出了不同位置处位移的 DQEM解与FEM解的比较,可以看到,在应用 DQEM时,6×6个节点实际上就可以达到精度要 求. 图6平面单元划分及布点 Fig.6Elementdivisionandnodescollocated 3.2几何非线性条件下层状土中桩一土系统问题的 DQEM解 运用上面的理论和方法,也可以求解几何非线 性条件下层状土中桩一土系统问题.其中,土被分为 3层:土层(1)的厚度h一0.4,土层(2)的厚度h.一 0.4,土层(3)的厚度h.===0.2(图8).将桩一土系统划 分为6单元(图8),其中每个单元分别布置12×12, 21×12,12×12,21×12,12×12,21×12个节点. :G(?.(2): {一?一}雹 :G'?,p(4: 一 卜t. iG(6,,p(8i :土层1:: :土层2:: {-土~层tt3i冀R 图8分层土的单元划分及布点 Fig.8Elementdivisionandnodescollocatedinlayeredsoil 并令G"一G一1000,G'===1,G一2,G ==G'===3,H一1,Hp一0.8,Rp一0.05,R..一2,Q一 第2期胡育佳等:求解几何非线性桩一土耦合系统的微分求积单元法?147? 1,'='.=0. 25,'='='一'.一0. 3.桩头 (Z—H===1)和桩底(Z—O.2)处方向和方向的 位移如图9所示.图1O和图11给出了桩一土系统的 方向的等位移线和等剪应力线. 图9桩头和桩底处的位移 Fig.9Displacementsatthetopandbottomofpile R 图1O桩一土系统方向等位移线 Fig.10Equi—displacementlinesintheZ-direction ofpilesoilsystem 图11桩一土系统的等剪应力线 Fig.11Equi-shearforcelinesofpile-soilsystem 4结论 本文建立了几何非线性条件下,空间轴对称耦 合桩一土系统的控制方程,并发展了求解空间轴对称 桩一土系统的DQEM,提出了在应用DQEM处理空 间轴对称问题单元之间的连接条件和边界条件的离 散化方法,并最终得到了桩一土系统非线性控制方程 的DQEM离散化形式,在此基础上利用Newton- Raphson迭代方法直接求解离散化的非线性方程 组,可得到结构在每个节点处的无量纲位移,并给出 了数值算利,与有限元单元的数值结果进行了比较. 可以看到,DQEM解和FEM解是非常一致的,说明 本文的方法是正确的.但是由于本文计算中所布置 的节点比较少,因此,计算工作量少.不但保留了 DQM计算的精度高,收敛性好的特点,同时也大大 拓宽了DQM的求解领域,使这种方法能够适合求 解具有各种间断性条件的问题.特别是由于DQM/ 能用很少的点得到令人满意的结果,因而它非常适 合岩土力学这种具有大范围计算区域的问题.还需 要说明的是,本文中提出的处理间断性条件的方法, 在数学上也是严谨的,它完全遵循求解具有间断性 条件的边值问题的原则和方法,而且编程方便,操纵 简单. [4] 参考文献 BellmamRE,CastiJ.Differentialquadratureand longtermintegrationEJ].JournalofMathematical AnalysisandApplications,1971,34:235—238. BellmanRE,KashefBG,CastiJ.Differentialquad— rature:Atechniquefortherapidsolutionofnonlin— earpartialdifferentialequations[J].JournalofCorn— putationalPhysics,1972,10:40—52. BertCW,MalikM.Differentialquadraturemethod incomputationalmechanics:areview[J].Applied MechanicsReviews,1996,49:1—28. 张纯,胡振东,仲政.集中载荷作用下梯度复合材料 梁的小波一微分求积法分析[J].固体力学, 2006,27(增刊):38—42.(ZhangChun,HuZhen— dong,ZhongZheng.Wavelet—baseddifferential quadraturemethodfortheanalysisofgradedcorn— positebeamunderconcentratedloading[J].Acta MechanicaSolinaSinica,2006(S.Issue):38—42.(in ?148?固体力学2008焦 Fsl [6] Chinese)) SunJian-An,ZhuZheng-You.Applicationofdiffer— entialquadraturemethodtosolveentryflowofvisco— elasticsecond—orderfluidEJ],InternationalJournal forNumericalMethodsinFluids,1999,30(8): 1109—1117. StrizAG,WangX,BertCW.Harmonicdifferential methodandapplicationstostructurecomponents[J]. [71 F81 AetaMechanica,1995,I11:85—94. LiuFL,LiewKM.Vibrationanalysisofdiscontinu— OUSMindlinplatesbydifferentialquadratureelement method[-J].VibrationandAcoustics,1999,121: 204—2O8. ChengChang-jun,ZhuYuan—yuan.TheoryofElas— ticity(RevisedEdition)FM].Shanghai:ShanghaiU— niversityPress,2005. DQEMFORSOLVINGPILE—SOILCOUPLINGSYSTEMSWITH GEoMETRICALNoNLINEARITY YujiaHuYuanyuanZhuChangjunCheng (Shn"g矗 ni[nstit"te0fAppliedMathematicsandMechanics,DepartmentofMechanics,ShanghaiUniversity,Shanghai,200072) (.DepartmentofComputerScienceandTechnology,ShanghaiNormalUniversity,Shanghai,200234) AbstractApile— soilsystemisregardedasanaxisymmetricspaceelasticbodyinwhichapilewithcir— cu1arcross—sectionisembedded.Anonlinearmathematicalmodelofthepile— soilsystemwithdiscontinuity conditionsisestablishedundertheconditionofgeometricalnonlinearity.Thedifferentialquadratureele— mentmethod(DQEM)isappliedtosolvetheproblem.Adiscretizationmethodispresentedtodealwith theiunctionconditions(includeddiscontinuityconditions)attheinterfacebetweenthepileandthesoilas wellastheboundaryconditionsintheapplicationoftheDQEMtosolvingnonlinearaxisymmetricprob lems.AsetofDQEMdiscretizationequationsareobtained.TheNewton-Raphsonmethodisusedtosolve thesystemofdiscretizednonlinearalgebraicequationsandthenodaldisplacementsareobtained,further thestressesandthestrainsofthesystemcanbeyielded.Twonumericalexamplesarepresented.Theob— tainedresultsarecomparedwiththoseobtainedbyFEMandtheyarecomparativelyaccordant.Duetofe wetnodesapplied,themethodpresentedinthispaperiswiththeadvantagesoflittleamountincomputa— tion,higherprecision,betterstabilityandconvergence,broaderapplicationandSOon,comparedtothe othercomparationtechniques.Atthesametime,Themethodfordealingwiththejunctionconditions presentedinthispaperisalsoageneralmethod,anditfollowsthetheoryandtheprinciplesolvingthe boundaryvalueproblemwithdiscontinuityconditions. Keywordspilesoilcouplingsystem,geometricalnonlinearity,discontinuitycondition,differential quadratureelementmethod(DQEM)