Sometimes it is hard to understand how numerical expressions are evaluated. I remember reading claims by Wolfram on how smart the Kernel is to evaluate expressions trees numerically by recognizing patterns, yet I don't see how it applies in very simple examples.
This question is about numerics, but to see if symbolics can help in an elegant way.
It is known that when working with finite precision, the function $\log(1+x)$ should have a special implementation for small $x$. That is why functions like
log1p exists in many libraries (on top of
log). For example: http://gnu.ist.utl.pt/software/gsl/manual/html_node/Elementary-Functions.html#Elementary-Functions
/*C code*/ log(1. + 1.e-15) == 1.11022e-15
/*C code*/ log1p(1.e-15) == 1.e-15
(The second version is more exact, the first is "wrong")
Log[1. + 1.*^-15] == 1.11022*10^-15
Mathematica doesn't have such
Log1P function. One can say, well, that is because it doesn't need to, because of the symbolic power. In fact one can know the answer.
N[Log[1 + 1/10^15], 100] = 9.999999999999995000000000000003333333333333330833333333333335333333333333331666666666666668095238095*10^-16
But this is not general, if I want to evaluate
Log[1 + x] and x has machine precision then I can't force to use something like
1+x will evaluate to a machine precision number.
These are my attempts:
x = 1.*^-15; Log[1 + x] Log[1. + x] N[Log[1. + N[x, 20]], 20] N[Log[1 + N[x, 20]], 20]
All evaluate to the wrong answer (
Finally I find this expression,
N[Log[1 + Rationalize[x, 0]], 20] 9.999999999999995000*10^-16
But really? is it that hard to get
log(1+x) for small x numerically?
Do I have to roll my own