程序无法运行,报错。请大神帮忙调试一下。
syms a b
syms x y z
syms z0 f
z0=405;
f=sqrt(3);
k1=1;k2=1;k3=1;k4=1;
K=[k1 0 0 0;0 k2 0 0;0 0 k3 0;0 0 0 k4];
S=[0,0,405,0]';
e1S=0:0.5:320;
e2S=0:0.5:210;
for i=1:length(e1S)
j=1:length(e2S)
a=e1S(i);
b=e2S(j);
for i=1:4
B1=[0,-f*b/2,0,1]';
B2=[b/2,0,0,1]';
B3=[0,f*b/2,0,1]';
B4=[-b/2,0,0,1]';
A1=[0,-f*a/2,0,1]';
A2=[a/2,0,0,1]';
A3=[0,f*a/2,0,1]';
A4=[-a/2,0,0,1]';
T=[1,0,0,x;0,1,0,y;0,0,1,z;0,0,0,1];
BBi=T*Bi+S;
for i=1:4
li=sqrt((BBi-Ai)'*(BBi-Ai));
D=[l1,0,0,0;0,l2,0,0;0,0,l3,0;0,0,0,l4];
E=[0,y+fa/2-fb/2,z;x-a/2+b/2,0,z;0,y-fa/2+fb/2,z;x+a/2-b/2,0,z];
J=inv(D)*E;
KC=J'*K*J;
K1=eig(KC'*KC);
m=[1;1;1;1];
M=inv(KC)*m;
K11=1/M(1,1);
K22=1/M(2,1);
K33=1/M(3,1);
K44=1/sqrt((M(1,1))^2+(M(2,1))^2+(M(3,1))^2);
% [V,D]=eig(KC);
% V1=V(:,1);
% V2=V(:,2);
% V3=V(:,3);
% P1=V1/sqrt(V1'*V1);
% P2=V2/sqrt(V2'*V2);
% P3=V3/sqrt(V3'*V3);
% P=[P1 P2 P3];
% Y=inv(P)*KC*P
% K11=Y(1,1)
% K22=Y(2,2)
% K33=Y(3,3)
Kmax(i,j)=double(sqrt(max(K1))); %最大刚度评价指标
Kmin(i,j)=double(sqrt(min(K1))); %最小刚度评价指标
k11(i,j)=double(K11);
k22(i,j)=double(K22);
k33(i,j)=double(K33);
k44(i,j)=double(K44);
K2(i,j)=K11+K22+K33+K44;
e1M(i,j)=a;
e2M(i,j)=b;
end
end
end
subplot(2,3,2)
mesh(e1M,e2M,Kmax)
colorbar
subplot(2,3,1)
mesh(e1M,e2M,K2)
colorbar
title('stiffness')
xlabel(' ','rotation',20)
ylabel(' ','rotation',-20)
zlabel(' ','rotation',90)
% Kmax
% Kmin
% V
% Y
```![图片说明](https://img-ask.csdn.net/upload/202006/20/1592656170_839724.png)