📄 plot_results_tri.m
字号:
function plot_results_tri ( in_data, SIGsys, resp, dof, NF)
% PLOT STATIC DISPLACEMENTS/STRESSES FROM 'CST' & 'CSQ' ELEMENTS (PLAIN FEM)
% =============================================================================
% INPUT:
% in_data - structure with input data: ND, EL LOAD_, CON
% SIGsys - element stresses vector
% resp - structure with displacements vector: resp.static.D
% dof - # of dof's
% NF - # of figure plot. For NF > 100 loads & constrains are plot !!!
% plot main stress =============================================================
% figure(2);
if 1
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 / 6); laby = (maxY / 6);
labx = min([labx laby]); laby = labx;
figure(NF);
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);
axis([(minX-labx) (maxX+labx) (minY-laby) (maxY+laby)]);
% plot elements -------------------------------------
subplot(121); plot(in_data.ND(:,2),in_data.ND(:,3),'b.','MarkerSize',4);
axis equal; axis off; hold on;
for i=1:size(in_data.EL)
yy = (MAX_SIG(i) - minSIGmax)*1/(maxSIGmax - minSIGmax); % white to black
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));
subplot(121);
fill([in_data.ND(node1,2) in_data.ND(node2,2) in_data.ND(node3,2) in_data.ND(node1,2)], ...
[in_data.ND(node1,3) in_data.ND(node2,3) in_data.ND(node3,3) in_data.ND(node1,3)],[yy yy yy]);
end;
if in_data.EL(i,2)==5 % "5" - CSQ element (quadrilateral)
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));
subplot(121);
fill([in_data.ND(node1,2) in_data.ND(node2,2) in_data.ND(node3,2) ...
in_data.ND(node4,2) in_data.ND(node1,2)],...
[in_data.ND(node1,3) in_data.ND(node2,3) in_data.ND(node3,3) in_data.ND(node4,3)...
in_data.ND(node1,3)],[yy yy yy]);
end;
end;
% plot lines/directions of MAX stresses -------------
for i=1:size(in_data.EL)
yy = (MAX_SIG(i) - minSIGmax)*1/(maxSIGmax - minSIGmax); % white to black
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));
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 % "5" - CSQ element (quadrilateral)
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)); % check if XX > YY
yi = abs(cos(angle_max)*yy*laby);
xi = abs(sin(angle_max)*yy*labx);
angle_max = angle_max*360/(2*pi); % radian to degrees
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-'; else cv='r-'; end; % blue -tension
% plot directions of main (MAX) stresses
subplot(121);
if vv<0 & angle_max>0 plot([Xi(2)-0.5*xi Xi(2)+0.5*xi],...
[Xi(3)+0.5*yi Xi(3)-0.5*yi],cv,'LineWidth',2);
end;
if vv<0 & angle_max<0 plot([Xi(2)+0.5*xi Xi(2)-0.5*xi],...
[Xi(3)+0.5*yi Xi(3)-0.5*yi],cv,'LineWidth',2);
end;
if vv>0 & angle_max>0 plot([Xi(2)+0.5*xi Xi(2)-0.5*xi],...
[Xi(3)+0.5*yi Xi(3)-0.5*yi],cv,'LineWidth',2);
end;
if vv>0 & angle_max<0 plot([Xi(2)-0.5*xi Xi(2)+0.5*xi],...
[Xi(3)+0.5*yi Xi(3)-0.5*yi],cv,'LineWidth',2);
end;
end;
title('Stress field - MAX (blue - tension)');
fg = ['max: ' num2str(maxSIGmax) ' ' 'min: ' num2str(minSIGmax)];
xlabel(fg);
% axes('Position',[.1 .92 .06 .07]); fill([0 1 1 0],[0 0 1 1],'w');
% set(gca,'XTick',[]); set(gca,'YTick',[]); fg=num2str(maxSIGmax);text(.1,.5,fg);
% plot lines/directions of MIN stresses -------------
subplot(122); plot(in_data.ND(:,2),in_data.ND(:,3),'b.','MarkerSize',4);
axis equal; axis off; hold on;
for i=1:size(in_data.EL)
yy = 1-(MIN_SIG(i) - minSIGmin)*1/(maxSIGmin - minSIGmin); % white to black
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));
subplot(122);
fill([in_data.ND(node1,2) in_data.ND(node2,2) in_data.ND(node3,2) in_data.ND(node1,2)], ...
[in_data.ND(node1,3) in_data.ND(node2,3) in_data.ND(node3,3) in_data.ND(node1,3)],[yy yy yy]);
end;
if in_data.EL(i,2)==5 % "5" - CST element (quadrilateral)
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));
subplot(122);
fill([in_data.ND(node1,2) in_data.ND(node2,2) in_data.ND(node3,2) in_data.ND(node4,2) in_data.ND(node1,2)], ...
[in_data.ND(node1,3) in_data.ND(node2,3) in_data.ND(node3,3) in_data.ND(node4,3) in_data.ND(node1,3)],[yy yy yy]);
end;
end;
for i=1:size(in_data.EL)
yy = 1-(MIN_SIG(i) - minSIGmin)*1/(maxSIGmin - minSIGmin); % white to black
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));
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 % "5" - CSQ element (quadrilateral)
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)); % compare XX and YY
yi = abs(cos(angle_max)*yy*laby);
xi = abs(sin(angle_max)*yy*labx);
angle_max = angle_max*360/(2*pi); % radian to degrees
if vv>0 txy=xi*laby/labx; xi=yi*labx/laby; yi=txy; end;
if MIN_SIG(i)>0 cv='b-'; else cv='r-'; end; % blue -tension
txy=xi*laby/labx; xi=yi*labx/laby; yi=-txy;
% plot directions of main (MIN) stresses
if vv<0 & angle_max>0 plot([Xi(2)-0.5*xi Xi(2)+0.5*xi],...
[Xi(3)+0.5*yi Xi(3)-0.5*yi],cv,'LineWidth',2);
end;
if vv<0 & angle_max<0 plot([Xi(2)+0.5*xi Xi(2)-0.5*xi],...
[Xi(3)+0.5*yi Xi(3)-0.5*yi],cv,'LineWidth',2);
end;
if vv>0 & angle_max>0 plot([Xi(2)+0.5*xi Xi(2)-0.5*xi],...
[Xi(3)+0.5*yi Xi(3)-0.5*yi],cv,'LineWidth',2);
end;
if vv>0 & angle_max<0 plot([Xi(2)-0.5*xi Xi(2)+0.5*xi],...
[Xi(3)+0.5*yi Xi(3)-0.5*yi],cv,'LineWidth',2);
end;
end;
title('Stress field - MIN (red - compression)');
fg = ['max: ' num2str(maxSIGmin) ' ' 'min: ' num2str(minSIGmin)];
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -