Solving a linear system of 1st order ODEs using the Forward Euler method and ode45

I need to solve a linear system of 1st order ODEs using the Forward Euler method and ode45.

The equations to solve:

dxdt = (-1/2)*x+y;

dydt = x+(1/3)*y;

this is what I have so far:

xf = 4; %final time
xspan = [0,xf];
h = .05; %stepsize
N = xf/h;
Y0 = [6,-4.2];
myx = zeros(N+1,1);
myY = zeros(N+1,2);
myY(1,:) = Y0;


% Integrate using Forward Euler
for k = 1:N
myY(k+1,:) = myY(k,:) +...
h*(jump(myx(k),myY(k,:)))';
myx(k+1) = myx(k) + h;
end

% Integrate using ode45
[x,y] = ode45(@jump,xspan,Y0);

plot(x,y(:,2),'ko',myx(1:skip:end),myY(1:skip:end,2),'r*')
title('Trajectory of an object');
xlabel('Time (s)'); ylabel('Position (m)');
legend('ode45','forward Euler')


function dxdt = jump(x,y)
dxdt = [-(1/2)*x+y;x+(1/3)*y];
end

Error I'm getting:

Unable to perform assignment because the size of the left side is 1-by-2 and the size of the

right side is 2-by-2.

Error (line 15)

myY(k+1,:) = myY(k,:) +...

I also need different initial positions for both methods, (6m, -4.2m) for Forward Euler and (6m, -3.8m) for the ode45 and I'm not sure how to implement these

Answer 1

[Alan Stevens]

You were confusing time and x position. Try

tf = 4; %final time  You are integrating wrt time not x
tspan = [0,tf];
h = .05; %stepsize
N = tf/h;
XY0 = [6,-4.2];   % [x0, y0]  for Euler
myXY = zeros(N+1,2);
myXY(1,:) = XY0;
tE = zeros(N+1,1); %Euler times

% Integrate using Forward Euler
for k = 1:N
    tE(k+1) = k*h; 
    myXY(k+1,:) = myXY(k,:) + h*jump(tE(k),myXY(k,:))';
end
xEuler = myXY(:,1);
yEuler = myXY(:,2);

XY0 = [6,-3.8];   % [x0, y0]  for ODE45

% Integrate using ode45
[t,XY] = ode45(@jump,tspan,XY0);
xODE45 = XY(:,1);
yODE45 = XY(:,2);

% Plot trajectories through time
plot(t,xODE45,'k',t,yODE45,'k--',tE,xEuler,'r',tE,yEuler,'r--'),grid
xlabel('time [s]'),ylabel('x and y positions')
legend('x ODE45','y ODE45','x Euler','y Euler')

% Plot trajectories in space
figure
plot(xODE45,yODE45,'ko',xEuler,yEuler,'r*'),grid
title('Trajectory of an object');
xlabel('x position (m)'); ylabel('y position (m)');
legend('ode45','forward Euler')

[Elle Rae]

I'm still getting the same error

Unable to perform assignment because the size of the left side is 1-by-2 and the size of the

right side is 2-by-2.

[Alan Stevens]

Sorry, I didn't include the modified jump function. Put the following on the end of my previous code:

function dxdt = jump(~,XY)
x = XY(1); y = XY(2);
dxdt = [-(1/2)*x+y;x+(1/3)*y];
end

[Elle Rae]

That worked! Thank you so much!