单元刚度矩阵等参元MATLAB编程Word文档下载推荐.docx

1、通用函数function B3,S3,K3 = ele_mat_tri3（xy3,mat）%生成平面3结点三角形单元的应变矩阵、应力矩阵和单元刚度矩阵的功能函数%*变量说明*%xy3-结点坐标数组%mat-材料参数矩阵（弹性模量，泊松比，壁厚）%B3-应变矩阵%S3-应力矩阵%K3-单元刚度矩阵%*xyh=1,xy3（1,1）,xy3（1,2）;1,xy3（2,1）,xy3（2,2）;1,xy3（3,1）,xy3（3,2）;A=0.5*det（xyh）;A=abs（A）;D=mat（1）/（1-mat（2）2）*1,mat（2）,0;mat（2）,1,0;0,0,（1-mat（2）/2;b=z

2、eros（1,3）;c=zeros（1,3）;for i=1:3 if i=1 j=2; m=3; elseif i=2 j=3; m=1; else j=1; m=2; end b（i）=xy3（j,2）-xy3（m,2）; c（i）=xy3（m,1）-xy3（j,1）;endB31=1/（2*A）*b（1）,0;0,c（1）;c（1）,b（1）;B32=1/（2*A）*b（2）,0;0,c（2）;c（2）,b（2）;B33=1/（2*A）*b（3）,0;0,c（3）;c（3）,b（3）;B3=B31,B32,B33;S3=D*B3;K3=A*mat（3）*B3*D*B3;主程序clear;

3、clc;%*输入结点坐标数组*xy3=0,0;5,1;1,4; mat=3e6,0.5,1.0;%*输入材料参数矩阵（弹性模量，泊松比，壁厚）*B3,S3,K3=ele_mat_tri3（xy3,mat）三、算例分析算例1：如图1所示三角形单元，结点坐标为1（0,0），2（5,2），3（1,4），弹性模量为200GPa，泊松比为0.35、厚度为0.5m。试求应变矩阵，应力矩阵和单元刚度矩阵。图1 算例1三角形单元解：根据如图1所示三角形单元及其几何和材料参数，编制主程序如下：5,2;mat=2e11,0.35,0.5;运行程序，得到应变矩阵B3如下：-0.1111 0.0000 0.2222

4、-0.2222 -0.0556 0.2778 得到应力矩阵S3（Pa）如下：-2.53E+10-1.77E+105.06E+10-4.43E+092.22E+10-8.86E+09-5.06E+101.77E+10-1.27E+106.33E+10-1.65E+10-8.23E+09-4.12E+091.65E+102.06E+10得到单元刚度矩阵K3（Pa）如下：2.91E+101.71E+10-2.12E+10-1.42E+10-7.91E+09-2.85E+095.48E+10-1.57E+104.43E+09-1.42E+09-5.92E+105.17E+10-8.55E+09-3.0

5、5E+102.42E+101.96E+102.28E+10-2.41E+103.84E+10-2.14E+108.33E+10算例2：如图2所示三角形单元，结点坐标为1（0,0），2（3,0），3（0,5），弹性模量为200GPa，泊松比为0.35、厚度为0.5m。图2 算例2三角形单元根据如图2所示三角形单元及其几何和材料参数，编制主程序如下：3,0;5,0;-0.3333 0.3333 -0.2000 0.2000 -7.60E+10-1.60E+107.60E+100.00E+001.60E+10-2.66E+10-4.56E+102.66E+104.56E+10-1.48E+10-2.

6、47E+102.47E+101.48E+101.06E+113.85E+10-9.50E+10-1.85E+10-1.11E+10-1.99E+106.51E+10-3.09E+10-3.42E+109.50E+101.99E+103.09E+101.85E+101.11E+103.42E+10算例3：如图3所示三角形单元，结点坐标为1（0,0），2（3,0），3（1.5,1.5），弹性模量为200GPa，泊松比为0.35、厚度为0.5m。图3 算例3三角形单元根据如图3所示三角形单元及其几何和材料参数，编制主程序如下：1.5,1.5*sqrt（3）; -0.1925 0.3849 -1.54E+103.07E+10-4.39E+108.77E+10-1.43E+102.85E+105.47E+101.92E+10-4.40E+107.12E+08-1.07E+103.25E+10-7.12E+084.11E+08-3.29E+10-1.92E+102.14E+106.58E+10实 验 二（30分）编写一个计算平面4结点四边形等参元的刚度矩阵的MATLAB函数文件K4 = ele_mat_quad4（xy4,mat），其中：输入变量xy4为结点坐标数组，mat为材料参数矩阵；输出变量K4为单元刚度矩阵。（要求给出3个不同算例进行验证，并绘