% przed użyciem funkcji należy ją zapisać do pliku o nazwie `rk4.m`

%%%%%%%%%%%%%%%%
%
% Metoda RK4 
%
%%%%%%%%%%%%%%%%

function x = rk4(f, x0, t0, T, N)
% Dla równania $x' = f(t,x), \quad t_0 < t \leq T$ z warunkiem początkowym $x(t_0) = x_0$  
% wyznacza x(T)
% metodą RK4 z krokiem $h = (T-t_0)/N$
x = x0; t = t0; h = (T-t0)/N;
for i = 1:N
	K1 = feval(f, t, x);
	K2 = feval(f, t + 0.5*h, x + 0.5*h*K1);
	K3 = feval(f, t + 0.5*h, x + 0.5*h*K2);
	K4 = feval(f, t + h, x + h*K3);
	x = x + (h/6)*(K1 + 2*(K2+K3) + K4); % nowy x to x(t+h)
	t = t+h;
end
end
