MATLAB Answers

RK4 Errors in MATLAB

2 ビュー (過去 30 日間)
Kyle Broder
Kyle Broder 2016 年 8 月 4 日
コメント済み: Torsten 2016 年 8 月 4 日

I'm trying to solve the following ODE using RK4 $$F = \frac{1}{50t^2}$$.

My code is the following

    function [t,y]=euler2(y0,dt,t0,tf)
    t=t0:dt:tf;
    y=zeros(1,length(t));
    y(:,1)=y0;
    F = @(t,r) 1/(50*Power(t,2));
    for i = 1:length(t)-1
    k1 = F(y(:,i),t(i));
    k2 = F(y(:,i)+0.5*dt*k1,t(i)+0.5*dt);
    k3 = F(y(:,i)+0.5*dt*k2,t(i)+0.5*dt);
    k4 = F(y(:,i)+dt*k3,t(i)+dt);
    y(:,i+1) = y(:,i) + dt*(k1+2.0*(k2+k3)+k4)/6.0;
    end
    end

When I run the code however, with y0 = 50 for example, or any other number, this yields just an array of zeros and eventually NaN.

Please help.

  0 件のコメント

サインイン to comment.

回答 (1 件)

Walter Roberson
Walter Roberson 2016 年 8 月 4 日
You use
F = @(t,r) 1/(50*Power(t,2));
Power appears to be undefined in the code you gave. Perhaps you are thinking of power
Your F function takes t (presumably time) as the first variable and r (unknown purpose) as its second column, and it ignores the second parameter. But you are invoking
k1 = F(y(:,i),t(i));
in which you appear to be passing y as the first parameter and a time as the second parameter. That k1 y(:,i) is going to be received into the parameter that F knows as t, and that k1 t(i) is going to be received into the parameter that F knows as r (where it is going to be ignored).
For the sake of the readers, you need to be consistent about whether it is your first parameter or the second parameter of F that is the time parameter.

  1 件のコメント

Torsten
Torsten 2016 年 8 月 4 日
... and t0 should be greater than 0.
Best wishes
Torsten.

サインイン to comment.

サインイン してこの質問に回答します。


Translated by