规则外形平面应力问题有限元分析.docVIP

规则外形的平面应力问题有限元分析 ㈠ 程序说明 本程序适用于规则外形受集中载荷作用的平面应力问题有限元分析，采用三节点三角形单元，不包括网格划分，因此前处理需要自行处理数据；根据最小位能原理进行求解，求解结果为节点位移，不包括后置处理。程序应用举例如下： ㈡ 程序举例 如图为矩形简支板，板厚为t=1m，板长为18m，板高为3m，受图示作用力，弹性模量，泊松比，容重。求位移？ ㈢ 程序求解过程 ⑴ 划分单元 划分单元，标出单元号码及节点，选取坐标。 ⑵ 输入数据 ①弹性模量：eo=2000000000； ②泊松比：co=0.167； ③梁的厚度：t=1m； ④节点总数：nnd=14； ⑤单元总数：nne=12； ⑥节点坐标：xynnd=[0 0;3 0;6 0;9 0;12 0;15 0;18 0;0 3;3 3;6 3;9 3;12 3;15 3;18 3]； ⑦单元节点编码：nonne=[1 2 9;2 3 10;3 4 11;4 5 12;5 6 13;6 7 14;1 9 8;2 10 9;3 11 10;4 12 11;5 13 12;6 14 13]； ⑧存在载荷的节点总数；w=7； ⑨已知位移的节点总数：d=2； ⑩载荷值：nwz=[0 -200000 8;0 -200000 9;0 -200000 10;0 -200000 11;0 -200000 12;0 -200000 13;0 -200000 14]； ⑩已知位移值：ndz=[0 0 1;0 0 7]。 ⑶ 程序框图及程序 ①程序框图 ②源程序代码 function [UV]=yxymainprogram() %有限元主程序，适用于规则外形受集中载荷作用的平面应力问题有限元分析。 %弹性模量eo %泊松比co %梁的厚度t %节点总数nnd %单元总数nne %节点坐标xynnd %单元节点编码nonne %存在载荷的节点总数w %已知位移的节点总数d %载荷值nwz %已知位移值ndz %载荷矩阵P %引入位移边界后的刚度矩阵KZ %引入位移边界后的载荷矩阵KP %总纲矩阵K %单元刚度矩阵KE %应变矩阵B %弹性矩阵D global UV; global KZ; global KP; global P; global w; global nnd; global nwz; global K; global B; global D; global A; global KE; global xynnd; global nonne; global eo; global co; global d; global ndz; clear; clc; [KZ,KP]=getKZP; UV=inv(KZ)*KP; function [KZ,KP]=getKZP(d) %引入位移边界后的矩阵 global UV; global KZ; global KP; global P; global w; global nnd; global nwz; global K; global B; global D; global A; global KE; global xynnd; global nonne; global eo; global co; global d; global ndz; d=input(d=); w=input(w=); nwz=input(nwz=); nne=input(nne=); nnd=input(nnd=); nonne=input(nonne=); xynnd=input(xynnd=); eo=input(eo=); co=input(co=); ndz=input(ndz=); [K]=getK(nne,nnd); [P]=getP(w); a=inf; for i=1:d ii=ndz(i,3); if ndz(i,1)==0 K(:,2*ii-1)=0; K(2*ii-1,:)=0; K(2*ii-1,2*ii-1)=1; P(2*ii-1,1)=0; else K(2*ii-1,2*ii-1)=a* K(2*ii-1,2*ii-1); P(2*ii-1,1)=a* K(2*ii-1,2*ii-1)*ndz(i,1); end if ndz(i,2)==0 K(:,2*ii)=0; K(2*ii,:)=0; K(2*ii,2*ii)=1; P(2*ii,1)=0; else K(2*ii,2*ii)=a*