声振论坛

 找回密码
 我要加入

QQ登录

只需一步,快速开始

查看: 1405|回复: 6

[编程技巧] subs替换的奇怪问题

[复制链接]
发表于 2009-5-20 22:27 | 显示全部楼层 |阅读模式

马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。

您需要 登录 才可以下载或查看,没有账号?我要加入

x
各位大侠帮忙看看这个奇怪的问题:
function yansuan
% L=[1 2 3]
% b=[1/4*pi 1/6*pi -1/4*pi]
% F= [1 2 1]
% m=[1 2 1]
% w=[1/2*pi -1/2*pi -1/2*pi]
% Mc=[0 0 0]
% n=length(L);
L=[0.4 0.2 0.3 0.4];
b=[1/2*pi 1/4*pi 1/6*pi -1/4*pi];
F=[2 1 1 2];
m=[2 1 1 1];
w=[-1/2*pi -1/2*pi  -1/2*pi  -1/2*pi]
Mc=[0 0 0 0]
n=length(L)
for i=1:n
              X(i)=L(i)*cos(b(i));   
              Y(i)=L(i)*sin(b(i));
         end
         p1=cumsum(X);
         q1=cumsum(Y);
for i=1:n
              varargout{i}=[p1(1,i),q1(1,i)];
         end
for i=1:n
             x(i+1)=varargout{i}(1);
             y(i+1)=varargout{i}(2);
         end
              x(1)=0;
              y(1)=0;

for i=1:n
    Fx(i)=F(i)*cos(w(i));
    Fy(i)=F(i)*sin(w(i));
    Fa(3*i-2)=Fx(i);
    Fa(3*i-1)=Fy(i);
    Fa(3*i)=Mc(i);
end
syms x1 y1 m1 x2 y2 m2 x3 y3 m3 x4 y4 m4 x5 y5 m5 x6 y6 m6
qq=[x1 y1 m1 x2 y2 m2 x3 y3 m3 x4 y4 m4 x5 y5 m5 x6 y6 m6];
q=qq(1:3*n);
MK(1,1)=q(1)-1/2*L(1)*cos(q(3));
MK(2,1)=q(2)-1/2*L(1)*sin(q(3));
for i=2:n
    MK(2*i-1,1)=q(3*i-2)-1/2*L(i)*cos(q(3*i))-q(3*i-5)-1/2*L(i-1)*cos(q(3*i-3));
    MK(2*i,1)=q(3*i-1)-1/2*L(i)*cos(q(3*i))-q(3*i-4)-1/2*L(i-1)*sin(q(3*i-3));
end
MK(2*n+1,1)=x(n+1)-q(3*n-2)-1/2*L(n)*cos(q(3*n));
MK(2*n+2,1)=y(n+1)-q(3*n-1)-1/2*L(n)*sin(q(3*n));
M=[MK]
Mq=jacobian(M,q);
syms t
Mt=diff(M,t);
for i=1:3*n
Kq(:,i)=diff(Mq(:,i),q(i));
end
Mtt=diff(M,t,2)
Kt=diff(Mq,t)
% vq=[0.2;0.3;0.5;0.4;1.2;1.5;1.3;2.3;0.4]
vq=[0.2;0.3;0.12;0.32;0.45;-0.5;0.312;-0.123;0.002;0.21;0215;0.225];
r=-Kq*(vq.*vq)-2*Kt*vq-Mtt;
F=[Fa.';r];
for i=1:n
J(i)=1/12*m(i)*L(i)^2;
z(3*i-2)=m(i);
z(3*i-1)=m(i);
z(3*i)=J(i);
end
Z=diag(z);
N=length(F)-length(z);
MZ=[Z Mq.';Mq zeros(N)];
AK=MZ\F;
aq=AK(1:3*n);
x2=AK(3*n+1:end);
% q2=[0.2;0.3;0.1;0.4;0.25;0.36;0.47;0.25;0.21]
q2=[0.21;0.2;0.36;0.32;0.51;0.02;0.321;0.65;0.25;0.21;0.32;0.5]
q=q2'
q1=mat2cell(q,1,ones(1,3*n))
% aq=double(subs(aq,{x1 y1 m1 x2 y2 m2 x3 y3 m3},q1))
% Q=double(subs(x2,{x1 y1 m1 x2 y2 m2 x3 y3 m3},q1))
aq=double(subs(aq,{x1 y1 m1 x2 y2 m2 x3 y3 m3 x4 y4 m4},q1))
Q=double(subs(x2,{x1 y1 m1 x2 y2 m2 x3 y3 m3 x4 y4 m4},q1))
clear all
那位大侠帮忙看看。主要被注视掉的部分以及与之相同的没有被注释掉的部分。两者其实是一样的,只是数组的长度不一样而已,一个是四维数组,一个是三维数组。可是为什么数组为三维的时候,可以算出结果,但数组为四维的时候,却是算不出结果来呢?好奇怪哦!恳求那位大侠帮帮忙,小弟感激不尽了。很着急着用,求求各位大侠啦!!
回复
分享到:

使用道具 举报

发表于 2009-5-20 22:31 | 显示全部楼层
有出错提示吗?有的话请给出
 楼主| 发表于 2009-5-20 23:00 | 显示全部楼层
??? Error using ==> maple
Error, (in evalm) evalm accepts only one argument, an expression

Error in ==> sym.subs at 230
         NEWf = sym(maple('evalm',strrep(NEWf,'MATRIX','array')));

Error in ==> yansuan at 86
aq=double(subs(aq,{x1 y1 m1 x2 y2 m2 x3 y3 m3 x4 y4 m4},q1))
这是出错信息!谢谢关注
发表于 2009-5-21 11:45 | 显示全部楼层
应该错在 x2=AK(3*n+1:end);这一行
此时改变了x2,已不是单纯的符号变量
 楼主| 发表于 2009-5-21 13:23 | 显示全部楼层
aq=AK(1:3*n);
x2=AK(3*n+1:end);
这两行为什么会出错呢?我只是取符号函数AK的前3*n行和后3*n行而已,为什么就会变成不再是单纯的符号变量呢?如不是单纯的符号变量那是什么呢?还有,当数组长度是3的时候是可以运行的。若只是单纯的x2=AK(3*n+1:end);出错,那aq也应该会有值呀。可是为什么每次都是运行到aq=double(subs(aq,{x1 y1 m1 x2 y2 m2 x3 y3 m3 x4 y4 m4},q1))这一行就会出现错误呢?这个问题我很困惑,还敬请能够详细解答,小弟在此在此拜谢!!
发表于 2009-5-21 13:31 | 显示全部楼层
看看
K>> whos x1 x2 m1 m2
  Name       Size              Bytes  Class    Attributes
  m1         1x1                 128  sym               
  m2         1x1                 128  sym               
  x1         1x1                 128  sym               
  x2       10x1             1957872  sym
 楼主| 发表于 2009-5-21 14:55 | 显示全部楼层
那可以怎么修改呀?
还请赐教!!
您需要登录后才可以回帖 登录 | 我要加入

本版积分规则

QQ|小黑屋|Archiver|手机版|联系我们|声振论坛

GMT+8, 2024-9-22 12:31 , Processed in 0.060325 second(s), 17 queries , Gzip On.

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

快速回复 返回顶部 返回列表