¿Cómo reducir el error numérico al calcular $\ln(x) - \ln(y)$ cuando $x \approx y$?

Question

$\newcommand{\Cond}{\operatorname{Cond}}$Habrá un gran error de cancelación al calcular $\ln(x)-\ln(y)$ para $x \approx y$.

Un método para resolver este problema es reescribir $\ln(x)-\ln(y)$ como $\log(\frac{x}{y})$ pero el número de condición para $f(u) = \ln(u)$ es infinito para $u = 1$ por lo que no es un método bien condicionado.

$$\Delta x = h$$ $$\Cond(f,x) = \frac{\frac{\Delta f}{f}}{\frac{\Delta x}{x}} = \frac{\frac{f(x+h)-f(x)}{h} \times \frac{h}{f(x)}}{\frac{h}{x}} = \frac{x f'(x)}{f(x)}$$ $$ \Cond(\ln(x), x) = \frac{1}{\ln(x)}$$ $$\Cond(\ln(x), 1) = \infty$$

¿Cómo podemos reducir el error al calcular $\ln(x) - \ln(y)$? ¿Ayuda reescribirlo como $\ln(1+\frac{x-y}{y})$?

Answer 1

Suponga que $x$, $y$ son exactos y $x\approx y$. Sea $z$ denote $(x-y)/y$ y sea $z^{*} = \text{fl}((x-y)/y)$ el valor computado de $z$. Como $x-y$ se puede calcular sin ningún error (por el lema de Sterbenz como se señaló en un comentario), $z^{*} = z\times(1+\alpha_1\delta_1)$, donde $|\delta_1|\leq \varepsilon$, $\varepsilon$ es la precisión de la máquina y $\alpha_1$ es la precisión de la división de punto flotante, usualmente $\alpha_1 = 1/2$.

Sea $\text{log1p}(x) = \text{ln}(1+x)$ y sea $\overset{\sim}{\text{log1p}}(x)$ el valor computado de $\text{log1p}(x)$.

Tenemos $\overset{\sim}{\text{log1p}}(z^{*}) = \text{log1p}(z^*)\times(1+\alpha_2\delta_2)$, donde $|\delta_2|\leq \varepsilon$, y $\alpha_2$ es la precisión de la función $\text{log1p}$.

Así: $$\overset{\sim}{\text{log1p}}(z^{*}) = \text{log1p}((1+\alpha_1\delta_1)z)\times(1+\alpha_2\delta_2) = \text{log1p}(z)\times(1+\text{cond}(\text{log1p}, z)\alpha_1\delta_1)\times(1+\alpha_2\delta_2) + O(\varepsilon^2) = \text{log1p}(z)\times(1+\text{cond}(\text{log1p}, z)\alpha_1\delta_1+\alpha_2\delta_2) + O(\varepsilon^2)$$ donde el número de condición de la función $\text{log1p}$ es $$\text{cond}(\text{log1p}, z) = \frac{z}{\text{ln}(z + 1)\times(z + 1)}$$ y para $z\approx 0$, $\text{cond}(\text{log1p}, z)= 1 - \frac{1}{2}z + O(z^2) \approx 1$. Esto muestra que el error relativo hacia adelante para $x\approx y$: $$\left|\frac{\overset{\sim}{\text{log1p}}(z^{*}) - \text{log1p}(z)}{\text{log1p}(z)}\right| \lesssim (\alpha_1 + \alpha_2)\varepsilon$$ es pequeño y muy cercano a la precisión de la función $\text{log1p}$.

Observa que no se puede usar $\text{ln}(1+z)$ en lugar de $\text{log1p}(z)$, debido al enorme número de condición de $\text{ln}(1+z)$ para $z\approx 0$.

Si $x$ o $y$ no son exactos, entonces también debemos tener en cuenta el número de condición de la función $\text{ln}(x)- \text{ln}(y)$. Como observaste, este número de condición es muy alto para $x\approx y$. En este caso es imposible obtener un algoritmo preciso para calcular esta función, y la única solución es usar aritmética de alta precisión.

Answer 2

Es correcto que el número de condición de la evaluación de la función en cualquier cero de la función es infinito. Esto es equivalente al hecho de que "ningún número está cerca de cero". Sin embargo, esto no es realmente un problema en la práctica. Por esta razón, Gautschi ha definido un "número de condición suavizado"

$$ \mathrm{cond}_{\mathrm{mol}}(f, x) := \frac{m(x)|f'(x)|}{m(f(x))} $$ donde el suavizador $m$ está dado por

$$ m(x) := \begin{cases} 1, & |x| < 1 \\ |x|, & \mathrm{otherwise} \end{cases} $$ Es fácil ver que $\mathrm{cond}_{\mathrm{mol}}(\ln, 1) = 1$, así que simplemente use $\log(x/y)$, o quizás $\mathrm{log1p}(x/y-1)$.

Aquí está el razonamiento completo de Gautschi sobre la introducción de esta definición ($y=f(x)$):

Si $x = 0$ y $y = 0$, se debe considerar una perturbación absoluta de $x$ y una perturbación relativa de $y, y viceversa si $x = 0$ y $y = f(x) = 0$. Los respectivos números de condición son entonces $ (\mathrm{cond} \, f)(x) = |f(x)/f(x)| $ y $ (\mathrm{cond} \, f)(x) = |xf(x)| $. Si $x = y = 0$, el número de condición apropiado es el de la Definición 2.3. Dado que estas modificaciones crean discontinuidades en el comportamiento del número de condición, es preferible, y de hecho bastante natural, utilizar error relativo solo si la cantidad en cuestión es mayor que 1 en módulo, y error absoluto de lo contrario.