function [back,dt,time]=TimeStepping(md,solution,dt,time)
%TIMESTEPPING - ...
%
%   Usage:
%      

if md.timestepping,
	
	time0=time-dt; %previous time. we might have to go back in time.

	back=0; %default is to not touch the time step.

	%Is thickness decreasing or increasing too fast? 
	dh=solution(end).thickness-solution(end-1).thickness;

	maxdeltah=max(abs(dh));
	if maxdeltah>md.DeltaH,
		disp(sprintf('%s%g','      thickness varying too fast: ',maxdeltah));
		disp(sprintf('%s%g%s%g%s','         decreasing time step from ',dt/md.yts,'(yr) to ',dt*md.timedec/md.yts,'(yr)'));
		back=1;
		dt=dt*md.timedec;
		time=time0+dt; %new time is previous time plus reduced time step
		return;
	end

	%Too slow?
	if maxdeltah<md.deltaH,
		disp(sprintf('%s%g','      thickness varying too slow: ',maxdeltah));
		disp(sprintf('%s%g%s%g%s','         increasing time step from ',dt/md.yts,'(yr) to ',dt*md.timeacc/md.yts,'(yr)'));
		back=0;
		dt=dt*md.timeacc;
		return;
	end
	
	
	%Is temperature decreasing or increasing too fast? 
	dT=solution(end).t_g-solution(end-1).t_g;

	maxdeltaT=max(abs(dT));
	if maxdeltaT>md.DeltaT,
		disp(sprintf('%s%g','      temperature varying too fast: ',maxdeltaT));
		disp(sprintf('%s%g%s%g%s','         decreasing time step from ',dt/md.yts,'(yr) to ',dt*md.timedec/md.yts,'(yr)'));
		back=1;
		dt=dt*md.timedec;
		time=time0+dt;
		return;
	end

	%Too slow?
	if maxdeltaT<md.deltaT,
		disp(sprintf('%s%g','      temperature varying too slow: ',maxdeltaT));
		disp(sprintf('%s%g%s%g%s','         increasing time step from ',dt/md.yts,'(yr) to ',dt*md.timeacc/md.yts,'(yr)'));
		back=0;
		dt=dt*md.timeacc;
		return;
	end

else
	%time stepping not activated.
	back=0;
	return;
end
	
