四边形八节点等参元matlab程序

悬臂钢梁，尺寸如图一所示；v0.3。h1，E2.1e11. 图一 悬臂钢梁 图二 单元划分与结点编号 MatlabMatlab 输出结果输出结果 附录Ⅰ附录Ⅰ 有限元 ANSYS 分析结果 采用 PLANE183PLANE183 单元（四边形八节点）单元（四边形八节点）单元得出的结构 Y 向最大位 移为-0.216E-04-0.216E-04。。 约等于 matlab 平面四边形八节点等参元结点 Y 向最 大位移-2.4024E-5-2.4024E-5。。 附录Ⅱ附录Ⅱ ------------------------------四边形八节点等参元四边形八节点等参元 matlabmatlab 计算程序计算程序-------------------------------------------------------- 主主程程序序 ************************************************************************************************************************************** ************************************************************************ 20122012 年年 本程序只能处理集中荷载作用下的情况本程序只能处理集中荷载作用下的情况 只输出了节点位移、单元中心点的应力只输出了节点位移、单元中心点的应力 ************************************************************************************************************************************** ****************************** 变量说明 Evh 弹性模量泊松比厚度 NPOINNELEMNVFIXNNODENFPOIN 总结点数,单元数,约束结点个数,单元节点数 ,受力结点数 COORDLNODS 结构节点整体坐标数组,单元定义数组, FPOINFORCEFIXED 结点力数组，总体荷载向量, 约束信息数组 HKDISP 总体刚度矩阵,结点位移向量 ****************************** clear all at short e FP1bjd.txt,rt;打开数据文件 读入控制数据 EfscanfFP1,f,1;弹性模量 vfscanfFP1,f,1; 泊松比 hfscanfFP1,f,1;厚度 NELEMfscanfFP1,d,1;单元数 NPOINfscanfFP1,d,1; 总结点数 NNODEfscanfFP1,d,1;单元节点数 NFPOINfscanfFP1,d,1;受力结点数 NVFIXfscanfFP1,d,1;约束结点个数 LNODSfscanfFP1,f,[NNODE,NELEM]; 单元定义单元结点号 （逆时针） COORDfscanfFP1,f,[2,NPOIN]; 结点号 x,y 坐标整体坐标下 FPOINfscanfFP1,f,[3,NFPOIN]; 节点力结点号、X 方向力向右正，Y 方向力向上正 FIXEDfscanfFP1,d,[3,NVFIX]; 约束信息数组n,3 n受约束节点数目, n,1约束点号 n,2与n,3分别为约束点 x 方向和 y 方向的约束情况,受约束为 1 否则为 0 ******************************************************************* ******************************************************************* 平面应力问题的求解 ******************************************************************* ******************************************************************* 刚度矩阵的生成 计算刚度矩阵，并对约束条件进行处理 Kezeros2*NNODE,2*NNODE; 单元刚度矩阵并清零 HKzeros2*NPOIN,2*NPOIN; 张成总刚矩阵并清零 调用子程序 生成单元刚度矩阵 for m1NELEMm 为单元号 KeKE,v,h,. COORDLNODSm,1,1,COORDLNODSm,1,2,. COORDLNODSm,3,1,COORDLNODSm,3,2,. COORDLNODSm,5,1,COORDLNODSm,5,2,. COORDLNODSm,7,1,COORDLNODSm,7,2;调用单元刚度矩阵 aLNODSm,;临时向量,用来记录当前单元的节点编号 对总刚度矩阵的处理 for j18 for k18 HKaj*2-1aj*2,ak*2-1ak*2HKaj*2-1aj*2,ak*2-1ak*2. Kej*2-1j*2,k*2-1k*2; end end end 对荷载向量进行处理 FORCEzeros2*NPOIN,1; 张成总荷载向量并清零 for i1NFPOIN b1FPOINi,1*2-1;b2FPOINi,1*2;FPIONi,1为作用点 FORCEb1FPOINi,2;FPIONi,2为 x 方向的节点力 FORCEb2FPOINi,3;FPIONi,3为 y 方向的节点力 end 将约束信息加入总刚,总荷载 for i1NVFIX if FIXEDi,21 c12*FIXEDi,1-1; HKc1,0;将一约束序号处的总刚列向量清 0 HK,c10;将一约束序号处的总刚行向量清 0 HKc1,c11;将行列交叉处的元素置为 1 FORCEc1