包括弯曲形变势能和中面应变势能8

1、基于Grbner基的圆薄板非线性大挠度分析Using Grbner bases to analyze the nonlinear deflection of a thin circular plate张迅炜，周建方，徐吉良 ZHANG Xun-wei, ZHOU Jian-fang, XU Ji-liang(河海大学 机电工程学院，常州 213022)摘 要：根据能量原理，用里茨法推导了等厚度轴对称圆薄板大挠度弯曲问题的变分列式，基于Grbner基法，应用计算代数系统Mathematica编制程序，分析了圆薄板受均布载荷时的最大挠度。计算结果表明，该方法求解过程易于实现程序化，计算精度较高，为

2、用变分法求解薄板大挠度问题提供了新的思路。关键词：等厚度圆薄板；非线性弯曲；里茨法；Grbner基；符号计算 中图分类号：TH114 文献标识码：A 文章编号：Doi:0 引言由于重量轻，易加工，又可承受较大的载荷，薄板结构在工程中使用广泛。为了准确地评估其工作性能，常常需要考察薄板在载荷作用下的力学行为，圆薄板非线性大挠度弯曲就是其中的一个重要问题。圆薄板大挠度弯曲问题有多种数学模型，其中卡门（Von Karman）方程组由于形式相对简单，又合理地考虑了非线性效应，因而得到广泛采纳。该方程组包括一对耦合的非线性偏微分方程，各国学者对其求解方法进行了大量研究15。Vincent1以载荷为摄动参

3、数，首先求解了均布载荷作用下的圆板大挠度弯曲问题，结果的误差较大。钱伟长2,3以中心位移为摄动参数求解了同样的问题，计算结果非常接近实验值，并且求解过程比Vincent的方法收敛更快。叶开沅等改进了摄动法的计算程序，提出了修正迭代法4，使得各阶摄动方程不再耦合，计算过程得到较大简化。上述方法均需要多次迭代求解微分方程组，以逐步提高近似解的精度，这要求工程设计人员有较高的数学分析和运算能力，不便于工程应用。变分法不直接求解卡门方程组，而是将大挠度弯曲问题转化为等价的弹性势能极值问题，通过选取适当的挠曲面函数，得到一个代数方程组，进而求出近似解，求解过程的数学分析较为简便。但是，当挠曲面函数待定系

4、数增多时，将形成多元高次非线性代数方程组。常用的数值方法6，如Newton-raphson迭代法，需要给出好的初始值才能保证非线性方程组的收敛。而解析法中的Grbner基法7，通过将一组非线性多项式方程化简为一个等价的三角化方程组，可求得封闭形式的解析解，不会产生增根或少根，能够最大限度地满足方程组的求解精度，克服了数值方法的不足。本文根据能量原理，用里茨法推出等厚度轴对称圆薄板非线性弯曲问题的变分列式，结合符号计算软件Mathematica，应用Grbner基法求解了非线性方程组，得到了圆薄板受均布载荷时的最大挠度。结果表明，采用Grbner基法求解大挠度问题，计算效率较高，且能求得高精度的

5、近似解。1 里茨法变分列式在轴对称分布载荷作用下，圆薄板大挠度弯曲时，板的形变势能包括弯曲形变势能和中面应变势能8： (1)(2) (3)其中， r为极坐标；w为板的中面挠度方程；为中面内点的径向位移；m为泊松比；弯曲刚度，式中：d为板厚，E为弹性模量。 周边固支时，板的位移边界条件及轴对称条件分别为， (4)， (5)应用里茨法求解时，要求假设的位移函数同时满足式（2）、（3），因此将中面内各点的径向与轴向位移分别取为 (6) (7)其中，、为互不依赖的待定系数；n为非负整数。根据最小势能原理，圆薄板轴对称弯曲问题的解应满足下列二组方程： (8) (9)其中，q为载荷集度；假设式（6）、（7

6、）中的参数值n，并代入式（8）、（9），得到一个关于与Cm的非线性方程组，它包括2n+1个方程。将该方程组的解代入式（6）、（7）即得到板中面各点的径向与轴向位移的近似解，据此可进一步求出应力、内力等其他未知量。随着参数n取值的增大，待定系数的数量相应增加，各级近似解将不断趋近于精确解。2 Grbner基法简介与求解步骤Grbner基法9是一种建立非线性多项式系统标准基的理论和方法。与求解线性方程组的高斯消去法思路相类似，Grbner基法将非线性系统转化为等价的三角形方程组，从而可用回代的方法求解。该方法在原非线性多项式系统所构成的多项式环内，通过对多项式变元的适当排序，求多项式系统中多项式对

7、的S-多项式，并进行彼此的消元，最终生成一个与原系统完全等价的标准基。Buchberger算法是构造多项式系统Grbner基的一种有效算法，在计算代数系统（如Mathematica等）中该算法及其改进版本均已实现，可直接以函数方式调用。本文用改进的Buchberger算法对式（8）、（9）组成的非线性方程组进行求解，步骤为：1) 假设n，得到关于与Cm的非线性方程组G；2) 将G中所有多项式赋于H； 3) 确定H中的变元的字典法排序，如；4) 用Grbner基函数将H化为最简Grbner标准基，进而得到与原多项式方程组等价的三角化方程组5) 求解其中一元多次方程得到C0，回代方程逐个求出其它变

8、元。3 算例与程序设等厚度轴对称圆板，半径为，厚度为，受均布载荷作用，周边固支，材料的弹性模量，泊松比，如图1所示。求圆板的最大挠度。图1 周边固支轴对称圆薄板由式（7），时圆板的中心挠度也就是最大挠度 (10) 设，式（6）、（7）分别表示为 (11) (12) 将式（11）、（12）代入式（1）（2）（3），得到形变势能表达式，再代入式（8）、（9）得到非线性方程组应用Grbner基法，按前述步骤求解并转换得到的等价三角化方程组为 (14)由式（14）的第一个方程式可求得C0的三个根。其中唯一的实根，据式（10），就是取时板的最大挠度再设，则 (11) (12) 类似地，求得最大挠度取不同

9、值时，求解Grbner基和最大挠度的Mathematica程序如下： (*设置参数n的值*)n=1;(*例题中的已知参数值*)Mu=3/10.;a=2;Delta=3/100;q=40000;(*求弯曲刚度，建立近似位移方程*)FR=2*1011*Delta3/(12*(1-Mu2);w=(1-Rho2/a2)2*Sumci*(1-Rho2/a2)i,i,0,n;uRho=(1-Rho/a)*Rho/a*SumAi*(Rho/a)i,i,0,n+1;(*求弯曲应变能*)VEpsilon1=Pi*FR*IntegrateRho*(Dw,Rho,2)2+1/Rho*(Dw,Rho)2+2*Mu*D

10、w,Rho,2*Dw,Rho,Rho,0,a;VEpsilon2=Pi*FR*12/Delta2*Integrate(DuRho,Rho+1/2*(Dw,Rho)2)2+uRho2/Rho2+2Mu*uRho/Rho(DuRho,Rho+(1/2*(Dw,Rho)2)*Rho),Rho,0,a;(*生成非线性方程组*)equ1=SimplifyTableDVEpsilon1+VEpsilon2,Ai,i,0,n+1; equ2Lef=SimplifyTableDVEpsilon1+VEpsilon2,ci,i,0,n;equ2Rig=SimplifyTableIntegrate2*Pi*q*(

11、1-Rho2/a2)(i+2)*Rho,Rho,0,a,i,0,n;equ2=equ2Lef-equ2Rig;equSys=Joinequ1,equ2(*设定字典序*)tab=ReverseJoinTableci,i,0,n,TableAi,i,0,n+1(*生成Grbner基*)equGro=GroebnerBasisequSys,tab(*逐个回代求解三角形方程组*)sol=SolveequGro1=0,c0;sol=c0-Selectc0/.sol,#ElementReals&1;Ifn=1,sol1=Solve(equGro2/.sol)=0,c1,sol1=c1-0;(*求最大挠度*

12、)wmax=(c0/.sol)+(c1/.sol1)/N为了检验Grbner基法以及Mathematica程序的求解精度，对算例中的q取不同值，求得最大挠度wmax，结果汇总于表1。表1 不同载荷集度q下的最大挠度比较q/ (kN/m2)wmax/ (m)本文解摄动解6n = 0n = 1100.0049890.0049800.004998200.0096310.0095720.009614300.013770.013610.01360400.017390.017120.01714600.023430.022880.02306800.028290.027470.02774计算结果表明：（1）应

13、用里兹法和Grbner基法求解等厚度圆薄板最大挠度，计算过程收敛速度快，时，与最大挠度摄动解的误差不超过2%。（2）随着参数n的增大，中面各点的位移表达式待定系数相应增加，时，最大误差减小为0.97%，最大挠度更加趋近于摄动解，且结果已能达到工程应用的精度要求。4 结语根据能量原理，本文结合里茨法和Grbner基法，应用符号计算软件Mathematica求解了等厚度圆薄板的大挠度问题。该方法原理清晰，求解步骤易于编程实现，无需应用复杂的数学分析方法和技巧，就可得到结果的解析符号表达式，或者代入参数求出高精度的数值结果，便于工程应用。本文的方法和程序可推广到其它载荷与约束条件下圆薄板的轴对称大挠度弯曲问题。参考文献：1 Vincent J. J. The bending of a thin circular plateJ. Phil.Mag., 1931, 12: 1851962 Chien Weizang. Large deflection of a circular plate under uniform pressureJ. Chinese J. Phys., 1947, 7(2):102-1133 钱伟长，叶开沅圆薄板大挠度问题J. 物理学报，1954, 10(3), 209-2384