% Function format for ODE solvers
  function f = fut(z)
% z contains the state vector and the instant of time
% z = [ u_1 ; u_2 ; ... ; u_N ; t]

  u = z(1:end-1,1) ; t = z(end,1) ;
  
% Harmonic oscillator (omega = 1)
  A = [0 1 ; -1 0] ; f = A*u;

 
