有限元平面矩形单元MATLAB程序设计精品资料.docx

《有限元平面矩形单元MATLAB程序设计精品资料.docx》由会员分享，可在线阅读，更多相关《有限元平面矩形单元MATLAB程序设计精品资料.docx（109页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、苏州科技学院本科生毕业设计（论文）有限元平面矩形单元MATLAB程序设计摘 要本论文主要研究内容是有限元平面矩形单元的基本原理和MATLAB软件的图形用户界面及函数编程的基本知识，并根据有限元平面矩形单元的计算单元刚度矩阵，结点位移等原理和方法，运用MATLAB 语言编写弹性力学平面矩形单元的计算程序，使得程序能够完成对不同荷载类型，不同结构单元的计算，同时设计了图形用户界面（GUI），使得程序能够在图形用户界面（GUI）上对应的输入框内输入相应数据，并能够计算出单元刚度矩阵，等效结点荷载，有限元上结点位移和单元应力；最后将集中荷载和均布荷载情况下的实例进行对MATLAB计算得出的位移和单元结

2、点应力与ANSYS计算值比较分析，验证其正确性和通用性，并总结了MATLAB软件在有限元的运用。关键词 MATLAB程序设计；有限元；矩形单元；刚度矩阵；单元应力The Rectangle Unit Plane of Finite Element of the MATLAB ProgrammingAbstractThis dissertation mainly studies the radical principle of the rectangle unit plane of the finite element and the basic knowledge of the Graphic

3、al User Interface and the function programming of the MATLAB software . And make use of MATLAB language to make up the GUI and the computation program of MATLAB according to the theory and method of the rectangle unit plane of the finite element to compute the element stiffness matrix , nodes displa

4、cement and so on about different types of loads and different elements . Imputing the datas into the GUI to figure out element stiffness matrix , equivalent nodal load , displacement of nodes and element stress . At last , comparing and analyzing the nodes displacement and nodes stress of element of

5、 MATLAB and ANSYS to make sure whether the programming is right and used universally , and making a summary the application of MATLAB software in finite element analysis . Keywords MATLAB programming ; finite element analysis ; rectangle unit of plane ;stiffness matrix ; element stress105目 录第1章 绪论11

6、.1研究背景11.1.1 有限元单元法概述11.1.2 MATLAB程序概述21.1.3 MATLAB应用于有限元程序设计的优点21.2 研究内容31.3 研究重点3第2章 矩形单元有限元原理42.1 位移函数42.2 几何矩阵52.3 应力矩阵62.4 单元刚度矩阵72.5 等效结点荷载82.6 整体等效结点荷载82.7 整体刚度矩阵82.8 结点位移的计算8第3章 MATLAB 图形用户界面（GUI）93.1 GUI界面的概述93.1.1 控件对象的描述103.1.2 部分常用控件通用属性113.2 GUI常用函数语句123.2.1 get和set语句123.2.2 if判断语句和for循

7、环语句133.3 M文件的一般结构13第4章 有限元平面矩形单元MATLAB程序设计154.1 图形用户界面设计154.2 设置控件的相关属性164.3 编写代码18第5章 应用实例315.1 9个单元16个结点计算335.2 6个单元13个结点计算395.4 均布荷载下9个相等单元44结 论51致 谢52参 考 文 献53附录A 译文54工程结构构件的有限元模型和剪应力分析54附录B 外文原文70第1章 绪论1.1研究背景1.1.1 有限元单元法概述 有限元法是伴随着电子计算机技术的进步而发展起来的一种工程数值计算方法，其适应性很强，可以解决各种各样的复杂工程问题。有限元法的基本思想就是将任

8、意一个连续体离散化为一组有限个、按一定方式互相连接且几何形状简单的单元的组合，单元本身又可以有不同形状，使复杂问题简化从而得到问题的解。有限元法作为数值分析方法的另一个重要特点是利用在每一个单元内假设的近似函数来分片地表示全求解区域上待求的未知场函数。单元内的近似函数通常由未知场函数，或这些函数与其导数在单元各个节点上的数值和其插值函数来表达。这样一来，未知场函数或未知场函数及其导数在各个节点上的数值就成为新的未知量（即自由度），从而使一个连续的无限自由度问题变成离散的有限自由度问题，一经求解出这些未知量，就可以通过插值函数计算出各个单元内场函数的近似值，从而得到整个求解域上的近似解。显然随着

9、单元数目的增加，解得近似程度将不断改进。如果单元是满足收敛要求的，近似解最后将收敛于精确解。有限元法的数学逻辑严谨，物理概念清晰，易于理解和掌握，应用范围广泛，能够灵活的处理和求解各种复杂问题，特别是它采用矩阵形式表达基本公式，便于运用计算机编程运算，这些优点赋予了有限单元法强大的生命力。一、结构物的离散化离散化是将待分析的结构物从集合上用线或面划分为有限个单元，即将结构物看出有限个单元构成的组合体。按结构物形状的不同和分析的要求，选取不同形式的单元，通常在单元的边界上设置结点，结点联接相邻的单元。二、单元分析所谓单元分析就是设法导出单元的结点位移和结点力之间的关系，即建立单元刚度矩阵。在分析

10、弹性力学平面问题时，每个平面单元内的任意点的位移需要按一定的函数关系用结点位移来表示，这种函数称为位移函数。选择的位移函数应保证解的收敛性，因此，建立合理的位移函数是单元分析点的关键。位移函数确定以后就可以利用弹性力学的基本方程推导出单元刚度矩阵。此外，还需要按静力等效原则将作用于每个单元上的外力简化到节点上，构成等效结点力。对于每一个单元进行上述分析以后，可建立起单元刚度方程。三、整体分析整体分析就是将各个单元组成结构整体进行分析。整体分析的目的在于导出整个结构结点位移与结点力之间的关系，即建立整个结构的刚度方程。整体分析的步骤为：首先按着一定的集成规则，将各个单元刚度矩阵集合成结构整体刚度

11、矩阵，并将单元等效结点荷载集合成整体等效节点荷载列阵；然后引入结构的位移边界条件，求解整体平衡方程，得出基本未知量结点位移列阵；最后计算各个单元的内力和变形1。1.1.2 MATLAB程序概述MATLAB是一种以矩阵为基本单位的高效数值计算语言，是一个集科学计算、图像处理、声音处理于一体的高度集成系统。在编程效率、程序可读性、可移植性和可扩充性上MATLAB远远优于其它的高级编程语言，而且其编程易学、直观，代码非常符合人们的思维习惯。另外MATLAB为用户提供了丰富的图形界面设计方法，使用户能够在利用其强大的数值计算功能的同时可设计出友好的图形界面，它受到了越来越多的用户的欢迎，成为了当今国际

12、上最流行的计算机辅助设计软件2-5。1.1.3 MATLAB应用于有限元程序设计的优点一、语言简洁紧凑、使用方便灵活、库函数极其丰富，几乎包括有限元编程中的所有基础程序，使开发者可以直接进行高、精、尖的研究。MATLAB是一个强大的数值计算软件，在数值计算方面，除了包括基本的数学函数、基本矩阵和数组运算函数以及多种插值函数之外，而且具有矩阵的求逆、LU分解、OR分解、矩阵指数等几乎所有的矩阵函数及矩阵分析函数。另外它还包括强大的稀疏矩阵的存储、初等变换、分解、特征值和奇异值的求法等功能，以及提供了系数阵为稀疏矩阵的线性方程组的各种解法。MATLAB提供的所有这些数值计算方面的功能，对有限元编程

13、中的数据存储、单元刚度矩阵生成、刚度矩阵集成、线性方程组的求解等方面大有益处，根本无需编程人员去编制有限元中基础程序，大大减少了工作量及提高了编程效率。二、可视化建模及图形功能强大。MATLAB可以给出数据的二维、三维、乃至四维等数据表现以及绘制一般科技绘图软件所能绘制的几乎所有图形，如曲线图、网格图、等直线图、表面图等，MATLAB这些功能为有限元后置处理数据可视化提供了充分的表现力度，这往往是有限元数据处理中最为困难的事。除此之外，MATLAB除了具有较强的应用程序应用界面编制功能，而且提供了专门的界面编制工具(GUI)，可以利用它编制理想的用户界面。三、含有多种学科的工具箱(ToolBo

14、x)以及程序代码的公开性，是MATLAB的一大特色。MATLAB提供的学科性工具箱可以使用户无需编写自己学科范围内的基础程序，太大节约编程时间。例如，MATLAB工具箱中的偏微分方程工具箱(PDE)就是利用MATLAB编制好的有限元基础程序库。四、程序可移植性好，MATLAB几乎可以在各种机型和操作系统上运行，所以它的可移植性和可扩充性MATLAB远远优于其它的高级编程语言6-10。1.2 研究内容本论文是根据有限元平面矩形单元的计算原理和方法，运用MATLAb语言的语法规则、程序结构和编程技巧，编写弹性力学平面矩形单元有限元计算程序和平面矩形单元有限元计算的图形用户界面；并通过在界面上输入不

15、同单元类型和荷载类型的数据，使得能够计算出单元的刚度矩阵、结点位移和单元应力。最后通过与其他程序软件（如ANSYS）的计算结果相比较分析，来验证本论文研究计算程序的相对精确性和通用性。1.3 研究重点1，深入了解有限元矩形单元计算刚度矩阵，结点位移等的基本原理；2，熟练掌握MATLAB语言程序的语法规则和编程技巧；3，编写图形用户界面和弹性力学平面矩形单元MATLAB的计算程序，正确计算出结点位移，刚度矩阵和单元应力。第2章 矩形单元有限元原理y矩形单元是以4个顶点为结点，每个结点有2个位移分量，共有8个结点位移分量。设图2.1所示位移一矩形单元，变长分别为2a和2b，四个结点以i、j、m、p

16、表示。为了简便，取单元形心为坐标原点，x轴、y轴与单元边界平行。ijmpox2b2a图2.1 矩形单元2.1 位移函数矩形单元的位移函数为u=1+2x+3y+4xyv=5+6x+7y+8xy (2-1)对于i、j、m、p4结点，将相应结点坐标带入上式，有ui=1a2b3+ab4uj=1+a2b3ab4um=1+a2+b3+ab4up=1a2+b3ab4 a vi=5a6b7+ab8vj=5+a6b7ab8vm=5+a6+b7+ab8vp=5a6+b7ab8 b由式（a）解出1、2、3、4；由式（b）解出5、6、7、8。将这8个常数带入式（2-1），可得u=Niui+Njuj+Nmum+Npup

17、v=Nivi+Njvj+Nmvm+Npvp 22其中Ni=141xa1ybNj=141+xa1ybNm=141+xa1+ybNp=141xa1+yb 23式（1-2）改写为矩阵形式fe=Ne 24式中fe=uvT N=Ni00NiNj00NjNm00NmNp00Np 25 e=uiviujvjumvmupvpT 26上述位移函数可以充分保证解答的收敛性。这是由于：式 （2-1）中包含1、2、3以及5、6、7，反映了刚体位移和常量应变；单元4条边界分别平行于x、y轴。在平行于x轴的边界上y=b，位移分量为坐标x的线性函数；在平行于y轴的边界上x=a，位移分量是y的线性函数，因此，可以保证两相邻单

18、元在公共边界上的位移是连续的1。2.2 几何矩阵 利用弹性力学平面问题点的几何方程和式（2-4），求得单元应变31e=B3881e 27 式中：仍然由x、y、xy组成的列阵；几何矩阵B为由结点位移求单元应变的转换矩阵，即 B=x00yyxN或改写成 B=BiBjBmBp 28 而 Bi=Nix00NiyNiyNix i，j，m，p 29上式中利用脚标i、j、m、p进行轮换，用以表示Bi、Bj、Bm、Bp。分别计算Niy，Nixi，j，m，p并带入式（2-8），可得B=14abby00axaxbyby00a+xa+xbyb+y00a+xa+xb+yb+y00axaxb+y 210 2.3 应力矩

19、阵利用弹性力学平面问题的物理方程=D，并将式（2-7）带入可得 31e=D33B3881e=S3881e 211 式中：仍然是由x、y、xy 组成的列阵；D为弹性矩阵，应力矩阵S为结点位移求应力的转换矩阵。S可写成 S=SiSjSmSp 212 而 Si=DBi i，j，m，p 113 D=E1210100012S=E4ab12byaxbyax12ax12bybya+xbya+x12a+x12by b+ya+xb+ya+x12a+x12b+yb+yaxb+yax12ax12b+y 214整理应力结果时，通常需要计算单元4个结点处的应力，将结点坐标值带入式（2-14）即可求出4个结点上的应力矩阵

20、，分别以Si、Sj、Sm、Sp表示。Si=E4ab122b2a2b2a1a1b2b02b001b00000002a02a1a0Sj=E4ab122b02b001b2b2a2b2a1a1b02a02a1a0000000Sm=E4ab1200000002a02a1a02b2a2b2a1a1b2b02b001bSp=E4ab1202a02a1a00000002b02b001b2b2a2b2a1a1b 215 2.4 单元刚度矩阵 矩形单元的单元刚度矩阵ke为88方阵，利用下式 ke=BTDBtdxdy 或 ke=BTStdxdy并带入弹性矩阵D、几何矩阵B、或几何矩阵B和应力矩阵S，可求出ke=Et

21、1213ba+16ab1+813ab+16ba13ba+112ab138对称13ba+16ab13816ab16ba16ba112ab1+81+816ab112ba1+813ab+16ba16ba16ab13813813ab+112ba16ba16ab13813813ab+112ba16ba112ab1+81+816ab112ba13ba+16ab1+813ab+16ba对称13ba+112ab13813816ab16ba13ba+16ab1+813ab+16ba 2162.5 等效结点荷载1，单元自重对于单元自重W，移置到每一个结点上都是W4，因此Re=W014014014014T 2172

22、，分布面力如果单元的任一边界上受有分布面力的作用，应将该分布面力的合力按静力等效原则向这条边界的两个结点上移置1。2.6 整体等效结点荷载将单元的等效结点力以“对号入座”的方式叠加，同时考虑结点上作用的集中反力，求得整体等效结点荷载列阵R。2.7 整体刚度矩阵建立了单元刚度矩阵ke之后，就可以利用刚度集成法将其集合成整体刚度矩阵K。具体做法是按照局部编号与与整体编号的对应关系，“对号入座”。2.8 结点位移的计算求解整体平衡方程K=R，考虑边界条件，引入支承条件后，可解方程得出未知结点位移，最后可得整体结点位移列阵1。第3章 Matlab 图形用户界面（GUI）3.1 GUI界面的概述图形用户

23、界面（GUI）是由窗口、光标、按键、菜单、文字说明等对象（Objects）构成的一个用户界面，用户通过一定的方法选择、激活这些对象，使计算机产生某种动作或变化，如实现计算、绘图等。创建MATLAB用户图形界面必须具有以下三类基本元素：组件、图形窗口、回应。实现一个GUI的过程包括两个基本任务：一是GUI的组件布局，另一个是GUI组件编程。这些任务都可以通过图形用户界面开发环境来完成。在Matlab中，GUIDE是一个组件布局工具集，能够生产所需的组件资源并保存在一个FIG文件中，如图3.1(a)所示；其次还可以生成一个包含GUI初始化和发布控制代码的M文件，如图3.1(b)所示，该文件作为回调

24、函数提供了一个框架11。GUI工具栏GUI对象Figure窗口(b)(a)图3.1 GUI编辑界面(a) FIG文件(b) M 文件3.1.1 控件对象的描述MATLAB中的控件大致可分为两种，一种为动作控件，鼠标点击这些控件时会产生相应的响应。一种为静态控件，是一种不产生响应的控件，如文本框等。（1） 按钮(Push Buttons) ：执行某种预定的功能或操作；图3.2 pushbutton回调函数在上图3.2语句中，pushbutton1_Callback(hObject, eventdata, handles)是用来执行pushbutton1按钮下的程序的调用语句。（2） 文本编辑器(

25、Editable Texts) ：用来使用键盘输入字符串的值，可以对编辑框中的内容进行编辑、删除和替换等操作；（3） 静态文本框(Static Texts) ：仅仅用于显示单行的说明文字；（4） 表格(uitable) ：用于数据的可视化；通过表格的属性设置窗口可以改变表格的行数，列数等各种参数值，如图3.3；也可以通过下面的编程来控制表格：column=str2num(get(handles.column_edit,String); %获取行数row=str2num(get(handles.row_edit,String); %获取列数num_elem=cell(column,row); %

26、设置num_elem的单元数为element9num_elem(:,:)=; %设置num_elem数组为空set(handles.uitable,Data,num_elem); %设置表格的句柄set(handles.uitable,ColumnEditable,true(1,column); %设置表格可编辑图3.3 表格属性编辑框（5） 面板(panel) ：用于把某些相关的交互控件组织在同一区域内，这样可以提高GUI界面的组织层次和易用性。3.1.2 部分常用控件通用属性每种控件都有一些可以设置的参数，用于表现控件的外形、功能及效果，即属性。属性由两部分组成：属性名和属性值，它们必须是

27、成对出现的。当用户创建好GUI界面需要的各个交互控件如图3.4，并调整好位置后，就需要设置这些控件的属性，双击控件会弹出该控件的属性设置窗口，如图3.5所示，用户可以通过控件设置窗口更加详细的了解和设置不同控件的属性。在一般的GUI界面设计中，各种控件的常用属性包括颜色设置、文字标签设置、回调函数设置等。图3.4 GUIDE 交互控件图3.5 组件属性设置窗口表3.1 部分常用控件的通用属性属性说明属性说明BackgroundColor背景色（空间颜色）Tag用户自定义的标记文字ForegroundColor前景色（标签文字颜色）Value控件当前取值FontName字体名称（标签文字）Sel

28、ected选中状态FontSize字号Visible可见状态FonfUnits字号单位Enable可用状态Sring标签文字Callback激活控件的回调函数Style交互控件种类Createfvn创建控件的回调函数3.2 GUI常用函数语句3.2.1 get和set语句get(handles, 属性名称)：获得特定属性名称对应的属性值；set(handles, 属性名称, 属性值)：设置句柄值为handles的对象属性。例如：set(handles.uitable3,ColumnEditable,true(1,2)设置句柄uitable3的第一列到第二列是可编辑的num=get(handle

29、s.num_edit,String)获得句柄num_edit内的字符串文字，并赋值给变量num。3.2.2 if判断语句和for循环语句1，If 结构的语法形式2如下：If logical_expressionstatementselseif logical_expressionstatementsend其中if elseif构成的各项分支里面，哪个的条件满足（逻辑表达式logical_expression的结果为真），就执行哪个分支后面紧跟的程序语句。2，for循环2用于知道循环次数的情况，其语法格式为：for index=start:increment:endstatementsendin

30、dex为循环变量，increment为增量，end用于判断循环是否应该终止。3.3 M文件的一般结构 如图3.6所示，MATLAB提供的一个函数M文件的全部内容，图中清晰显示了一般的M文件包括的各部分结构。Callback程序编辑区Callback程序编辑区图3.6 M文件编辑器界面Pushbutton为GUI最常使用也是最简单的对象，当用户单击压下pushbutton时，MATLAB即会立即依据其对应的Callback程序来执行操作。第4章 有限元平面矩形单元Matlab程序设计根据任务书上要求该程序具有如下功能及要求：1、 给定矩形单元，编写程序计算平面矩形单元的单元刚度矩阵；2、 给定荷

31、载作用下，计算结点位移；3、 计算有限单元上的单元应力，单元结点应力；4、 编写matlab语言的图形用户界面，使得程序的执行能够在图形用户界面下运行。4.1 图形用户界面设计在布局编辑器中布置控件，并使用几何位置排列工具对控件的位置进行调整如图4.1。1建立10个静态文本，用于显示程序标题，标注相应控件的提示和标注等内容；2建立7个可编辑文本框，用于输入各参数值；3建立5个表格框，用于输入单元的信息和显示计算位移，单元应力结果；4建立4个按钮，用来定义表格大小，计算，清除数据和结束程序。5界面设计如下图4.1所示，保存为untitled.fig文件。图4.1 图形用户界面草图4.2 设置控件

32、的相关属性控件的标识Tag是对于各控件的识别，每个控件创建时都会由开发环境自动产生一个一个标识，但在程序设计中，为了编辑、记忆和维护的方便，一般会为控件设置新的标识。1， 设置第1个按钮的“Tag”为“shuru_pushbutton，“String”为“输入单元、结点信息”，用来定义uitable的大小，如图4.2所示；图4.2 按钮shuru_pushbutton属性设置2， 设置第2个按钮的“Tag”为“jisuanxianshi_pushbutton”，“String”为“计算并显示结果”，用来计算刚度矩阵，结点位移，节点荷载，单元应力等并显示计算结果；3， 设置第3个按钮的“Tag”

33、为“clear_pushbutton”，“String”为“清除数据”，用来清除编辑文本框内的数据；4， 设置第4个按钮的“Tag”为“quit_pushbutton”，“String”为“退出”，用来结束程序；5， 设置10个静态文本的“String”，分别为图显示的内容，并分别定义其字体大小；6， 设置7个可编辑文本框的“Tag”分别为“e_edit”、“u_edit”、“t_edit”、“element_edit”、“node_edit”、“restriction_edit”、“num_edit”，用来输入计算相关参数；7， 设置3个表格框（uitable）的属性，如图4.3所示，用来输

34、入单元结点的信息。图4.3 表格行数、列数的设置属性4.3 编写代码开始读取数据计算单元刚度矩阵形成整体刚度矩阵计算有效结点荷载引入支承条件解平衡方程得出结点位移计算单元应力结束图4.4 程序流程图程序标识符说明element 单元数目 e 弹性模量node 结点数目 u 泊松比res 约束结点数 t 单元厚度k 单元刚度矩阵 num 结点荷载数K 整体刚度矩阵 sigma 单元应力f 结点位移 PM 等效结点荷载完成程序中变量的赋值等工作后，打开untitled.m文件，系统自动生成M文件程，按流程图4.4编写序代码如下：function varargout = untitled(varar

35、gin)gui_Singleton = 1;gui_State = struct(gui_Name, mfilename, . gui_Singleton, gui_Singleton, . gui_OpeningFcn, untitled_OpeningFcn, . gui_OutputFcn, untitled_OutputFcn, . gui_LayoutFcn, , . gui_Callback, );if nargin & ischar(varargin1) gui_State.gui_Callback = str2func(varargin1);end if nargout var

36、argout1:nargout = gui_mainfcn(gui_State, varargin:);else gui_mainfcn(gui_State, varargin:);end%定义可编辑文本框为空及更新句柄结构：function untitled_OpeningFcn(hObject, eventdata, handles, varargin)set(handles.t_edit,String,);set(handles.e_edit,String,);set(handles.u_edit,String,);set(handles.element_edit,String,);se

37、t(handles.node_edit,String,);set(handles.num_edit,String,);set(handles.restriction_edit,String,);handles.output = hObject;guidata(hObject, handles);%按钮shuru_pushbutton定义uitable的行数和列数：function shuru_pushbutton_Callback(hObject, eventdata, handles)element=str2num(get(handles.element_edit,String);res=s

38、tr2num(get(handles.restriction_edit,String);num=str2num(get(handles.num_edit,String);num_elem2=cell(element,9);num_elem2(:,:)=;num_elem3=cell(res,2);num_elem3(:,:)=;num_elem4=cell(num,3);num_elem4(:,:)=;set(handles.uitable2,Data,num_elem2);set(handles.uitable2,ColumnEditable,true(1,9);set(handles.ui

39、table3,Data,num_elem3);set(handles.uitable3,ColumnEditable,true(1,2);set(handles.uitable4,Data,num_elem4);set(handles.uitable4,ColumnEditable,true(1,3);按钮jisuanxianshi_pushbutton的代码编写：function jisuanxianshi_pushbutton_Callback(hObject, eventdata, handles)%计算单元刚度矩阵及整体刚度矩阵，并使计算结果分别显示到Excel：e=str2num(g

40、et(handles.e_edit,String);u=str2num(get(handles.u_edit,String);t=str2num(get(handles.t_edit,String);element=str2num(get(handles.element_edit,String);node=str2num(get(handles.node_edit,String);K=zeros(2*node,2*node);data2=get(handles.uitable2,data);for n=1:element; i=str2num(data2n,2); j=str2num(data

41、2n,3); m=str2num(data2n,4); p=str2num(data2n,5); a=str2num(data2n,6); b=str2num(data2n,7); c1=1/3*b/a+(1-u)/6*a/b; c2=1/3*a/b+(1-u)/6*b/a; d1=(1+u)/8; d2=(1-3*u)/8; e1=-1/3*b/a+(1-u)/12*a/b; e2=-1/3*a/b+(1-u)/12*a/b; f1=-1/6*b/a-(1-u)/12*a/b; f2=-1/6*a/b-(1-u)/12*a/b; g1=1/6*b/a-(1-u)/6*a/b; g2=1/6*

42、a/b-(1-u)/6*a/b; k1=c1,d1,e1,-d2,f1,-d1,g1,d2; d1,c2,d2,g2,-d1,f2,-d2,e2; e1,d2,c1,-d1,g1,-d2,f1,d1; -d2,g2,-d1,c2,d2,e2,d1,f2; f1,-d1,g1,d2,c1,d1,e1,-d2; -d1,f2,-d2,e2,d1,c2,d2,g2; g1,-d2,f1,d1,e1,d2,c1,-d1; d2,e2,d1,f2,-d2,g2,-d1,c2; k=e*t/(1-u2)*k1 xlswrite(刚度矩阵.xlsx,k,n) K(2*i-1,2*i-1)=K(2*i-1,2*i-1)+k(1,1); K(2