Suma por pares

En análisis numérico , la suma por pares , también llamada suma en cascada , es una técnica para sumar una secuencia de números de punto flotante de precisión finita que reduce sustancialmente el error de redondeo acumulado en comparación con la acumulación ingenua de la suma en secuencia. ^[1] Aunque existen otras técnicas, como la suma de Kahan , que normalmente tienen errores de redondeo aún más pequeños, la suma por pares es casi tan buena (difiere solo en un factor logarítmico) y tiene un costo computacional mucho menor: se puede implementar de manera que tienen casi el mismo costo (y exactamente el mismo número de operaciones aritméticas) que la suma ingenua.

En particular, la suma por pares de una secuencia de n números x _n funciona dividiendo recursivamente la secuencia en dos mitades, sumando cada mitad y sumando las dos sumas: un algoritmo de divide y vencerás . Sus errores de redondeo en el peor de los casos crecen asintóticamente como máximo como O (ε log  n ), donde ε es la precisión de la máquina (suponiendo un número de condición fijo , como se analiza a continuación). ^[1] En comparación, la técnica ingenua de acumular la suma en secuencia (sumando cada x _i uno a la vez para i  = 1,...,  n ) tiene errores de redondeo que, en el peor de los casos, crecen como O (ε n ). ^[1] La suma de Kahan tiene un error en el peor de los casos de aproximadamente O (ε), independiente de n , pero requiere varias veces más operaciones aritméticas. ^[1] Si los errores de redondeo son aleatorios y, en particular, tienen signos aleatorios, entonces forman un paseo aleatorio y el crecimiento del error se reduce a un promedio de para la suma por pares. ^[2] ${\displaystyle O(\varepsilon {\sqrt {\log n}})}$

En muchos algoritmos de transformada rápida de Fourier (FFT) se encuentra una estructura recursiva de suma muy similar , y es responsable de la misma acumulación lenta de redondeo de esas FFT. ^{[2] [3]}

el algoritmo

En pseudocódigo , el algoritmo de suma por pares para una matriz x de longitud n ≥ 0 se puede escribir:

s = por pares ( x [1... n ]) si  n ≤ N  caso base: suma ingenua para una matriz suficientemente pequeña  s = 0 para  i = 1 to n  s = s + x [ i ] de lo contrario  divide y vencerás: recursivamente sumar dos mitades de la matriz  m = piso ( n / 2) s = por pares ( x [1... m ]) + por pares ( x [ m +1... n ]) fin si

Para algunos N suficientemente pequeños , este algoritmo cambia a una suma ingenua basada en bucles como caso base , cuyo límite de error es O(Nε). ^[4] La suma completa tiene un error en el peor de los casos que crece asintóticamente como O (ε log  n ) para n grande , para un número de condición dado (ver más abajo).

En un algoritmo de este tipo (como en los algoritmos de divide y vencerás en general ^[5] ), es deseable utilizar un caso base más grande para amortizar la sobrecarga de la recursividad. Si N  = 1, entonces hay aproximadamente una llamada de subrutina recursiva para cada entrada, pero de manera más general hay una llamada recursiva para (aproximadamente) cada N /2 entradas si la recursión se detiene exactamente en  n  =  N. Al hacer que N sea lo suficientemente grande, la sobrecarga de la recursividad puede hacerse insignificante (precisamente esta técnica de un caso base grande para la suma recursiva es empleada por implementaciones FFT de alto rendimiento ^[3] ).

Independientemente de N , se realizan exactamente n −1 adiciones en total, lo mismo que para la suma ingenua, por lo que si la sobrecarga de recursividad se hace insignificante, entonces la suma por pares tiene esencialmente el mismo costo computacional que para la suma ingenua.

Una variación de esta idea es dividir la suma en b bloques en cada etapa recursiva, sumando cada bloque de forma recursiva y luego sumando los resultados, lo que sus proponentes denominaron algoritmo de "superbloque". ^[6] El algoritmo por pares anterior corresponde a b  = 2 para cada etapa excepto para la última etapa que es  b  =  N.

Exactitud

Supongamos que uno está sumando n valores x _i , para i  = 1, ...,  n . La suma exacta es:

${\displaystyle S_{n}=\sum _{i=1}^{n}x_{i}}$

(calculado con precisión infinita).

Con la suma por pares para un caso base N  = 1, se obtiene , donde el error está limitado arriba por: ^[1] ${\displaystyle S_{n}+E_{n}}$ ${\displaystyle E_{n}}$

${\displaystyle |E_{n}|\leq {\frac {\varepsilon \log _{2}n}{1-\varepsilon \log _{2}n}}\sum _{i=1}^{n}|x_{i}|}$

donde ε es la precisión de la máquina de la aritmética que se emplea (por ejemplo, ε ≈ 10 ⁻¹⁶ para coma flotante estándar de doble precisión ). Generalmente, la cantidad de interés es el error relativo , que por lo tanto está limitado por: ${\displaystyle |E_{n}|/|S_{n}|}$

${\displaystyle {\frac {|E_{n}|}{|S_{n}|}}\leq {\frac {\varepsilon \log _{2}n}{1-\varepsilon \log _{2}n}}\left({\frac {\sum _{i=1}^{n}|x_{i}|}{\left|\sum _{i=1}^{n}x_{i}\right|}}\right).}$

En la expresión del límite de error relativo, la fracción (Σ| x _i |/|Σ x _i |) es el número de condición del problema de suma. Esencialmente, el número de condición representa la sensibilidad intrínseca del problema de suma a los errores, independientemente de cómo se calcule. ^[7] El límite de error relativo de cada método de suma ( estable hacia atrás ) mediante un algoritmo fijo con precisión fija (es decir, no aquellos que usan aritmética de precisión arbitraria , ni algoritmos cuyos requisitos de memoria y tiempo cambian según los datos), es proporcional a este número de condición. ^[1] Un problema de suma mal condicionado es aquel en el que esta relación es grande y, en este caso, incluso la suma por pares puede tener un error relativo grande. Por ejemplo, si los sumandos x _i son números aleatorios no correlacionados con media cero, la suma es un paseo aleatorio y el número de condición crecerá proporcionalmente a . Por otro lado, para entradas aleatorias con valores distintos de cero, el número de condición asíntotas de una constante finita es . Si todas las entradas no son negativas , entonces el número de condición es 1. ${\displaystyle {\sqrt {n}}}$ ${\displaystyle n\to \infty }$

Tenga en cuenta que el denominador es efectivamente 1 en la práctica, ya que es mucho menor que 1 hasta que n se vuelve del orden 2 ^1/ε , que es aproximadamente 10 ^{10 15} en doble precisión. ${\displaystyle 1-\varepsilon \log _{2}n}$ ${\displaystyle \varepsilon \log _{2}n}$

En comparación, el error relativo límite para la suma ingenua (simplemente sumando los números en secuencia, redondeando en cada paso) crece al multiplicarse por el número de condición. ^[1] En la práctica, es mucho más probable que los errores de redondeo tengan un signo aleatorio, con media cero, de modo que formen un paseo aleatorio; en este caso, la suma ingenua tiene un error relativo cuadrático medio que crece a medida que la suma por pares tiene un error que crece en promedio. ^[2] ${\displaystyle O(\varepsilon n)}$ ${\displaystyle O(\varepsilon {\sqrt {n}})}$ ${\displaystyle O(\varepsilon {\sqrt {\log n}})}$

Implementaciones de software

La suma por pares es el algoritmo de suma predeterminado en NumPy ^[8] y el lenguaje de computación técnica Julia , ^[9] donde en ambos casos se encontró que tiene una velocidad comparable a la suma ingenua (gracias al uso de un caso base grande).

Otras implementaciones de software incluyen la biblioteca HPCsharp ^[10] para el lenguaje C# y la suma de bibliotecas estándar ^[11] en D.

Referencias

1. Higham, Nicholas J. (1993), "La precisión de la suma en coma flotante", SIAM Journal on Scientific Computing , 14 (4): 783–799, CiteSeerX  10.1.1.43.3535 , doi :10.1137/0914050
2. Manual de métodos analítico-computacionales de Manfred Tasche y Hansmartin Zeuner en matemáticas aplicadas Boca Raton, FL: CRC Press, 2000).
3. SG Johnson y M. Frigo, "Implementación de FFT en la práctica, en transformadas rápidas de Fourier , editado por C. Sidney Burrus (2008).
4. Higham, Nicolás (2002). Precisión y estabilidad de algoritmos numéricos (2 ed) . SIAM. págs. 81–82.
5. Radu Rugina y Martin Rinard, "Desarrollo de recursión para programas de divide y vencerás", en Lenguajes y compiladores para computación paralela , capítulo 3, págs. Apuntes de conferencias sobre informática vol. 2017 (Berlín: Springer, 2001).
6. Anthony M. Castaldo, R. Clint Whaley y Anthony T. Chronopoulos, "Reducción del error de punto flotante en el producto escalar utilizando la familia de algoritmos de superbloques", SIAM J. Sci. Computadora. , vol. 32, págs. 1156-1174 (2008).
7. LN Trefethen y D. Bau, Álgebra lineal numérica (SIAM: Filadelfia, 1997).
8. ENH: implementar la suma por pares, github.com/numpy/numpy pull request #3685 (septiembre de 2013).
9. RFC: use suma por pares para suma, cumsum y cumprod, solicitud de extracción github.com/JuliaLang/julia n.º 4039 (agosto de 2013).
10. https://github.com/DragonSpit/HPCsharp Paquete HPCsharp nuget de algoritmos C# de alto rendimiento
11. "std.algorithm.iteration - Lenguaje de programación D". dlang.org . Consultado el 23 de abril de 2021 .