📄 plot_results_tri.m
字号:
function plot_results_tri ( in_data, SIGsys, resp, dof, NF, showPS)
%=============================================================
warning off
maxX = max(in_data.ND(:,2)); minX = min(in_data.ND(:,2));
maxY = max(in_data.ND(:,3)); minY = min(in_data.ND(:,3));
labx = (maxX / 10); laby = (maxY / 10);
labx = max([labx laby]); laby = labx;
if 1
if in_data.EL(1,2)==4
MAX_SIG = (SIGsys(1:3:length(SIGsys))+SIGsys(2:3:length(SIGsys))).*0.5 +...
sqrt((SIGsys(1:3:length(SIGsys))-SIGsys(2:3:length(SIGsys))).^2 +...
SIGsys(3:3:length(SIGsys)).^2.*4);
MIN_SIG = (SIGsys(1:3:length(SIGsys))+SIGsys(2:3:length(SIGsys))).*0.5 -...
sqrt((SIGsys(1:3:length(SIGsys))-SIGsys(2:3:length(SIGsys))).^2 +...
SIGsys(3:3:length(SIGsys)).^2.*4);
maxSIGmax = max(MAX_SIG); minSIGmax = min(MAX_SIG);
maxSIGmin = max(MIN_SIG); minSIGmin = min(MIN_SIG);
SIG_nd = zeros(size(in_data.ND,1),1); SIG_nd_nr = zeros(size(in_data.ND,1),1);
end
if in_data.EL(1,2)==5 | in_data.EL(1,2)==51
MAX_SIG = (SIGsys(1:3:end,:)+SIGsys(2:3:end,:)).*0.5 +...
sqrt((SIGsys(1:3:end,:)-SIGsys(2:3:end,:)).^2 + SIGsys(3:3:end,:).^2.*4);
MIN_SIG = (SIGsys(1:3:end,:)+SIGsys(2:3:end,:)).*0.5 -...
sqrt((SIGsys(1:3:end,:)-SIGsys(2:3:end,:)).^2 + SIGsys(3:3:end,:).^2.*4);
maxSIGmax = max(MAX_SIG); minSIGmax = min(MAX_SIG);
maxSIGmin = max(MIN_SIG); minSIGmin = min(MIN_SIG);
SIG_nd = zeros(size(in_data.ND,1),1); SIG_nd_nr = zeros(size(in_data.ND,1),1);
end
for i=1:size(in_data.EL)
if in_data.EL(i,2)==4
node1 = find(in_data.ND(:,1)==in_data.EL(i,3));
node2 = find(in_data.ND(:,1)==in_data.EL(i,4));
node3 = find(in_data.ND(:,1)==in_data.EL(i,5));
SIG_nd(node1,1:3) = SIG_nd(node1,1) + MAX_SIG(i)';
SIG_nd(node2,1:3) = SIG_nd(node2,1) + MAX_SIG(i)';
SIG_nd(node3,1:3) = SIG_nd(node3,1) + MAX_SIG(i)';
SIG_nd_nr(node1,1:3) = SIG_nd_nr(node1,1) + [1];
SIG_nd_nr(node2,1:3) = SIG_nd_nr(node2,1) + [1];
SIG_nd_nr(node3,1:3) = SIG_nd_nr(node3,1) + [1];
end
if in_data.EL(i,2)==5
node1 = find(in_data.ND(:,1)==in_data.EL(i,3));
node2 = find(in_data.ND(:,1)==in_data.EL(i,4));
node3 = find(in_data.ND(:,1)==in_data.EL(i,5));
node4 = find(in_data.ND(:,1)==in_data.EL(i,6));
SIG_nd(node1,1:3) = SIG_nd(node1,1) + MAX_SIG(i,2)';
SIG_nd(node2,1:3) = SIG_nd(node2,1) + MAX_SIG(i,3)';
SIG_nd(node3,1:3) = SIG_nd(node3,1) + MAX_SIG(i,4)';
SIG_nd(node4,1:3) = SIG_nd(node4,1) + MAX_SIG(i,5)';
SIG_nd_nr(node1,1:3) = SIG_nd_nr(node1,1) + [1];
SIG_nd_nr(node2,1:3) = SIG_nd_nr(node2,1) + [1];
SIG_nd_nr(node3,1:3) = SIG_nd_nr(node3,1) + [1];
SIG_nd_nr(node4,1:3) = SIG_nd_nr(node4,1) + [1];
end
end
SIG_nd_max = SIG_nd./SIG_nd_nr;
figure(NF); hold off;
% ****************************************************************************************
for i=1:size(in_data.EL)
if in_data.EL(i,2)==4
lp=subplot(121);
node1 = find(in_data.ND(:,1)==in_data.EL(i,3));
node2 = find(in_data.ND(:,1)==in_data.EL(i,4));
node3 = find(in_data.ND(:,1)==in_data.EL(i,5));
trisurf([1 2 3],in_data.ND([node1 node2 node3],2),...
in_data.ND([node1 node2 node3],3),SIG_nd_max([node1 node2 node3],1));
hold on;
end
if in_data.EL(i,2)==5 | in_data.EL(i,2)==51
lp=subplot(121);
node1 = find(in_data.ND(:,1)==in_data.EL(i,3));
node2 = find(in_data.ND(:,1)==in_data.EL(i,4));
node3 = find(in_data.ND(:,1)==in_data.EL(i,5));
node4 = find(in_data.ND(:,1)==in_data.EL(i,6));
trisurf([1 2 3 4],in_data.ND([node1 node2 node3 node4],2),...
in_data.ND([node1 node2 node3 node4],3), SIG_nd_max([node1 node2 node3 node4],1) );
hold on;
end
end
% ****************************************************************************************
% ****************************************************************************************
shading interp; axis equal; colormap(hsv); lighting phong;
axis off; colorbar('vert'); title('Stress field - MAX (tension)'); view(2); hold on
if isfield(in_data.mater,'ex')
axis([min(in_data.ND(:,2)) max(in_data.ND(:,2)) ...
min(in_data.ND(:,3)) max(in_data.ND(:,3)) ...
min(in_data.mater.rb_) max(in_data.mater.rb) ]);
caxis([ min(in_data.mater.rb_) max(in_data.mater.rb) ]);
end
SIG_nd = zeros(size(in_data.ND,1),1); SIG_nd_nr = zeros(size(in_data.ND,1),1);
for i=1:size(in_data.EL)
if in_data.EL(i,2)==4 % "4" - CST element (triangle)
node1 = find(in_data.ND(:,1)==in_data.EL(i,3));
node2 = find(in_data.ND(:,1)==in_data.EL(i,4));
node3 = find(in_data.ND(:,1)==in_data.EL(i,5));
SIG_nd(node1,1:3) = SIG_nd(node1,1) + MIN_SIG(i)';
SIG_nd(node2,1:3) = SIG_nd(node2,1) + MIN_SIG(i)';
SIG_nd(node3,1:3) = SIG_nd(node3,1) + MIN_SIG(i)';
SIG_nd_nr(node1,1:3) = SIG_nd_nr(node1,1) + [1];
SIG_nd_nr(node2,1:3) = SIG_nd_nr(node2,1) + [1];
SIG_nd_nr(node3,1:3) = SIG_nd_nr(node3,1) + [1];
end
if in_data.EL(i,2)==5 | in_data.EL(i,2)==51
node1 = find(in_data.ND(:,1)==in_data.EL(i,3));
node2 = find(in_data.ND(:,1)==in_data.EL(i,4));
node3 = find(in_data.ND(:,1)==in_data.EL(i,5));
node4 = find(in_data.ND(:,1)==in_data.EL(i,6));
SIG_nd(node1,1:3) = SIG_nd(node1,1) + MIN_SIG(i,2)';
SIG_nd(node2,1:3) = SIG_nd(node2,1) + MIN_SIG(i,3)';
SIG_nd(node3,1:3) = SIG_nd(node3,1) + MIN_SIG(i,4)';
SIG_nd(node4,1:3) = SIG_nd(node4,1) + MIN_SIG(i,5)';
SIG_nd_nr(node1,1:3) = SIG_nd_nr(node1,1) + [1];
SIG_nd_nr(node2,1:3) = SIG_nd_nr(node2,1) + [1];
SIG_nd_nr(node3,1:3) = SIG_nd_nr(node3,1) + [1];
SIG_nd_nr(node4,1:3) = SIG_nd_nr(node4,1) + [1];
end
end
SIG_nd_min = SIG_nd./SIG_nd_nr;
figure(NF); hold on;
% ****************************************************************************************
for i=1:size(in_data.EL)
if in_data.EL(i,2)==4
rp=subplot(122);
node1 = find(in_data.ND(:,1)==in_data.EL(i,3));
node2 = find(in_data.ND(:,1)==in_data.EL(i,4));
node3 = find(in_data.ND(:,1)==in_data.EL(i,5));
trisurf([1 2 3],in_data.ND([node1 node2 node3],2),...
in_data.ND([node1 node2 node3],3),SIG_nd_min([node1 node2 node3],1));
hold on;
end
if in_data.EL(i,2)==5 | in_data.EL(i,2)==51
rp=subplot(122);
node1 = find(in_data.ND(:,1)==in_data.EL(i,3));
node2 = find(in_data.ND(:,1)==in_data.EL(i,4));
node3 = find(in_data.ND(:,1)==in_data.EL(i,5));
node4 = find(in_data.ND(:,1)==in_data.EL(i,6));
trisurf([1 2 3 4],in_data.ND([node1 node2 node3 node4],2),...
in_data.ND([node1 node2 node3 node4],3), SIG_nd_min([node1 node2 node3 node4],1) );
hold on;
end
end
% ****************************************************************************************
shading interp; axis equal; colormap(hsv); axis off;
colorbar('vert'); title('Stress field - MIN (compression) '); view(2); hold on
set(NF,'name',[' Stress Field. MAX = ' num2str(max(SIG_nd_max(:,1))) ...
' MIN = ' num2str(min(SIG_nd_min(:,1)))],'NumberTitle','off','position',[ 50 50 800 420 ]);
if isfield(in_data.mater,'ex')
axis([min(in_data.ND(:,2)) max(in_data.ND(:,2)) min(in_data.ND(:,3)) max(in_data.ND(:,3)) ...
min(in_data.mater.rb_) max(in_data.mater.rb) ]);
caxis('manual');
caxis([ min(in_data.mater.rb_) max(in_data.mater.rb) ]);
end
figure(NF); hold on;
for i=1:size(in_data.EL)
if in_data.EL(i,2)==4
yy = (MAX_SIG(i) - minSIGmax)*1/(maxSIGmax - minSIGmax);
node1 = find(in_data.ND(:,1)==in_data.EL(i,3));
node2 = find(in_data.ND(:,1)==in_data.EL(i,4));
node3 = find(in_data.ND(:,1)==in_data.EL(i,5));
Xi = [1 1 1; in_data.ND(node1,2) in_data.ND(node2,2) in_data.ND(node3,2); ...
in_data.ND(node1,3) in_data.ND(node2,3) in_data.ND(node3,3)]*[1/3; 1/3; 1/3];
end;
if in_data.EL(i,2)==5 | in_data.EL(i,2)==51
yy = (MAX_SIG(i,1) - minSIGmax(1))*1/(maxSIGmax(1) - minSIGmax(1));
node1 = find(in_data.ND(:,1)==in_data.EL(i,3));
node2 = find(in_data.ND(:,1)==in_data.EL(i,4));
node3 = find(in_data.ND(:,1)==in_data.EL(i,5));
node4 = find(in_data.ND(:,1)==in_data.EL(i,6));
xi = 0; niu=0;
N1 = 1/4*(1-xi)*(1-niu); % shape function of quadrilateral
N2 = 1/4*(1+xi)*(1-niu); % shape function of quadrilateral
N3 = 1/4*(1+xi)*(1+niu); % shape function of quadrilateral
N4 = 1/4*(1-xi)*(1+niu); % shape function of quadrilateral
Xi = [1; in_data.ND(node1,2)*N1 + in_data.ND(node2,2)*N2 + in_data.ND(node3,2)*N3 + in_data.ND(node4,2)*N4;...
in_data.ND(node1,3)*N1 + in_data.ND(node2,3)*N2 + in_data.ND(node3,3)*N3 + in_data.ND(node4,3)*N4];
end;
angle_max = 0.5*atan(2*(SIGsys(i*3))/(SIGsys(i*3-2)-SIGsys(i*3-1)));
vv = (SIGsys(i*3-2)) - (SIGsys(i*3-1));
yi = abs(cos(angle_max)*yy*laby);
xi = abs(sin(angle_max)*yy*labx);
angle_max = angle_max*360/(2*pi);
ANm(i) = angle_max;
if vv>0 txy=xi*laby/labx; xi=yi*labx/laby; yi=txy; end;
if MAX_SIG(i)>0 cv='b-'; cv='k-'; else cv='r-'; cv='k-'; end;
% ****************************************************************************************
if showPS(1)==1
figure(NF); subplot(lp); hold on
if vv<0 & angle_max>0 plot3([Xi(2)-0.5*xi Xi(2)+0.5*xi],...
[Xi(3)+0.5*yi Xi(3)-0.5*yi],[max(max(SIG_nd_max))+100 max(max(SIG_nd_max))+100],cv,'LineWidth',1);
end;
if vv<0 & angle_max<0
plot3([Xi(2)+0.5*xi Xi(2)-0.5*xi],...
[Xi(3)+0.5*yi Xi(3)-0.5*yi],[max(max(SIG_nd_max))+100 max(max(SIG_nd_max))+100],cv,'LineWidth',1);
end;
if vv>0 & angle_max>0 plot3([Xi(2)+0.5*xi Xi(2)-0.5*xi],...
[Xi(3)+0.5*yi Xi(3)-0.5*yi],[max(max(SIG_nd_max))+100 max(max(SIG_nd_max))+100],cv,'LineWidth',1);
end;
if vv>0 & angle_max<0 plot3([Xi(2)-0.5*xi Xi(2)+0.5*xi],...
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -