Euler’s Method Or ode45 for solving ODE for control systems

Question

The dominant approach for solving ODE in control systems books is ode45 since the majority of these books use Matlab. I'm not acquainted with how the ode45 works but lately I started reading about Euler's method in this book Numerical Methods for Engineers. If the step size is very small, then the results are satisfactory. For simulation, one can actually set the step size to be very small value. I've used ode45 in here for regulation and tracking problems. I faced some difficulties for using ode45 for tracking problem since the step size is not fixed. Now for the same experiment, I've used the Euler's method with step size 0.001 sec. The results are amazing and so friendly in comparison with ode45. This is a snapshot from the result

And this is the code

clear all;
clc;

dt = 0.001;
 t = 0;

% initial values of the system
 a = 0; % angular displacement
da = 0; % angular velocity 

% PID tuning
Kp = 50;
Kd = 18.0;
Ki = 0.08;

error  = 0;

%System Parameters:
m = 0.5;       % mass (Kg)
d = 0.0023e-6; % viscous friction coefficient
L = 1;         % arm length (m)
I = 1/3*m*L^2; % inertia seen at the rotation axis. (Kg.m^2)
g = 9.81;      % acceleration due to gravity m/s^2

% Generate Desired Trajectory 
y = 0:dt:(3*pi)/2;
AngDes = y; % Ang: angle , Des: desired
AngDesPrev = 0;

for i = 1:numel(y)

   % get the first derviative of the desired angle using Euler method.
   dAngDes  =  ( AngDes(i) - AngDesPrev )/ dt;
   AngDesPrev = AngDes(i);

   % torque input
   u = Kp*( AngDes(i) - a ) + Kd*( dAngDes - da ) + Ki*error;
   % accumulated error
   error = error + ( AngDes(i) - a );

   %store the erro
   E(i) = ( a - AngDes(i) );
   T(i) = t;


   dda = 1/I*(u - d*da - m*g*L*sin(a));

   % get the function and its first dervative
   da = da + dda*dt;
    a = a  +  da*dt;

   %store data for furhter investigation 
   A(i) = a;
   dA(i) = da;

   t = t + dt; 

end

plot(T, AngDes, 'b', T, A, 'g', 'LineWidth', 1.0)
h = legend('$\theta_{d}(t)$', '$\theta(t)$');
set(h, 'Interpreter','LaTex')

My question is why ode45 is preferred in many control books assuming the step size is very small.

Answer 1

Euler method is horrible for general ODE timestepping because it is globally only first order accurate. This means that to decrease the error by a factor of 1/2, you'll need to decrease the timestep by half. Effectively, this also means that you'll need twice as many timesteps to improve your accuracy at the same final timestep. Furthermore, the error tends to increase quite dramatically as you take more and more timesteps. This is especially true for more complex systems.

The system you posted in your code appears to be very well behaved; no extreme oscillations or sharp changes in the gradient of the solution. It doesn't surprise me that euler and runge-kutta yield nearly the same solution in this case. I caution you that this is generally not the case for most systems. I highly discourage the use of forward euler method for general purpose ode timestepping.

The Runge-Kutta (ode45) method is more popular because it is fourth order accurate globally (and fifth order locally, hence the numbers 45 in the name). This means that if you decrease your timestep by 1/2, your error decreases by a factor of 1/16. The extra accuracy you gain from the 4th order accuracy of runge-kutta far outweighs the extra computational cost. This is especially true if you need adaptive timestepping, as @fibonatic suggests, because you may not need as many refinement steps if you use a higher order accurate method like runge-kutta.