Mathematica代数系统与岩土工程理论演绎推导

第24卷增1岩石力学与工程学报Vol.24Supp.12005年8月ChineseJournalofRockMechanicsandEngineeringAug.，2005收稿日期：2005–04–02；修回日期：2005–05–30基金项目：国家自然科学基金资助项目(50278051)；上海市重点学科建设项目作者简介：罗仁安(1948–)，男，博士，1982年毕业于武汉建材工业学院力学工程系，现任教授，主要从事应用数学与力学方面的教学与研究工作。E-mail：[e-mail]。Mathematica代数系统与岩土工程理论演绎推导罗仁安1，3，宋喜艳1，3，黄理兴2，程昌钧1，3(1.上海大学，上海200444；2.中国科学院武汉岩土力学研究所，湖北武汉430071；3.上海大学上海市应用数学和力学研究所，上海200072)摘要：研究使用第二代计算机代数系统Mathematica对岩土工程中某些理论或公式进行解析解的演绎推导，以嵌岩桩基为例，将其简化为粘弹性半空间桩土相互作用数学模型，计算机经过符号运算、数值计算、绘图和编程四大功能运算，演绎推导结果与Novak的研究结果和实际工程相吻合。关键词：岩土工程；计算机代数系统；解析解；嵌岩桩；桩土相互作用中图分类号：TU45文献标识码：A文章编号：1000–6915(2005)增1–4688–06MATHEMATICAALGEBRASYSTEMANDDEDUCINGTHEORYOFGEOTECHNICALENGINEERINGLUORen-an1，3，SONGXi-yan1，3，HUANGLi-xing2，CHENGChang-jun1，3(1.ShanghaiUniversity，Shanghai200444，China；2.InstituteofRockandSoilMechanics，ChineseAcademyofSciences，Wuhan430071，China；3.ShanghaiInstituteofAppliedMathematicsandMechanics，ShanghaiUniversity，Shanghai200072，China)Abstract：ThesecondalgebrasystemMathematicaisusedtodeduceanalyticalsolutioninsometheoriesofrockengineeringinthispaper.Themethod，forexample，takesthesoilandthepileasawholesystem(visco-elastichalf-spaceembeddedpile).Usingthemethodofthecontinuummechanics，thenonlinearvisco-elasticpile-soilmodeisestablished.AssumingthatthepileisaEulerbeamwithcircularcross-section，thevisco-elasticmaterialandthesoilarecomposedoftwoparts：anouterinfinitelinearvisco-elasticzoneandaninnernonlinearvisco-elasticzone(hollowcylinder).Inthethree-dimensionalcylindricalcoordinates，theboundaryvalueproblemandjoiningconditionarepresented.ThecomputerwiththesystemMathematicagoesthroughthefourmainprocessesofsymbolicoperation，numericalcalculation，plottingandprogramming.TheresultsofthismethodagreewellwithNovak′sresults.Keywords：geotechnicalengineering；computeralgebrasystem；analyticalsolution；embeddedpile；interactionofpile-soil1引言计算机自问世以来已被成功地用于科学计算、文字和管理工作。迄今为止，计算机除了编制程序进行了一些工程数值计算外，还不能直接进行理论研究及其推导工作。因此，人们一直盼望有一个可以进行理论演绎推导的计算机系统。现在，一种被称为“计算机代数”的工具在符号计算中已经被广泛应用。Mathematica正在应用于各研究领域。文[1]用Mathematica对考虑集中载荷和带阻尼的梁进行动第24卷增1罗仁安等.Mathematica代数系统与岩土工程理论演绎推导•4689•力分析。文[2]用Mathematica模拟机械系统从而进行定性分析，解决了刚性机械和各类电路中的许多问题。本文企图运用计算机代数方法研究岩土工程理论，例如，将嵌岩桩简化为粘弹性半空间桩土相互作用的数学模型，应用Mathematica软件进行进行数学公式演绎和解析理论推导，分析嵌岩桩的非线性力学行为，并把研究结果与工程实际相比较。2Mathematica简述Mathematica第2代计算机代数系统，是一个符号计算系统，是美国WolframResearch公司开发的数学软件[3]。它可以用于解决各种领域中涉及复杂的符号运算和数值推导问题，代替许多人工演绎和推导工作，这种思维和解题工具的革新将对理论研究领域和工程领域产生深远的影响。Mathematica包括符号运算、数值计算、绘图和编程四大功能，更重要的是它把这些功能有机地结合一个系统。在使用本系统时，研究人员可以根据需要，在符号演算、图形绘制和数值计算之间灵活方便的转换，使一些非常复杂的问题变得简单。表1列出几种常用的计算机代数系统。其中，Mathematica以符号运算为主，建立了输入各种数学符号和函数的专用模板，输入数学公式和各种操作命令变得极为简便直观，而且允许用户自制模板，在输入方面大大优于Matlab和Mathcad。Mathematica的绘图功能比Matlab和Mathcad使用简单，且具有相当强的绘图功能。另外，它带有扩展的绘图软件包，使得绘图功能更加完善。在数值计算方面，Mathematica注重保证计算精度。此外，Mathematica吸取了不同类型软件的优点：具有BASIC语言简单易学的交互式操作方式；有Macsyma，Maple，Reduce和SMP的符号计算功能；有APL和LISP表1常用计算机代数系统Table1CommonComputerAlgebraSystems计算机代数系统功能介绍Mathematica功能齐全、图形功能强，可用C语言和FORTRAN语言输出结果、可用于微机Maple功能齐全、文本处理功能灵活、可用于微机Reduce代换功能丰富、用LISP语言，无绘图功能Derive小巧灵活、有绘图功能，但接口不佳，适宜于解决简单问题Macsyma规模昀大、功能昀齐全；但机种要求高人工智能及列表处理能力；有像C语言与PASCAL语言一样的结构化程序设计等等。3用Mathematica推导岩土理论举例3.1嵌岩桩数理模型在粘弹性半空间，假设嵌岩桩桩身由线性粘弹性材料组成，其半径为0R；桩周土由内层土域和外层土域构成粘弹性柱体，0R和iR分别为内层土域(空心圆柱体)的内外半径；同时iR也是外层土域的内半径，其外半径为无限大[4]。在空间柱()zr，，θ坐标系中，其中坐标原点与桩底重合，oz轴沿桩基的轴线，桩基的长度H。如图1所示。图1数理模型Fig.1Physicalmathematicalmodel3.2内层土域动力方程根据弹性理论[5]，内层土柱坐标下的动力方程为⎪⎪⎭⎪⎪⎬⎫∂∂=∂∂+∂∂+∂∂+∂∂=∂∂+∂∂−∂∂+∗∗∗∗∗∗∗∗22n22nnnn22n22nnnn21)2(12)2(tuzuGrwGΔrGtuzuGwrGrΔGiizirirzθθρθλρθλ(1)式中：iiruuθ，分别为内层土在θ，r方向的位移；∗∗nnλ，G为复拉梅常数；θθ∂∂+∂∂=iirurrurrΔ1)(1，⎥⎦⎤⎢⎣⎡∂∂−∂∂=θθirizurrurrw1)(121。当土为各向同性时，)21/(2)i21(nnnnnnvvGG−=+=∗∗λξ，。求解式(1)时要引进势函数iiΨΦ，，则iiruuθ，可简化为tiitirirΨrrΦuuωωθiie1e⎟⎟⎠⎞⎜⎜⎝⎛∂∂+∂∂==(2a)tiitiirΨΦruuωωθθθiie1e⎟⎟⎠⎞⎜⎜⎝⎛∂∂−∂∂==(2b)对势函数iiΨΦ，进行分离变量可得z土桩y外域土xN内域土R0RioxoH•4690•岩石力学与工程学报2005年)()()(zZΘrRΦiθ=(3a))()()(111zZΘrRΨiθ=(3b)将表达式(3a)，(3b)代入式(1)进行求解可得到势函数iiΨΦ，的通解：++=)sin()][()([2m1m1hzAqrIBqrKAΦi)sin()][cos(32θmAhzB)]cos(3θmB+(4a)++=)sin()][()([5m4m4hzAsrIBsrKAΨi)sin()][cos(65θmAhzB)]cos(6θmB+(4b)式中：mI(·)和mK(·)分别为一阶和二阶的修正贝塞尔函数。3.3外层土域动力方程外层土柱坐标下的动力方程为⎪⎪⎭⎪⎪⎬⎫∂∂=∂∂+∂∂−∂∂+∂∂=∂∂+∂∂−∂∂+22n22nnnn22n22nnnn121)2(12)2(tuZuGrrGrΔrGtuZuGrGrΔGoozororzθθρωλρθωλ(5)式中：oru和ouθ分别是内层土在θ，r方向的位移；nn212Gννλ−=，)i21(nξ+=GG为外层土的复剪切模量；ν为外层土的泊松比；=Δ+∂∂)(1orrurrθθ∂∂our1；⎥⎦⎤⎢⎣⎡∂∂−∂∂=θθorozurrurrw1)(121。同理引入势函数ooΨΦ，，并求得其通解如下：++=)sin()][()([211hzCqrIDqrKCΦmmo)sin()][cos(32θmChzD)]cos(3θmD+(6a)++=)sin()][()([544hzCsrIDsrKCΨmmo)sin()][cos(65θmChzD)cos(6θmD+(6b)3.4Mathematica求解势函数3.4.1计算机演绎推导势函数在推导过程中，势函数求导后有一个很长很繁琐的方程，因此可以用Simplify[·]来化简。Simplify[·]能将表达式化简成含项数昀少的昀简形式，化简后的表达式可以立即得到某些信息参数。在表达式极其烦琐的情况下，还可以用FullSimplify[·]来化简。FullSimplify[·]与Simplify[·]的功能相似，可对某些特殊函数进行化简，其功能更强，但它对内存的要求也更高，化简的时间往往要比后者多很多。上述工作假如用手工进行推导，不仅需要很耐心、很细心，要有很好的数学基础，而且涉及到对贝塞尔函数、求导等概念熟悉掌握，即使完成了也可能会出错，耗费很长的时间是勿庸置疑的。本文以内层土的求解为例，演绎计算机推导势函数方法[6]。求解时先输入势函数的通解，Mathematica自身定义了很多函数，贝塞尔函数也是其中之一，它不仅定义了一阶和二阶贝塞尔函数，还定义了一阶和二阶的修正贝塞尔函数。程序中的BesselK[m，qr]就代表二阶的修正贝塞尔函数mK(·)，在程序运行过程中，计算机将BesselK[m，qr]按照二阶修正贝塞尔函数的特性进行处理。3.4.2内层土势函数演绎推导求解内层土势函数待定系数时，用到如下边界条件和连接条件。桩土连接条件。假设内层土域与桩之间不脱离，当0Rr=时，有)(sincospzuuuuiirx=−=θθθ(7a)0cossin=+=θθθii