\SweaveOpts{engine=R, keep.source=TRUE}
\SweaveOpts{eps=FALSE, pdf=TRUE, width=9, height=6, strip.white=true}
\setkeys{Gin}{width=\textwidth}

<<preliminaries, echo=FALSE, results=hide>>=
options(SweaveHooks= list(fig=function() par(mar=c(5.1, 4.1, 1.1, 2.1))),
        width = 75,
        digits = 7, # <-- here, keep R's default!
        prompt = "> ",
        continue="  ")
try(Mlibrary(Rmpfr))
stopifnot(require("Rmpfr"))
@
\subsection{Why $\log(1+x)$ is not good enough, but $\mathrm{log1p}(x)$ is}
\begin{frame}[fragile]{Why $\log(1+x)$ is not good enough, but $\mathrm{log1p}(x)$ is}
  $1+x$ cannot be numerically accurate when $|x| \ll 1$.
  In double precision (53 bits $\approx 16$ digits) accuracy,
  $1+x$ ``sees'' only 2--3 digits of $x$ when $x = 10^{-14}$,
<<ex-in-R>>=
u <- 1 + (e <- 4e-13/9)  ## then   u - 1 == e  mathematically:
rbind(`u-1` = u - 1, e)
@

And the consequence for $\log(1 + x)$,
<<log1p-err, eval=FALSE>>=
curve(abs(1 - log(1+x) / log1p(x)), 1e-17, .2,  log = 'xy', main = "..", ..)
sfsmisc :: eaxis(1); eaxis(2)
<<log1p-err-do, fig=TRUE, height=3, echo=FALSE>>=
curve(abs(1 - log(1+x)/log1p(x)), 1e-17, .2,  log='xy',
      n=1024, ylab='', main="| relative error of log(1+x) |",
      lwd = .25, col="red", axes=FALSE, frame=TRUE)
sfsmisc :: eaxis(1); eaxis(2)
@
\end{frame}

\begin{frame}{Why $\mathrm{log1p}(x)$ beats $\log(1+x)$}
  \begin{center}
    \includegraphics[width=\textwidth]{log1p-ex-log1p-err-do}
  \end{center}
  Solution: Expand $\log(1+x)$ around $x=0$.
  Well known
  \[  \log(1+x) = x - x^2/2 + x^3/3 \pm \ldots = \sum_{n=1}^\infty
  (-1)^{n+1} \frac{x^n}{n},
  \]
  for $|x| < 1$.

  \medskip

  Fast version of this expansion: typically used in
  % the  implementation of
  \Rfun{log1p}.
\end{frame}

%%% Local Variables:
%%% TeX-command-default: "LaTeX PDF"
%%% TeX-master: "Maechler-2014.tex"
%%% End: