[已答复] 错误使用 + 矩阵维度必须一致

[复制链接]
lzs_j4nVy 发表于 2021-1-10 10:36:47
求助,编写的矩阵微分方程组,然后在运行的时候报错:错误使用  +  矩阵维度必须一致。之后我将式子中求导的项去掉了,然后仍然要求维度必须一致,请问这该如何解决。程序如下:
opt=odeset;opt.RelTol=1e-8;
Tmh=0;Tmp=0;Mh=4100;Ihx=4032;Ihy=3461;Ihz=4980;Ihxy=31;Ihxz=637;Ihyz=26;Mp=1500;Ipx=117;Ipy=861;
Ipz=754;Ipxy=0.615;Ipxz=152;Ipyz=14.5;
g=9.81;xh=0.237;yh=0.013;zh=0.404;xp=0.5;yp=0;zp=0;xhp=0;yhp=0;zhp=0.96;xch=0;ych=0;zch=0.9;u=12*sin((pi/3)*t);
v=0;w=12*sin((pi/3)*t);
p=12*sin((pi/3)*t);q=12*sin((pi/3)*t);r=12*sin((pi/3)*t);theta=12*sin((pi/3)*t);varphi=12*sin((pi/3)*t);Ch=100;Cp=90;Kp=100;
A=p*cos(x(1))*cos(x(2))+q*sin(x(1))*cos(x(2))-r*sin(x(2));
B=p*cos(x(1))*sin(x(2))+q*sin(x(1))*sin(x(2))+r*cos(x(2));
E=-p*sin(x(1))+q*cos(x(1));
F=p*cos(x(1))+q*sin(x(1));
e=-diff(p)*sin(x(1))+diff(q)*cos(x(1));
f=diff(p)*cos(x(1))+diff(q)*sin(x(1));
Tcwh=Mh.*xh.*(v+xch.*r-zch.*p+xh.*r-zh.*f)+Mh.*yh.*(u-r.*ych+q.*zch-r.*yh+zh.*e)+....
    Mh.*zh.*F.*(u-r.*ych+q.*zch-r.*yh+zh.*E)+Mh.*zh.*E.*(v+r.*xch-p.*zch+r.*xh-zh.*F)-Mh.*xh.*F.*(w-q.*xch+p.*ych-E.*xh+yh.*F)...
    -Mh.*yh.*E.*(w-q.*xch+p.*ych-E.*xh+yh.*F)...
    +Mp.*((diff(u)-diff(r).*ych+diff(q).*zch-diff(r).*yhp+zhp.*e).*(-yhp-cos(x(2)).*yp)+(diff(v)+xch.*diff(r)-zch.*diff(p)+...
    xhp.*diff(r)-zhp.*f+b.*xp-a.*zp).*(xhp+xp.*cos(x(2))+sin(x(2)).*zp)+(diff(w)-diff(q).*xch+diff(p).*ych-e.*xhp+yhp.*f-xp.*e+...
    yp^2.*sin(x(2)).*a)-((u-r.*ych+q.*zch-r.*yhp+zhp.*E-yp.*B+zp.*E).*(-(zhp+zp).*F-yp.*sin(x(2)).*E)+(v+r.*xch-p.*zch+r.*xhp-zhp.*F+...  
    xp.*B-zp.*A).*(-zhp.*E+xp.*sin(x(2)).*E-zp.*cos(x(2)).*E)+(w-q.*xch+p.*ych-E.*xhp+yhp.*F-xp.*E+yp.*A)*(xhp*F+yhp*E+xp*F+yp*cos(x(2))*E)))...
    +Ihxz*f+Ihyz*e+Ihz*diff(r)-Ihx*E*F+Ihy*E*F+Ihxy*F^2+Ihxy*E^2-Ihxz*r*E+Ihyz*r*F-Ipx*sin(x(2))*a+Ipz*cos(x(2))*b-...
    Ipxy*sin(x(2))*e-Ipxz*sin(x(2))^2*f+Ipxz*cos(x(2))*a+Ipyz*cos(x(2))*e-Ipx*cos(x(2))*E*A+Ipy*E*F-Ipz*A*E*sin(x(2))-...
    Ipxy*cos(x(2))*E^2+Ipxy*A*F-Ipxz*cos(x(2))*B*E-Ipxz*sin(x(2))*A*E+Ipyz*B*F-Ipyz*sin(x(2))*E^2;

289014845476818569.jpg
503381164817086127.jpg

4 条回复


20141303 发表于 2021-1-10 10:51:23
代码不全,t未定义

lzs_j4nVy 发表于 2021-1-10 10:54:34
20141303 发表于 2021-1-10 10:51
代码不全,t未定义
  1. function dx=shuangzhou(t,x,Tmh,Tmp,Mh,Ihx,Ihy,Ihz,Ihxy,Ihxz,Ihyz,Mp,Ipx,Ipy,...
  2.     Ipz,Ipxy,Ipxz,Ipyz,g,xh,yh,zh,xp,yp,zp,xhp,yhp,zhp,xch,ych,zch,u,v,w,...
  3.     p,q,r,theta,varphi,Ch,Cp,Kp)
  4. A=p*cos(x(1))*cos(x(2))+q*sin(x(1))*cos(x(2))-r*sin(x(2));
  5. B=p*cos(x(1))*sin(x(2))+q*sin(x(1))*sin(x(2))+r*cos(x(2));
  6. E=-p*sin(x(1))+q*cos(x(1));
  7. F=p*cos(x(1))+q*sin(x(1));
  8. e=-diff(p)*sin(x(1))+diff(q)*cos(x(1));
  9. f=diff(p)*cos(x(1))+diff(q)*sin(x(1));
  10. Ib=Mh*(xh^2-yh^2)+Ihz+Mp*((yhp+cos(x(2))*yp)^2+(xhp+cos(x(2))*xp+sin(x(2))*zp)^2+(sin(x(2))*yp)^2)...
  11.     +Ipx*sin(x(2))^2+Ipz*cos(x(2))^2-2*Ipxz*sin(x(2))*cos(x(2));
  12. Ibc=Mp*(-zp*(yhp+yp*cos(x(2)))+xp*yp*sin(x(2))-Ipxy*sin(x(2))+Ipyz*cos(x(2)));
  13. Ic=Mp*(-zp*yhp+yp*(xp*sin(x(2))-zp*cos(x(2))))+Ipyz*cos(x(2))-Ipxy*sin(x(2));
  14. Chb=0;
  15. Chbc=Mp*(-2*yp*yhp*sin(x(2))+2*(zp*cos(x(2))-xp*sin(x(2)))*(xhp+xp*cos(x(2))+zp*sin(x(2))))...
  16.     +2*(Ipx-Ipz)*sin(x(2))*cos(x(2))+2*Ipxz*(sin(x(2))^2-cos(x(2))^2);
  17. Ch1=Ch+Chb+Chbc*x(2);
  18. Chc=Mp*(yp*(yhp+yp*cos(x(2)))*A+yp*sin(x(2))*(u-r*ych+q*zch-r*yhp+zhp*E-yp*B+zp*E)+(xhp+xp*cos(x(2))...
  19.     +sin(x(2))*zp)*(xp*A+zp*B)-(xp*sin(x(2))-zp*cos(x(2)))*(v-r*xch-p*zch+r*xhp-zhp*F+xp*B-zp*A)...
  20.     +yp^2*sin(x(2))*B-yp*cos(x(2))*(w-q*xch+p*ych-E*xhp+yhp*F-xp*E+yp*A)+zp*(zhp*F+yp*sin(x(2))*E+zp*F)...
  21.     +xp*(xhp*F+yhp*E+xp*F+yp*cos(x(2))*E))+((Ipx-Ipz)*(sin(x(2))^2-cos(x(2))^2)-4*Ipxz*sin(x(2))*cos(x(2)))*F...
  22.     +2*(Ipx-Ipz)*r*sin(x(2))*cos(x(2))-2*E*(Ipxy*cos(x(2))+Ipyz*sin(x(2)))+Ipy*F+2*Ipxz*r*(sin(x(2))^2-cos(x(2))^2);
  23. Chcc=Mp*yp*(zp*sin(x(2))+xp*cos(x(2)))-Ipxy*cos(x(2))-Ipyz*sin(x(2));
  24. Ch2=Chc+Chcc*x(2);
  25. Cpb=Mp*((-(zhp+zp)*F-yp*sin(x(2))*E)*zp-xp*(F*(xhp+xp)+F*(yhp+yp*cos(x(2))))-(A*yp*(yhp+yp*cos(x(2)))+(xhp+xp*cos(x(2))...
  26.     +sin(x(2))*zp)*(xp*A+zp*B)+yp^2*sin(x(2))*B+yp*sin(x(2))*(u-r*ych+q*zch-r*yhp+zhp*E-yp*B+zp*E)-(xp*sin(x(2))-zp*cos(x(2)))...
  27.     *(v+r*xch-p*zch+r*xhp-zhp*F+xp*B-zp*A)-yp*cos(x(2))*(w-q*xch+p*ych-E*xhp+yhp*F-xp*E+yp*A)))+(-Ipy*F+2*E*(Ipxy*cos(x(2))...
  28.     +Ipyz*sin(x(2)))-(Ipx-Ipz)*sin(x(2))*B+(Ipz-Ipx)*cos(x(2))*A-4*Ipxz*sin(x(2))*cos(x(2))*F+2*Ipxz*r*(sin(x(2))^2-cos(x(2))^2));
  29. Cpcc=-(Mp*(-yp*yhp*sin(x(2))+(zp*sin(x(2))+xp*cos(x(2)))*(zp*cos(x(2))-xp*sin(x(2))))+(Ipx-Ipz)*sin(x(2))*cos(x(2))+...
  30.     Ipxz*(sin(x(2))^2-cos(x(2))^2));
  31. Cp1=Cpb+Cpcc*x(1);
  32. Cpc=2*zp*yp*A;
  33. Cpbc=Mp*(zp*sin(x(2))+xp*cos(x(2)));
  34. Cp2=Cp+Cpc+Cpbc*x(1);
  35. Tgh=Mh.*g.*xh.*...
  36.     (sin(x(1)).*sin(theta)+cos(x(1)).*sin(varphi).*cos(theta))...
  37.     +Mh.*g.*yh.*(cos(x(1)).*sin(theta)-sin(x(1)).*sin(varphi).*cos(theta))...
  38.     +Mp.*g.*(0.5.*sin((varphi)+(theta)).*(xhp.*cos(x(1))+xp.*cos(x(2)).*cos(x(1))-yp.*sin(x(1))+zp.*sin(x(2)).*cos(x(1))-yhp.*sin(x(1)))+...
  39.     sin(theta).*sin(x(1)).*(xhp+xp.*cos(x(2))+zp.*sin(x(2)))+cos(x(1)).*(yh+yp).*sin(theta).*cos(x(1))+0.5.*(xhp+xp.*cos(x(2))+...
  40.     zp.*sin(x(2))).*cos(x(1)).* sin((varphi)-(theta))-0.5.*(yhp+yp).*sin((varphi)-(theta)).*sin(x(1)));
  41. Tcwh=Mh.*xh.*(v+xch.*r-zch.*p+xh.*r-zh.*f)+Mh.*yh.*(u-r.*ych+q.*zch-r.*yh+zh.*e)+....
  42.     Mh.*zh.*F.*(u-r.*ych+q.*zch-r.*yh+zh.*E)+Mh.*zh.*E.*(v+r.*xch-p.*zch+r.*xh-zh.*F)-Mh.*xh.*F.*(w-q.*xch+p.*ych-E.*xh+yh.*F)...
  43.     -Mh.*yh.*E.*(w-q.*xch+p.*ych-E.*xh+yh.*F)...
  44.     +Mp.*((diff(u)-diff(r).*ych+diff(q).*zch-diff(r).*yhp+zhp.*e).*(-yhp-cos(x(2)).*yp)+(diff(v)+xch.*diff(r)-zch.*diff(p)+...
  45.     xhp.*diff(r)-zhp.*f+b.*xp-a.*zp).*(xhp+xp.*cos(x(2))+sin(x(2)).*zp)+(diff(w)-diff(q).*xch+diff(p).*ych-e.*xhp+yhp.*f-xp.*e+...
  46.     yp^2.*sin(x(2)).*a)-((u-r.*ych+q.*zch-r.*yhp+zhp.*E-yp.*B+zp.*E).*(-(zhp+zp).*F-yp.*sin(x(2)).*E)+(v+r.*xch-p.*zch+r.*xhp-zhp.*F+...  
  47.     xp.*B-zp.*A).*(-zhp.*E+xp.*sin(x(2)).*E-zp.*cos(x(2)).*E)+(w-q.*xch+p.*ych-E.*xhp+yhp.*F-xp.*E+yp.*A)*(xhp*F+yhp*E+xp*F+yp*cos(x(2))*E)))...
  48.     +Ihxz*f+Ihyz*e+Ihz*diff(r)-Ihx*E*F+Ihy*E*F+Ihxy*F^2+Ihxy*E^2-Ihxz*r*E+Ihyz*r*F-Ipx*sin(x(2))*a+Ipz*cos(x(2))*b-...
  49.     Ipxy*sin(x(2))*e-Ipxz*sin(x(2))^2*f+Ipxz*cos(x(2))*a+Ipyz*cos(x(2))*e-Ipx*cos(x(2))*E*A+Ipy*E*F-Ipz*A*E*sin(x(2))-...
  50.     Ipxy*cos(x(2))*E^2+Ipxy*A*F-Ipxz*cos(x(2))*B*E-Ipxz*sin(x(2))*A*E+Ipyz*B*F-Ipyz*sin(x(2))*E^2;
  51. Th=Tmh-Tgh-Tcwh;
  52. Tgp=Mp*g*(xp*cos(x(1))*sin(x(2))*sin(theta)-zp*cos(x(1))*sin(x(2))*sin(theta)-xp*cos(x(2))*cos(varphi)*sin(theta)-...
  53.     zp*sin(x(2))*cos(varphi)*cos(theta)-xp*sin(x(1))*cos(varphi)*sin(theta)*sin(x(2))+zp*sin(x(1))*sin(varphi)*...
  54.     cos(theta)*cos(x(2)));
  55. Tcwp=Mp*(zp*(diff(u)-diff(r)*ych+diff(q)*zch+zh*e-yp*b-zp*sin(x(1))*diff(p))-xp*(diff(w)-diff(q)*xch+diff(p)*ych-...
  56.     e*xhp+yhp*f-xp*e+yp*a)-(-(u-r*ych+q*zch-r*yhp+zhp*E-yp*B+zp*E)*yp*A+(v+r*xch-p*zch+r*xhp-zhp*F+...  
  57.     xp*B-zp*A)*(xp*A+zp*B)-(w-q*xch+p*ych-E*xhp+yhp*F-xp*E+yp*A)*yp*B))+Ipy*e+Ipxy*a+Ipyz*b+Ipx*A*B-Ipz*A*B+Ipxy*E*B...
  58.     +Ipxz*(B^2-A^2)-Ipyz*E*B;
  59. Tp=Tmp-Tgp-Tcwp;
  60. K=[0, 0; 0, Kp];
  61. M=[Ib, Ibc; Ibc, Ic];
  62. C=[Ch1, Ch2; Cp1, Cp2];
  63. F=[Th; Tp];
  64. dx=[x(3:4);inv(M)/(1/(F-C*x(3:4)-K*x(1:2)))];
  65.    


复制代码

这是函数的代码。
下面是用ode45求解的代码
  1. opt=odeset;opt.RelTol=1e-8;
  2. Tmh=0;Tmp=0;Mh=4100;Ihx=4032;Ihy=3461;Ihz=4980;Ihxy=31;Ihxz=637;Ihyz=26;Mp=1500;Ipx=117;Ipy=861;
  3. Ipz=754;Ipxy=0.615;Ipxz=152;Ipyz=14.5;
  4. g=9.81;xh=0.237;yh=0.013;zh=0.404;xp=0.5;yp=0;zp=0;xhp=0;yhp=0;zhp=0.96;xch=0;ych=0;zch=0.9;u=12*sin((pi/3)*t);
  5. v=0;w=12*sin((pi/3)*t);
  6. p=12*sin((pi/3)*t);q=12*sin((pi/3)*t);r=12*sin((pi/3)*t);theta=12*sin((pi/3)*t);varphi=12*sin((pi/3)*t);Ch=100;Cp=90;Kp=100;
  7. [t,x]=ode45(@shuangzhou,[0,0.5],zeros(4,1),opt,Tmh,Tmp,Mh,Ihx,Ihy,Ihz,Ihxy,Ihxz,Ihyz,Mp,Ipx,Ipy,...
  8.     Ipz,Ipxy,Ipxz,Ipyz,g,xh,yh,zh,xp,yp,zp,xhp,yhp,zhp,xch,ych,zch,u,v,w,...
  9.     p,q,r,theta,varphi,Ch,Cp,Kp);
  10. plot(t,x(:,1:2)),figure;
  11. Mh.*xh.*(diff(v)+xch.*diff(r)-zch.*diff(p)+xh.*diff(r)-zh.*f)+Mh.*yh.*(diff(u)-diff(r).*ych+diff(q).*zch-diff(r).*yh+zh.*e)
复制代码

lzs_j4nVy 发表于 2021-1-10 10:56:40
lzs_j4nVy 发表于 2021-1-10 10:54
这是函数的代码。
下面是用ode45求解的代码

希望您能帮忙看一下哪里错了,昨天和其他人讨论,他们说,我的那个里面虽然有的是常数,但是有的是向量。但是即便是向量的话也只是一维的向量,这怎么会出现维度问题呢?

20141303 发表于 2021-1-10 15:33:57
用ode45求解的代码中第三行的t未定义
您需要登录后才可以回帖 登录 | 注册

本版积分规则

相关帖子
热门教程
站长推荐
快速回复 返回顶部 返回列表