急急急大伙帮我看看准双曲面小轮的程序 老出问题

飞龙在天88866

syms a t pi  rt f g j q e01 xb1 m1 x1  p  io1   s1  n    r1  s  b c d
pi=3.1415926
a=10*pi/180
rt= 96.52 ;       %sh(407)
f=(15+21.4723)*pi/180   ;     %sh(396)刀转角
g=350.111*pi/180  ;      %sh(417)刀倾角
q=96.444*pi/180    ;       %sh(431)
e01=40.2217   ;    %sh(430)
xb1=9.0661   ;       %sh(429)
m1= 357*pi/180 ;   %sh(402)齿陪安装角
x1=-50.8138    ;     %sh(425)%ffffffff安装距
    io1=  3.7051;        %sh(436)
    s1=  102.902 ;          %sh(4;16)
r1=[];
for p=0:0.1:0.8;
  for t= 0:pi/30:2*pi;
syms s real
rr1=[(rt+s*sin(a))*cos(t);((rt+s*sin(a)))*sin(t);-s*cos(a);1];
nt=[-cos(a)*cos(t);-cos(a)*sin(t);-sin(a);0];
mbt=[cos(f) 0 sin(f) 0;0 1 0 0;-sin(f) 0 cos(f) 0;0 0 0 1];
mcb=[-sin(g) -cos(g) 0  s1;cos(g) -sin(g) 0 0;0 0 1 0;0 0 0 1];
moc=[cos(q) sin(q) 0 0;-sin(q) cos(q) 0 0;0 0 1 0;0 0 0 1];
ro=moc*mcb*mbt*rr1;
no=moc*mcb*mbt*nt;
wc=-io1.*[0;0;1;1];%%%%%%%%%%%%%%%%%
wp=-[cos(m1);0;sin(m1);1];
ooon=[0;-e01;xb1;1];%%%%%%%%%%%%%%%%%%
wc1=wc(1:3,:);
wp1=wp(1:3,:);
no1=no(1:3,:);
ooon1=ooon(1:3,:);
ro1=ro(1:3,:);
no11=no1.';
s2=no11*(cross((wc1-wp1),ro1)+cross(ooon1,wp1));
s=solve(s2,s);   
s=vpa(s,5);
mnm=[1 0 0 0;0 1 0 e01;0 0 1 -xb1;0 0 0 1];
mqn=[cos(m1) 0 sin(m1) -x1;0 1 0 0;-sin(m1) 0 cos(m1) 0;0 0 0 1];
mpq=[1 0 0 0;0 cos(p) -sin(p) 0;0 sin(p) cos(p) 0;0 0 0 1];
rrr1=[((rt+s*sin(a)))*cos(t) ;((rt+s*sin(a)))*sin(t); -s*cos(a); 1];
r_1=mpq*mqn*mnm*moc*mcb*mbt*rrr1;
r_1=vpa(r_1,5);
r_1=double(r_1);
r1=[r1 r_1]
   end
end
  r1=r1';
  r1(:,4)=[];
xlswrite('1234.xls',r1,'Sheet1','A1');

q364228550

系统还原一下看看。。。