Un poco de historia sobre Python
Python es un gran lenguaje. Intenté varios idiomas antes: Pascal en la escuela; C, C con clases, C ++ en la universidad. Los últimos dos (tres) infundieron una fuerte aversión a la programación: en lugar de resolver el problema, te metas con asignaciones y destructores (palabras aterradoras del pasado), piensas en términos de primitivas de bajo nivel. Mi opinión es que C no es adecuado para resolver problemas educativos y científicos (en cualquier caso, en el campo de las matemáticas). Estoy seguro de que se opondrán a mí, pero no estoy tratando de imponer nada a nadie, simplemente estoy expresando mi opinión.
Python fue una vez una revelación. Por primera vez en mi vida, escribí varios niveles de abstracción superiores a los habituales en C. La distancia entre la tarea y el código se ha reducido como nunca antes.
Probablemente habría hecho esto toda mi vida en Python si no hubiera tenido que implementar de repente las pruebas estadísticas NIST. Parecería que la tarea es muy simple: hay una matriz de longitud de varios (> = 10) megabytes, hay un conjunto de pruebas que deben aplicarse a esta matriz.
¿Para qué sirve un buen numpy?
En python, hay un buen paquete numpy para trabajar con matrices, que es esencialmente una interfaz de alto nivel para bibliotecas rápidas como LAPACK. Parece que todo el conjunto de caballeros para la informática científica está disponible (Sympy, Numpy, Scipy, Matplotlib), ¿qué más se puede pedir?
Todos los que han lidiado con numpy saben que él es bueno, pero no en todo. También debemos tratar de asegurarnos de que las operaciones estén vectorizadas (como en R), es decir. en cierto sentido, uniforme para toda la matriz. Si de repente su problema por alguna razón no encaja en este paradigma, entonces tiene problemas.
¿De qué tipo de tareas no vectorizadas estoy hablando? Sí, al menos tome el mismo NIST: calcule la longitud del registro de desplazamiento lineal utilizando el algoritmo Berlekamp-Messi; calcule la longitud de la subsecuencia más larga de unidades consecutivas y así sucesivamente. No excluyo la posibilidad de que haya algún tipo de solución ingeniosa que reduzca el problema a uno vectorizado.
Astucia?Como un ejemplo del mismo NIST: era necesario calcular el número de secuencias de "conmutación", donde por "conmutación" me refiero a cambiar las unidades consecutivas (... 1111 ...) a ceros consecutivos (... 000 ... ) y viceversa. Para hacer esto, puede tomar la secuencia original sin el último elemento (x [: -1]) y restarle la secuencia desplazada por 1 (x [2:]), y luego calcular el número de elementos distintos de cero en la nueva secuencia resultante. Todos juntos se verán así:
np.count_nonzero(x[:-1] - x[1:])
Puede parecer un entrenamiento entretenido para la mente, pero esencialmente está sucediendo algo antinatural aquí, un truco que no será obvio para el lector después de un corto tiempo. Sin mencionar que esto todavía es lento: nadie ha cancelado la asignación de memoria.
Las operaciones no vectorizadas en Python son mucho tiempo. ¿Cómo lidiar con ellos si numpy no es suficiente? Por lo general, ofrecen varias soluciones:
- Numba JIT. Si ella trabajó como se describe en la página de título de Numba, entonces creo que valdría la pena terminar la historia. Había olvidado un poco en el pasado que lo que había salido mal con ella; La conclusión es que la aceleración no fue tan impresionante como esperaba, por desgracia.
- Cython OK, levanten la mano aquellos que creen que Cython es una solución realmente hermosa y elegante que no destruye el significado y el espíritu de Python. No lo creo si llegas a Cython, ya puedes dejar de engañarte y cambiar a algo menos sofisticado como C ++ y C.
- Reescribe los cuellos de botella en C y sácalos de tu amado Python. De acuerdo, pero ¿qué pasa si tengo todo el programa, se trata de cálculos y cuellos de botella? Xi no ofrece! No estoy hablando del hecho de que en esta solución necesitas saber no solo uno, sino dos idiomas: Python y C.
¡Aquí viene la JULIA!
Habiendo sido contemplado con las soluciones existentes y no haber encontrado (incapaz de programar) nada bueno, decidí reescribirlo con algo "más rápido". De hecho, si escribe "trillador de números" en el siglo XXI con soporte normal para matrices, operaciones vectorizadas "listas para usar", etc. etc., entonces su elección no es muy densa:
- Fortran Y no te rías, ¿cuál de nosotros no sacó BLAS / LAPACK de nuestros idiomas favoritos? Fortran es un lenguaje realmente bueno (¡mejor SI!) Para la informática CIENTÍFICA. No lo tomé por la razón de que desde la época de Fortran se han inventado y agregado muchas cosas a los lenguajes de programación; Esperaba algo más moderno.
- R sufre los mismos problemas que Python (vectorización).
- Matlab : probablemente sí, desafortunadamente no tengo dinero para verificar.
- Julia : el caballo está oscuro, despegará, no despegará (y era natural hasta la versión estable 1.0)
Algunas ventajas obvias de Julia
- Parece Python, al menos el mismo "de alto nivel", con la capacidad, si es necesario, de reducirse en números.
- Sin problemas con las asignaciones de memoria y similares.
- Potente sistema de tipos. Los tipos se prescriben opcionalmente y muy dosificados. Se puede escribir un programa sin especificar un solo tipo, y si lo hace, podrán hacerlo rápidamente. Pero hay matices.
- Es fácil escribir tipos personalizados que serán tan rápidos como los tipos integrados.
- Despacho múltiple. Puedes hablar de esto durante horas, pero en mi opinión, esto es lo mejor que tiene Julia, es la base del diseño de todos los programas y, en general, la base de la filosofía del lenguaje.
Gracias al envío múltiple, muchas cosas se escriben de manera muy uniforme.
Múltiples ejemplos de despacho rand() # 0 1 rand(10) # 10 0 1 rand(3,5) # 3 5 .... using Distributions d = Normal() # 0, 1 rand(d) # rand(d, 10) # 10 ...
Es decir, el primer argumento puede ser cualquier distribución (unidimensional) de Distribuciones, pero la llamada a la función seguirá siendo literalmente la misma. Y no solo distribución (es posible transmitir el RNG en sí mismo como un objeto MersenneTwister y mucho más). Otro ejemplo (en mi opinión, ilustrativo) es la navegación en DataFrames sin loc / iloc.
- 6. Las matrices son nativas, incorporadas. Vectorización fuera de la caja.
Contras para guardar silencio sobre cuál sería un crimen!
- Nuevo lenguaje Por un lado, por supuesto, un signo menos. En algo, quizás inmaduro. Por otro lado, tuvo en cuenta el rastrillo de muchos idiomas pasados, se encuentra "sobre los hombros de gigantes", ha absorbido muchas cosas interesantes y nuevas.
- Inmediatamente, es poco probable que los programas rápidos escriban. Hay algunas cosas no obvias con las que es muy fácil lidiar: pisa el rastrillo, pide ayuda a la comunidad, da un paso nuevamente ... etc. Estos son principalmente inestabilidad de tipo, inestabilidad de tipo y variables globales. En general, por lo que puedo decir por mí mismo, un programador que quiere escribir rápidamente en Julia pasa por varias etapas: a) escribe en Python. Esto es genial, y lo es, pero a veces habrá problemas de rendimiento. b) escribe como en C: prescribiendo tipos de forma manual siempre que sea posible. c) comprende gradualmente dónde es necesario (muy medido) prescribir tipos y dónde interfieren.
- Ecosistema Algunos paquetes son crudos, en el sentido de que en algún lugar constantemente algo se cae; algunos son maduros, pero inconsistentes entre sí (por ejemplo, es necesaria la conversión a otros tipos). Por un lado, esto es obviamente malo; Por otro lado, Julia tiene muchas ideas interesantes y audaces solo porque "estamos sobre los hombros de los gigantes": hemos acumulado una tremenda experiencia de pisar un rastrillo y "cómo no hacerlo", y esto es tenido en cuenta (parcialmente) por los desarrolladores de paquetes.
- Numeración de matrices de 1. Ja, bromeando, ¡esto es, por supuesto, una ventaja! No, bueno, en serio, ¿qué sucede? ¿Con qué índice comienza la numeración? Te acostumbras en 5 minutos. Nadie se queja de que los enteros se llaman enteros, no enteros. Todos estos son asuntos de gustos. Y sí, tome al menos el mismo Cormen de acuerdo con los algoritmos: todo está numerado a partir de uno allí, y a veces es inconveniente rehacerlo en Python, por el contrario.
Bueno, ¿por qué velocidad?
Es hora de recordar por qué todo comenzó.
Nota: Olvidé Python con seguridad, así que escriba sus mejoras en los comentarios, intentaré ejecutarlas en mi computadora portátil y ver el tiempo de ejecución.
Entonces, dos pequeños ejemplos con microbenchmarks.
Algo vectorizadoProblema: se suministra un vector X de 0 y 1 a la entrada de la función. Es necesario convertirlo en un vector de 1 y (-1) (1 -> 1, 0 -> -1) y calcular cuántos coeficientes de la transformada de Fourier de este vector El valor absoluto excede el límite del límite.
def process_fft(x, boundary): return np.count_nonzero(np.abs(np.fft.fft(2*x-1)) > boundary) %timeit process_fft(x, 2000) 84.1 ms ± 335 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)
function process_fft(x, boundary) count(t -> t > boundary, abs.(fft(2*x.-1))) end @benchmark process_fft(x, 2000) BenchmarkTools.Trial: memory estimate: 106.81 MiB allocs estimate: 61 -------------- minimum time: 80.233 ms (4.75% GC) median time: 80.746 ms (4.70% GC) mean time: 85.000 ms (8.67% GC) maximum time: 205.831 ms (52.27% GC) -------------- samples: 59 evals/sample: 1
No veremos nada sorprendente aquí; de todos modos, no lo consideran ellos mismos, sino que lo pasan a programas fortran bien optimizados.
Algo no vectorizableSe alimenta a la entrada una matriz de 0 y 1. Encuentre la longitud de la subsecuencia más larga de unidades consecutivas.
def longest(x): maxLen = 0 currLen = 0
function longest(x) maxLen = 0 currLen = 0 # This will count the number of ones in the block for bit in x if bit == 1 currLen += 1 maxLen = max(maxLen, currLen) else maxLen = max(maxLen, currLen) currLen = 0 end end return max(maxLen, currLen) end @benchmark longest(x) BenchmarkTools.Trial: memory estimate: 0 bytes allocs estimate: 0 -------------- minimum time: 9.094 ms (0.00% GC) median time: 9.248 ms (0.00% GC) mean time: 9.319 ms (0.00% GC) maximum time: 12.177 ms (0.00% GC) -------------- samples: 536 evals/sample: 1
La diferencia es obvia a simple vista. Consejos para "terminar" y / o vectorizar la versión numpy son bienvenidos. También quiero señalar que los programas son casi idénticos. Por ejemplo, no registré un solo tipo en Julia (compárelo con el anterior); de todos modos, todo funciona rápidamente.
También quiero señalar que las versiones presentadas no se incluyeron en el programa final, pero se optimizaron aún más; aquí se dan a modo de ejemplo y sin complicaciones innecesarias (reenvío de memoria adicional en Julia para hacer el rfft in situ, etc.).
¿Qué salió al final?
No mostraré el código final para Python y para Julia: esto es un secreto (al menos por ahora). En el momento en que me detuve para terminar la versión de Python, funcionó todas las pruebas NIST en una matriz de 16 megabytes de caracteres binarios en ~ 50 segundos. En este momento, en Julia, todas las pruebas se ejecutan en el mismo volumen ~ 20 segundos. Bien puede ser que soy un tonto, y había espacio para optimizaciones, etc. etc. Pero este soy yo, como soy, con todas sus ventajas y desventajas, y en mi opinión, no se debe considerar la velocidad esférica de los programas en el vacío en los lenguajes de programación, sino cómo puedo programarlo personalmente (sin errores realmente graves).
¿Por qué escribí todo esto?
La gente
aquí se interesó; Decidí escribir cuando llegue el momento. En el futuro, pienso analizar Julia más a fondo, con un ejemplo de resolución de algunos problemas típicos. Por qué Más idiomas, buenos y diferentes, y deja que todos los que quieran encontrar algo que le convenga personalmente.