2

我对 Matlab 很陌生。我有一个脚本,在 Matlab 上使用ode45并显示摆动弹簧在通过 3-D 空间时的运动。arrow.m该程序几乎正在做我想做的事。现在,钻石的密度正在显示有效的弹簧速度(除了当ode45将其作为个人最喜欢的数字样本),并且速度几乎可以通过函数的步长(至少以我的计算机运行代码的速度)准确计算。我想要做的是让我在代码中注释掉的位置矢量只显示在质量的瞬时位置,即曲线的端点,而不是钻石出现的每个点。我四处寻找帮助,但我尝试的一切似乎都会导致错误。如果有人能指出我正确的方向,那将不胜感激。请尝试运行程序,你会明白我的意思,也可以玩弄函数的参数。

function elasticPendulum(totalTime)
time=1;
totalTime=20
X=1
Vx=0
Y=0
Vy=0
Z=1.1
Vz=0
w=3;
g=9;
l=1;
while (time<totalTime)
    tspan=[0,time];
    x0=[X,Vx,Y,Vy,Z,Vz];
    [t,x]=ode45(@F,tspan,x0);
    plot3(x(:,1),x(:,3),x(:,5),'dk'), axis([-2 2 -2 2 -10 2]);
    grid on;
    o=[0, 0, 0];
    Xeq=[0, 0, -(g/(w^2)+l)];
    arrow(o,Xeq,'Length',5);
    %{ 
    Xt=[x(:,1) x(:,3) x(:,5)]
    arrow(o,Xt,'Length',5);
    %}
    time=time+.1*(((x(2))^2+(x(4))^2+(x(6))^2)^(1/2))
    pause(.1);
end

    function xprime=F(t,x)
    r=sqrt(x(1)^2+x(3)^2+x(5)^2);
    xprime=[x(2);-w*(r-l)/r*x(1);x(4);-w*(r-l)/r*x(3);x(6);-w*(r-l)/r*x(5)-g];
    end

end
4

1 回答 1

1

我认为你可以简单地通过设置Xt来完成你想要的,这样每次迭代只绘制最后一个向量:

Xt = [x(end,1) x(end,3) x(end,5)];

此外,您提到“除了 ode45 采用个人最喜欢的数字样本时”,事情正在发挥作用。您可以ode45通过更改简单地告诉使用固定步长输出tspan

tspan = 0:dt:time;

在哪里dt = 1/time(或者你可以使用linspace)。这与使用固定步长求解器不同——阅读我对这个问题的回答,或许可以对此有所了解。

我认为您没有正确更新动画。您更新时间,但不更新初始条件。因此,您在循环的每次迭代中通过相同的路径进行积分while。您是否有理由在每次迭代中调整结束时间?你是想找到振荡周期还是什么?

此外,您的动画方法相当粗糙/效率低下。您应该了解可以通过Matlab 的 ODE 求解器设置的输出选项。odeset特别是'OutputFcn'。您甚至可以使用和/或查看某些内置输出函数的代码——例如,在命令窗口中键入edit odeplot或查看。edit odephas3这是我为您的案例创建的一个简单的输出函数:

function elasticPendulum
totalTime = 20;
stepsPerSec = 10;
X = 1; Vx = 0;
Y = 0; Vy = 0;
Z = 1.1; Vz = 0;
w = 3; g = 9; l = 1;

opts = odeset('OutputFcn',@(t,x,flag)arrowplot(t,x,flag,w,l,g));
tspan = linspace(0,totalTime,totalTime*stepsPerSec);
x0 = [X,Vx,Y,Vy,Z,Vz];
ode45(@(t,x)F(t,x,w,l,g),tspan,x0,opts);

function xprime=F(t,x,w,l,g)    %#ok<INUSL>
r=sqrt(x(1)^2+x(3)^2+x(5)^2);
xprime=[x(2);-w*(r-l)/r*x(1);x(4);-w*(r-l)/r*x(3);x(6);-w*(r-l)/r*x(5)-g];

function status=arrowplot(t,x,flag,w,l,g)   %#ok<INUSL>
persistent h; % Handle for moving arrow
status = 0;
o = [0, 0, 0];
switch flag
    case 'init'
        axis([-2 2 -2 2 -10 2]); grid on; hold on;
        Xeq = [0, 0, -(g/(w^2)+l)];
        arrow(o,Xeq,'Length',5);
        plot3(x(1,:),x(3,:),x(5,:),'dk');  % Initial position
        Xt = [x(1,end) x(3,end) x(5,end)];
        h = arrow(o,Xt,'Length',5);        % Initial arrow, get handle
    case 'done'
        hold off; status = 1;
    otherwise
        plot3(x(1,:),x(3,:),x(5,:),'dk');        % Plot new positions
        Xt = [x(1,end) x(3,end) x(5,end)];
        arrow(h,'Start',o,'Stop',Xt,'Length',5); % Update arrow
        pause(0.1);
end

plot3每次迭代的调用(中的otherwise情况arrowplot)仍然效率低下,因为它添加了新的高级绘图对象,占用了更多内存。更好/更快的方法是从第一次调用plot3和使用中获取句柄setget添加新数据点。请参阅代码odeplot以及odephas3如何执行此操作(它有点高级)。

请注意我如何通过匿名函数而不是通过创建子函数来传递参数。这有点风格问题。

于 2014-03-02T18:33:22.103 回答