clear all;
close all;
Delta_t=0.2;
t=0:Delta_t:8;
x=zeros(length(t),2);
x(1,:)=[1.0 1.0];
for n=1:length(t)-1
k1=f(t(n),x(n,:))’;
k2=f(t(n)+Delta_t,x(n,:)+Delta_t*k1)’;
x(n+1,:)=x(n,:)+Delta_t*(k1/2.0+k2/2.0);
end
[tmat,xmat]=ode23(@f,[0 8],[1 1]);
hold off
plot(t,x(:,2),’bo-‘,’linewidth’,2,’markersize’,10)
hold on
plot(tmat,xmat(:,2),’k-‘,’linewidth’,2)
xlabel(‘t’);
ylabel(‘x_1(t)’);
function dxdt=f(t,x)
dxdt=[-6*x(1)-3*x(2) ;5*x(1)+2*x(2)];
end