lunes, 24 de octubre de 2016

Conjetura de Collatz


Ya se trató esta conjetura en este blog, pero desde el punto de vista de su experimentación en un Taller de Matemáticas

http://hojaynumeros.blogspot.com.es/2011/05/la-conjetura-de-collatz-en-un-taller-de.html

En esta ocasión se buscarán rutinas y funciones que nos ayuden a comprobar la conjetura para muchos números, así como encontrar sus cúspides y órbitas. Se ha escrito mucho sobre esta conjetura, por lo que aquí se desarrollará sólo ese aspecto.

Planteamiento

Para quienes no conozcan esta conjetura recordaremos su planteamiento:

Se toma un número entero positivo N cualquiera, por ejemplo el 13, y se le aplica la siguiente operación, a la que llamaremos función COLL(N):
  • Si el número es par, se divide entre 2.
  • Si el número es impar, se multiplica por 3 y se le suma 1.
En el caso del 13, como es impar, se le aplicará la segunda, y quedará COLL(13)=13*3+1=40.

La idea de la conjetura es que sigamos aplicando esta operación a todos los resultados que obtengamos, En nuestro caso sería COLL(40)=20 (por ser par), COLL(20)=10, COLL(10)=5, COLL(5)=3*5+1=16, COLL(16)=8, COLL(8)=4, COLL(4)=2, COLL(2)=1, y a partir del 1 se entra en el ciclo {4, 2, 1}

La conjetura afirma que este final en el 1 y el ciclo posterior ocurre para cualquier otro entero positivo. Sea cual sea el comienzo, se llegará al número 1. Todas las sucesiones construidas así terminarán en el ciclo 4, 2, 1.

Lo vemos con otro ejemplo:

COLL(6)=3, COLL(3)=10, COLL(10)=5, COLL(5)=16, COLL(16)=8, COLL(8)=4, COLL(4)=2, COLL(2)=1

Insistimos en que existen muchas publicaciones sobre esta conjetura y aquí sólo nos limitaremos a pequeñas comprobaciones. La más sencilla es mediante celdas en la hoja de cálculo.

Comprobación con celdas de hoja de cálculo

Escribimos el entero positivo inicial (lo podemos nombrar como semilla) en una celda, sea por ejemplo, la D4. En la celda inferior D5 escribimos SI(RESIDUO(D4;2)=0;D4/2;3*D4+1), que divide entre 2 el valor de la D4 si es par (porque entonces se verifica RESIDUO(D4;2)=0) y lo multiplica por 3 añadiendo 1 si es impar. En la imagen vemos el resultado para el número 132:



Al ser par el 132 se ha dividido entre 2. Ahora lo único que tenemos que hacer es rellenar esa fórmula hacia abajo y parar cuando aparezca un 1:

(Troceamos la imagen porque aparecen muchos números)

La conjetura ha sido comprobada hasta números muy grandes, por lo que puedes tener la seguridad de que llegarás siempre al valor 1. Al conjunto de números que se recorren hasta llegar a ese valor le podemos llamar órbita del número dado, que aquí son los 29 números que aparecen en la imagen. Al número mayor que hayamos alcanzado en la órbita le llamaremos cúspide. En este ejemplo la cúspide es el mismo 132.

De esta forma tan simple podemos comprobar la conjetura dentro del alcance de la herramienta que usamos. Si la órbita tiene una longitud grande este procedimiento puede alargarse. Por ello acudiremos ahora a la definición de funciones:

Funciones sobre la conjetura de Collatz

Ya hemos presentado COLL(N). Sería bueno introducir su versión en VBA para poder construir sobre ella otras funciones más complicadas. Su código es muy sencillo:

Public Function coll(n)
If n / 2 = n \ 2 Then coll = n / 2 Else coll = 3 * n + 1
End Function

No necesita grandes explicaciones. La condición n/2=n\2 equivale a indicar que n es par, ya que entonces el resultado de la división n/2 es idéntico al de la división entera n\2. El resto se entiende bien. Con esta función podemos reproducir las órbitas en columna de forma idéntica a como procedimos en el primer ejemplo.

En PARI el código es similar:

coll(n)=if(n/2==n\2,n/2,3*n+1)

En la imagen vemos el resultado de pedir coll(132)



Función orbicoll

A cualquier número entero le podemos asignar una cadena que contenga todos los números por los que “pasa” hasta llegar al 1, es decir, su órbita. En VBA podía ser esta:

Public Function orbicoll(n)
Dim b
Dim s$ ‘Cadena (string) para recoger los resultados
b = n: s$ = Str$(b) ‘La cadena comienza con el número inicial (semilla)
While b <> 1 ‘Se trabaja hasta llegar al 1
b = coll(b) ‘ En cada paso se aplica la función COLL
s$ = s$ + Str$(b) ‘Se incorpora el resultado al string
Wend
orbicoll = s$
End Function

En la imagen puedes comprobar la creación de la órbita del número 132:


Después puedes usar la prestación de convertir texto en columnas (antes debes copiar la órbita en otra celda mediante copiar valores) y crear un gráfico que abarque toda la órbita:



Puedes observar que arranca en 132 y termina en 1. Los tramos ascendentes representan números impares, y los descendentes a los pares.
Puedes seguir este proceso con cualquier otro número entero y seguir la evolución de su órbita.

En la imagen aparece la órbita del número 127, más compleja y con una cúspide cercana a 4500:



Función orbicoll en PARI

Es fácil la traducción de esta función a PARI:

coll(n)=if(n/2==n\2,n/2,3*n+1)
orbicoll(n)=my(b=n,s=Str(n));while(b<>1,b=coll(b);s=concat(concat(s," "),Str(b)));s

En la imagen se ha pedido la órbita de 127:



Estudia al código siguiente para la función lorbicoll, que devuelve el número de elementos de una órbita:

Public Function lorbicoll(n)
Dim a, b
a = 1: b = n
While b <> 1
b = coll(b)
a = a + 1
Wend
lorbicoll = a
End Function

Con ella podemos comprobar lo que ya sabemos, que 132 tiene una órbita de 29 elementos.

Con la versión PARI puedes abordar casos con números mayores:

coll(n)=if(n/2==n\2,n/2,3*n+1)
lorbicoll(n)=my(b=n,a=1);while(b<>1,b=coll(b);a+=1);a

Te proponemos comprobar que el número  871 es el que posee la órbita de más longitud entre los de tres cifras, y que contiene 179 elementos. De los de cuatro cifras el de órbita de más longitud es el número 6171, con 262 elementos.

Función cuspicoll

Del mismo modo que construimos la órbita de un número entero positivo, podemos encontrar su cúspide. El procedimiento será similar, pero, en lugar de añadir resultados a un string, tomaremos nota en cada paso del máximo valor que ha aparecido:

Public Function cuspicoll(n)
Dim a, b
a = n: b = n
While b <> 1
b = coll(b)
If b > a Then a = b
Wend
cuspicoll = a
End Function

La instrucción clave es If b > a Then a = b, que convierte a en el nuevo máximo si aparece un elemento b mayor que los precedentes.

Con las funciones definidas podemos construir un esquema en el que se analice el comportamiento de la conjetura de Collatz para una semilla dada:


Finales previsibles

Algunos tipos de números presentan una órbita bastante previsible. Por ejemplo:

Potencias de 2: a partir de ellos se entra en una ruta descendente y previsible que finaliza en el 1.

Números tipo (2^n-1)/3: Desembocan en una potencia de 2, por lo que también inician una ruta directa,  y esto ocurre una potencia sí y otra no, porque 2^n es del tipo 3k+1 o 3k+2 y al multiplicar por Si es 3k+1 es candidato a que (2^n-1)/3 sea entero, y si es del tipo 3k+2, la siguiente potencia será 6k+4, o sea del tipo 3k+1

Otros tienen recorrido corto, como el 6, el 10 o el 20.

Es normal que pensemos en que muchas órbitas pasarán por ellos, y existan pares de órbitas  que pasan ambas por el mismo punto de entrada, más o menos primario.

Podemos intentar ver si dos números presentan alguna coincidencia en sus órbitas, porque entonces compartirán final. No es difícil programar una función que nos devuelva un punto de coincidencia en las órbitas de dos números. En primer lugar necesitamos una función que nos indique si el número n pertenece a la órbita del número m. Puede ser esta:

Public Function enlacoll(m, n) As Boolean
Dim e As Boolean
Dim p
If m = n Then
e = True
Else
e = False
p = m  ‘p recorrerá la órbita de m
While Not e And p <> 1
p = coll(p)
If p = n Then e = True  ‘Si p es igual a n, sí pertenece
Wend
End If
enlacoll = e
End Function

Con esta función enlacoll(m,n) podemos saber si n pertenece a la órbita de m. Se puede organizar un esquema de cálculo:


En la imagen se ha verificado que 52 pertenece a la órbita de 57.

Coincidencia en las órbitas

Con la anterior función ya estamos preparados para encontrar la primera coincidencia entre dos órbitas. Si no existe una coincidencia anterior, se nos devolverá un cero. El código de la función es:

Public Function coincicoll(m, n)
Dim c, q
q = m
c = 0
While q > 1 And c = 0 ‘q recorre la órbita de m
If enlacoll(n, q) Then c = q ‘si q pertenece a la órbita de n, hay coincidencia
q = coll(q)
Wend
coincicoll = c
End Function

Con ella también podemos construir otro esquema de cálculo. En la imagen se comprueba que 125 y 126 comparten una subórbita de 95 elementos, que comienzan en 364.


Con estas ideas puedes construir otras muchas funciones, o emprender otras búsquedas. Sólo se ha pretendido en esta entrada dar ideas para comprobaciones y experimentaciones sobre la conjetura, ya de por sí bastante estudiada en otros aspectos.

jueves, 13 de octubre de 2016

Hipotenuseando


En mis exploraciones por la página OEIS (Enciclopedia On-line de sucesiones de números enteros, http://oeis.org/?language=spanish), me encontré con la sucesión http://oeis.org/A104863

10, 30, 31, 43, 53, 68, 86, 109, 138, 175, 222, 282, 358, 455, 578, 735, 935, 1189,…

En ella, a partir de los valores a(1)=10 y a(2)=30, se van formando los siguientes como la parte entera de la raíz cuadrada de la suma de cuadrados de los dos anteriores.

Así, el tercer elemento es igual a

ENTERO(RAIZ(10^2+30^2))=ENTERO(RAIZ(1000))=ENTERO(31,62)=31.

Prueba a justificar que el siguiente es 43.

Esta definición equivale a que cada término es la hipotenusa truncada correspondiente a los términos anteriores tomados como catetos. Por eso le hemos llamado a esta sucesión la de “hipotenusear”. Su expresión recurrente sería:


Esta sucesión presenta algunas características que le hacen merecer esta entrada:

(1) El uso de la parte entera obliga a renunciar a las fórmulas teóricas. De hecho, en términos avanzados de la sucesión, el cumplimiento del Teorema de Pitágoras es deplorable. Observa este ejemplo, en el que a(n-2)=8348, a(n-1)= 10618 y a(n)=13506. Si elevamos al cuadrado obtenemos:

13506^2=182412036
8348^2+10618^2=182431028

Restamos y nos queda un error de 18992 unidades. Así que las hipotenusas que obtengamos con esta recurrencia no son tales, aunque seguiremos llamándolas así.

(2) Como también ocurre en las recurrencias lineales, el cociente entre dos términos consecutivos se va estabilizando y tiende a un límite. Esto nos permite “razonar en el límite”, aunque sepamos que es una técnica aproximada, que sólo nos valdrá para explicar (y no demostrar) algunas conjeturas que se han afirmado para esta sucesión y otras similares.

(3) La sucesión presentada es sólo un caso particular de toda una familia en la que podemos fijar a(1) y a(2) como deseemos. La mayoría de las propiedades se mantendrán. Vemos la primera:

Conjetura: El límite del cociente a(n+1)/a(n) es la raíz cuadrada del número áureo.
La llamamos conjetura por causa de la parte entera, que nos impide mejores razonamientos. Esta cuestión, de manejarnos entre aproximaciones y conjeturas, es uno de los objetivos de esta entrada.

Podemos comprobar lo anterior con hoja de cálculo. Escribimos dos catetos uno debajo de otro, como 2 y 5, y después, en columna rellenamos la fórmula

=ENTERO(RAIZ(CATETO1^2+CATETO2^2)).

No es difícil de organizar. Después, en la columna de la derecha calculamos los cocientes entre dos términos consecutivos. Algo así:



Al llegar al término 14 ya se adivina el valor deseado. Si seguimos bajando, la aproximación mejora mucho



Podemos razonarlo en el límite. Llamamos k al cociente a(n+1)/a(n). Por tanto, en la expresión de a(n) podemos escribir:


O bien, pasando a(n-1) al primer miembro,


Elevando al cuadrado y agrupando, tenemos que k se debe aproximar a la solución de la ecuación k4 – k2 – 1 = 0, una bicuadrada cuya solución es el límite sugerido, la raíz cuadrada del número áureo.

Incluimos las cuatro soluciones tal como las da WolframAlpha:



Elegimos la real positiva, y, efectivamente, resulta 1,27201964951407…

Manteniendo el razonamiento en el límite, si a(n-1) y a(n) se comportan como cateto e hipotenusa respectivamente con esa razón dada, el otro cateto, a(n-2), se podrá aproximar (también en el límite) de esta forma:

Plantéate como ejercicio demostrar el último paso. Recuerda que F-1=1/F
En esta sucesión a(n) tiende en el límite a a(n-2)*F

Lo hemos demostrado en el párrafo anterior. También lo podemos razonar mediante la idea de que si el cociente entre dos términos consecutivos se aproxima a la raíz del número áureo, el correspondiente a a(n) y a(n-2) será dicho número F.

Por tanto, en el límite, cada tres términos consecutivos forman un triángulo rectángulo cuyos lados son proporcionales a (1, 1,272019…, 1,618033…) y cuyo ángulo menor es de 38,17º.

Un ejercicio: ¿Cuál es, en el límite, el cociente entre el área del triángulo (a(n+1), a(n), a(n-1)) y el correspondiente a (a(n), a(n-1), a(n-2))?

Para quienes conozcáis el lenguaje PARI, con una línea de código similar a esta podéis estudiar la sucesión hasta términos más avanzados:

a=1;b=7;for(i=1,30,c=truncate(sqrt(a^2+b^2));a=b;b=c;print1(c,", "))

Podéis estudiar los cocientes añadiendo el código adecuado.

Conjetura: A partir de un término mínimo, a(n) se diferencia de a(n-2)+a(n-4) en a lo sumo una unidad.

Esta conjetura está publicada en la página OEIS citada para el caso a(1)=10 y a(2)=30, en el que la diferencia se estabiliza en 1. Su verificación no depende de los términos iniciales, salvo, quizás, el tope inferior de 1. Por ejemplo, lo comprobaremos con hoja de cálculo y los términos iniciales a(1)=4 y a(2)=7:



En este caso vemos que a(n) tiende a coincidir con a(n-2)+a(n-4)-1

En el límite se puede justificar usando todos los cocientes presentados más arriba:

a(n-2)+a(n-4) = a(n)/F+a(n)/F^2 = a(n)*(F+1)/F^2 = a(n)

Así que en el límite la coincidencia es exacta: a(n-2)+a(n-4) = a(n), y la unidad como error aparece por los truncamientos.

Puedes cambiar la función ENTERO por la de REDONDEAR. Así lo hacen las sucesiones A104803,  A104804, A104805 y A104806, con resultados similares.

lunes, 3 de octubre de 2016

¿Variaciones o combinaciones?


En el mes de julio pasado descubrí que el número 1716 equivale a un número de variaciones sin repetición y también de combinaciones sin repetición, ambas con el mismo índice superior. En efecto, 1716 = 13*12*11 = V(13,3), pero también equivale a C(13,6) o C(13,7), ya que


Se ha producido la feliz casualidad de que 13*12*11 = 6*5*4*3*2*1, y por eso se ha podido simplificar.

En esta propiedad el verdadero protagonista, a efectos de construcción de algoritmos, es el número 13, que es el que participa en ambas fórmulas, de combinaciones y variaciones. No existen muchos índices que cumplan esto. Los primeros son 8, 13, 27, 124, 725 y 5046, si imponemos la condición razonable de que el índice inferior sea mayor que 1, para evitar trivialidades.

Búsqueda “ingenua”

Para encontrar estos índices superiores podemos acudir a la definición que hemos insinuado: “Números n para los que existen dos índices k y h tales que V(n,k) = C(n,h)”. Si disponemos de las funciones C y V, basta recorrer índices para cada candidato y parar cuando se dé una coincidencia. Podemos acotar la búsqueda eligiendo para h el intervalo (2, n/2), por cuestión de simetría en los números combinatorios. Por otra parte, es claro que k ha de ser menor que h, para que tenga lugar la igualdad. En Basic de hojas de cálculo podía usarse algo así:

Public Function vari(n, k)
Dim v, i
v = 1
For i = 0 To k - 1: v = v * (n - i): Next i
vari = v
End Function

Public Function combi(n, k)
Dim v, w, i
v = 1: w = 1
For i = 0 To k - 1: v = v * (n - i): w = w * (i + 1): Next i
combi = v / w
End Function

Código para cada valor de n (aquí representado por la variable i):
a = 0 ‘variable para parar la búsqueda
k = 2
While k <= i / 2 And a = 0 ‘ se recorren los valores de k
b = combi(i, k) ‘ se encuentra el número de combinaciones b
h = 2
While h < k And a = 0 ‘se recorren los valores de h
c = vari(i, h) ‘ se encuentra el número de variaciones c
If b = c Then a = 1: m = k: n = h ‘en caso de igualdad, se para y toma nota
h = h + 1
Wend
k = k + 1
Wend

La salida será el valor m del índice k y el n del índice h. En forma de tabla estos serían los primeros valores:



Los comprobamos (salvo el 13 que ya se ha visto):

V(8,2)=8*7 = 56, C(8,3)=(8*7*6)/(3*2*1)=56

V(27,3)=27*26*25=17550, C(27,4)=(27*26*25*24)/(4*3*2*1) =27*26*25=17550

V(124,4)=124*123*122*121=225150024 = C(124,5)

Algoritmo con recursividad

En este primer intento estamos realizando más operaciones de lo debido. No es necesario calcular C(n,k) y V(n,h) en cada paso. Es mejor generar cada intento recursivamente a partir del anterior:

a = 0
k = 2
b = i
While k <= i / 2 And a = 0
b = b * (i - k + 1) / k ‘recursividad para combinaciones
h = 2
c = i
While h < k And a = 0
c = c * (i - h + 1) ‘recursividad para variaciones
If b = c Then a = 1: m = k: n = h
h = h + 1
Wend
k = k + 1
Wend

Con hoja de cálculo se llega pronto al desbordamiento de decimales. Deberemos cambiar a PARI:

for(i=3,1000,a=0;j=2;m=i;while(j<=i/2+1&&a==0,m=m*(i-j+1)/j;k=2;n=i;while(k<j&&a==0&&n<=m,n*=i-k+1;if(m==n,a=1);k+=1);j+=1);if(a==1,print(i,", ",j,", ",k,", ",m)))

Es poco legible. Contiene las mismas ideas desarrolladas con VBA, pero escritas de forma excesivamente compacta.

El resultado te será familiar, y aparece un nuevo índice, el 725:



Para seguir avanzando se requiere ya mucha paciencia, porque los cálculos se van haciendo lentos y complejos. El siguiente índice superior en aparecer es el 5046:



El valor de V(5046,7) = C(5046,8) da idea de cómo se va complicando esto. Sin embargo, nos conduce al hecho de que en

V(5046,7)=5046*5045*5044*5043*5042*5041*5040,

el último factor es el factorial de 7, lo que permite la simplificación que da lugar a la igualdad entre variaciones y combinaciones.

Casos particulares

El último ejemplo nos da una idea de la naturaleza de algunos de los índices superiores con la propiedad buscada. Es fácil entender que todo índice del tipo n!+(n-1) da lugar a un número m en el que el número de variaciones de n-1 elementos coincide con el de combinaciones de n elementos. Así ha ocurrido con muchas de las soluciones presentadas. En la siguiente tabla se han destacado en rojo las que ya conocíamos. Nos hemos detenido en el último factorial que Excel puede expresar de forma entera:


Por tanto, el número de índices adecuados es infinito, y crece a ritmo de factorial.
Los casos que faltan, como el 13, provienen de la casualidad de que un producto de números consecutivos equivalga a un factorial, que es lo que ocurre con 10*9*8 = 6! ¿Existirán más casos? Mediante una búsqueda manual descubrimos: 6*5*4=5!, lo que nos da de nuevo el candidato 8, pero esta vez con la expresión C(8,5) y la solución 56 ya vista. De ella podemos extraer la coincidencia 10*9*8*7=7!, que nos llevaría al índice 13.

Este estudio es un ejemplo más de una forma clásica de abordar problemas:

  •  Usar un algoritmo sencillo, que no hace uso de propiedades especiales. Es un buen método para comenzar, pero suele ser largo y poco interesante. Así ha sido nuestro procedimiento “ingenuo”.
  •  Perfeccionar el algoritmo a fin de conseguir mayor velocidad de búsqueda. En este caso se ha logrado con la recursividad.
  •  Acudir a la teoría o el razonamiento. Hemos descubierto así que existen infinitos casos con la fórmula n!+n-1, más unos cuantos casos aislados. Así le hemos quitado el misterio a la cuestión planteada.