function y = logdq(x)
% function y = logdq(x)
% Accurately evaluates the function
% 
% y = log(1+x)/x   by   y = log(1+x)/[(1+x)-1]
%
% Follows Question 1.9 in Demmel, pg 25.
% 
d = 1+x;
y = ones(size(x));
% 
% Vectorized version. Initialize y to all 1's;
% then change those that correspond to a d ~= 1.
% 
j = d ~= 1;
if ~isempty(j)
  y(j) = log(d(j)) ./ (d(j)-1);
end