function u = Timestepping(u,T)

Globals1D;

u_temp = zeros(N_p, K);

% time-step choice for stable scheme
CLF = 0.75; dt = (CLF / abs(a_advec)) * min(abs(x(1,:) - x(2,:)));
% dt = 0.25 * dt; % kvuli omezeni chyby v dusledku casove diskretizace
N_steps = ceil(T / dt); dt = T / N_steps;

current_time = 0;
for timestep = 1:N_steps
	for i = 1:5 % 5 stages of low-storage Runge-Kutta method
		time_RK = current_time + rk4c(i) * dt;
		u_rhs = ODEs_RHS(u,time_RK);
		u_temp = rk4a(i) * u_temp + dt * u_rhs; % rk4a(1) = 0, proto bude u_temp = dt * u_rhs pri kroku i = 1
		u = u + rk4b(i) * u_temp;
	end
	current_time = current_time + dt;
    
    if (rem(timestep,100) == 0)
        fprintf(2,'%i of total %i timesteps processed.\n',timestep,N_steps);
    end

end

end