1 Introducción a la programación con Sage

1 Introducción a la programación con Sage

1.1 Tipos de datos en Sage

1.1.1 Introducción a Sage

Página de Sage.
Formas de trabajar con Sage: consola, emacs y navegador.
Interacción con el navegador: SageMath
- Cuadros de texto: admite textos en HTML y Mardown
  - se inician con %html o %md
  - se procesan con mayúscula+Enter.
- Cuadros de código:
  - se inician pulsando en la línea azul
  - se evalúa con mayúscula+Enter

1.1.2 Operaciones

Se puede usar Sage como calculadora con las siguientes:
- operaciones aritméticas: +, -, *, /, // (cociente entero), % (resto) y ^ (potencia).
- operaciones booleanas: and, or, not
- relaciones: < , <= , > , >= , <> , != , ==
La precedencia y asociatividad son las habituales.
Las funciones se llaman como es habitual. Por ejemplo, max(2,3), sin(pi/2).

1.1.3 Variables

Definición de variables: La forma de asignar valores a las variables es

variable = expresión

Presentación de valores: El valor de las variables se muestra con la función print. Por ejemplo,

print variable

Los nombres de las variables deben empezar por una letra o un guión bajo (_), pueden contener números, y son distintos si se usan mayúsculas o minúsculas.
Liberación de variables: Con función del se libera una variable, y a partir de ese punto el nombre de la variable deja de estar deﬁnido.

del variable

1.1.4 Tipos de datos

Ejemplos de tipos de datos:
- Booleanos: True o False .
- Enteros: por ejemplo, 1, 10, -30.
- Racionales: por ejemplo, 1/2, -3/4
- Números de coma ﬂotante: por ejemplo, 1.25, -1.5e6
- Expresiones simbólicas: por ejemplo, pi/4, (1+sqrt(2))/2
Sage usa tipado dinámico: todos los datos tienen necesariamente un tipo.
n(x) es una aproximación del valor de x. Por ejemplo,

sage:  n(sqrt(2))
1.41421356237310

Métodos específicos de cada tipo: Cada tipo de datos tiene sus propios métodos y se llaman escribiendo primero la variable que contiene el dato, después un punto (.), y después el método. Por ejemplo,

sage: n = 12
sage: n.factor()
2^2 * 3

Se puede completar el nombre del método usando el tabulador.
Los números de Python se definen con int o float.

1.1.5 Secuencias de datos

1.1.5.1 Cadenas de caracteres

Las cadenas de caracteres se escriben entre comillas simple, dobles o triples.
x + y es la concatenación de las cadenas x e y. Por ejemplo.

sage: print ’Hola’ + ’tu’
Holatu

x * n es la cadena obtenida repitiendo n veces la cadena x. Por ejemplo,

sage: print ’Hola’*2
HolaHola

len(x) es la longitud de la cadena x. Por ejemplo,

sage: print len(cadena1), len(cadena2), len(cadena3), len(cadena4)
4  8  12  15

x[j] es el j-ésimo elemento (empezando en 0) de la cadena x. Por ejemplo,

sage: cadena = ’Si miras al abismo...’
sage: print cadena[0], cadena[1], cadena[2], cadena[3]
S i   m

x[j:k] es la subcadena de x desde el índice j (inclusive) hasta el k (exclusive). Por ejemplo,

sage: cadena = ’Si miras al abismo...’
sage: print cadena[3:10]
miras a

1.1.5.2 Tuplas

Las tuplas se crean escribiendo los elementos entre paréntesis y separados por comas.
x[j] es el j-ésimo elemento (empezando en 0) de la tupla x. Por ejemplo,

sage: ej_tupla = (1.5,'a',7,'d')
sage: ej_tupla[2]
7

x[j:k] es la subtupla de x desde el índice j (inclusive) hasta el k (exclusive). Por ejemplo,

sage: ej_tupla = (1.5,'a',7,'d',9,6)
sage: ej_tupla[1:3]
('a', 7)

x + y es la concatenación de las tuplas x e y. Por ejemplo.

sage: (1,3) + (2,5,7)
(1, 3, 2, 5, 7)

x * n es la tupla obtenida repitiendo n veces la tupla x. Por ejemplo,

sage: (5,2)*3
(5, 2, 5, 2, 5, 2)

Las tuplas se pueden desempaquetar con patrones. Por ejemplo,

sage: a,b,c = (2,5,3)
sage: (c,b,a)
(3, 5, 2)

Las tuplas son immutables, es decir, que no se puede quitar ni añadir elementos a una tupla una vez ha sido creada, ni siquiera sustituir un elemento por otro.

1.1.5.3 Listas

Las listas se crean escribiendo los elementos entre corchetes y separados por comas.
x[j] es el j-ésimo elemento (empezando en 0) de la lista x. Por ejemplo,

sage: ej_lista = [1.5,'a',7,'d']
sage: ej_lista[2]
7

x[j:k] es la sublista de x desde el índice j (inclusive) hasta el k (exclusive). Por ejemplo,

sage: ej_lista = [1.5,'a',7,'d',9,6]
sage: ej_lista[1:3]
['a', 7]

x + y es la concatenación de las tuplas x e y. Por ejemplo.

sage: [1,3] + [2,5,7]
[1, 3, 2, 5, 7]

x * n es la tupla obtenida repitiendo n veces la tupla x. Por ejemplo,

sage: [5,2]*3
[5, 2, 5, 2, 5, 2]

Se puede modificar el valor de un elemento de una lista. Por ejemplo,

sage: ej_lista = [1,'a',7,'d',9,6]
sage: ej_lista[2] = 8
sage: ej_lista
[1, 'a', 8, 'd', 9, 6]

del x[j] elimina el elemento j-ésimo de la lista x. Por ejemplo,

sage: ej_lista = [1,'a',7,'d',9,6]
sage: del ej_lista[2]
sage: ej_lista
[1, 'a', 'd', 9, 6]

xs.insert(j,y) inserta en la posición j de la lista xs el elemento y. Por ejemplo,

sage: ej_lista = [1,'a',7,'d',9,6]
sage: ej_lista.insert(2,'e')
sage: ej_lista
[1, 'a', 'e', 7, 'd', 9, 6]

Las listas son contenedores dinámicos de datos: se pueden quitar y añadir elementos en cualquier momento.
srange(j,k,d) es la lista de los números entre j (inclusive) y k (exclusive), pero contando de d en d elementos. Por ejemplo,

sage: srange(3,9,2)
[3, 5, 7]

srange(k) es equivalente a srange(0,k,1). Por ejemplo,

sage: srange(9)
[0, 1, 2, 3, 4, 5, 6, 7, 8]

srange(j,k) es equivalente a srange(j,k,1). Por ejemplo,

sage: srange(3,9)
[3, 4, 5, 6, 7, 8]

1.1.5.4 Pertenencia

x in ys se verifica si x pertenece al contenedor (cadena, tupla o lista) ys. Por ejemplo,

sage: 'a' in "Hoy"
False
sage: 'o' in "Hoy"
True
sage: 7 in (3,2,5)
False
sage: 7 in (3,7,5)
True
sage: 7 in [3,2,5]
False
sage: 7 in [3,7,5]
True

1.1.5.5 Conversiones entre secuencias

tuple(x) es la tupla correspondiente a la cadena o lista x. Por ejemplo,

sage: tuple("abc")
('a', 'b', 'c')
sage: tuple([1,2,3])
(1, 2, 3)

list(x) es la lista correspondiente a la cadena o tupla x. Por ejemplo,

sage: list("abc")
['a', 'b', 'c']
sage: list((1,2,3))
[1, 2, 3]

c.join(xs) es la cadena obtenida uniendo los elementos de xs con la cadena c entre ellos. Por ejemplo,

sage: ej_lista = ['a','b','c']
sage: ''.join(ej_lista)
'abc'
sage: '-'.join(ej_lista)
'a-b-c'

1.1.6 Mostrar información por pantalla

print var1, var2, ... escribe los valores de las variables var1, var2, .... Por ejemplo,

sage: a = 3
sage: b = 7
sage: print a, b
3 7

print cadena valores escribe la cadena sustituyendo los códigos por sus valores. Los códigos usuales son
- %d: número entero
- %f: número de coma ﬂotante, con decimales
- %.3f: número de coma ﬂotante, con 3 decimales exactamente
- %s: cadena de caracteres (o cualquier dato que no sea entero ni de coma ﬂotante)

sage: x = "Esto"
sage: print "%s es una cadena de %d elementos" %(x,len(x))
Esto es una cadena de 4 elementos
sage: a=2 ; b=7
sage: print "La media de %d y %d es %f" %(a,b,(a+b)/2)
La media de 2 y 7 es 4.500000
sage: print "La media de %d y %d es %.1f" %(a,b,(a+b)/2)
La media de 2 y 7 es 4.5

1.2 Control del flujo del programa

1.2.1 Alternativas y bucles

1.2.1.1 Alternativas lógicas

La sintaxis de los boques if ... elseif ... else es

if condicion1:
  bloque de instrucciones 1
elif condicion2:
  bloque de instrucciones 2
...
else:
  bloque de instrucciones final

Ejemplos de if

sage: x = 6
sage: if x%2 == 0:
....:       print '%d es par' %x
....: 
6 es par
sage: x = 7
sage: if x%2 == 0:
....:       print '%d es par' %x
....: else:
....:       print '%d es impar' %x
....: 
7 es impar

La condición puede ser cualquier expresión que al evaluar devuelva un booleano, y no se pueden hacer asignaciones.
La condición puede ser un número, en cuyo caso 0 se interpreta como False , y cualquier otra cosa como True.
La condición puede ser un contenedor (lista, tupla, cadena...), en cuyo caso un contenedor vacío se interpreta como False, y cualquier otra cosa como True.

sage: xs = [5]
sage: if xs:
....:       print "la lista es no vacía"
....: else:
....:       print "la lista es vacía"
....: 
la lista es no vacía
sage: xs = []
sage: if xs:
....:       print "la lista es no vacía"
....: else:
....:       print "la lista es vacía"
....: 
la lista es vacía

1.2.1.2 Bucle for

La sintaxis del bucle for es

for elemento in lista:
    instrucciones

Ejemplos de for

sage: suma = 0
sage: for x in [3,2,5]:
....:   suma = suma + x
....: 
sage: suma
10

1.2.1.3 Bucles while

La sintaxis de los bucles while es

while condicion:
    instrucciones

Ejemplo de bucle while

sage: i = 0
sage: while i < 7:
....:   i = i+1
....: 
sage: i
7

1.2.2 Bloques anidados

Podemos anidar bloques if...else, for y while de cualquier forma posible.
Ejemplo if dentro de for para calcular la suma de los inversos de los primos menores que k.

sage: k = 100
sage: suma = 0
sage: for j in range(k):
....:   if is_prime(j):
....:      suma = suma + 1/j
....: 
sage: print 'La suma de los inversos de los primos menores que  %d es  %.3f' %(k, suma)
La suma de los inversos de los primos menores que  100 es  1.803

Ejemplo de for dentro de for para escribir un triángulo numérico

sage: for j in srange(10):
....:   for k in srange(j):
....:     print k, 
....:   print
....: 
0
0 1
0 1 2
0 1 2 3
0 1 2 3 4
0 1 2 3 4 5
0 1 2 3 4 5 6
0 1 2 3 4 5 6 7
0 1 2 3 4 5 6 7 8

1.2.2.1 Tablas de verdad

Ejemplo de la tabla de verdad de $(p \wedge q) \vee (p \wedge \neg r)$

sage: for p in [True,False]:
....:    for q in [True,False]:
....:       for r in [True,False]:
....:          print 'p = %s, q = %s, r = %s => e = %s' %(p,q,r,e)
....: 
p = True, q = True, r = True => e = e
p = True, q = True, r = False => e = e
p = True, q = False, r = True => e = e
p = True, q = False, r = False => e = e
p = False, q = True, r = True => e = e
p = False, q = True, r = False => e = e
p = False, q = False, r = True => e = e
p = False, q = False, r = False => e = e

1.2.2.2 Buscar un número con ciertas propiedades

Ejemplo: Buscar un número N primo relativo a un número k dado y tal que k/3 ≤ N ≤ 2k/3.
Búsqueda por fuerza bruta (saliendo del bucle con break):

sage: k = 196
sage: tercio = ceil(k/3)
sage: dosTercios = floor(2*k/3)
sage: for x in range(tercio, dosTercios):
....:    if gcd(x,k) == 1:
....:       break
....: 
sage: print x
67

Búsqueda aleatoria usando randint (tal que randint(a,b) es un número aleatorio entre a y b)

sage: k = 196
sage: tercio = ceil(k/3)
sage: dosTercios = floor(2*k/3)
sage: x = randint(tercio,dosTercios)
sage: while gcd(x,k) != 1:
....:    x = randint(tercio,dosTercios)
....: 
sage: print x
69

1.2.3 Deﬁnición de funciones

La sintaxis para deﬁnir una función en Python es:

def funcion(argumento1,argumento2,...):
   documentación
   instrucciones

Ejemplo de definición: Definir el procedimiento informacion tal que informacion(n) escriba información del número n, indicando si es primo, par o potencia de dos. Por ejemplo,

sage: informacion(2)
Información sobre 2
    es primo
    es par
    es potencia de dos
sage: informacion(7)
Información sobre 7
    es primo
sage: informacion(8)
Información sobre 8
    es par
    es potencia de dos

Solución (Informacion_del_numero.sage):

def informacion(n):
    print 'Información sobre %d'%n
    if is_prime(n):
        print '    es primo'
    if n % 2 == 0:
        print '    es par'
    if is_power_of_two(n):
        print '    es potencia de dos'

Las funciones tienen la posibilidad de devolver un valor, usando la palabra clave return.
Ejemplo de función con valor: Definir la función radio tal que radio(x,y) es el radio del círculo con centro en el origen y que pasa por el punto (x,y). Por ejemplo,

sage: radio(3,4)
5
sage: r1 = radio(1,1); r2 = radio(5,0)
sage: r2-r1
-sqrt(2) + 5

Solución (Radio.sage):

def radio(x,y):
    return sqrt(x^2 + y^2)

1.2.3.1 Documentación

Se le puede añadir documentación a las definiciones. Por ejemplo,

def radio(x,y):
    '''el radio del círculo con centro en el origen y que pasa por el
    punto (x,y).'''
    return sqrt(x^2 + y^2)

La documentación se puede consultar terminando el nombre de la función con ?.

1.2.3.2 ¿Cuándo usar una función?

Como regla general, siempre que, resolviendo un problema, podamos identiﬁcar una subtarea claramente delimitada, debemos escribirla como función.
Ejemplo: Definir la función sumaPrimos tal que sumaPrimos(n) es la suma de los inversos de los números primos menores que n. Por ejemplo,

sage: sumaPrimos(3)
1/2
sage: n(sumaPrimos(3))
0.500000000000000
sage: sumaPrimos(50)
1021729465586766997/614889782588491410
sage: n(sumaPrimos(50))
1.66164651701580

Solución (Suma_inversos_primos.sage):

def sumaPrimos(k):
  suma = 0
  for j in srange(k):
    if is_prime(j):
      suma = suma + 1/j
  return suma

# 2ª definición (usando prime_range)
# ==================================

def sumaPrimos2(k):
  suma = 0
  for j in prime_range(k):
      suma = suma + 1/j
  return suma

1.2.3.3 Tratar con argumentos erróneos

Si algún fragmento de código produce un error, lanza una excepción.
Si una función espera datos de cierto tipo y recibe datos de otro tipo distinto, lanza un TypeError. Por ejemplo

sage: range('a')
TypeError: range() integer end argument expected, got str.

Si realiza una división por cero lanza un ZeroDivisionError. Por ejemplo,

sage: 1/0
ZeroDivisionError: Rational division by zero

Los errores se propagan desde la línea que produjo el error a través de la pila de ejecución. Por ejemplo,

def f(a):
    return 1/a

def g(a,b):
    return a+f(b)

# Ejemplo de la propagación del error
#    sage: g(1,0)
#    .../codigo/<string> in g(a, b)
#    .../codigo/<string> in f(a)
#    ZeroDivisionError: Rational division by zero

Si en algún momento sentimos la necesidad de asegurarnos de que los datos de entrada de nuestras funciones veriﬁcan ciertas restricciones, utilizaremos la sintaxis raise Exception.
Ejemplo, definir la función factorial tal que factorial(n) es el factorial de n, si n es un entero positivo y da un error si n es negativo. Por ejemplo,

sage: factorial(3)
6
sage: factorial(-3)
Exception: El argumento de factorial debe ser positivo

Solución:

def factorial(n):
    if n <= 0:
        raise Exception, 'El argumento de factorial debe ser positivo'
    acumulador = 1
    k = n
    while k != 0:
        acumulador *= k
        k -= 1
    return acumulador

1.3 Un poco de programación funcional

1.3.1 Funciones recursivas

Ejercicio. [Factorial recursivo.] Definir, por recursión, la función factorial tal que factorial(n) es el factorial de n. Por ejemplo,

sage: factorial(5)
120
sage: factorial(50)
30414093201713378043612608166064768844377641568960512000000000000

Solución:

def factorial(n):
    if n == 0:
        return 1
    else:
        return n*factorial(n-1)

Nota. Si en la definición falta el caso base, termina después de 1.000 llamadas recursivas. Por ejemplo, con la definición

def factorial(n):
    return n*factorial(n-1)

Se tiene

sage: factorial(5)
RuntimeError: maximum recursion depth exceeded

1.3.2 Funciones de orden superior

Las funciones son ciudadanos de primera clase.
Ejemplo de funciones como argumento:

sage: funciones = [is_prime, is_prime_power, is_power_of_two, is_square]
sage: k = 5
sage: for f in funciones:
....:     print f(k)
....: 
True
True
False
False

Ejemplo de función de orden superior: [Iteración hasta que se cumpla una condición] Definir la función hasta tal que hasta(f,p,x) itera la función f, comenzando por el valor x, hasta que se verifique la condición p y devuelve el elemento final de la iteración. Por ejemplo,

sage: hasta(lambda x: x*2, lambda x: x > 1000, 1)
1024

Solución:

def hasta(f,p,x):
    t = x
    while not p(t):
        t = f(t)
    return t

Ejemplo de uso de hasta [El algoritmo de Herón de cálculo de raíces cuadradas con hasta] Definir, usando hasta, la función raiz tal que raiz(a,e) es la raíz cuadrada de a, con un error menor que e calculada mediante al algoritmo de Herón. Por ejemplo,

sage: raiz(49,0.1)
7.00140647524394

Solución:

# 1ª definición (con funciones auxiliares):
def raiz(a,e):
    def aceptable(t): 
        return abs(t**2-a) < e
    def mejora(t):  
        return (t+a/t)/2
    return n(hasta(mejora,aceptable,a))

# 2ª definición (con funciones anónimas):
def raiz2(a,e):
    return n(hasta(lambda t: (t+a/t)/2,
                   lambda t: abs(t**2-a) < e,
                   a))

Ejemplo de aplicación de hasta. [Cálculo de ceros por el método de bisección]
El método de la bisección para calcular un cero de una función en el intervalo [a,b] se basa en el teorema de Bolzano:

"Si f(x) es una función continua en el intervalo [a, b], y si, además, en los extremos del intervalo la función f(x) toma valores de signo opuesto (f(a) * f(b) < 0), entonces existe al menos un valor c en (a, b) para el que f(c) = 0".

La idea es tomar el punto medio del intervalo c = (a+b)/2 y considerar los siguientes casos:

Si f(c) es aproximadamente igual a 0, hemos encontrado una aproximación del punto que anula f en el intervalo con un error aceptable.
Si f(c) tiene signo distinto de f(a), repetir el proceso en el intervalo [a,c].
Si no, repetir el proceso en el intervalo [c,b].

Definir la función ceroBiseccion tal que ceroBiseccion(f,a,b,e) es una aproximación del punto del intervalo [a,b] en el que se anula la función f, con un error menor que e, aplicando el método de la bisección (se supone que f(a)*f(b)<0). Por ejemplo,

sage: ceroBiseccion(cos,0,2,0.0001)
1.57080078125000
sage: ceroBiseccion(lambda x: x^2-49,0,10,0.1)
7.00195312500000

Solución:

load Hasta.sage

# 1ª solución (por recursión):
def ceroBiseccion(f,a,b,e):
    def aceptable(x): 
        return abs(f(x)) < e 
    def aux(c,d): 
        m = (c+d)/2
        if aceptable(m):    
            return m
        elif f(c) * f(m)  < 0:  
            return aux(c,m)
        else: 
            return aux(m,d)
    return n(aux(a,b))


# 2ª solución (con until):
def ceroBiseccion2(f,a,b,e): 
    def aceptable(t):
        (x,y) = t
        return abs(f(x)) < e
    def mejora(t): 
        (x,y) = t
        m = (x+y)/2
        if f(x) * f(m) < 0:  
            return (x,m)
        else:
            return (m,y)
    return n(hasta(mejora,aceptable,(a,b))[0])

Ejemplo [Reducción (o plegado) de listas]
Definir la función reducir tal que reducir(f,xs) es el resultado de aplicar la operación f entre los elementos de la lista xs; es decir, si xs es [x1,x2,...,xn], entoncesreducir(f,xs)esf(...(f(f(x1,x2) ,x3),...,xn))`. Por ejemplo,

sage: def suma (x,y):
....:    return x+y
....: 
sage: reducir(suma,[2,5,3])
10
sage: reduce(lambda x,y: x+y, [2,5,3])
10
sage: reducir(operator.add,[2,5,3])
10
sage: reducir(operator.mul,[2,5,3])
30
sage: reducir(operator.concat,[[3,5],[2],[4,6,7]])
[3, 5, 2, 4, 6, 7]

Solución:

def reducir(f,xs):
    r = xs[0]
    for x in xs[1:]:
        r = f(r,x)
    return r

En Python existe la función reduce que es equivalente a reducir.

1.4 Secuencias y flujos de datos

1.4.1 Transformar los elementos de una lista

Una lista de comprensión es una lista de la forma [f(x) for x in xs]. Por ejemplo,

sage: [x^2 for x in range(10)]
[0, 1, 4, 9, 16, 25, 36, 49, 64, 81]
sage: [is_prime(x) for x in range(10)]
[False, False, True, True, False, True, False, True, False, False]
sage: [len(x) for x in ['En', 'un', 'lugar', 'de', 'la', 'Mancha']]
[2, 2, 5, 2, 2, 6]
sage: [x.upper() for x in ['En', 'un', 'lugar', 'de', 'la', 'Mancha']]
['EN', 'UN', 'LUGAR', 'DE', 'LA', 'MANCHA']

1.4.2 Filtrar una lista

Lista de comprensión con condición: La lista [f(x) for x in xs if p(x)] es la lista obtenida aplicando f a los elementos de xs que cumplen la propiedad p. Por ejemplo,

sage: [x for x in srange(10) if is_prime(x)]
[2, 3, 5, 7]
sage: [x for x in ['En', 'un', 'lugar', 'de', 'la', 'Mancha'] if len(x)==2]
['En', 'un', 'de', 'la']

1.4.3 Todo a la vez

Se puede usar varios generadores y combinar transformaciones y filtrados. Por ejemplo,

sage: [[x,y] for x in [2,3,5,6] if is_even(x) for y in [4,7,6,9] if is_odd(y)]
[[2, 7], [2, 9], [6, 7], [6, 9]]

1.5 Funciones que trabajan sobre listas

Funciones:
- max(lista) es el máximo de la lista.
- min(lista) es el mínimo de la lista
- sum(lista) es la suma de los elementos de la lista.
- all(lista) se verifica si todos los elementos de la lista son ciertos.
- any(lista) se verifica si algún elemento de la lista es cierto.
Ejemplo 1: ¿Cuántos divisores de 21.000 son cuadrados perfectos?

sage: len([x for x in divisors(21000) if is_square(x)])
4

Ejemplo 2: ¿Hay algún número k entre 20 y 100 tal que $2^k+1$ es primo?

sage: any([is_prime(x) for x in srange(20,101)])
True
sage: [x for x in srange(20,101) if is_prime(x)]
[23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97]

Ejempl 3: Encuentra el tercer cuadrado perfecto n=k^2 tal que n+1 es primo.

sage: [n^2 for n in srange(1,1000) if is_prime(n+1)][2]
16

1.6 Evaluación perezosa

Los generadores son objetos que generan los elementos de la lista uno por uno, según se vayan solicitando desde otras partes del programa.
A esta técnica se le denomina evaluación perezosa , porque no se hacen cálculos hasta que no son necesarios.
xsrange es la versión perezosa de srange. Una vez creado el generador, podemos pedir los números uno por uno, usando el método next(), o recorrerlos todos, usando un bucle for.

sage: xlista = xsrange(100)
sage: print xlista.next()
0
sage: print xlista.next()
1
sage: acumulador = 0
sage: for j in xsrange(101):
....:     acumulador += j
....: 
sage: acumulador
5050

Definición de generadores: (f(x) for x in xs if p(x)). Por ejemplo,

sage: cuadrados = (x^2 for x in xsrange(10))
sage: sum(cuadrados)
285

Ejemplos de ganancia de eficiencia mediante evaluación perezosa.: En los 4 siguientes ejercicios se muestra la ganancia.
Ejercicio 1. Definir la función primo tal que primo(n,x,y) es el n-ésimo primo entre x e y. Por ejemplo,

sage: primo(2,10,100)
17

Solución:

def primo(n,x,y):
    ls = srange(x,y)
    primos = [t for t in ls if is_prime(t)]
    return primos[n]

Ejercicio 2. Definir la función tiempo tal que tiempo(f,n,x,y) es el tiempo para calcular f(n,x,y). Por ejemplo,

sage: tiempo(primo,2,10,100)
0.015407085418701172

Solución:

def tiempo(f,n,x,y):
    inicio = time.time()
    f(n,x,y)
    final = time.time()
    return final - inicio

Ejercicio 3. Definir, mediante evaluación perezosa, la función primo tal que primo(n,x,y) es el n-ésimo primo entre x e y. Por ejemplo,

sage: primoP(2,10,100)
17

Solución:

# 1ª definición
def primoP(n,x,y):
    ls = xsrange(x,y)
    primos = (t for t in ls if is_prime(t))
    for k in srange(n):
        primos.next()
    return primos.next()

# 2ª definición
def primoP2(n,x,y):
    k = n
    for j in xsrange(x,y):
        if is_prime(j):
            k -= 1
            if k == -1:
                break
    return j

Ejercicio 4. Comparar el tiempo necesario para calcular el segundo primo entre 10.000 y 100.000 con las distintas definiciones.
Solución: La comparación es

sage: tiempo(primo,2,10000,100000)
8.730193138122559

sage: tiempo(primoP,2,10000,100000)
0.00787806510925293

sage: tiempo(primoP2,2,10000,100000)
0.006257057189941406

1.7 Combinación de generadores y flujos

Nota En los siguientes ejercicios se vuelve a mostrar la ganancia de eficiencia de la evaluación perezosa junto con la simplicidad de los flujos.
Ejercicio 1. Definir la función tiene_especiales tal que tiene_especiales(x,y) se verifica si existe algún k entre x e y tal que 2*k+1 es primo.

def tiene_especiales(x,y):
    return any([is_prime(2*k+1) for k in srange(x,y)])

Ejercicio 2. Definir la función tiempo tal que tiempo(f,x,y) es el tiempo para calcular f(x,y).

def tiempo(f,x,y):
    inicio = time.time()
    f(x,y)
    final = time.time()
    return final - inicio

Ejercicio 3. Definir, con evaluación perezosa, la función tiene_especiales tal que tiene_especiales(x,y) se verifica si existe algún k entre x e y tal que 2*k+1 es primo.

def tiene_especialesP(x,y):
    return any(is_prime(2*k+1) for k in xsrange(x,y))

Ejercicio 4. Comparar el tiempo necesario para comprobar si existe algún k entre 10.000 y 100.00 tal que 2*k+1 es primo, con las distintas definiciones.

sage: tiempo(tiene_especiales,10000,100000)
9.206285953521729

sage: tiempo(tiene_especialesP,10000,100000)
0.0020439624786376953

1.8 Ejercicios

Ejercicio [Raices de la ecuación de segundo grado] Definir la función raices de forma que raices(a,b,c) devuelve la lista de las raices reales de la ecuación $ax^2 + bx + c = 0$ . Por ejemplo,

sage: raices(1,3,2)
[-1, -2]
sage: raices(1,-2,1)
[1]
sage: raices(1,2,3)
Exception: No tiene raices reales

Solución:

def raices(a,b,c):
    d = b^2-4*a*c
    e = sqrt(d) 
    if d > 0:
        return [(-b+e)/(2*a), (-b-e)/(2*a)]
    elif d == 0: 
        return [(-b+e)/(2*a)]
    else:
        raise Exception, "No tiene raices reales"

Ejercicio. [Restos de dividir los números de una lista por uno dado]. Definir, usando un bucle for, el procedimiento restos tal que restos(xs,n) escribe los restos de los elementos de xs entre n. Por ejemplo,

sage: restos([1,3,7,13,21,31,43,57,73,91],5)
1 3 2 3 1 1 3 2 3 1

Solución:

def restos(xs,n):
    for x in xs:
        print x % n,

Ejercicio. [Elementos de una lista que son múltiplos de uno dado]. Definir, usando un bucle for con un bloque if anidado, el procedimiento multiplos tal que multiplos(xs,n) escribe los elementos de xs que son múltiplos de n. Por ejemplo,

sage: multiplos([1,3,7,13,21,31,43,57,73,91],3)
3 21 57

Solución:

def multiplos(xs,n):
    for x in xs:
        if x % n == 0:
            print x,

Ejercicio. [Divisores de los números de una lista]. Definir, usando un bloque for con otro bloque for anidado, el procedimiento divisores tal que divisores(xs) escribe cada número de n junto a sus divisores. Por ejemplo,

sage: divisores([1,3,7,13,21,31,43,57,73,91])
Divisores de 1: 1
Divisores de 3: 1 3
Divisores de 7: 1 7
Divisores de 13: 1 13
Divisores de 21: 1 3 7 21
Divisores de 31: 1 31
Divisores de 43: 1 43
Divisores de 57: 1 3 19 57
Divisores de 73: 1 73
Divisores de 91: 1 7 13 91

Solución:

def divisores(xs):
    for x in xs:
        print "Divisores de %d:"%x,
        for k in srange(1,x+1):
            if x % k == 0:
                print k,
        print

Ejercicio. [Suma de inversos hasta superar la cota]. Definir, usando un bucle while, el procedimiento sumaInversos tal que sumaInversos(n) suma los inversos de los números naturales $1/k$ hasta que la suma sea mayor que n y escribe la suma y el número $k$ en el que la suma se detiene. Por ejemplo,

sage: sumaInversos(1)
La suma es 1.500000 y el número de sumandos es 2
sage: sumaInversos(2)
La suma es 2.083333 y el número de sumandos es 4
sage: sumaInversos(3)
La suma es 3.019877 y el número de sumandos es 11
sage: sumaInversos(10)
La suma es 10.000043 y el número de sumandos es 12367

Solución

def sumaInversos(n):
    k = 0
    suma = 0
    while suma <= n:
        k = 1 + k
        suma = suma + 1/k
    print "La suma es %f y el número de sumandos es %d"%(suma,k)

Ejercicio. [El algoritmo de Herón de cálculo de raíces cuadradas.] El algoritmo de Herón de cálculo de la raíz cuadrada de un número real positivo $a$ es comienza con una aproximación burda,
$t_0 = a$
y después en cada iteración del algoritmo la sustituye por otra aproximación mejor
$t_{n+1} = \frac{1}{2}\left(t_n + \frac{a}{t_n}\right)$
La sucesión $t_n$ converge a la raíz cuadrada de $a$ .

Definir la función raiz tal que raiz(a,e) es la raíz cuadrada de n, con un error menor que e calculada mediante al algoritmo de Herón. Por ejemplo.

sage: raiz(49,0.1)
7.00140647524394

Solución:

def raiz(a,e):
    t = a
    while (abs(t**2-a) >= e):
        t = (t+a/t)/2
    return n(t)

Ejercicio. [Órbitas de Collatz] El sucesor de Collatz de un número n es n/2, si n es par y 3n+1, en caso contrario. La órbita de un número se obtiene escribiendo los sucesivos sucesores de Collatz hasta llegar a 1. Por ejemplo, la órbita de Collatz de 7 es 7, 22, 11, 34, 17, 52, 26, 13, 40, 20, 10, 5, 16, 8, 4, 2, 1.

Definir la función collatz tal que collatz(n) es la órbita de Collatz de n. Por ejemplo,

sage: collatz(7)
[7, 22, 11, 34, 17, 52, 26, 13, 40, 20, 10, 5, 16, 8, 4, 2, 1]

Solución:

def collatz(n):
    s = [n]
    while (n != 1):
        if is_even(n):
            n = n//2
        else:
            n = 3*n+1
        s = s + [n]
    return s

Ejercicio. Definir la función longitud_collatz tal que longitud_collatz(n) es la longitud de la órbita de Collatz de n. Por ejemplo,

longitud_collatz(7)  ==  17

Solución:

# 1ª definición (usando collatz)
def longitud_collatz(n):
    return len(collatz(n))

# 2ª definición (sin usar collatz)
def longitud_collatz_2(n):
    x = 1
    while (n != 1):
        if is_even(n):
            n = n//2
        else:
            n = 3*n+1
        x = x + 1
    return x

Ejercicio. Definir la función collatz_mayor tal que collatz_mayor(n) es el menor número cuya órbita de Collatz tiene más de n elementos. Por ejemplo,

collatz_mayor(16)  ==  7
collatz_mayor(50)  ==  27

Solución:

def collatz_mayor(n):
    x = 1
    while (longitud_collatz_2(x) <= n):
        x = x+1
    return x

Introducción a la programación con Sage