矩阵法求解压杆稳定的有限差分欧拉方程

胡章咏

(黄冈师范学院 机电工程学院，湖北 黄州 438000)



矩阵法求解压杆稳定的有限差分欧拉方程

胡章咏

(黄冈师范学院 机电工程学院，湖北 黄州 438000)

将压杆稳定的欧拉方程用差分的形式表达后，根据压杆稳定的实际要求，通过矩阵的方法求解该差分方程的系数矩阵，最后求出简支梁等截面直杆的压杆稳定临界应力。将求出的结果与解析解进行对比并进行误差分析，找出产生误差的原因和减小误差的措施，并指出该方法在解决变截面杆的实际意义。

压杆稳定；有限差分；欧拉方程；矩阵

承压杆件的稳定性问题一直是材料力学当中的三大问题之一(材料力学研究杆件的三大问题为：杆件的强度、刚度和稳定性)，相比较于杆件的强度和刚度问题，杆件的稳定性在计算和模拟方面有着一定的难度，特别是对于变截面杆件的稳定性问题，其临界应力的计算目前用解析法还是比较困难。

目前，解决压杆稳定问题的主要方法有：针对直杆定截面的压杆稳定问题主要用欧拉方程来求解，也有用能量法和图乘法来进行求解[1]；针对变截面杆，求解的方法有很多，主要有能量法、数值法等[2-3]。

1　问题的引入

如图1所示，细长压杆两端为球铰支座，轴线为直线，压力P与轴线重合。当压力达到临界值时，压杆将由直线平衡形态转变为曲线平衡形态。选取坐标系如图1所示，距离原点为x的任意截面的挠度为v，弯矩M的绝对值为Pv，对于微小的弯曲变形挠曲线的近似方程为：

图1　简支梁力学模型

(1)

引入记号：

(2)

于是式(1)可以写成[1]：

(3)

式(2)的解析解法在参考文献[4]中有详细讨论，本文不再赘述。

2　建立差分方程并确定求解方法

式(2)的差分离散化问题做如下讨论：

首先将梁等分成等份(如图2所示)，

图2　简支梁差分离散模型

v在x=xi处的一阶和二阶倒数可以表达为差分公式[5]：

(4)

(5)

将式(5)代入式(3)中,有差分方程：

(6)

不妨设将细杆分成n等分，x对应的坐标分别为[x0,x1,x2,…,xn-1,xn]，v对应xi的坐标分别为：[v0,v1,v2,…,vn-1,vn]，且令：

(7)

则式(6)可以分别离散成为下面的方程组：

(8)

又因为根据杆件的边界条件：x0=0时，v0=0；xn=l时，vn=0；

所以式(8)变成了以下的一个n元一次其次方程组：

(9)

只有当[v0,v1,v2,…,vn-1,vn]有异于0的解时，才是失稳的情况，因此就要求式(8)系数行列式等于0，即：

(10)

将式(10)解得的a的值反代入式(7)和式(2)中，求出P：

(11)

但由于式(10)是一个一元高次方程，由此可以得出a的解时多个数值，根据该问题的实际意义，即在多个与a值对应的P值中，使杆件保持微小弯曲的最小压力，才是临界应力Pcr。

3　计算结果及误差讨论

表1　计算结果数据表

从表1和图3可以看出，当n≤10之前，Pcr的收敛速度很快，误差可以控制在1.0%左右，但是随着n值的增大，其收敛速度明显放慢，而且还会出现数值上的波动。其主要原因是因为式(10)展开后得如下所示方程：

(12)

通过式(12)知：当n增大时，对a值的求解精度越来越高，但是当n>15时会出现误差反而略有增大的现象(如当n=15时误差为0.37%，但是当n=16时误差反而增大到了0.39%)，其主要的原因是因为当n增大时，式(12)是一个高阶的多项式方程，所以当多项式的阶数增大时会出现高阶震荡现象。不过从计算的结果上来看，当n=15时，0.37%的误差对于工程上的精度书完全能够满足，因此不建议求解n>15时的结果。由计算的数值说明该计算方法具有很高的精度，因此该方法具有很好的实际应用价值。

图3　n值与对应的a值和Pcr关系图

4　对于轴线为直线的变截面压杆的讨论

以变截面圆形压杆为例略作讨论，设其截面的直径函数为：d=φ(x)，所以在xi处的惯性矩为：

(13)

将式(13)代入式(2)中，则

(14)

式(6)则变为如下：

(15)

式(7)变为：

(16)

将式(16)代入式(10)中，式(10)变为如下：

(17)

由于a1,a2,…,an-1,an都是变量，所以需要将式(16)后面包含P的式子直接代入式(17)，然后直接求出Pcr的解即可。

用矩阵法求解压杆稳定的欧拉方程有一定的精确性，为解决变截面压杆的稳定性问题提供了一种数值求解的方法，并且为求解压杆失稳后的大挠度杆件变形问题提供了一个求解的方向。

[1]宣海洋.弹性介质上等截面压杆稳定分析[J].山西建筑，2010,36(9)：45-46.

[2]洪振德. 变截面压杆稳定临界力能量计算方法[J].江苏建筑,2011,(3):28-30.

[3]熊容，苏培东. 特定边界条件下压杆稳定问题的解析解与数值解[J].四川建筑,2012,32(2):114-116.

[4]刘鸿文.材料力学(第三版·下册)[M].北京:高等教育出版社，1992.

[5]刘鸿文.高等材料力学[M].北京:高等教育出版社，1985:237.

责任编辑王菊平

The matrix method for solving the finite differential form Euler equation of column stability

HU Zhang-yong

(College of Mechanical & Electrical Engineering, Huanggang Normal University, Huangzhou 438000, Hubei, China)

Euler equation of column stability can be expressed in the form of finite differential. According to the requirements of column stability, coefficient matrices the differential equations are solved with the method of matrices calculation. Subsequently, the critical stress of column stability of constant section of the simply supported beams is obtained. The results are compared with the analytical results. Also presented in the paper is an error analysis including the causes of the error and measures for diminishing it. The method proposed has positive value for the variable section beams.

column stability; finite differential form; Euler equation; matrices

TB12

A

1003-8078(2016)03-0076-04

2016-02-27

10.3969/j.issn.1003-8078.2016.03.19

胡章咏，男，湖北麻城人，讲师，主要研究方向为力学和材料成型工艺与数值模拟等。