So I'm doing a project for school, where i have to solve the x''-x'-x=sin(t) equation with Runge-Kutta.
intervalmin=0;
intervalmax=1;
numnodes=5;
f=@(t,x) [x(2); x(2)+x(1)+sin(t)];
inival=0;
h=0.1;
t=zeros(1,numnodes);
t(1)=intervalmin;
x=zeros(1,numnodes);
X=zeros(1,numnodes);
x(1)=inival;
X(1)=inival;
for i = 2:numnodes
k1 = f(t(i-1),x(i-1));
k2 = f(t(i-1)+h/2,x(i-1)+(h/2)*k1);
k3 = f(t(i-1)+h,x(i-1)-h*k1+2*h*k2);
x(i)=(x(i-1)+h/6)*(k1+4*k2+k3);
X(i)=(X(i-1)+h/6)*((f(t(i-1),X(i-1)))+4*(f(t(i-1)+h/2,X(i-1)+(h/2)*k1))+(f(t(i-1)+h,X(i-1)-h*k1+2*h*k2)));
t(i) = t(i-1)+h;
end
figure
plot(t,x,'.-b',t,X,'-.r');
The code is this so far, but it has multiple errors like:
Index exceeds the number of array elements (1).
Error in Untitled5>@(t,x)[t(2);t(2)+t(1)+sin(x)]
Error in Untitled5 (line 14) k1=f(t(i-1),x(i-1));
I'm pretty new to matlab, and just can't seem to find a way to make it work.

 Respuesta aceptada

darova
darova el 23 de Mzo. de 2020

0 votos

Here is the mistake
SOlution

6 comentarios

Dávid Bak
Dávid Bak el 23 de Mzo. de 2020
Thank you, but after these corrections, I get another error: Unable to perform assignment because the left and right sides have a different number of elements.
Error in masodikprojekt_elsocsoport (line 18)
x(i)=(x(i-1)+h/6)*(k1+4*k2+k3);
darova
darova el 23 de Mzo. de 2020
Should be changed everywhere
Dávid Bak
Dávid Bak el 24 de Mzo. de 2020
intervalmin=0;
intervalmax=1;
h=0.1;
g=0.01;
numnodes=ceil(((intervalmax-intervalmin)/h)+1);
numnodes2=ceil(((intervalmax-intervalmin)/g)+1);
inival=0;
%0.1
t=zeros(1,numnodes);
x=zeros(2,numnodes);
t(1)=intervalmin;
x(:,1)=inival;
f=@(t,x) [x(2); x(2)+x(1)+sin(t)];
for i = 2:numnodes
t(i) = t(i-1)+h;
k1 = f(t(i-1), x(:,i-1));
k2 = f(t(i-1)+h/2, x(:,i-1) + (h/2)*k1);
k3 = f(t(i-1)+h/2, x(:,i-1) + (h/2)*k2);
x(:,i)=x(:,i-1)+(h/6)*(k1+4*k2+k3);
end
%0.01
o = zeros(1,numnodes2);
z = zeros(2,numnodes2);
o(1) = intervalmin;
z(:,1) = inival;
m=@(o,z) [z(2); z(2)+z(1)+sin(o)];
for i = 2:numnodes
o(i) = o(i-1)+g;
k1 = m(o(i-1), z(:,i-1));
k2 = m(o(i-1)+h/2, z(:,i-1) + (g/2)*k1);
k3 = m(o(i-1)+h/2, z(:,i-1) + (g/2)*k2);
z(:,i)=z(:,i-1)+(g/6)*(k1+4*k2+k3);
end
[TT,XX] = ode45(f,[intmin intmax],x(:,1));
figure
plot(t,x(2,:),'.-',o,z(2,:),'.-',TT,XX(:,2),'.-')
hold on
syms x(t)
dz = diff(x);
ode = diff(x,t,2) == (x(2)+x(1)+sin(t));
cond1 = x(0) == 2;
cond2 = dz(0) == 0;
dsolve(ode,cond1,cond2)
tt = linspace(0,1);
yy = 4 - 2./(tan(tt/2) + 1) - tt;
legend('h=0.1','h=0.01','ode45')
ax=gca;
ax.XAxisLocation='origin';
ax.YAxisLocation='origin';
grid on;
Whit this code, can you tell me if the red curve is supposed to be that small? I can't really work out the math part of this, so that's why i'm wondering.
darova
darova el 24 de Mzo. de 2020
Please use code button
darova
darova el 24 de Mzo. de 2020
Mistake
Dávid Bak
Dávid Bak el 24 de Mzo. de 2020
Oh, now everything's alright, thank you very much, you helped a lot!

Iniciar sesión para comentar.

Más respuestas (0)

Categorías

Más información sobre Mathematics en Centro de ayuda y File Exchange.

Productos

Versión

R2019a

Etiquetas

Preguntada:

el 23 de Mzo. de 2020

Comentada:

el 24 de Mzo. de 2020

Community Treasure Hunt

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

Start Hunting!

Translated by