stressconcentrationforholedefectbymatlab.rar_firstv54_孔应力_带圆孔缺陷平

%用三角形单元求解中心有洞的平板中心边缘受分布载荷时应力集中问题 %------------------------------------ % 输入控制参数 %------------------------------------ clear %清除workplace残留数据 nel=288; % 单元数 nnel=3; % 每个单元的节点数 ndof=2; % 每个节点的自由度数 nnode=168; % 系统总节点数 sdof=nnode*ndof; % 系统总自由度数 edof=nnel*ndof; % 每个单元的自由度数 E=1; % 杨氏弹性模量 u=0; % 泊松比 Q=1; % 外力集度 t=1; % 薄板厚度 %--------------------------------------------- % 从程序读入坐标数值 % zb(i,j) i：节点号　j：x,y值 %--------------------------------------------- %从程序读入坐标数值 jdzb=zeros(nnode,2); jdzb=getzb(nnode); printzb; %--------------------------------------------------------- % 从程序读入每个单元对应的节点号（逆时针排列） % dyjd(i,j) i:节点号 j:对应的单元号 %--------------------------------------------------------- dyjd=getdyjd(nel); printdyjd; %------------------------------------- % 输入边界条件 %------------------------------------- bcdof=getbcdof(); % 约束的自由度 bcval=zeros(1,length(bcdof)); % 对应的值 %---------------------------------------------- % 初始化矩阵和矢量 %---------------------------------------------- lhz=zeros(sdof,1); % 载荷矢量,矩阵荷载 hz=zeros(168,2); kk=zeros(sdof,sdof); % 系统刚度矩阵 disp=zeros(sdof,1); % 系统位移矢量 index=zeros(edof,1); % 每个单元所包含的自由度 %---------------------------- % 载荷矢量 %---------------------------- lhz=getlhz(Q,nnode,sdof,jdzb); printhz; %----------------------------------------------------------------- % 单元刚度矩阵计算及其组合 %----------------------------------------------------------------- for iel=1:nel % 对所有单元数的循环 %------------------------------------------------------ % 计算单元 刚度矩阵 %-------------------------------- k=zeros(edof,edof); % 初始化单元刚度矩阵。 mlwjgetk;%得到单元刚度矩阵。 %-------------------------------- getkz; % 合成系统刚度矩阵 end %----------------------------- % 加载边界条件 %----------------------------- [kk,lhz]=bjtj(kk,lhz,bcdof,bcval); %---------------------------- % 求解 %---------------------------- disp=kk\lhz; printdisp;% 输出节点位移 %---------------------------- % 求y应力 %---------------------------- for m =1:nel ss=getss(u,E,m,dyjd,jdzb); stress=getstress( ss,disp,dyjd,m); printstress; end