Convert higher order to first order system

Question

kingsley 2017-5-17

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/340631-convert-higher-order-to-first-order-system

回答： Torsten 2017-5-17

I'm trying to solve the higher order ode by using RK4 method. Here is the code I have so far.

function [y] = rk4_high_ode(a,t0,n,h,f)
%f = @(t,y) 2*y(2)-y(3)+2*y(4);
   % Do I need to define y(2) derivative first?
   F=@(t,y)[y(2:end);f]; % Convert the higher order to the 1st order system
   t(1)=t0;
  for i=1:n
      % update time
       t(i+1)=t(i)+h;
      k1=F(t(i)      ,y(i)         );
      k2=F(t(i)+0.5*h,y(i)+0.5*h*k1);
      k3=F(t(i)+0.5*h,y(i)+0.5*h*k2);
      k2=F(t(i)+h    ,y(i)+h*k1    );
      y(i+1)=y(i)+h/6*(k1+2*k2+2*k3+k4);
  end
  end

And this is the test program:

clear
% test y=sin(t)
% y^(4) = 2*y'-y"+2*y^(3)
t0 = 0.1; n = 100; h = 1e-2;
a = [sin(t0) cos(t0) -sin(t0) -cos(t0)]';
f = @(t,y) 2*y(2)-y(3)+2*y(4);
y = rk4_high_ode(a,t0,n,h,f);
ye = sin(t0+n*h);
error2 = abs(y-ye)