%
% All vectors are assumed to be *column* vectors.
%
- alpha = -(g' * p)/(p' * Q * p);
+ denom = (p' * Q * p);
+
+ % denom is non-negative, since it's a Q-norm. No need to abs() it.
+ if (denom < eps)
+ % Catch divide-by-zeros. If denom is effectively zero, set it to
+ % something tiny instead. This trick is also used in the PCGM.
+ denom = eps;
+ end
+
+ alpha = -(g' * p)/denom;
end