[HELP] A Classical Numerical Computing Question

f(x) = (exp(x)-1)/x; g(x) = (exp(x)-1)/log(exp(x))

Analytically, f(x) = g(x).

When x is approaching to 0, both f(x) and g(x) are approaching to 1. However, g(x) works better than f(x).

% Compute y against x
for k = 1:15
    x(k) = 10^(-k);
    f(k) =(exp(x(k))-1)/x(k); 
    De(k) = log(exp(x(k)));
    g(k)= (exp(x(k))-1)/De(k);
end

% Plot y
plot(1:15,f,'r',1:15,g,'b');

f(x) actually diverges when x approaches to 0.But shouldn't them be the same??

Answer 1

[Matt Fig]

No, they shouldn't be the same at the fringes. This is an example of why we often have to look for more stable ways of doing in floating point arithmetic what is analytically simple. Look how f oscillates:

f = @(x) (exp(x)-1)./x;
g = @(x) (exp(x)-1)./log(exp(x));
x = 0:1e-13:1e-7;  % Try with x = 0:2e-14:1e-7;
ax(1) = subplot(1,2,1);
plot(x,f(x),'r')
title('(exp(x)-1)./x')
ax(2) = subplot(1,2,2);
plot(x,g(x),'b');
title('(exp(x)-1)./log(exp(x))')
L = get(gca,{'xlim','ylim'});
axis(ax(1),[L{:}])

[Xiao Tang]

That's awesome... But could you please explain why there are so many oscillates in f?

[Matt Fig]

I think what we have is a case of near perfect cancellation of errors. Take a look:

x = 1e-12:1e-12:1e-9;  % Double values
X = vpa(1/10^12:1/10^12:1/10^9,80);  % Symbolic values
syms Z
D = abs(X - log(exp(x)));
E = abs(X-x);
plot(D)
figure
plot(E)  % Very little difference
F = abs((exp(x)-1) - subs(exp(Z)-1,X)); 
figure
plot(F)  % Notice similarity to D!
G = F-D;
max(double(G))  % Cancellation of errors.

[Xiao Tang]

Thank you so much!