Automágica: durante 2017 estoy trabajando bastante en Automágica, mi software para editar libros: Más información - Posts relacionados

La soledad de los números primos

La soledad de los números primos

En las últimas semanas leí la novela La soledad de los números primos. Compré el libro atraído por el título y por ser la ópera prima de un escritor de mi edad, Paolo Giordano, italiano.

El libro habla sobre la soledad, y con saltos de varios años va contando la historia de sus dos personajes.

En una clase de primer curso Mattia había estudiado que entre los números primos hay algunos aún más especiales. Los matemáticos los llaman números primos gemelos: son parejas de números primos que están juntos, o mejor dicho, casi juntos, pues entre ellos media siempre un número par que los impide tocarse de verdad. Números como el 11 y el 13, el 17 y el 19, o el 41 y el 43. Mattia pensaba que Alice y él eran así, dos primos gemelos, solos y perdidos, juntos pero no lo bastante para tocarse de verdad.

Dos notas: It's not about maths. No me gustó el final.

Una tercera: tengo la película para ver en italiano y creo que va a estar buena.


Los 6 magníficos: divertimentos matemáticos

Cuando estaba en los primeros años de ingeniería, en la sala de consultas de Materias Básicas (el departamento que reúne las matemáticas, las físicas, las químicas, etc...) vi una calcamonía que me llamó la atención. Planteaba una igualdad entre cinco números: pi, i, e, 0 y 1.

Ayer, durante una reunión en la misma sala, volví a encontrarla. Seguía en el mismo rincón de un pizarrón, con una punta ajada. No resistí la tentación de fotografiarla y quedármela para siempre. Se las comparto:

Imagen778

Una versión, más desprolija, pero en dónde se lee mejor:

Imagen779

Una primera pregunta es: ¿Es verdad esa igualdad? Resulta que sí. Lo que me llamó la atención por años es conocido como Identidad de Euler.

Pero lo que realmente desvela es que la calcamonía se titulaba "Los seis magníficos" y remata "¿y el que falta?".¿ Alguien tiene la respuesta?


Conjuntos en Python

Cuando estaba preparando mi charla introductoria a Python, le pasé mis slides a un amigo para que me diga su opinión y una de las cosas que me dijo fue

Nunca senti que set sea algo nativo de python, le daria mas importancia a los diccionarios, aunque tal vez tu publico este mas interesado en sets, no lo se.

Me sorprendió el comentario. Para mi, set es un tipo de dato muy útil y poderoso. En este post voy a intentar hacer una apología de set, el tipo de dato que incorpora Python para representar la noción matemática de conjunto.

Presentación

Un objeto set es una colección sin orden de objetos hasheables. Puede contener objetos de todos los tipos inmutables de Python, pero no los contenedores mutables como listas. También puede contener objetos de clases definidas por el usuario (Los objetos instancias de clases definidas por el usuario son por defecto hasheables.).

Los conjuntos se pueden crear, por ejemplo, a partir de una lista. Podemos quitar elementos (al azar o uno en concreto) o agregarlos:

>>> heladera = ['huevo', 'huevo', 'queso', 'leche', 'pera', 'pera', 'pera']

>>> alimentos = set(heladera)

>>> alimentos

set(['queso', 'leche', 'huevo', 'pera'])

>>> alimentos.pop()

'queso'

>>> alimentos.remove('leche')

>>> alimentos

set(['huevo', 'pera'])

>>> alimentos.add('empanada')

>>> alimentos

set(['empanada', 'huevo', 'pera'])

Prueba de pertenencia

Otra función muy común y útil es probar la pertenencia de objetos al conjunto:

>>> 'empanada' in alimentos

True

>>> 'leche' in alimentos

False

>>> 'leche' not in alimentos

True

Iterar sobre conjuntos

Podemos interar sobre conjuntos de la misma forma que lo hacemos sobre listas:

>>> for a in alimentos:

...     "Debo comer " + a

...

'Debo comer empanada'

'Debo comer huevo'

'Debo comer pera'

Operaciones sobre conjuntos

Los conjuntos en Python soportan las operaciones típicas de conjuntos: restas, intersección, unión y diferencia simétrica (los elementos que están en uno de los conjuntos, pero no en ambos). Repasar operaciones con conjuntos.

>>> frutas = set(['banana', 'naranja', 'pera'])

>>> frutas - alimentos

set(['banana', 'naranja'])

>>> alimentos - frutas

set(['huevo', 'empanada'])

>>> frutas & alimentos

set(['pera'])

>>> frutas | alimentos

set(['huevo', 'empanada', 'pera', 'banana', 'naranja'])

>>> frutas ^ alimentos

set(['huevo', 'empanada', 'banana', 'naranja'])

También podemos preguntar sin un conjunto es subconjunto de otro. En los ejemplos se utiliza set(), el conjunto vacío:

>>> alimentos < set()

False

>>> set() < alimentos

True

>>> set() > alimentos

False

>>> alimentos <= alimentos

True

El problema de las dos comisiones

Este problema está basado en un caso real y lo escuché en la charla Escribí menos código, pensá como un (buen) matemático de Gustavo Carmona (FCEYN - UBA) bio y Matías A Graña (FCEyN - UBA) bio.

Se tienen dos archivos de textos con una lista de e-mails en cada uno. Cada archivo tiene los mails de los funcionarios de una comisión; los archivos no están bien depurados, por lo que pueden contener direcciones repetidas; hay funcionarios trabajando en las dos comisiones. Luego de una reunión en la que trabajaron ambas comisiones, se generó un material que se necesita enviar a todos los participantes. ¿Cómo obtener la lista de destinatarios?

archivo1

dir1@mail.com

dir2@mail.com

dir3@mail.com

dir1@mail.com
archivo2

dir21@mail.com

dir23@mail.com

dir3@mail.com

dir1@mail.com

Queremos una tercera lista que tenga la primera más la segunda, pero que no estén repetidos.

El enfoque tradicional que utilizaría un programador para resolver este problema mediante bucles es:

def unionlarga(l1, l2):

    l3 = []

    for x in l1:

        if not x in l3:

            l3.append(x)

    for x in l2:

        if not x in l3:

            l3.append(x)

    return l3

Supongamos que ya tenemos los archivos leídos y almacenamos el contenidos en listas:

>>> unionlarga(archivo1, archivo2)

['dir1@mail.com', 'dir2@mail.com', 'dir3@mail.com',

'dir21@mail.com', 'dir23@mail.com']

La solución es genérica para cualquier lenguaje. Sin embargo, puede lograr una mejor solución utilizando sets en Python:

>>> list(set(archivo1) | set(archivo2))

['dir21@mail.com', 'dir3@mail.com', 'dir2@mail.com',

'dir1@mail.com', 'dir23@mail.com']

Convertimos ambas listas a conjuntos (con lo que se eliminan los repetidos dentro de las listas), realizamos la unión de ambos conjuntos (con lo que se eliminan los repetidos entre listas y finalmente se convierte el resultado en una nueva lista.

Más sobre conjuntos en la referencia del lenguaje.


All you need is love!

Desde el etorno gráfico de tu GNU/Linux abrí una terminal (gnome-terminal, Konsole, Eterm u otra). Con tu editor preferido (vi, emacs, gedit u otro) creá el archivo love con este contenido:

plot a*(1-sin(t))

a = a - 2

pause 1

if(a == 0) a = 16

reread

Desde la terminal ejecutá el comando gnuplot y allí:

gnuplot> set polar

        dummy variable is t for curves

gnuplot> set xrange[-25:25]; set yrange[-35:5]

gnuplot> a = 16

gnuplot> set title "All you need is love"

gnuplot> load 'love'

Done! All you need is love!.. and gnuplot! :)

Gnuplot es un potente Software Libre que permite graficar funciones matemáticas tanto en 2 como en 3 dimenciones.

a*(1-sin(t)) es una función matemática en coordenadas polares que se denomina cardioide por su forma de corazón.

Este pequeño ejemplo ilustra como utilizar esta herramienta para crear animaciones en base a la parametrización de funciones matemáticas que nos permitan observar comportamientos según cambian los valores de los parámetros.


Fermat vs. Pythagoras mejorado

Logré bajar la complejidad de mi solución de un orden cúbico a un orden cuadrático, la diferencia en cuanto a tiempo de ejecución requerido es sorprendente:

comparacion de 3 soluciones

Pero todavía no e suficiente, con un N igual a 10000 la respuesta se obtiene luego de 10 minutos (en mi computadora de 700 Mhz).

En esta solución ya no uso la función pythagoras(x,y,z). Solo hay 2 ciclos for anidados y el valor de z se calcula como z = sqrt( pow(x,2) + pow(y,2) ). Si z es entero, se tiene un triplete válido y solo falta verificar si x,y,z son primos relativos.

def is_int(i):
return i == int(i)

def solv106f(N):

miembros = []

tripletes_rprim = 0

for x in range(1,N-1):

    for y in range(x+1,N):

        z = sqrt( pow(x,2) + pow(y,2) )

        if (is_int(z) and z &lt;= N):

            for a in [x,y,z]:

                if a not in miembros:

                    miembros.append(a)

            if (rprim(x,y,z)):

                tripletes_rprim +=1



no_miembros = N - len(miembros)

return N, tripletes_rprim, no_miembros


Fermat vs. Pythagoras

Esta es mi solución al problema número 6 de una larga lista de problemas de programación. A pesar de que, como estudiante de Ingeniería, pasé los últimos 3 años y medio estudiando distintas formas de las matemáticas, no conocí el Último Teorema de Fermat hasta que mi amigo Joel, estudiante de Filosofía, me lo comentó.

El problema en realidad no tiene mucho que ver con Fermat, pero sí con Pitágoras: dado un valor N, hay que encontrar cuántos tripletes x,y,z satisfacen x² + y² = z² tales que

    <li>x,y,z sean menores o iguales a N</li>
    
    <li>x,y,z sean primos relativos (es decir que 1 es su único divisor común)</li>
    
    <li>x &lt; y &lt; z</li>
    

    Además hay que encontrar cuántos números mayores que 0 y menores o iguales a N no pertenecen a ninguno de los tripletes (no solo a los tripletes de primos relativos).

    El enunciado original aquí.

    Una posibilidad sería explorar exasutivamente todas las combinaciones de números x,y,z menores o iguales a N, y para cada uno verificar si cumplen las características pedidas. Pero esto sería de altísima complejidad, sería de orden cúbico. Hay que pensar una mejor solución.

    Empecemos con una función que dado 3 valores x,y,z me diga si satisfacen el Teorema de Pitágoras:

    from math import pow
    
    
    
    def pythagoras(x,y,z):
    
        return (pow(x,2)+pow(y,2) == pow(z,2))
    
    
    
    >>> pythagoras(3,4,5)
    
    True
    
    >>> pythagoras(1,1,1)
    
    False
    
    >>> pythagoras()
    
    True
    
    

    También necesito una función, que dado 3 valores x,y,z me diga sin son primos relativos entre sí:

    def mcd(a,b): # algoritmo de euclides
    
        if (a == b): return a
    
        elif (a < b): a,b = b,a # a es el más grande
    
        while(b != 0):
    
            r = a % b
    
            a = b
    
            b = r
    
            return a
    
    
    
    def rprim(x,y,z):
    
        return (mcd(mcd(x,y),z) == 1)
    
    

    La solución que explora todas las combinaciones de 3 números constaría de 3 ciclos for anidados y en el curerpo del más profundo se haría las preguntas pertitenes a la resolución del problema: ¿ese triplete satisface el teorema de pitagoras? ¿son primos relativos? ¿x<y<z?

    for(x=1;x=<N;x++)
    
        for(y=1;y=<N;y++)
    
            for(z=1;z=<N;z++)
    
                # ¿...?
    
    

    Una de las condiciones es que x<y<z, esto permite una redefinición tal que NO se exploren todas las combinaciones de 3 números menores que N sino solo las combinaciones que cumplen con esta restricción:

    for(x=1;x=<N-2;x++)
    
        for(y=x+1;y=<N-1;y++)
    
            for(z=y+1;z=<N;z++)
    
                # ¿..?
    
    

    ¿Esta solución es buena o no es mucho mejor a una en la que todos los ciclos for vayan de 1 a N y por cada uno se pregunte si x,y,z?

    En el archivo pythagoras.py están definidas dos funciones: solv106s(N), basada en la primera aproximación y solv106(N) que espera ser mejor.

    Este programita sirve para medir el tiempo de las dos ejecuciones y compararlos.

    import timing
    
    
    
    def test(N=1000):
    
        timing.start()
    
        solv106s(N)
    
        timing.finish()
    
        print "la primer solución tardo " + str(timing.seconds()) + " segundos"
    
        timing.start()
    
        solv106(N)
    
        timing.finish()
    
        print "la segunda solución tardo " + str(timing.seconds()) + " segundos"
    
    

    Hago una prueba con un valor relativamente grande, N=1000 el valor por defecto:

    >>> test()
    
    la primer solución tardo  1879 segundos
    
    la segunda solución tardo 1561 segundos
    
    

    31 minutos contra 26 minutos.. no me parece una diferencia muy grande. Nota: mi máquina es un AMD K6 II de 700 Mhz.

    Esta gráfica compara las dos soluciones luego de haberlas probado con diferentes valores de N, obtenido unos puntos e interpolar una curva representativa:

    comparacion

    mmm, hay que pensar una mejor solución.


    3n+1

    ACM (Association for Computing Machinery) organiza unas competencias de programación llamadas ICPC (International Collegiate Programming Contest). En esta dirección hay muchos problemas de las competencias: http://acm.uva.es/problemset/.

    Este es el enunciado del primer problema de la guía: http://acm.uva.es/p/v1/100.html

    El fin de semana hice esta solución:

    (nota: los lenguajes soportados por las competencias son C y Java, mi solución está escrita en Python por que me resulta más ágil. De todas formas importa más el algoritmo que el lenguaje.)

    def odd(n):
    
        return not (n%2 == 0)
    
    
    def tnu(n=1):
    
        print n
    
        if (n == 1): return
    
        elif (odd(n)): n = 3*n+1
    
        else: n = n/2
    
        tnu(n) # la recursividad por la cola es evitable
    
    
    >>> tnu(22)
    
    22
    
    11
    
    34
    
    17
    
    52
    
    26
    
    13
    
    40
    
    20
    
    10
    
    5
    
    16
    
    8
    
    4
    
    2
    
    1
    
    

    La longitud del ciclo para un n dado es la cantidad de valores que la función tnu() develve (incluido el 1 final)

    El problema pide leer un archivo con pares i,j y para cada par devolver la mayor longitud de ciclo para todos los enteres enter i y j incluyendo a i y a j.

    Solución Exaustiva

    Una solución exaustiva podría consistir en calcular para cada par i,j todos las longitudes de ciclo asociadas y devolver la mayor.

    Redefino tnu para que en lugar de imprimir los sucesivos valores de n devulva una lista con los valores, el segundo argumento de la función es una lista vacía.

    def tnu(n,l):
    
        #print n
    
        l.append(n)
    
        if (n == 1): return l
    
        elif (odd(n)): n = 3*n+1
    
        else: n = n/2
    
        return tnu(n,l) # la recursividad por la cola es evitable
    
    

    Suponiendo que ya se leyeron los pares del archivo correspondiente y se tienen en una lista llamada pairs:

    def solv1(pairs):
    
        for i,j in pairs:
    
            maxc = 0
    
            for n in range(i,j+1):
    
                maxc = max(len(tnu(n,[])),maxc)
    
            print i, j, maxc
    
    
    >>> pairs = [(1,10),(88,123),(7,1),(55,455),(10000,10003)]
    
    >>> solv1(pairs)
    
    1 10 20
    
    88 123 119
    
    7 1 0
    
    55 455 144
    
    10000 10003 180
    
    >>> pairs = [(1,10),(100,200),(201,210),(900,1000)]
    
    >>> solv1(pairs)
    
    1 10 20
    
    100 200 125
    
    201 210 89
    
    900 1000 174
    
    

    Bien, esta solución ANDA, pero no es óptima. De hecho probablemente sea la solución menos eficiente :-) la de la fuerza bruta.

    Solución Dinámica

    La idea es almacenar los resultados ya calculados para no tener que hacerlos una y otra vez. Si se pide la máxima longitud de ciclo (mlc) enter 100 y 200 y obtengo 125 y luego se pide la máxima longitud de ciclo enter 1 y 300, voy a calcular la mlc para el par 1,100, luego lo haré para 200,300 y luego compararé los resultados con 125 quedándome con el mayor de los 3.

    Pero.. que pasa si ya he calculado la mlc para los valores entre 2 y 99?

    def solv1din(pairs,dic):
    
        for i,j in pairs:
    
            maxc = f_maxc(i,j,dic)
    
            print i, j, maxc
    
            dic[i,j]=maxc
    
    
    def f_maxc(i,j,dic):
    
        wide = 0
    
        for a,b in dic.keys():
    
            if a == i and b == j: return dic[a,b]
    
            if (i<=a and b<=j and b-a > wide):
    
                f,g = a,b
    
                wide = b-a
    
        if (wide == 0): #ningún rango del diccionario estaba dentro de i,j
    
            maxc = 0
    
            for n in range(i,j+1):
    
                maxc = max(len(tnu(n,[])),maxc)
    
            return maxc
    
        else:
    
            return max(f_maxc(i,f,dic),dic[f,g],f_maxc(g,j,dic))
    
    
    >>> solv1din(pairs,{})
    
    1 10 20
    
    100 200 125
    
    201 210 89
    
    900 1000 174
    
    

    Una versión no recursiva de tnu podría ser algo como esto:

    def tnu(n):
    
        l = []
    
        while(n != 1):
    
            l.append(n)
    
            if (odd(n)): n = 3*n+1
    
            else: n = n/2
    
        l.append(1)
    
        return l
    
    

    Download

    • Una implementación completa que resuelve el problema: 3n+1.py
    • Un archivo de entrada de ejemplo: input100.txt

    Ejemplo de uso

    juanjo@sarge:~/problemas$ ./3n+1.py input100.txt
    
    1 10 20
    
    100 200 125
    
    201 210 89
    
    900 1000 174
    
    88 123 119
    
    7 1 0
    
    55 455 144
    
    10000 10003 180
    
    


    Quiero obtener el mínimo cuadrado perfecto mayor a n

    Programando mi generador de constelaciones llegué a este problema:

    Obtener el mínimo cuadrado perfecto mayor a n (con n

    ¿Cómo hubiese sido la forma iterativa de hacerlo?

    
    >>> def min_perf_sqrt(n):
    
    ...    i = s = 1
    
    ...    while (s 
    
    Este es mi one-liner hecho en python con un toque de programación funcional:
    
    
    
    min(filter(lambda i: i > n, map(lambda j: j*j, range(13))))
    
    

    Mmm probablemente la versión iterativa sea más eficiente.. programadores? comentarios?


    gnuplot

    El año pasado escribí una guía sobre este programa para el proyecto Gleducar, en la introducción decía algo que cuando lo leí reien me gustó:

    Los autores iniciales de gnuplot son Thomas Williams y Colin Kelly, quienes decidieron crear un programa que les permitiera visualizar las ecuaciones matemáticas de las clases de electromagnetismo y ecuaciones diferenciales. Ellos eran estudiantes igual que ustedes, se encontraron con un problema (ya que no contaban con una herramienta como esta) y le dieron una solución. Gracias a su generosidad, hoy está disponible para todos.

    Lo encontré en la web con el formato de media wiki (el sistema que usan en Wikipedia) -- muy lindo, gracias a quien lo haya pasado!! les dejo el link por si quieren leerlo http://wiki.gleducar.org.ar/wiki/index.php/Manual_GNUPlot


    tartaglia.py

    La primer vez que escribí un código para el triángulo de tartaglia lo hice en Prolog y se me ocurrió que podría hacerlo en Python para practivar, lo hice de 4 formas distintas: tartaglia.py