基于matlab的有限元法分析平面应力应变问题刘刚

合集下载

有限元法MATLAB作业

例题：如图所示的受力均匀分布载荷作用的薄平板结构，将平板离散成两个线性三角元，如图所示，假定2210,0.3,0.025,3000/m E GPa v t m w KN ====。

求（1）该结构的整体刚度矩阵。

（2）节点2和节点3的水平位移和垂直位移。

（3）节点1和节点4的支反力。

（4）每个单元的应力。

（5）每个单元的主应力和主应力方向角。

图2. 用双线性三角元离散化的薄木板解：（1）离散化域我们将平板分为两个单元，4个节点，如图2所示，分布载荷的总作用力平均分给节点2和节点3，由于结构是薄平板，所以假定其属于平面应力情况。

MATLAB 中采用的计算单位是KN 和m 。

表1给出了该题的单元连通性。

（2）单元刚度矩阵通过调用MATLAB 的LinearTriangleElementStiffness 函数，得到两个单元矩阵k 1和k 2，每个矩阵都是66⨯的。

k1=gangdujuzhen(@LinearTriangleElementStiffness, 210000000, 0.3, 0.025, 0, 0, 0.5, 0.25, 0, 0.25, 1)的源程序：function k1= gangdujuzhen( f,E,NU,t,xi,yi,xj,yj,xm,ym,p ) %UNTITLED4 此处显示有关此函数的摘要% 此处显示详细说明k1=gangdujuzhen(@LinearTriangleElementStiffness, 210000000, 0.3, 0.025, 0, 0, 0.5, 0.25, 0, 0.25, 1)图1. 薄平板结构k1= f(E,NU,t,xi,yi,xj,yj,xm,ym,p);>>k1=gangdujuzhen(@LinearTriangleElementStiffness, 210000000, 0.3, 0.025, 0, 0, 0.5, 0.25, 0, 0.25, 1)k1 =1.0e+06 *2.0192 0 0 -1.0096 -2.0192 1.00960 5.7692 -0.8654 0 0.8654 -5.76920 -0.8654 1.4423 0 -1.4423 0.8654-1.0096 0 0 0.5048 1.0096 -0.5048-2.0192 0.8654 -1.4423 1.0096 3.4615 -1.87501.0096 -5.7692 0.8654 -0.5048 -1.8750 6.2740k2=gangdujuzhen2(@LinearTriangleElementStiffness, 210000000, 0.3, 0.025, 0, 0, 0.5, 0, 0.5, 0.25, 1)的源程序：function k2= gangdujuzhen2( f,E,NU,t,xi,yi,xj,yj,xm,ym,p )%UNTITLED3 此处显示有关此函数的摘要% 此处显示详细说明k2=gangdujuzhen2(@LinearTriangleElementStiffness, 210000000, 0.3, 0.025, 0, 0, 0.5, 0, 0.5, 0.25, 1)k2= f(E,NU,t,xi,yi,xj,yj,xm,ym,p);>>k2=gangdujuzhen2(@LinearTriangleElementStiffness, 210000000, 0.3, 0.025, 0, 0, 0.5, 0, 0.5, 0.25, 1)k2 =1.0e+06 *1.4423 0 -1.4423 0.8654 0 -0.86540 0.5048 1.0096 -0.5048 -1.0096 0-1.4423 1.0096 3.4615 -1.8750 -2.0192 0.86540.8654 -0.5048 -1.8750 6.2740 1.0096 -5.76920 -1.0096 -2.0192 1.0096 2.0192 0-0.8654 0 0.8654 -5.7692 0 5.7692（3）集成整体刚度矩阵⨯的，因此为了得到整体刚度矩阵K，我们由于结构有4个节点，所以刚度矩阵是88⨯的零矩阵。

(完整版)有限元大作业matlab---课程设计例子

有限元大作业程序设计学校：天津大学院系：建筑工程与力学学院专业：01级工程力学姓名：刘秀学号：\\\\\\\\\\\指导老师：连续体平面问题的有限元程序分析[题目]：如图所示的正方形薄板四周受均匀载荷的作用，该结构在边界上受正向分布压力，m kNp 1=，同时在沿对角线y 轴上受一对集中压力，载荷为2KN ，若取板厚1=t ，泊松比0=v 。

[分析过程]：由于连续平板的对称性，只需要取其在第一象限的四分之一部分参加分析，然后人为作出一些辅助线将平板“分割”成若干部分，再为每个部分选择分析单元。

采用将此模型化分为4个全等的直角三角型单元。

利用其对称性，四分之一部分的边界约束，载荷可等效如图所示。

[程序原理及实现]：用FORTRAN程序的实现。

由节点信息文件NODE.IN和单元信息文件ELEMENT.IN，经过计算分析后输出一个一般性的文件DATA.OUT。

模型基本信息由文件为BASIC.IN生成。

该程序的特点如下：问题类型：可用于计算弹性力学平面问题和平面应变问题单元类型：采用常应变三角形单元位移模式：用用线性位移模式载荷类型：节点载荷，非节点载荷应先换算为等效节点载荷材料性质：弹性体由单一的均匀材料组成约束方式：为“0”位移固定约束，为保证无刚体位移，弹性体至少应有对三个自由度的独立约束方程求解：针对半带宽刚度方程的Gauss消元法输入文件：由手工生成节点信息文件NODE.IN，和单元信息文件ELEMENT.IN结果文件：输出一般的结果文件DATA.OUT程序的原理如框图：（1）主要变量：ID：问题类型码，ID＝1时为平面应力问题，ID=2时为平面应变问题N_NODE：节点个数N_LOAD：节点载荷个数N_DOF：自由度，N_DOF=N_NODE*2（平面问题）N_ELE：单元个数N_BAND：矩阵半带宽N_BC：有约束的节点个数PE：弹性模量PR：泊松比PT：厚度LJK_ELE(I,3)：单元节点编号数组，LJK_ELE(I,1)，LJK_ELE(I,2)，LJK_ELE(I,3)分别放单元I的三个节点的整体编号X(N_NODE), Y(N_NODE)：节点坐标数组，X(I),Y(I)分别存放节点I的x，y 坐标值P_LJK(N_BC,3)：节点载荷数组，P_LJK(I,1)表示第I个作用有节点载荷的节点的编号，P_LJK(I,2)，P_LJK(I,3)分别为该节点沿x,y方向的节点载荷数值AK(N_DOF,N_BAND)：整体刚度矩阵AKE(6,6)：单元刚度矩阵BB(3,6)：位移……应变转换矩阵（三节点单元的几何矩阵）DD(3,3)：弹性矩阵SS(3,6)；应力矩阵RESULT_N(N_NOF)：节点载荷数组，存放节点载荷向量，解方程后该矩阵存放节点位移DISP_E(6):：单元的节点位移向量STS_ELE(N_ELE,3)：单元的应力分量STS_ND(N_NODE,3)：节点的应力分量（2）子程序说明：READ_IN：读入数据BAND_K：形成半带宽的整体刚度矩阵FORM_KE：计算单元刚度矩阵FORM_P：计算节点载荷CAL_AREA：计算单元面积DO_BC：处理边界条件CLA_DD：计算单元弹性矩阵SOLVE：计算节点位移CLA_BB：计算单元位移……应变关系矩阵CAL_STS：计算单元和节点应力（3）文件管理：源程序文件：chengxu.for程序需读入的数据文件：BASIC.IN，NODE.IN，ELEMENT.IN（需要手工生成）程序输出的数据文件：DATA.OUT（4）数据文件格式：需读入的模型基本信息文件BASIC.IN的格式如下表需读入的节点信息文件NODE.IN的格式如下表需读入的单元信息文件ELEMENT.IN的格式如下表输出结果文件DATA.OUT格式如下表[算例原始数据和程序分析]：（1）模型基本信息文件BASIC.IN的数据为1,4,6,5,31.,0.,1.1,1,0,2,1,0,4,1,1,5,0,1,6,0,11,-0.5,-1.5,3.,-1.,-1,6,-0.5,-0.5（2）手工准备的节点信息文件NODE.IN的数据为1 0.0 2.02 0.0 1.03 1.0 1.04 0. 0.5 1.0 0.6 2.0 0.（3）手工准备的单元信息文件ELEMENT.IN的数据为1 2 3 3 0 0 0 0 1 1 1 1 0 12 4 5 5 0 0 0 0 1 1 1 1 0 25 3 2 2 0 0 0 0 1 1 1 1 0 33 5 6 6 0 0 0 0 1 1 1 1 04 （4）源程序文件chengxu.for为：PROGRAM FEM2DDIMENSION IJK_ELE(500,3),X(500),Y(500),IJK_U(50,3),P_IJK(50,3),&RESULT_N(500),AK(500,100)D IMENSION STS_ELE(500,3),STS_ND(500,3)OPEN(4,FILE='BASIC.IN')OPEN(5,FILE='NODE.IN')OPEN(6,FILE='ELEMENT.IN')OPEN(8,FILE='DATA.OUT')OPEN(9,FILE='FOR_POST.DAT')READ(4,*)ID,N_ELE,N_NODE,N_BC,N_LOADIF(ID.EQ.1)WRITE(8,20)IF(ID.EQ.2)WRITE(8,25)20 FORMAT(/5X,'=========PLANE STRESS PROBLEM========')25 FORMAT(/5X,'=========PLANE STRAIN PROBLEM========')CALL READ_IN(ID,N_ELE,N_NODE,N_BC,N_BAND,N_LOAD,PE,PR,PT, & IJK_ELE,X,Y,IJK_U,P_IJK)CALL BAND_K(N_DOF,N_BAND,N_ELE,IE,N_NODE,& IJK_ELE,X,Y,PE,PR,PT,AK)CALL FORM_P(N_ELE,N_NODE,N_LOAD,N_DOF,IJK_ELE,X,Y,P_IJK, & RESULT_N)CALL DO_BC(N_BC,N_BAND,N_DOF,IJK_U,AK,RESULT_N)CALL SOLVE(N_NODE,N_DOF,N_BAND,AK,RESULT_N)CALL CAL_STS(N_ELE,N_NODE,N_DOF,PE,PR,IJK_ELE,X,Y,RESULT_N, & STS_ELE,STS_ND)c to putout a data fileWRITE(9,70)REAL(N_NODE),REAL(N_ELE)70 FORMAT(2f9.4)WRITE(9,71)(X(I),Y(I),RESULT_N(2*I-1),RESULT_N(2*I),& STS_ND(I,1),STS_ND(I,2),STS_ND(I,3),I=1,N_NODE)71 FORMA T(7F9.4)WRITE(9,72)(REAL(IJK_ELE(I,1)),REAL(IJK_ELE(I,2)),&REAL(IJK_ELE(I,3)),REAL(IJK_ELE(I,3)),&STS_ELE(I,1),STS_ELE(I,2),STS_ELE(I,3),I=1, N_ELE)72 FORMAT(7f9.4)cCLOSE(4)CLOSE(5)CLOSE(6)CLOSE(8)CLOSE(9)E NDcc to get the original data in order to model the problemSUBROUTINE READ_IN(ID,N_ELE,N_NODE,N_BC,N_BAND,N_LOAD,PE,PR, &PT,IJK_ELE,X,Y,IJK_U,P_IJK)DIMENSION IJK_ELE(500,3),X(N_NODE),Y(N_NODE),IJK_U(N_BC,3), & P_IJK(N_LOAD,3),NE_ANSYS(N_ELE,14)REAL ND_ANSYS(N_NODE,3)READ(4,*)PE,PR,PTREAD(4,*)((IJK_U(I,J),J=1,3),I=1,N_BC)READ(4,*)((P_IJK(I,J),J=1,3),I=1,N_LOAD)READ(5,*)((ND_ANSYS(I,J),J=1,3),I=1,N_NODE)READ(6,*)((NE_ANSYS(I,J),J=1,14),I=1,N_ELE)DO 10 I=1,N_NODEX(I)=ND_ANSYS(I,2)Y(I)=ND_ANSYS(I,3)10 CONTINUEDO 11 I=1,N_ELEDO 11 J=1,3IJK_ELE(I,J)=NE_ANSYS(I,J)11 CONTINUEN_BAND=0DO 20 IE=1,N_ELEDO 20 I=1,3DO 20 J=1,3IW=IABS(IJK_ELE(IE,I)-IJK_ELE(IE,J))IF(N_BAND.LT.IW)N_BAND=IW20 CONTINUEN_BAND=(N_BAND+1)*2IF(ID.EQ.1) THENELSEPE=PE/(1.0-PR*PR)PR=PR/(1.0-PR)END IFR ETURNENDcC to form the stiffness matrix of elementSUBROUTINE FORM_KE(IE,N_NODE,N_ELE,IJK_ELE,X,Y,PE,PR,PT,AKE) DIMENSION IJK_ELE(500,3),X(N_NODE),Y(N_NODE),BB(3,6),DD(3,3), & AKE(6,6), SS(6,6)CALL CAL_DD(PE,PR,DD)CALL CAL_BB(IE,N_NODE,N_ELE,IJK_ELE,X,Y,AE,BB)DO 10 I=1,3DO 10 J=1,6SS(I,J)=0.0DO 10 K=1,310 SS(I,J)=SS(I,J)+DD(I,K)*BB(K,J)DO 20 I=1,6DO 20 J=1,6AKE(I,J)=0.0DO 20 K=1,320 AKE(I,J)=AKE(I,J)+SS(K,I)*BB(K,J)*AE*PTRETURNENDcc to form banded global stiffness matrixSUBROUTINE BAND_K(N_DOF,N_BAND,N_ELE,IE,N_NODE,IJK_ELE,X,Y,PE, & PR,PT,AK)DIMENSIONIJK_ELE(500,3),X(N_NODE),Y(N_NODE),AKE(6,6),AK(500,100)N_DOF=2*N_NODEDO 40 I=1,N_DOFDO 40 J=1,N_BAND40 AK(I,J)=0DO 50 IE=1,N_ELECALL FORM_KE(IE,N_NODE,N_ELE,IJK_ELE,X,Y,PE,PR,PT,AKE)DO 50 I=1,3DO 50 II=1,2IH=2*(I-1)+IIIDH=2*(IJK_ELE(IE,I)-1)+IIDO 50 J=1,3DO 50 JJ=1,2IL=2*(J-1)+JJIZL=2*(IJK_ELE(IE,J)-1)+JJIDL=IZL-IDH+1IF(IDL.LE.0) THENELSEAK(IDH,IDL)=AK(IDH,IDL)+AKE(IH,IL)END IF50 CONTINUERETURNENDcc to calculate the area of elementSUBROUTINE CAL_AREA(IE,N_NODE,IJK_ELE,X,Y,AE)DIMENSION IJK_ELE(500,3),X(N_NODE),Y(N_NODE)I=IJK_ELE(IE,1)J=IJK_ELE(IE,2)K=IJK_ELE(IE,3)XIJ=X(J)-X(I)YIJ=Y(J)-Y(I)XIK=X(K)-X(I)YIK=Y(K)-Y(I)AE=(XIJ*YIK-XIK*YIJ)/2.0RETURNENDcc to calculate the elastic matrix of elementSUBROUTINE CAL_DD(PE,PR,DD)DIMENSION DD(3,3)DO 10 I=1,3DO 10 J=1,310 DD(I,J)=0.0DD(1,1)=PE/(1.0-PR*PR)DD(1,2)=PE*PR/(1.0-PR*PR)DD(2,1)=DD(1,2)DD(2,2)=DD(1,1)DD(3,3)=PE/((1.0+PR)*2.0)RETURNENDcc to calculate the strain-displacement matrix of elementSUBROUTINE CAL_BB(IE,N_NODE,N_ELE,IJK_ELE,X,Y,AE,BB) DIMENSION IJK_ELE(500,3),X(N_NODE),Y(N_NODE),BB(3,6)I=IJK_ELE(IE,1)J=IJK_ELE(IE,2)K=IJK_ELE(IE,3)DO 10 II=1,3DO 10 JJ=1,310 BB(II,JJ)=0.0BB(1,1)=Y(J)-Y(K)BB(1,3)=Y(K)-Y(I)BB(1,5)=Y(I)-Y(J)BB(2,2)=X(K)-X(J)BB(2,4)=X(I)-X(K)BB(2,6)=X(J)-X(I)BB(3,1)=BB(2,2)BB(3,2)=BB(1,1)BB(3,3)=BB(2,4)BB(3,4)=BB(1,3)BB(3,5)=BB(2,6)BB(3,6)=BB(1,5)CALL CAL_AREA(IE,N_NODE,IJK_ELE,X,Y,AE)DO 20 I1=1,3DO 20 J1=1,620 BB(I1,J1)=BB(I1,J1)/(2.0*AE)RETURNENDcc to form the global load matrixSUBROUTINE FORM_P(N_ELE,N_NODE,N_LOAD,N_DOF,IJK_ELE,X,Y,P_IJK, & RESULT_N)DIMENSION IJK_ELE(500,3),X(N_NODE),Y(N_NODE),P_IJK(N_LOAD,3), & RESULT_N(N_DOF)DO 10 I=1,N_DOF10 RESULT_N(I)=0.0DO 20 I=1,N_LOADII=P_IJK(I,1)RESULT_N(2*II-1)=P_IJK(I,2)20 RESULT_N(2*II)=P_IJK(I,3)RETURNENDcc to deal with BC(u) (here only for fixed displacement) using "1-0" method SUBROUTINE DO_BC(N_BC,N_BAND,N_DOF,IJK_U,AK,RESULT_N) DIMENSION RESULT_N(N_DOF),IJK_U(N_BC,3),AK(500,100)DO 30 I=1,N_BCIR=IJK_U(I,1)DO 30 J=2,3IF(IJK_U(I,J).EQ.0)THENELSEII=2*IR+J-3AK(II,1)=1.0RESULT_N(II)=0.0DO 10 JJ=2,N_BAND10 AK(II,JJ)=0.0DO 20 JJ=2,II20 AK(II-JJ+1,JJ)=0.0END IF30 CONTINUERETURNENDcc to solve the banded FEM equation by GAUSS eliminationSUBROUTINE SOLVE(N_NODE,N_DOF,N_BAND,AK,RESULT_N) DIMENSION RESULT_N(N_DOF),AK(500,100)DO 20 K=1,N_DOF-1IF(N_DOF.GT.K+N_BAND-1)IM=K+N_BAND-1IF(N_DOF.LE.K+N_BAND-1)IM=N_DOFDO 20 I=K+1,IML=I-K+1C=AK(K,L)/AK(K,1)IW=N_BAND-L+1DO 10 J=1,IWM=J+I-K10 AK(I,J)=AK(I,J)-C*AK(K,M)20 RESULT_N(I)=RESULT_N(I)-C*RESULT_N(K)RESULT_N(N_DOF)=RESULT_N(N_DOF)/AK(N_DOF,1)DO 40 I1=1,N_DOF-1I=N_DOF-I1IF(N_BAND.GT.N_DOF-I-1)JQ=N_DOF-I+1IF(N_BAND.LE.N_DOF-I-1)JQ=N_BANDDO 30 J=2,JQK=J+I-130 RESULT_N(I)=RESULT_N(I)-AK(I,J)*RESULT_N(K)40 RESULT_N(I)=RESULT_N(I)/AK(I,1)WRITE(8,50)50 FORMAT(/12X,'* * * * * RESULTS BY FEM2D * * * * *',//8X,&'--DISPLACEMENT OF NODE--'//5X,'NODE NO',8X,'X-DISP',8X,'Y-DISP') DO 60 I=1,N_NODE60 WRITE(8,70) I,RESULT_N(2*I-1),RESULT_N(2*I)70 FORMAT(8X,I5,7X,2E15.6)RETURNENDcc calculate the stress components of element and nodeSUBROUTINECAL_STS(N_ELE,N_NODE,N_DOF,PE,PR,IJK_ELE,X,Y,RESULT_N, &STS_ELE,STS_ND)DIMENSION IJK_ELE(500,3),X(N_NODE),Y(N_NODE),DD(3,3),BB(3,6), &SS(3,6),RESULT_N(N_DOF),DISP_E(6)DIMENSION STS_ELE(500,3),STS_ND(500,3)WRITE(8,10)10 FORMAT(//8X,'--STRESSES OF ELEMENT--')CALL CAL_DD(PE,PR,DD)DO 50 IE=1,N_ELECALL CAL_BB(IE,N_NODE,N_ELE,IJK_ELE,X,Y,AE,BB)DO 20 I=1,3DO 20 J=1,6SS(I,J)=0.0DO 20 K=1,320 SS(I,J)=SS(I,J)+DD(I,K)*BB(K,J)DO 30 I=1,3DO 30 J=1,2IH=2*(I-1)+JIW=2*(IJK_ELE(IE,I)-1)+J30 DISP_E(IH)=RESULT_N(IW)STX=0STY=0TXY=0DO 40 J=1,6STX=STX+SS(1,J)*DISP_E(J)STY=STY+SS(2,J)*DISP_E(J)40 TXY=TXY+SS(3,J)*DISP_E(J)STS_ELE(IE,1)=STXSTS_ELE(IE,2)=STYSTS_ELE(IE,3)=TXY50 WRITE(8,60)IE,STX,STY,TXY60 FORMAT(1X,'ELEMENT NO.=',I5/18X,'STX=',E12.6,5X,'STY=',&E12.6,2X,'TXY=',E12.6)c the following part is to calculate stress components of nodeWRITE(8,55)55 FORMAT(//8X,'--STRESSES OF NODE--')DO 90 I=1,N_NODEA=0.B=0.C=0.II=0DO 70 K=1,N_ELEDO 70 J=1,3IF(IJK_ELE(K,J).EQ.I) THENII=II+1A=A+STS_ELE(K,1)B=B+STS_ELE(K,2)C=C+STS_ELE(K,3)END IF70 CONTINUESTS_ND(I,1)=A/IISTS_ND(I,2)=B/IISTS_ND(I,3)=C/IIWRITE(8,75)I,STS_ND(I,1),STS_ND(I,2),STS_ND(I,3)75 FORMAT(1X,'NODE NO.=',I5/18X,'STX=',E12.6,5X,'STY=',&E12.6,2X,'TXY=',E12.6)90 CONTINUERETURNENDc FEM2D programm end[算例结果]：chengxu.for所输出的数据文件DATA.OUT数据内容如下：=========PLANE STRESS PROBLEM========* * * * * RESULTS BY FEM2D * * * * *--DISPLACEMENT OF NODE--NODE NO X-DISP Y-DISP1 .000000E+00 -.525275E+012 .000000E+00 -.225275E+013 -.108791E+01 -.137363E+014 .000000E+00 .000000E+005 -.824176E+00 .000000E+006 -.182418E+01 .000000E+00--STRESSES OF ELEMENT--ELEMENT NO.= 1STX=-.108791E+01 STY=-.300000E+01 TXY= .439560E+00ELEMENT NO.= 2STX=-.824176E+00 STY=-.225275E+01 TXY= .000000E+00ELEMENT NO.= 3STX=-.108791E+01 STY=-.137363E+01 TXY= .307692E+00ELEMENT NO.= 4STX=-.100000E+01 STY=-.137363E+01 TXY=-.131868E+00--STRESSES OF NODE--NODE NO.= 1STX=-.108791E+01 STY=-.300000E+01 TXY= .439560E+00NODE NO.= 2STX=-.100000E+01 STY=-.220879E+01 TXY= .249084E+00NODE NO.= 3STX=-.105861E+01 STY=-.191575E+01 TXY= .205128E+00NODE NO.= 4STX=-.824176E+00 STY=-.225275E+01 TXY= .000000E+00NODE NO.= 5STX=-.970696E+00 STY=-.166667E+01 TXY= .586081E-01NODE NO.= 6STX=-.100000E+01 STY=-.137363E+01 TXY=-.131868E+00[结论与体会]:通过本次的课程设计，我对有限元的概念有了更加深刻的理解，同时也弥补了平时学习是疏忽的地方，充实了有限元知识。

基于matlab的平面应力问题研究

基于matlab的平面应力问题研究近年来，随着科学技术的发展，研究者不断对平面应力问题进行研究，并采用不同的方法来模拟和解决这一问题。

其中，借助 MATLAB 术可以高效地完成这项研究。

本文首先简要介绍了MATLAB的基本知识，然后根据现有研究成果，总结出了应用MATLAB解决平面应力问题的典型方法，并结合实际应用，进一步分析其优缺点。

此外，本文还提出了在计算机程序中避免平面应力问题的策略以及可行性研究。

MATLAB是计算机软件，适用于数值计算、科学计算、可视化和交互式编程等多种方面。

它的实用性应用和强大的数据处理能力，使得MATLAB成为当前应用最广泛的计算软件。

目前，MATLAB于平面应力问题的研究，主要有以下两种方法：第一种是使用 MATLAB非线性方程求解方法，即数值求解法。

要解决一个复杂的平面应力问题，需要先建立一个非线性方程组，然后使用 MATLAB行求解。

这种方法的优点是：计算快速，求解精度高，并且可以解决大规模、非线性问题。

但是，它的缺点是：由于非线性方程的求解非常复杂，需要花费较长的时间才能获得准确的结果。

第二种研究方法是使用MATLAB的有限元分析方法。

有限元分析是一种采用数学建模技术来对物体受力过程详细分析的方法。

它可以用来分析平面应力问题，从而更好地了解平面应力的变化趋势，为后续的控制策略提供理论依据。

优点是：可以在计算机上快速求解非线性方程，以及通过二维或三维数值模拟，更准确地分析问题；缺点是：求解精度较低，计算时间较长，技术要求较高。

在计算机程序开发过程中，为了避免出现平面应力问题，应坚持以下步骤：首先，要持续关注平面应力的变化趋势，并分析各种可能影响应力变化的因素，以便找到平面应力的增加的原因；其次，要确定加强应力控制的策略，以避免发生平面应力大小的变化；最后，结合实际情况，应控制计算机程序的开发进度，确保软件系统的性能达到预期目标。

本文分析了MATLAB技术解决平面应力问题的典型方法，以及计算机程序中避免平面应力问题的策略，对现有的研究工作提供了一定的参考。

基于matlab的平面应力问题研究

基于matlab的平面应力问题研究
把握平面应力问题的发展趋势，不仅可以更好地解决工程问题，而且可以更精准地预测结构安全性。

Matlab作为一种计算数学软件，可以模拟平面应力问题，并有效解决相关问题，因此基于Matlab技术的平面应力分析研究受到社会各界的高度重视。

一、Matlab技术在平面应力分析中的应用
（1）Matlab技术可以求解平面应力分析问题，实现非线性状态的力学分析，可以模拟结构的动态变化，以最大程度地测量结构的稳定性。

（2）Matlab技术可以模拟复杂的平面应力分析问题，预测应力的变化情况，确定各种不同环境平面应力的影响因素，以确保结构的安全性。

（3）Matlab技术可以在非线性应力状态下进行模拟，实现对材料强度及变形性能的准确评估，以分析材料对结构安全性的影响，从而更好地保证安全可靠。

二、基于Matlab的平面应力分析的研究方法
（1）首先，建立数学模型，包括应力场模型、变形模型等，以用于模拟应力场分布，研究基于Matlab的平面应力分析问题；
（2）其次，运用Matlab对建立的数学模型进行对应的数值求解，将平面应力分析问题转化为适当的线性或非线性方程；
（3）最后，运用Matlab进行数据分析，并提出有效的解决方案，实现平面应力分析问题的有效解决。

三、结论
Matlab技术在平面应力分析中的应用使其能够实现非线性状态的力学分析，借助Matlab对复杂的平面应力分析问题进行有效模拟，研究分析实现了结构安全性的预测。

Matlab技术为平面应力分析提供了有效的计算方法，使得分析更准确，有助于提高工程的安全性及可靠性。

在未来，基于Matlab 的平面应力分析研究将不断发展，以更好地解决工程应用中的问题。

基于matlab的平面应力问题研究

基于matlab的平面应力问题研究近年来，随着研究的深入及科技的迅速发展，应力学问题受到了越来越多的关注。

应力学问题是工程和物理科学研究的基本内容，它主要涉及金属、机械、材料、结构力学等多个领域。

应力学问题可以从结构力学分析、材料性能预测及有关地震、水力、气动力学等几方面考虑，其中研究应力分布最为重要。

由于应力学问题涉及到多个领域，要研究其解析解和数值模拟解是一项极具挑战性的工作，针对应力分布的研究，主要是利用应力分析的计算机软件，如ANSYS、Comsol等进行模拟研究。

本研究采用Matlab作为计算工具，利用其强大的编程功能以及各种数值算法，构建应力分布模拟系统，对平面应力问题进行研究分析。

首先通过详细的理论推导，导出二维平面应力的基本的计算公式。

然后，利用Matlab的编程，构建一个应力分布模拟系统，根据一维变形规律，分析平面应力的分布情况，采用数值的方法研究不同数据的变化规律。

在此基础上，结合经典的应力模型，探讨不同材料的应力分布情况，并对相关的细节参数进行深入研究，为进行模型模拟提供更加准确的参数。

既然研究计算机模拟，则需要考虑实际情况下的问题，本研究中采用Matlab作为编程语言，首先利用Matlab的视觉技术、绘图工具和图表化工具，构建出一个多维的动态模拟系统，实时展现出平面应力随着外界影响变化的情况，以便分析应力分布中出现的问题。

其次，利用Matlab的数值计算工具，对计算出来的数据进行处理，最后可以求出应力随着外界因素变化的参数分析结果，从而为实际工程提供设计依据。

本文结合Matlab软件与其他相关软件，通过数值模拟的方式，研究了平面应力问题的分布规律，从而确定应力的大小及受力是如何变化的，从而为实际工程设计提供参考。

本研究所得出的结果可以提供一种新的应力分析方法，为现有应力分析手段的改进提供一种有效性的检验方法，能够更加准确、快捷的对平面应力问题进行分析。

综上所述，本文基于Matlab软件，利用其强大的编程功能和数值算法，构建应力分布模拟系统，研究了平面应力问题的分布情况，并对不同材料的应力分布情况作了详细的分析，为工程设计提供了参考。

根据MATLAB的有限元法分析平面应力应变问答刘刚

姓名：刘刚学号：15平面应力应变分析有限元法Abstruct:本文通过对平面应力/应变问题的简要理论阐述，使读者对要分析的问题有大致的印象，然后结合两个实例，通过MATLAB软件的计算，将有限元分析平面应力/应变问题的过程形象的展示给读者，让人一目了然，快速了解有限元解决这类问题的方法和步骤！一.基本理论有限元法的基本思路和基本原则以结构力学中的位移法为基础，把复杂的结构或连续体看成有限个单元的组合，各单元彼此在节点出连接而组成整体。

把连续体分成有限个单元和节点，称为离散化。

先对单元进行特性分析，然后根据节点处的平衡和协调条件建立方程，综合后做整体分析。

这样一分一合，先离散再综合的过程，就是把复杂结构或连续体的计算问题转化简单单元分析与综合问题。

因此，一般的有限揭发包括三个主要步骤：离散化单元分析整体分析。

二.用到的函数1. LinearTriangleElementStiffness(E,NU,t,xi,yi,xj,yj,xm,ym,p)2.LinearBarAssemble(K k I f)3.LinearBarElementForces(k u)4.LinearBarElementStresses(k u A)5.LinearTriangleElementArea(E NU t) 三.实例例1.考虑如图所示的受均布载荷作用的薄平板结构。

将平板离散化成两个线性三角元，假定E=200GPa ，v=0.3，t=0.025m,w=3000kN/m.1.离散化2.写出单元刚度矩阵通过matlab 的LinearTriangleElementStiffness 函数，得到两个单元刚度矩阵1k 和2k ，每个矩阵都是6 6的。

>> E=210e6 E =210000000>> k1=LinearTriangleElementStiffness(E,NU,t,0,0,0.5,0.25,0,0.25,1) k1 =1.0e+006 *Columns 1 through 52.0192 0 0 -1.0096 -2.01920 5.7692 -0.8654 0 0.86540 -0.8654 1.4423 0 -1.4423-1.0096 0 0 0.5048 1.0096 -2.0192 0.8654 -1.4423 1.0096 3.46151.0096 -5.7692 0.8654 -0.5048 -1.8750 Column 61.0096-5.76920.8654-0.5048-1.87506.2740>> NU=0.3NU =0.3000>> t=0.025t =0.0250>> k2=LinearTriangleElementStiffness(E,NU,t,0,0,0.5,0,0.5,0.25,1) k2 =1.0e+006 *Columns 1 through 51.4423 0 -1.4423 0.8654 00 0.5048 1.0096 -0.5048 -1.0096-1.4423 1.0096 3.4615 -1.8750 -2.01920.8654 -0.5048 -1.8750 6.2740 1.00960 -1.0096 -2.0192 1.0096 2.0192-0.8654 0 0.8654 -5.7692 0 Column 6-0.86540.8654-5.76925.76923.集成整体刚度矩阵8*8零矩阵K =0 0 0 0 0 0 0 00 0 0 0 0 0 0 00 0 0 0 0 0 0 00 0 0 0 0 0 0 00 0 0 0 0 0 0 00 0 0 0 0 0 0 00 0 0 0 0 0 0 00 0 0 0 0 0 0 0>> K=LinearTriangleAssemble(K,k1,1,3,4)K =1.0e+006 *Columns 1 through 52.0192 0 0 0 00 5.7692 0 0 -0.86540 0 0 0 00 0 0 0 00 -0.8654 0 0 1.4423-1.0096 0 0 0 0 -2.0192 0.8654 0 0 -1.44231.0096 -5.7692 0 0 0.8654 Columns 6 through 8-1.0096 -2.0192 1.00960 0.8654 -5.76920 0 0* *0 0 00 -1.4423 0.86540.5048 1.0096 -0.50481.0096 3.4615 -1.8750-0.5048 -1.8750 6.2740>> K=LinearTriangleAssemble(K,k1,1,2,3)K =1.0e+007 *0.4038 0 0 -0.1010 -0.2019 0 -0.2019 0.10100 1.1538 -0.0865 0 0 -0.5769 0.0865 -0.57690 -0.0865 0.1442 0 -0.1442 0.0865 0 0-0.1010 0 0 0.0505 0.1010 -0.0505 0 0-0.2019 0 -0.1442 0.1010 0.4904 -0.1875 -0.1442 0.08650 -0.5769 0.0865 -0.0505 -0.1875 0.6779 0.1010-0.0505-0.2019 0.0865 0 0 -0.1442 0.1010 0.3462 -0.18750.1010 -0.5769 0 0 0.0865 -0.0505 -0.1875 0.62744.引入边界条件.用上一步得到的整体刚度矩阵，可以得到该结构的方程组如下形式⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡4y 4X 3y 3X 2y 2X 1y 1X 4y 4X 3y 3X 2y 2X 1y 1X 6F F F F F F F F U U U U U U U U 6.2740 1.8750- 0.5048- 0.8654 0 0 5.7692- 1.0096 1.8750- 3.4615 1.0096 1.4423- 0 0 0.8654 2.0192- 0.5048- 1.0096 6.2740 0 5.7692- 0.8654 01.8750- 0.8654 1.4423- 0 3.4615 1.00962.0192- 1.8750- 0 0 0 5.7692- 1.0096 6.2740 1.8750- 0.5048- 0.8654 0 0 0.8654 2.0192- 1.8750-3.4615 1.0096 1.4423- 5.7692- 0.8654 0 1.8750- 0.5048- 1.0096 6.2740 0 1.0096 2.0192- 1.8750- 0 0.8654 1.4423- 03.461510 本题的边界条件：04411====y x y x U U U U0,375.9,0,375.93322====y x y x F F F F将边界条件带入，得到：⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡4y 4X 1y 1X 3y 3X 2y 2X 6F F 09.37509.375F F 0 0 U U U U 0 0 6.2740 1.8750- 0.5048- 0.8654 0 0 5.7692- 1.0096 1.8750- 3.4615 1.0096 1.4423- 0 0 0.8654 2.0192- 0.5048- 1.0096 6.2740 0 5.7692- 0.8654 01.8750- 0.8654 1.4423- 0 3.4615 1.00962.0192- 1.8750- 0 0 0 5.7692- 1.0096 6.2740 1.8750- 0.5048- 0.8654 0 0 0.8654 2.0192- 1.8750-3.4615 1.0096 1.4423- 5.7692- 0.8654 0 1.8750- 0.5048- 1.0096 6.2740 0 1.0096 2.0192- 1.8750- 0 0.8654 1.4423- 03.4615105.解方程分解上述方程组，提取总体刚度矩阵K 的第3-6行的第3-6列作为子矩阵⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎣⎡⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡ 09.37509.375 U U U U 6.2740 0 5.7692- 0.8654 0 3.4615 1.0096 2.0192- 5.7692- 1.0096 6.2740 1.8750- 0.8654 2.0192- 1.8750-3.4615103y 3X 2y 2X 6 Matlab 命令>> k=K(3:6,3:6) k =1.0e+006 *3.4615 -1.8750 -2.0192 0.8654 -1.8750 6.2740 1.0096 -5.7692 -2.0192 1.0096 3.4615 0 0.8654 -5.7692 0 6.2740>> f=[9.375;0;9.375;0] f =9.3750 0 9.3750 0>> u=k\fu =1.0e-005 *0.71110.11150.65310.0045现在可以清楚的看出，节点2的水平位移和垂直位移分别是0.7111m和0.1115m。

基于有限元法的物体受力变形和应力研究

基于有限元法的物体受力变形和应力研究作者：王重来源：《科技视界》 2014年第17期王重（上海飞机设计研究院环控氧气系统设计研究部，中国上海 201210）【摘要】本文以工程中常见的弯曲悬臂梁作为研究对象，基于有限元分析方法对悬臂梁进行计算模型建立，并求解出悬臂梁在受力后的变形轮廓和应力分布，为其设计和优化提供了验证和指导。

结果表明，该悬臂梁在受力时中部位置应力较大，为该工况下的薄弱区域，需要进行局部加强，以防止在使用中发生失效。

【关键词】弯曲悬臂梁；有限元；刚度矩阵；变形；应力；应变0 引言物体受力变形和应力分析是工程中的常见问题，它在受力零件设计过程中是不可或缺的重要工作，例如在飞机承力梁、高压导管支架等的设计阶段对其将来在飞行中承受载荷后会出现的变形量、变形后轮廓进行估计以及完成强度校核。

这类问题的实质是经典弹性问题，它们的数学模型一般都是一组具有相应边界条件和初值的微分方程，这些微分方程组的解析解能够向我们展示出精确且完整的系统行为，也就是我们所需要的分析结果[1-2]。

但由于这些微分方程组的复杂性，我们又往往无法得到它们的解析解，这时我们就需要利用数值方法来求出近似解，这时在系统中各“节点”的数值解近似于解析解。

因此我们在使用数值方法进行求解前需要对“节点”和“单元”进行合理划分和定义，也就是“离散化”[2-3]。

在此过程之后我们再使用数值解法对问题进行求解。

有限元法是工程中常用的一种数值解法，它使用积分方法来建立系统的代数方程组，用一个连续的函数来近似描述每个单元的解，正因为每个单元的边界是连续的，因此整个系统的解可以由每一个单个的解“组装”起来[2-3]。

不难理解，当我们所划分的单元趋近于无穷多时，使用有限元法所得的解会趋近于精确解。

本文将基于有限元法的基本原理，对具体受力物体进行变形和应力的分析和研究，求解出物体的形变轮廓和相应的应力分布。

1 物体受力工况图1所示为工程中常见的弯曲悬臂梁。

基于matlab的平面应力问题研究

基于matlab的平面应力问题研究近年来，随着物理模拟技术和计算机技术的发展，应用Matlab 从事物理模拟及仿真的研究已经成为重要的学科。

Matlab综合多种数学运算功能及图像处理技术，能够快速、有效地解决各种复杂物理模拟及仿真问题，因此，在电力、能源及制造等多种领域得到了广泛的应用。

本文以Matlab为工具，从物理学的角度，结合一些物理模型，研究平面应力问题。

平面应力问题是几何力学理论中最基本的问题之一，它是研究物体内部力学系统状态的重要数学方法。

根据力学原理，平面应力问题可分为三种，分别是平面压强问题、拉力问题和推力问题，它们的计算主要着重于求解物体内部的力学系统的各种分量的模型。

首先，基于Matlab应用程序，利用拉格朗日乘子法和有限元法，构建平面应力问题的数学模型，并对模型进行参数估计，以求出最优解。

然后，利用Matlab中内置的物理仿真引擎对模型进行数值仿真。

最后，通过Matlab软件对模型进行可视化和调试，实现更为直观的应力分布及力学分析。

Matlab可以解决复杂的物理模拟问题，因此，基于Matlab的平面应力分析在工程实践中广泛使用。

近年来，许多工程师和数学家均使用Matlab研究平面应力问题，在结构力学、机械设计及结构抗震等方面取得重大进展，丰富和拓展了Matlab应用领域，提高了结构性能，为结构抗震和可靠性分析提供了有力的技术支持。

除此之外，本文还介绍了Matlab的其它应用，如图像处理、声音识别等，可用于平面应力问题的分析和计算。

基于Matlab的图像处理技术可以模拟和模拟任何形状的结构以及物体的动态运动，为平面应力问题提供了一种新的解决方法，大大提高了计算的精确性和可靠性。

另外，基于Matlab的声音处理技术也可用于研究平面应力问题，通过捕捉声音信号，可以使这些信号可视化，从而更好地理解平面应力问题产生的原因。

本文介绍了基于Matlab的平面应力分析的基本方法及其在工程实践中的重要性，同时介绍了Matlab的其它应用，如图像处理和声音处理等。

基于matlab的平面应力问题研究

基于matlab的平面应力问题研究随着科学技术的飞速发展，微机和计算机在工程计算和分析领域的应用也逐步发展为一种重要的工具。

而Matlab是一种流行的仿真工具，它利用简单的编程语言和图形化界面，可以方便快捷地实现对复杂工程问题的仿真和分析。

因此，基于Matlab工具的研究主题在近几年得到了广泛的重视和研究。

本文旨在利用Matlab工具，探讨解决平面结构中应力和应变问题。

我们首先详细介绍了Matlab的基本安装、操作和控制，然后重点讨论了Matlab中弹性力学理论和在平面结构中应力应变分析的理论基础。

此后，本文提出了一系列关于应力分析的问题以及求解这些问题所需要的编程步骤和程序。

在求解过程中，本文采用了Matlab程序来辅助完成，加快了计算的速度，准确的计算了平面结构中的应力和应变。

最后，本文给出了几个实例，让读者可以对实际应用进行熟悉，并进一步验证了所提出的方法的正确性和可行性。

Matlab作为一种功能强大的工具，在解决复杂工程问题中有着广泛的应用。

因此，在平面结构中的应力应变分析中，Matlab也可以得到广泛的应用，以有效提高计算效率和提高分析精度。

本文就是基于Matlab工具对平面结构中应力问题的研究。

通过介绍Matlab的基本操作，讨论弹性力学理论和应力分析的基本步骤，以及给出典型实例，本文针对平面结构中应力问题的研究给出了实用性的指导意见。

未来，随着计算机技术和仿真技术的不断发展，Matlab工具分析复杂工程问题的作用会越来越大。

因此，本文对基于Matlab工具的平面应力问题研究具有指导作用，旨在为 s 今后研究相关问题提供有价值的参考。

基于matlab的平面应力问题研究

基于matlab的平面应力问题研究近年来，随着机械工程的发展，解决平面应力问题已经成为一个重要的任务。

近年来，MATLAB一直是机械工程领域中最有用的工具之一，因为它提供了一个高效的方法来求解各种复杂问题。

本文拟讨论基于MATLAB的平面应力问题研究，主要包括：MATLAB的技术原理，各种对平面应力问题的模拟，以及对研究成果的分析与评价。

首先，让我们来看看MATLAB的技术原理。

MATLAB是一种特殊的编程语言，它可以用来解决数学和工程方面的有关问题，包括数值计算、计算机图形学、数据处理和可视化分析等。

MATLAB还提供了强大的矩阵运算系统，可以用来处理复杂的数学运算。

因此，MATLAB 可以用来分析复杂的工程问题，包括平面应力问题。

其次，我们来看看如何用MATLAB模拟平面应力问题。

首先，我们需要利用MATLAB的几何计算模块，根据实际应用场景，构造出网格模型；然后，我们可以使用MATLAB的特定函数，通过给定的尺寸和特性，计算出材料的应力水平；最后，根据计算结果，利用MATLAB 的可视化功能，绘制出应力分布图，将二三维应力空间中的每一点的应力情况进行有效可视化。

最后，我们要从实际应用上进行分析评价，以了解MATLAB对于解决平面应力问题的有效性。

首先，MATLAB在模拟结果上非常准确，因为它可以利用精确的数学模型，高精度的计算；其次，MATLAB可以很好地实现平面应力的可视化，使研究人员可以更加容易地分析研究结果；最后，MATLAB也可以用来解决其他复杂的工程问题，比如热传导问题、力学分析问题等，通过对不同问题的模拟，可以更深入地理解平面应力问题的特性。

综上所述，MATLAB是一种功能强大的编程语言，可以用来解决复杂的工程应用问题，包括平面应力问题。

它的几何计算模块、矩阵运算系统和可视化分析功能为解决平面应力问题提供了有效的解决办法，其准确性和可视化效果也被证明是非常可行的。

因此，基于MATLAB的平面应力问题研究绝对是机械工程领域中一个具有重要意义的方向。

基于MATLAB有限元法结构分析

移是大势所趋 , 坚信建立与社会主义市场经济体制相适应的政府行政管理制度, 就是要加强对交通造价行业协会等社会中介组织的培育与发展。只有统一思想观念 , 不等不靠, 做到洞悉先机、未雨绸缪 , 争取更多的政府主管部门最初的帮助和扶持 , 保证持续发展壮大的能力, 交通造价行业协会就必将在市场经济的大潮中大展宏图 , 完成社会发展赋予的历史使命。
2Hale Waihona Puke 实例分析下面通过一个实例验证 MATLAB 的有限元求 [ 5] 解全过程。如图 1 所示的连续梁 , 已知 E = 210 GPa , I = 50 10 m 。求节点 2 、节点 3 和节点
-6 4
BeamE lem enM t om entD icgram ( ) 函数分别绘出单元的剪力图和弯矩图。本文仅给出单元 ! 的剪力图 ( 图 2)和弯矩图 ( 图 3) 。 2 . 3 M ATLAB 实现过程 ( 程序代码 ) clear E = 210e6 ; I= 50e- 6 ; L 1= 3; L 2= 3; L 3= 4; k1= BeamE le m entS tiffness( E, I , L1) ; k2= BeamE le m entS tiffness( E, I , L2) ; k3= BeamE le m entS tiffness( E, I , L3) ; K = zeros( 8 , 8) ; K = BeamA ssem ble( K, k1, 1 , 2) ; K = BeamA ssem ble( K, k2, 2 , 3) ; K = BeamA ssem ble( K, k3, 3 , 4) k= [ K ( 4 , 4) K( 4 , 6) K ( 4 , 8) ; K ( 6 , 4) K( 6 , 6) K (6 , 8) ; K ( 8 , 4) K( 8 , 6) K ( 8 , 8) ] f= [ 7 . 5 ; - 15 ; 15] u= k \ f U = [ 0; 0 ; 0 ; u( 1) ; 0 ; u( 2) ; 0; u( 3) ] F= K* U u1= [ U ( 1) ; U ( 2) ; U ( 3) ; U ( 4) ] ; u2= [ U ( 3) ; U ( 4) ; U ( 5) ; U ( 6) ] ; u3= [ U ( 5) ; U ( 6) ; U ( 7) ; U ( 8) ] ; f1= B eamE lem entForces( k1, u1) ; f1= f1- [ - 15 ; - 7 . 5 ; - 15 ; 7 . 5] ( 下转第 136 页 )

基于MATLAB的有限元法分析平面应力应变问题--刘刚

把连续体分成有限个单元和节点，称为离散化。

先对单元进行特性分析，然后根据节点处的平衡和协调条件建立方程，综合后做整体分析。

这样一分一合，先离散再综合的过程，就是把复杂结构或连续体的计算问题转化简单单元分析与综合问题。

因此，一般的有限揭发包括三个主要步骤：离散化单元分析整体分析。

二.用到的函数1. LinearTriangleElementStiffness(E,NU,t,xi,yi,xj,yj,xm,ym,p)2.LinearBarAssemble(K k I f)3.LinearBarElementForces(k u)4.LinearBarElementStresses(k u A)5.LinearTriangleElementArea(E NU t)三.实例例1.考虑如图所示的受均布载荷作用的薄平板结构。

>> E=210e6 E =210000000>> k1=LinearTriangleElementStiffness(E,NU,t,0,0,0.5,0.25,0,0.25,1) k1 =1.0e+006 *Columns 1 through 52.0192 0 0 -1.0096 -2.01920 5.7692 -0.8654 0 0.86540 -0.8654 1.4423 0 -1.4423-1.0096 0 0 0.5048 1.0096 -2.0192 0.8654 -1.4423 1.0096 3.46151.0096 -5.7692 0.8654 -0.5048 -1.8750 Column 61.0096-5.76920.8654-0.5048-1.87506.2740>> NU=0.3NU =0.3000>> t=0.025t =0.0250>> k2=LinearTriangleElementStiffness(E,NU,t,0,0,0.5,0,0.5,0.25,1) k2 =1.0e+006 *Columns 1 through 51.4423 0 -1.4423 0.8654 00 0.5048 1.0096 -0.5048 -1.0096-1.4423 1.0096 3.4615 -1.8750 -2.01920.8654 -0.5048 -1.8750 6.2740 1.00960 -1.0096 -2.0192 1.0096 2.0192-0.8654 0 0.8654 -5.7692 0 Column 6-0.86540.8654-5.76925.76923.集成整体刚度矩阵8*8零矩阵K =0 0 0 0 0 0 0 00 0 0 0 0 0 0 00 0 0 0 0 0 0 00 0 0 0 0 0 0 00 0 0 0 0 0 0 00 0 0 0 0 0 0 00 0 0 0 0 0 0 00 0 0 0 0 0 0 0 >> K=LinearTriangleAssemble(K,k1,1,3,4)K =1.0e+006 *Columns 1 through 52.0192 0 0 0 00 5.7692 0 0 -0.86540 0 0 0 00 0 0 0 00 -0.8654 0 0 1.4423-1.0096 0 0 0 0 -2.0192 0.8654 0 0 -1.44231.0096 -5.7692 0 0 0.8654Columns 6 through 8-1.0096 -2.0192 1.00960 0.8654 -5.76920 0 00 0 00 -1.4423 0.86540.5048 1.0096 -0.50481.0096 3.4615 -1.8750-0.5048 -1.8750 6.2740>> K=LinearTriangleAssemble(K,k1,1,2,3)K =1.0e+007 *0.4038 0 0 -0.1010 -0.2019 0 -0.2019 0.10100 1.1538 -0.0865 0 0 -0.5769 0.0865 -0.57690 -0.0865 0.1442 0 -0.1442 0.0865 0 0-0.1010 0 0 0.0505 0.1010 -0.0505 0-0.2019 0 -0.1442 0.1010 0.4904 -0.1875 -0.1442 0.08650 -0.5769 0.0865 -0.0505 -0.1875 0.6779 0.1010 -0.0505-0.2019 0.0865 0 0 -0.1442 0.1010 0.3462 -0.18750.1010 -0.5769 0 0 0.0865 -0.0505 -0.1875 0.62744.引入边界条件.用上一步得到的整体刚度矩阵，可以得到该结构的方程组如下形式⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡4y 4X 3y 3X 2y 2X 1y 1X 4y 4X 3y 3X 2y 2X 1y 1X 6F F F F F F F F U U U U U U U U 6.2740 1.8750- 0.5048- 0.8654 0 0 5.7692- 1.0096 1.8750- 3.4615 1.0096 1.4423- 0 0 0.8654 2.0192- 0.5048- 1.0096 6.2740 0 5.7692- 0.8654 01.8750- 0.8654 1.4423- 0 3.4615 1.00962.0192- 1.8750- 0 0 0 5.7692- 1.0096 6.2740 1.8750- 0.5048- 0.8654 0 0 0.8654 2.0192- 1.8750-3.4615 1.0096 1.4423- 5.7692- 0.8654 0 1.8750- 0.5048- 1.0096 6.2740 0 1.0096 2.0192- 1.8750- 0 0.8654 1.4423- 03.461510 本题的边界条件：04411====y x y x U U U U0,375.9,0,375.93322====y x y x F F F F将边界条件带入，得到：⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡4y 4X 1y 1X 3y 3X 2y 2X 6F F 09.37509.375F F 0 0 U U U U 0 0 6.2740 1.8750- 0.5048- 0.8654 0 0 5.7692- 1.0096 1.8750- 3.4615 1.0096 1.4423- 0 0 0.8654 2.0192- 0.5048- 1.0096 6.2740 0 5.7692- 0.8654 0 1.8750- 0.8654 1.4423- 0 3.4615 1.0096 2.0192- 1.8750- 0 0 0 5.7692- 1.0096 6.2740 1.8750- 0.5048- 0.8654 0 0 0.8654 2.0192- 1.8750- 3.4615 1.0096 1.4423- 5.7692- 0.8654 0 1.8750- 0.5048- 1.0096 6.2740 0 1.0096 2.0192- 1.8750- 0 0.8654 1.4423- 03.4615105.解方程分解上述方程组，提取总体刚度矩阵K 的第3-6行的第3-6列作为子矩阵⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎣⎡⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡ 09.37509.375 U U U U 6.2740 0 5.7692- 0.8654 0 3.4615 1.0096 2.0192- 5.7692- 1.0096 6.2740 1.8750- 0.8654 2.0192- 1.8750- 3.4615103y 3X2y 2X 6 Matlab 命令>> k=K(3:6,3:6) k =1.0e+006 *3.4615 -1.8750 -2.0192 0.8654 -1.8750 6.2740 1.0096 -5.7692 -2.0192 1.0096 3.4615 0 0.8654 -5.7692 0 6.2740>> f=[9.375;0;9.375;0] f =9.37509.3750>> u=k\fu =1.0e-005 *0.71110.11150.65310.0045现在可以清楚的看出，节点2的水平位移和垂直位移分别是0.7111m和0.1115m。

基于matlab的平面应力问题研究

基于matlab的平面应力问题研究凡是浮现在我们身边的物体，都会受到一些形式的外力，这些外力也就是叫做应力。

应力可以分为压力、拉力和转矩，这些应力可以影响物体的状态和受力情况。

因此，应力的研究在工程物理和机械结构力学中是一个重要的研究课题。

而最常见的应力影响是平面应力问题，因为它涉及到应力的分布特性和力学性质，因此历来受到广泛关注。

基于Matlab的平面应力问题研究，主要是通过Matlab计算出来的数据来研究平面应力问题，探究平面应力问题的机理，从而更好的理解整个应力的分布情况。

首先，针对于平面应力问题，我们要明确其定义，平面应力指的是通过平面上受力系统中的单元把均匀的外力分布到每个单元上的应力，即将外力以单元的形式表示出来，能够更好的描述平面应力问题。

其次，我们可以使用Matlab软件，针对应力分布中的参数进行建模，也就是使用Matlab创建无源模型，其中涉及到了材料特征参数、应力分布参数等多种参数，然后根据这些参数，进行有限元算法计算，从而得出完整的平面应力分布。

最后，为了更好的理解应力分布特性和力学性质，我们可以利用Matlab进行可视化处理，利用图像生成工具，可以直观的查看应力的分布情况，以及统计其分布的特性，由此可以进一步分析和探究平面应力的规律，为应力的结构设计提供了可靠的依据。

从以上讨论可以清晰的看出，使用Matlab的平面应力问题研究可以更好的帮助我们理解应力的分布特性和本质性质，从而更好的进行结构设计。

总之，基于Matlab的平面应力研究可以为有关应力分布特性和力学性质的研究提供有力的支持，是一项重要的研究课题。

Matlab 有着易于使用的界面和完整的图形语言，从而可以轻松的分析和探究应力的分布特性和力学特性，是对结构设计起到不可或缺的作用。

基于matlab的平面应力问题研究

基于matlab的平面应力问题研究近年来，随着科技的进步，机械制造领域中逐渐复杂的计算问题越来越受到重视。

其中，最为重要的是平面应力问题的研究，它的研究关乎更大的机械制造发展的未来。

在研究平面应力问题的过程中，matlab是一款广泛使用的强大的编程工具。

它的出现极大地改变了机械分析的流程，使得工程师能够进行更有效率的运算分析。

本文将介绍matlab在平面应力问题研究中的应用，以期能够给大家带来一些有用的认识。

基于matlab的平面应力问题研究，主要是通过matlab编程来实现的。

平面应力问题的研究，基本等同于对于板材形状和材料性质的复杂运算分析。

有了matlab的帮助，用户可以使用matlab的计算和可视化工具，对于需要分析的平面应力系统进行结构参数分析、材料性质分析和应力分析等复杂运算分析，从而可以获得完整的分析结果。

matlab的可视化技术，不仅可以更好地反映出实际应力问题的变化规律，而且可以帮助用户更加直观地了解这些变化规律。

此外，matlab还提供了一系列实用的功能，可以帮助用户更有效地研究平面应力问题。

例如，matlab的数据分析功能可以帮助用户更有效地提取有用的信息。

matlab的图形化工具，可以帮助用户更好地可视化应力系统的变化规律。

此外，matlab的优化算法还可以帮助用户更有效地求解出一系列最优的解决方案。

最后，matlab在复杂运算分析领域发挥了巨大作用，在平面应力问题研究中也发挥了重要作用。

为了更加有效地使用matlab，用户可以利用matlab提供的各种实用功能，以及matlab编程语言提供的有用脚本，开发出更有效率的平面应力问题研究程序。

从而有助于机械分析中的进一步研究。

综上所述，matlab在平面应力问题研究中的应用，使得平面应力问题的研究变得更加有效率，而且能够给大家带来更多有用的认识。

未来，matlab将会在机械分析中发挥极大的作用，以期能够为机械制造发展带来更大的帮助。

matlab平面应力问题例题

一、介绍Matlab是一种用于数学计算和工程科学的高级技术计算语言和交互式环境。

它在工程，科学和经济领域广泛应用，包括线性代数，统计学，积分和微分方程等数学计算。

在工程科学领域，Matlab也被广泛应用于解决力学问题，其中包括平面应力问题。

二、问题描述平面应力问题是指在平面内受力体系中，应力分量只存在于两个相互垂直的方向上，而第三个方向不受力。

在工程实践中，平面应力问题经常出现在工程结构设计中。

建筑结构、飞机机身等工程结构在受力时可能产生平面应力问题。

解决平面应力问题需要确定应力分布和位移场。

三、解题步骤1. 确定边界条件：首先需要明确物体的几何形状和边界条件。

边界条件包括受力情况和位移情况。

2. 应力分析：利用弹性力学理论，可以推导出平面应力问题的一般解。

在Matlab中可以使用有限元分析方法来求解应力分布。

3. 位移场分析：根据边界条件和应力分布，可以求解出平面应力问题的位移场。

位移场分析可以帮助工程师了解结构的稳定性和变形情况。

4. 结果分析：最后需要对得到的结果进行分析，可以评估结构的安全性和稳定性。

四、例题分析接下来我们以一个简单的例题来演示如何使用Matlab解决平面应力问题。

考虑一个长方形梁，在两端分别受到水平方向的拉力和压力。

假设梁的材料为钢，弹性模量为200 GPa，泊松比为0.3。

梁的几何参数为宽度b=100mm，高度h=200mm，长度L=1000mm。

求解在梁内部的应力分布和位移场。

其中，拉力和压力的大小分别为P=10kN和P=5kN。

在Matlab中，我们可以按照以下步骤进行求解。

1. 确定边界条件。

将拉力作用于梁的左端，压力作用于梁的右端。

2. 进行应力分析。

利用有限元分析方法，可以求解出梁内部的应力分布。

3. 进行位移场分析。

根据边界条件和应力分布，可以求解出梁的位移场。

4. 对结果进行分析。

可以评估梁的稳定性和变形情况。

五、结论通过Matlab的求解，我们可以得到梁内部的应力分布和位移场。

MATLAB有限元分析与应用可编辑全文

%
modulus of elasticity E, cross-sectional
%
area A, and length L. The size of the
%
element stiffness matrix is 2 x 2.
y = [E*A/L -E*A/L ; -E*A/L E*A/L];
2019/11/28
y = k * u/A;
2019/11/28
18
§3-2 线性杆元
3、实例计算分析应用
如图所示二线性杆元结构，假定E=210MPa，A=0.003m^2，P=10kN，节点3的右位移为0.002m。
求：系统的整体刚度矩阵；节点2的位移；节点1、3的支反力；每个杆件的应力
解：
步骤1：离散化域
%LinearBarElementStresses This function returns the element nodal
%
stress vector given the element stiffness
%
matrix k, the element nodal displacement
%
vector u, and the cross-sectional area A.
? ?
?
? ?
?
10??
630000 ????0.002?? ?? F3 ??
线性杆元也是总体和局部坐标一致的一维有限单元，用线性函数描述
每个线性杆元有两个节点（node）
? EA
单刚矩阵为：k
?
? ?
L
???
EA L
?
EA L
? ? ?

基于matlab的平面应力问题研究

基于matlab的平面应力问题研究平面应力问题是指具有给定的空间位置和尺寸的刚体所受到的平面外的平面力作用。

本文介绍了一种求解平面应力问题的常用方法——matlab分离变量法，并以它为例说明了应用matlab软件求解平面应力问题的一般步骤。

一、确定微分方程根据已知条件，列出应力微分方程。

(1)式中， Y为拉力，τ为应力。

(2)式中，θ为力的作用点到截面的距离; x为力作用点到质点的距离;y为应力; M为材料的弹性模量; a为材料的泊松比;ε为弹性系数。

(3)式中， M为结构的重量;G为物体的重力。

(4)式中， H为结构自身的重量。

(5)式中， F 为结构对支座的力。

(6)式中， B为结构对基础的力。

(7)式中， E为结构对地面的压力。

(8)式中， F′为结构对外部荷载的力矩。

(9)式中， E′为结构对外部荷载的力矩。

(10)式中， E′′为结构对外部荷载的力矩。

(11)式中， P′为结构的外加荷载; P为地面上物体的压力。

(12)式中， E′为结构对物体的力矩; F为结构对基础的力矩; a为结构的惯性矩; E为结构对地面的压力。

二、求解微分方程(1)分别在微分方程的两边取和，代入数值得到微分方程的近似解。

(2)计算系数即式中， X和M分别为求解的一对多边形元素，γ为导数。

(3)令，代入上式得到近似的积分关系式。

(4)代入解得到求解的积分表达式。

(5)设，代入式求解即可得到近似解。

(6)因此，由式求解。

(7)令，代入式得到积分表达式。

(8)令，代入上式得到近似的积分关系式。

(9)令，代入式得到近似的积分关系式。

(10)令，代入式得到积分表达式。

(11)因此，由式求解。

(12)因此，由式求解。

(13)若应力计算是线性的，则应力计算是线性的。

(14)若应力计算是非线性的，则应力计算是非线性的。

(15)三、整理数据进行回归分析整理出的应力分布曲线，并对应力分布曲线进行回归分析。

四、最后结论通过应用matlab软件对有限元方法求解平面应力问题的应用，实现了对平面应力问题的求解，将有限元方法应用于求解平面应力问题，能够使得复杂问题更容易被解决。

基于matlab的平面应力问题研究

基于matlab的平面应力问题研究随着社会发展，基于计算机的工程计算正变得越来越重要。

工程计算被用于多个不同的领域，包括结构设计、材料性能评估、计算流体动力学、力学计算等。

这些领域当中，应力问题在结构设计及分析领域被广泛用于，如何有效、准确地解决应力问题成为必要解决难题。

MATLAB是一个专业的工程计算和数据可视化软件，已经被广泛应用于工程计算。

MATLAB拥有完整的数学计算能力及全面高效的编程辅助功能，能够方便地解决复杂的问题。

因此，MATLAB在应力问题的计算中起到了重要的作用。

基于MATLAB的平面应力问题研究是必要的科学探索。

平面应力问题的研究可以分为两个部分，一部分是计算部分，另一部分是可视化部分。

计算部分可以使用MATLAB提供的函数和程序库来实现，如求解应力问题所需要的偏微分方程、求解矩阵方程等。

此外，可以使用MATLAB编程实现有限元法的求解，有效的强化平面应力问题的计算能力。

可视化部分可以使用MATLAB提供的可视化工具，如绘制应力分布图，绘制变形和应力曲线等。

在实际的计算中，需要考虑多种因素的影响。

例如，在计算应力分布时要考虑材料的弹性模量、温度、应变等变化；考虑屈服条件及受力条件；考虑材料的局部结构及外部环境情况等。

MATLAB可以提供数学模型解析以及复杂的编程能力，能够很好的处理复杂的工程计算问题，使用MATLAB可以更加有效地求解应力问题。

本文针对基于MATLAB的平面应力问题进行研究，主要讨论了MATLAB在应力问题计算中的应用。

本文分析了MATLAB在解决应力问题中编程能力及可视化工具的作用；此外，还进行了复杂因素的考虑和讨论，指出了MATLAB在处理复杂应力问题中具有强大的计算能力。

最后，本文提出了使用MATLAB改进应力问题求解及计算效率的建议。

由于MATLAB提供的计算及可视化能力，它受到了工程学界的广泛应用，并可以帮助我们更好的理解应力问题，在应力分析与设计领域有重要作用。

《有限元方法与MATLAB程序设计》-课件--4-平面问题

x j ym xm y j ( y j ym )x (xm x j ) y
ai bi x ci y
L1
L2
1 2A
ai
a
j
bi bj
ci 1
c
j
x
L3
am bm cm y
1 1
x
xi
1 xj
1 xm
L1 L2
y yi y j ym L3
0 Ni
假设位移模式 u(x, y) 1 2x 3 y
m (xm , ym )
v(x, y) 4 5x 6 y
结点位移
j(xj, yj )
ui u(xi , yi ) 1 2 xi 3 yi
vi v(xi , yi ) 4 5xi 6 yi o
x
u j u(x j , y j ) 1 2 x j 3 y j
3
bi 0
3
3
3
3
3
2 Ai (x, y) (ai bi x ci y) ai bi x ci y 2A
i 1
i 1
i 1
i 1
i 1
i 1
3
ci 0
i 1
3
ai 2A
i 1
8
面积坐标
1x y 2Ai (x, y) 1 x j y j
1 xm ym
Ai (xi , yi ) 2 A Ai (x j , y j ) 0 Ai (xm , ym ) 0
x
Ni y
Ni y
Ni x
N j x
N j y
N j y
N j x
Nm x
Nm y
ui
Nm y
Nm x

1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。

把连续体分成有限个单元和节点，称为离散化。

先对单元进行特性分析，然后根据节点处的平衡和协调条件建立方程，综合后做整体分析。

这样一分一合，先离散再综合的过程，就是把复杂结构或连续体的计算问题转化简单单元分析与综合问题。

因此，一般的有限揭发包括三个主要步骤：离散化单元分析整体分析。

二.用到的函数1. LinearTriangleElementStiffness(E,NU,t,xi,yi,xj,yj,xm,ym,p)(K k I f)(k u)(k u A)(E NU t)三.实例例1.考虑如图所示的受均布载荷作用的薄平板结构。

将平板离散化成两个线性三角元，假定E=200GPa，v=，t=0.025m,w=3000kN/m.1.离散化2.写出单元刚度矩阵通过matlab 的LinearTriangleElementStiffness 函数，得到两个单元刚度矩阵1k 和2k ，每个矩阵都是6 6的。

>> E=210e6 E =>> k1=LinearTriangleElementStiffness(E,NU,t,0,0,,,0,,1) k1 = +006 *Columns 1 through 50 0 0 0 0 0 0 0Column 6 >> NU= NU = >> t= t =>> k2=LinearTriangleElementStiffness(E,NU,t,0,0,,0,,,1)k2 =+006 *Columns 1 through 50 00 0Column 63.集成整体刚度矩阵 8*8零矩阵K =0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 >> K=LinearTriangleAssemble(K,k1,1,3,4)K =+006 *Columns 1 through 50 0 0 00 0 00 0 0 0 00 0 0 0 0 0 0 00 0 0 0 0 0 0 0Columns 6 through 80 0 0 0 0 0 0>> K=LinearTriangleAssemble(K,k1,1,2,3) K =+007 *0 0 0 0 0 00 0 0 0 0 0 0 0 0 00 0 0 04.引入边界条件.用上一步得到的整体刚度矩阵，可以得到该结构的方程组如下形式⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡4y 4X 3y 3X 2y 2X 1y 1X 4y 4X 3y 3X 2y 2X 1y 1X 6F F F F F F F F U U U U U U U U 6.2740 1.8750- 0.5048- 0.8654 0 0 5.7692- 1.0096 1.8750- 3.4615 1.0096 1.4423- 0 0 0.8654 2.0192- 0.5048- 1.0096 6.2740 0 5.7692- 0.8654 01.8750- 0.8654 1.4423- 0 3.4615 1.00962.0192- 1.8750- 0 0 0 5.7692- 1.0096 6.2740 1.8750- 0.5048- 0.8654 0 0 0.8654 2.0192- 1.8750-3.4615 1.0096 1.4423- 5.7692- 0.8654 0 1.8750- 0.5048- 1.0096 6.2740 0 1.0096 2.0192- 1.8750- 0 0.8654 1.4423- 03.461510 本题的边界条件：04411====y x y x U U U U0,375.9,0,375.93322====y x y x F F F F将边界条件带入，得到：⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡4y 4X 1y 1X 3y 3X 2y 2X 6F F 09.37509.375F F 0 0 U U U U 0 0 6.2740 1.8750- 0.5048- 0.8654 0 0 5.7692- 1.0096 1.8750- 3.4615 1.0096 1.4423- 0 0 0.8654 2.0192- 0.5048- 1.0096 6.2740 0 5.7692- 0.8654 01.8750- 0.8654 1.4423- 0 3.4615 1.00962.0192- 1.8750- 0 0 0 5.7692- 1.0096 6.2740 1.8750- 0.5048- 0.8654 0 0 0.8654 2.0192- 1.8750-3.4615 1.0096 1.4423- 5.7692- 0.8654 0 1.8750- 0.5048- 1.0096 6.2740 0 1.0096 2.0192- 1.8750- 0 0.8654 1.4423- 03.4615105.解方程分解上述方程组，提取总体刚度矩阵K 的第3-6行的第3-6列作为子矩阵⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎣⎡⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡ 09.37509.375 U U U U 6.2740 0 5.7692- 0.8654 0 3.4615 1.0096 2.0192- 5.7692- 1.0096 6.2740 1.8750- 0.8654 2.0192- 1.8750-3.4615103y 3X 2y 2X 6 Matlab 命令>> k=K(3:6,3:6) k =+006 *0 0>> f=[;0;;0] f = 0 0 >> u=k\f u =*现在可以清楚的看出，节点2的水平位移和垂直位移分别是0.7111m 和0.1115m 。

节点3的水平位移和垂直位移分别是0.6531m 和0.0045m 。

6.后处理用matlab 命令求出节点1和节点4的支反力以及每个单元的应力。

首先建立总体节点位移矢量U ，U=[0;0;u;0;0] U = *0 0 0 0 >> F=K*U F =由以上知，节点1的水平反力和垂直反力分别是（指向左边）和（作用力方向向下），节点4的水平反力和垂直反力分别是（指向左边）和（作用力方向向下）.满足力平衡条件。

接着，建立单元节点位移矢量21u 和u ，然后调用matlab 命令LinearTriangleElementStresses 计算单元应力sigma1和sigma2>> u1=[U(1);U(2);U(5);U(6);U(7);U(8)]u1 =*>> u2=[U(1);U(2);U(3);U(4);U(5);U(6)]u2 =*>> sigma1=LinearTriangleElementStresses(E,NU,,0,0,,,0,,1,u1) sigma1 =+003 *>> sigma2=LinearTriangleElementStresses(E,NU,,0,0,,0,,,1,u2) sigma2 =+003 *由以上可知，单元1的应力)(0144.3拉应力MPa x =σ，)(9043.0y 拉应力MPa =σ，)(0072.0y 正值MPa x =τ 。

单元2的应力是)(9856.2拉应力MPa x =σ)(0036.0y 压应力MPa =σ)(0072.0y 负值MPa x =τ。

显然，在x方向的应力（拉应力）接近于正确的值3MPa （拉应力）。

接着调用LinearTriangleElementStresses 函数计算每个单元的主应力和主应力方向角。

>> s1= LinearTriangleElementPStresses(sigma1) s1 = +003 *>> s2= LinearTriangleElementPStresses(sigma2) s2 = +003 *)(0144.31拉应力MPa =σ，)(9043.02拉应力MPa =σ主应力方向角ο2.0=p θ )(M 9856.21拉应力Pa =σ，)(0036.02压应力MPa =σ，ο1.0-=p θ例 2.考虑如图所示的由均匀分布载荷和集中载荷作用的薄平板结构。

将平板离散化成12个线性三角单元，如图4所示。

假定E=210GPa，v=，t=0.025m，w=100kN/m和P=。

1.离散化2.写出单元刚度矩阵>> E=201e6;>> NU=;>> t=;>> k1= LinearTriangleElementStiffness(E,NU,t,0,,,,,,1); >> k2= LinearTriangleElementStiffness(E,NU,t,0,,0,,,,1); >> k3= LinearTriangleElementStiffness(E,NU,t,,,,,,,1); >> k4= LinearTriangleElementStiffness(E,NU,t,,,0,,,,1); >> k5= LinearTriangleElementStiffness(E,NU,t,0,,,,,,1); >> k6= LinearTriangleElementStiffness(E,NU,t,0,,0,0,,,1); >> k7= LinearTriangleElementStiffness(E,NU,t,,,,,,0,1); >> k8= LinearTriangleElementStiffness(E,NU,t,,,0,0,,0,1); >> k9= LinearTriangleElementStiffness(E,NU,t,025,,,0,,,1); >> k10= LinearTriangleElementStiffness(E,NU,t,,,,,,,1); >> k11= LinearTriangleElementStiffness(E,NU,t,,0,,0,,,1); >> k12= LinearTriangleElementStiffness(E,NU,t,,,,0,,,1)k1 =+006 *3.集成整体刚度矩阵：>>K=zero(22,22);>>K=LinearTriangleAssemble(K,k1,1,3,2);>>K=LinearTriangleAssemble(K,k2,1,4,3);>>K=LinearTriangleAssemble(K,k3,3,5,2);>>K=LinearTriangleAssemble(K,k4,3,4,5);>>K=LinearTriangleAssemble(K,k5,4,6,5);>>K=LinearTriangleAssemble(K,k6,4,7,6);>>K=LinearTriangleAssemble(K,k7,5,6,8);>>K=LinearTriangleAssemble(K,k8,6,7,8);>>K=LinearTriangleAssemble(K,k9,5,8,9);>>K=LinearTriangleAssemble(K,k10,5,9,10);>>K=LinearTriangleAssemble(K,k11,8,11,9);>>K=LinearTriangleAssemble(K,k12,9,11,10)运行得+008 *Columns 1 through 70 00 0 00 00 0 00 0 0 0 0 00 0 0 0 0 00 0 0 0 0 00 0 0 0 0 00 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0Columns 8 through 140 0 0 0 0 00 0 0 0 0 00 0 0 0 00 0 0 0 00 0 0 00 0 0 00 00 0 00 00 00 00 0 0 0 00 0 0 0 00 0 0 0 00 0 0 0 00 0 0 0 0 0 0 0 0 0 0 0 0 0Columns 15 through 210 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 00 0 0 0 00 0 0 0 00 0 0 0 00 0 0 0 00 00 00 00 0Column 224.引入边界条件：U1x= U1y= U4x= U4y=U7x=U7y=0F2x= F2y= F3x= F3y=F6x=F6y=F8x= F8y= F9x= F9y=F10x=F10y= F11x= F11y= 0F5x= 0,F5y=5.解方程：>>k=[K(3:6,3:6),K(3:6,9:12),K(3:6,15:22);K(9:12,3:6),K(9:12,9:12),K(9:12,15:22) ;K(15:22,3:6),K(15:22,9:12) ,K(15:22,15:22)];K =+008 *Columns 1 through 80 00 00 0 00 0 00 0 00 0 00 0 0 0 0 00 0 0 0 0 00 0 0 0 0 00 0 0 0 0 00 0 0 0 0 0 0 00 0 0 0 0 0 0 00 0 0 0 0 0 0 00 0 0 0 0 0 0 00 0 0 0 0 0 0 00 0 0 0 0 0 0 00 0 0 0 0 0 0 00 0 0 0 0 0 0 0Columns 9 through 160 0 0 0 0 0 0 00 0 0 0 0 0 0 00 0 0 0 0 00 0 0 0 0 00 0 0 0 0 00 0 0 0 0 00 0 00 0 00 00 00 00 00 0 0 00 0 0 00 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0Columns 17 through 220 0 0 0 0 00 0 0 0 0 00 0 0 0 0 00 0 0 0 0 00 0 0 0 0 00 0 0 0 0 00 0 0 0 0 00 0 0 0 0 00 00 00 0 0 0 0 00 0 0 0 0 00 0 0 0 0 00 0 0 0 0 00 00 0>>f=[0;0;0;0;0;;0;0;0;0;0;0;0;0;0;0];>>u=k\fu=*[ ]T6.后处理>>U=[0;0;u(1:4);0;0;u(5:8);0;0;u(9:16)];>>F=K*UF=[ 0 0 ]T三.小结通过这次练习，对有限元结合MATLAB软件解决平面应力/应变问题的方法和过程逐渐清晰，特别是在应用MATLAB软件的过程中学到了很多东西。