function [fnew, wnew, Awnew, gnew, D, stepsize] = restricted_linesearch(loss, grad, A, b, sigma, w, p, g, f, Aw, Ap, init_step, epsilon)
	succeed = 0;

% 	if (~exist('init_step'))
% 		init_step = 1; %Default is Newton, so start from unit step size
% 	end

	delta = sigma*(g'*p);

	stepsize = init_step;
	if (delta> 0)
		stepsize = 0;
	end

	while (stepsize > 1e-4)
		wnew = w + stepsize * p;
		Awnew = Aw + stepsize * Ap;
		[fnew, ~, gnew, D] = loss_fun_and_grad(loss, grad, wnew, A, b, epsilon, Awnew);
		if (fnew < f + stepsize * delta)
			succeed = 1;
			break;
		end
		stepsize = stepsize / 2;
	end

	if (succeed == 0)
		stepsize = 0;
		gnew = g;
		wnew = w;
		D = 0;
		Awnew = Aw;
		fnew = f;
		return;
	end
end
