Solving 2nd degree ODE with Euler method in MATLAB

I am trying to solve the equation below;

$$\ddot{x}= -x + sin(t)$$

by the initial conditions;

$$x(0) = 0 \\ \dot{x}(0)= 0$$

my MATLAB code is as follows:

clc
clear all

%define time step and tend

dt = 0.01;
tf = 1.0;

%create arrays for tvec, xvec and uvec

tvec = 0:dt:tf;
xvec = zeros(length(tvec),1);
xvec(1) = 0;
uvec = zeros(length(tvec),1);
uvec(1) = 0;

%run euler
for i = 1:length(tvec)-1
xvec(i+1) = xvec(i) + dt*uvec(i);
uvec(i+1) = uvec(i) + dt*(-xvec(i)+sin(i+1));
end


but this code gives me the solution at $x(1.0)$ wrong. Where am I being failed ?

-

You may want to try this:

tf = 5;
Nt = 150;
dt = tf/Nt;

t = 0:dt:tf;

x0 = 0;
u0 = 0;

x = zeros(Nt+1,1);
u = x;

x(1) = x0;
u(1) = u0;

for i = 1:Nt;

u(i+1) = u(i) + dt*(-x(i)+sin(t(i)));
x(i+1) = x(i) + dt*u(i);

end

plot(t,x,t,u,'k');
legend('x(t)','u(t)');


You use the RHS of the equation at the time $t_i$ because you are using an explicit or forward Euler's scheme. An alternative for improving consistence would be computing $x_{i+1}$ with the information of $u$ at the time $i+1$, i.e:

x(i+1) = x(i) + dt*u(i+1);


You can check that there are some differences but don't expect big changes. You can derive a matrix-form system of equations for a backward or implicit scheme in order to solve the unknowns of: $g(x_{n+1},u_{n+1},t_{n+1}) = f(x_n,u_n)$.

Cheers!

-
thank you Dmoreno! I was placing the counter in the $sin$ function as you mentioned.. –  Yirmidokuz May 21 at 22:25
@Dmoreno: You for loop should be from 1 to Nt. Fixing that may require changing other things. –  horchler May 21 at 23:40
Of course, the counter is wrong, but note that this doesn't affect much since I initialized the state vector to zero. Thanks for noting this! On the other hand, you're welcome @Yirmidokuz. Glad this helps! –  Dmoreno May 22 at 7:49
Edited my comment. Fixed now. –  Dmoreno May 22 at 7:50
1. Use sin(dt*(i)) or sin(dt*(i+1)), not sin(i+1).

2. This is a first order method, so don't expect to be very close.

-
how should I decide to choose sin(dt*(i)) or sin(dt*(i+1)) ? –  Yirmidokuz May 21 at 21:55
I am already multiplying sin(i+1) with $\text{dt}$, as $...+dt*(-xvec(i)+sin(i+1))$ –  Yirmidokuz May 21 at 21:58
Note that you are now evaluating $\sin{i}$ where $i$ is a counter. You must compute $\sin{t}$ at $t = t_i$ or, in Matlab notation, $\text{sin(t(i))}$. You might evaluate it at $i+1$ as well. –  Dmoreno May 21 at 22:12
@Yirmidokuz - As noted by Dmoreno, you should evaluate the sine at t(i), not at i. As to accuracy - it doesn't make any big difference whether you use sin(t(i)) or sin(t(i+1)). The method is first order in either case and the stability of the method is not affected by this choice. –  Hans Engler May 21 at 23:14