有限元算例分析..docVIP

一、平面3节点三角形单元分析的算例 F=10000N作用，材料常数为：弹性模量、泊松比，板的厚度为，试按平面应力问题计算各个节点位移及支座反力。 解： 结构的离散化与编号 对该结构进行离散，单元编号及节点编号如图4-20(b)所示，即有二个3节点三角形单元。载荷F按静力等效原则向节点1、节点2移置等效。 节点位移列阵： 节点外载列阵： 约束的支反力列阵： 其中和分别为节点3和节点4的两个方向的支反力。 各个单元的描述 当两个单元取图示中的局部编码(i,j,m)时，其单元刚度矩阵完全相同，即 = 建立整体刚度方程 按单元的位移自由度所对应的位置进行组装可以得到整体刚度矩阵，该组装过程可以写成 具体写出单元刚度矩阵的各个子块在总刚度矩阵中的对应位置如下 代入整体刚度方程Kq=P中，有 边界条件的处理及刚度方程求解 该问题的位移边界条件为，将其代入上式中，划去已知节点位移对应的第5行至第8行（列），有 由上式可求出节点位移如下 支反力的计算 将所求得的节点位移式代入总刚度方程中，可求得支反力如下 二、MATLAB—平面3节点三角形单元分析的算例(Triangle2D3Node) 3节点三角形单元，单元编号及节点编号如图4-20(b)所示。 计算各单元的刚度矩阵(以国际单位) 首先在MATLAB环境下，输入弹性模量E、泊松比NU，薄板厚度为t，平面应力问题性质指示参数ID，然后针对单元1和单元2，分别调用两次函数Triangle2D3Node_Stiffness，就可以得到单元的刚度矩阵k1(6×6)和k2(6×6)。 建立整体刚度方程 由于该结构共有4个节点，则总共的自由度数为8，因此，结构总的刚度矩阵为KK(8×8)，先对KK清零，然后两次调用函数Triangle2D3Node_Assembly进行刚度矩阵的组装。 KK = zeros(8,8); KK=Triangle2D3Node_Assembly(KK,k1,2,3,4); KK=Triangle2D3Node_Assembly(KK,k2,3,2,1) KK = 1.0e+006 * 0.6563 0.3750 -0.3750 -0.1875 -0.2813 -0.1875 0 0 0.3750 1.2188 -0.1875 -1.1250 -0.1875 -0.0938 0 0 -0.3750 -0.1875 0.6563 0 0 0.3750 -0.2813 -0.1875 -0.1875 -1.1250 0 1.2188 0.3750 0 -0.1875 -0.0938 -0.2813 -0.1875 0 0.3750 0.6563 0 -0.3750 -0.1875 -0.1875 -0.0938 0.3750 0 0 1.2188 -0.1875 -1.1250 0 0 -0.2813 -0.1875 -0.3750 -0.1875 0.6563 0.3750 0 0 -0.1875 -0.0938 -0.1875 -1.1250 0.3750 1.2188 边界条件的处理及刚度方程求解 由图4-20(b)可以看出，节点3和节点4的两个方向的位移将为零，即。因此，将针对节点1和节点2的位移进行求解，节点1和节点2的位移将对应KK矩阵中的前4行和前4列，则需从KK(8×8)中提出，置给k，然后生成对应的载荷列阵p，再采用高斯消去法进行求解，注意：MATLAB中的反斜线符号“\”就是采用高斯消去法。 k=KK(1:4,1:4) k = 1.0e+006 * 0.6563 0.3750 -0.3750 -0.1875 0.3750 1.2188 -0.1875 -1.1250 -0.3750 -0.1875 0.6563 0 -0.1875 -1.1250 0 1.2188 p=[0;-5000;0;-5000]; u=k\p u = 0.0188 -0.0899 -0.0150 -0.0842 由此可以看出，所