运行matlab代码表明没有给出足够的输入,并且我们的RHS函数有问题



此代码的目标是用数字模拟三个三维物体。所有的物体都有恒定的质量。我们将确定物体的位置随时间的变化。然而,当我们运行代码时,RHS函数出现故障;没有给出足够的输入";我们不知道为什么。

    % initial conditions and parameters
    clc; 
    clear all; 
    M1=1;   %mass of body 1
    M2=2;   %mass of body 2 
    M3=0.5; %mass of body 3 

    % the position of the three bodies in the initial conditions 
    P1_0=[0,0,0];
    P2_0=[1,0,0];
    P3_0=[0,0,1];
    P1=P1_0;
    P2=P2_0;
    P3=P3_0;
    % the initial speed of the three bodies
    a=0.05;
    V1_0=[a,a,a];
    V2_0=[-a,a,a];
    V3_0=[-a,-a,-a];
    V1=V1_0;
    V2=V2_0;
    V3=V3_0;
    X0=[P1_0,V1_0,P2_0,V2_0,P3_0,V3_0]'; % the initial positions and speeds of the three bodies in 
    one column vector.
    [t, X]= ode45(@(t,X) RHS(t,X), [0 3], X0);

    h=0.1;
    K1=(h*RHS(X0));
    K2=h*(RHS(X0+0.5*K1));
    K3=h*(RHS(X0-K1+(2*K2)));
    X1=X0+(1/6)*(K1+4*K2+K3);
    function Y=RHS(t,X)
    V1=X(4:6);
    V2=X(10:12);
    V3=X(16:18);
    P1=X(1:3);
    P2=X(7:9);
    P3=X(13:15);
    M1=1;
    M2=2;
    M3=0.5;
    Y=[V1;Force(P1,M1,P2,M2,P3,M3)/M1;V2;Force(P2,M2,P1,M1,P3,M3)/M2;V3;Force(P3,M3,P1,M1,P2,M2)/M3];
    end 

    function F=Force(P,M,Pa,Ma,Pb,Mb)
    F=(-M*Ma*(P-Pa))/(norm(P-Pa)^3)-(M*Mb*(P-Pb))/(norm(P-Pb)^3);
    end 

Thomas:

我刚刚在我的端上运行了您的代码,它完美地工作,直到并包括对ode45的调用。然而,我确实得到了与你在问题中声称的相同的错误。该错误是由K1=(h*RHS(X0))行触发的,我认为这是一个块的开始,您正试图在该块中设置Runge-Kutta方法。出现错误的原因是,在这一行中,您只给RHS函数一个输入,然而,您已经将RHS定义为接受两个输入(tX(:

    function Y=RHS(t,X)
    V1=X(4:6);
    V2=X(10:12);
    V3=X(16:18);
    P1=X(1:3);
    P2=X(7:9);
    P3=X(13:15);
    M1=1;
    M2=2;
    M3=0.5;
    Y=[V1;Force(P1,M1,P2,M2,P3,M3)/M1;V2;Force(P2,M2,P1,M1,P3,M3)/M2;V3;Force(P3,M3,P1,M1,P2,M2)/M3];
    end 

由于RHS函数根本不使用t(ODE的自治系统(,以下将避免错误:

    h=0.1;
    K1=(h*RHS([],X0));
    K2=h*(RHS([],X0+0.5*K1));
    K3=h*(RHS([],X0-K1+(2*K2)));
    X1=X0+(1/6)*(K1+4*K2+K3);

MATLAB只是在寻找额外的输入,仅此而已。注意:输入一个空向量[]没有什么特别的,实际上这个";"修复";可以使用任何输入(包括string或复数(。诀窍是占用占位符输入。

另一个变体不是用伪参数进行扩展,而是减少到

function Y=RHS(X)

然后RK3步骤可以保持原来的状态,而ode45呼叫需要适应

[t, X]= ode45(@(t,X) RHS(X), [0 3], X0);

相关内容

  • 没有找到相关文章

最新更新