Need assistance with incorrect output in plot.

1 次查看（过去 30 天）
Sean Sweeney
Sean Sweeney 2021-3-27
评论： Sean Sweeney 2021-3-29
Hi all. I am working on a class project that takes a set of four first order equations and integrates them into a plot. The plot should appear as this below but my graph is totally different and I am not sure what is causing the issue. Can anyone point out where I am going wrong? I am including my full code below. Thanks in advance
%% Creating custom R-K Integrator
clear
clc
dt=0.0001 %0.01 sec was used initially,
x_vect=[0,0,0,0]; %Initial Conditions
t_span=[0:dt:2];
npoints=length(t_span);
x_dot = @deriv
t0=0;
t(1)=t0;
x_answer(1)=x_vect(1,4);% Initialize answer vector.
for i=1:npoints
t(i+1) = t(i)+dt; %update t, this value will not be used until next iteration.
x_vect = integ(x_vect,x_dot,t(i),dt);
x_answer(i+1)=x_vect(1,4);
t0=t(i);
end
plot(t,x_answer)
function out=integ(x_vect,x_dot,t,dt) %fourth order runge kutta
%out=x_vect+x_dot*dt;
k1=x_dot(t,x_vect)*dt;
k2=x_dot(t+dt/2,x_vect+k1/2)*dt;
k3=x_dot(t+dt/2,x_vect+k2/2)*dt;
k4=x_dot(t+dt,x_vect+k3)*dt;
out=x_vect+1/6*(k1+2*k2+2*k3+k4);
end
function out=deriv(t,y_init)
out(1) = y_init(2);
out(2) = y_init(3);
out(3) = y_init(4);
out(4) = -206*y_init(4)-1341.4*y_init(3)-848.4*y_init(2)+60000;
end
  5 个评论
显示 3更早的评论
darova
darova 2021-3-29
@Sean Sweeney
x_dot according to my reacher is supposed to represent the derivitive.

请先登录，再进行评论。

请先登录，再回答此问题。

回答（1 个）

darova
darova 2021-3-29
Try to change this part
x_answer = zeros(npoints,4); % preallocate answer vector
x_answer(1,:) = x_vect(1,4); % Initialize answer vector.
And this part inside for loop
x_answer(i+1,:) = x_vect
  1 个评论
Sean Sweeney
Sean Sweeney 2021-3-29
x_dot = @deriv
t0=0;
t(1)=t0;
x_answer = zeros(npoints,4);
x_answer(1,:) = x_vect(1,4);
for i=1:npoints
t(i+1) = t(i)+dt; %update t, this value will not be used until next iteration.
x_vect = integ(x_vect,x_dot,t(i),dt);
x_answer(i+1,:) = x_vect
t0=t(i);
end
The bold is what I updated from the original code, but it just ran forever until I manually ended the script. Was this the change you were suggesting?

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 Programming 的更多信息

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by