弹性力学平面程序设计

1、第九章弹性力学平面程序设计第一节概述一、序编制说明（1）本程序用于计算弹性力学平面应力和应变两类问题。这两类问题用类型码ISP 区别：ISP=1 平面应力问题ISP=2 平面应变问题2）本程序适用于弹性力学平面问题在给定载荷，如结点荷载、自重作用下的静力计算问题。 若有其它非结点载荷作用，则应事先换算成等效结点荷载。3） 所计算的平面弹性体至少应由三个独立的支承，以保证弹性体的几何性质的不变性。支承方式可以在任何一结点的水平支承和竖直支承。4）弹性体由单一的弹性材料组成。5） 本程序计算方法为有限元位移法。组集整体刚度矩阵时采用直接刚度法和后处理法。6） 由于整体刚度矩阵 K 的非零元素分布在

2、以主对角线为中心的斜带形区域，为带形矩阵，根据带形矩阵的特点和对称性，采用半带存储方式。7） 解方程求解结点位移采用等带宽半带存储高斯消元法。二、计算模型及计算方法1计算模型、杆端力及杆端位移采用三角形单元进行分析。由于每个结点有两个位移分量，因此每个单元共有六个位移分量。 以这六个结点位移分量作为基本未知量。单元中结点力和结点位移列阵为FeUiViUjVjUmVmTeuiviujvjumTvm2单元刚度矩阵的建立在每个单元上，结点力用结点位移表示，可建立单元刚度方程Fekee确定单元刚度矩阵采用步骤：1）应用式（ 8-9）或式（ 8-10）求弹性矩阵 D ；2）应用式（ 8-26）求应变转换

3、矩阵B；（ 3）依下式求单元刚度矩阵keBTDBA t3整体刚度矩阵的建立同杆件结构位移法一样，由单元刚度矩阵按直接刚度法集成。4支承条件的引入首先将结点力向量F 换成结点载荷向量P，结构的整体刚度方程变为K=P根据支承条件，把K 和 P 中与支承有杆的元素进行修改。5方程的求解及应力计算采用等带宽半带存贮高斯消元法解方程，可求得结构的结点位移，并存于P 中。为求解单元应力，应从结构位移中取出单元的结点位移，进而求解单元应力。计算式为式（8-27），即eDBe主应力和主平面角的计算公式可参照图9-1 中的应力圆得出。平面应力：AST(xy) / 2(xy)22应力圆半径：RST4xy最大主应力

4、： STMA=AST RST最小主应力： STMI=AST RST主平面角：180arccos(xy)ySTMI9057.29578arctg (xy)STMIy第二节弹性力学平面问题计算框图与程序一、序标识符说明平面应力问题和平面应变问题静力分析程序名PASSAP（Plane Stress/Strain Problem AnalysisProgram）。其主要标识符说明如下：TL （20）算例标题。实型数组，输入参数。NJ结点总数。整形变量，输入参数。NE 单元数组。整形变量，输入参数。NZ 约束结点数。整形变量，输入参数。NPJ结点荷载数。整形变量，输入参数。NPJ0存放结点荷载数NPJ

5、信息，以保证当NPJ0=0时，亦可以按规定定义数组。IPS问题类型码。ND 半带宽。整形变量。E弹性模量。实型变量，输入参数。PR |泊松比。实型变量，输入参数。T弹形体厚度。实型变量，输入参数。V容量。实型变量，输入参数。LND （ NE，3）单元编号组数。LND （ NE ， 1）、LND （ NE，2）、 LND （ NE ， 3）分别为单元 NE 的三个结点编号。整型数组，输入参数。X （ NJ）， Y （ NJ ）结点坐标数组。X （ I）、 Y （ I）分别为I 号结点的x坐标、 y 坐标。实型数组，输入参数。JN（NZ ，3）支承结点数组。JZ（ I， 1）为第 I 个支承的节点

6、编号。JZ（ I，2）、JZ（ I，3）分别是第 I 个支承结点在u，v 位移方向的约束信息。该信息为1 时，表示有约束；为零时，则无约束。实型数组，输入参数。PJ（ NPJ， 3）结点载荷数组。PJ（ I， 1）表示第I 个支承的结点编号。PJ（ I， 2）、 PJ（I ， 3）分别表示该结点沿x， y 方向作用的结点荷载数值。实型数组，输入参数。N自由度数AKE （ 6，6）单元刚度矩阵B（ 3， 6）位移 -应变转换矩阵（几何矩阵）D（ 3，3）弹性矩阵S（ 3， 6）应力矩阵， S=DBP（ N ）结点荷载数组，存放结点荷载列向量；解方程后，存结点位移。DE(6) 单元e 的结点位移。

7、实型数组。ST(3)单元赢利。实型数组。READ 子程序名，求单元面积。MKE 子程序名，形成单元刚度矩阵。MA 子程序名，求单元面积。MD 子程序名，计算弹性矩阵。MB 子程序名，计算位移-应变转换矩阵。MF 子程序名，求结点荷载。RKR 子程序名，引入支承条件。SLOV 子程序名，求结点位移并输出。MADE 子程序名，计算应力并输出。二、框图弹性力学平面问题计算程序总框图如图9-2 所示。开始标题及数组说明子程序 READ读入数据并打印子程序 MKE形成整体刚度矩阵子程序 MA子程序 MD形成结点荷载向量子程序 MB引进支承条件子程序 MF子程序 RKR解方程并打印结点位移子程序 SLOV计算并打印应力子程序 MADE结束NJ K+ND-1IM=NK1=K+1I =K1, IML = I-K+1C = KS(K, L) / KS(K, 1)IW = NW-L+1J =1,IWM = J+I-K消元轮码循环否IM = K+ND-1行码循环修改系数列码循环AK(I,J) =AK(I,J)+C*AK(K,M)修改AKP(I) = P(I)+C*P(K)P(NJ