jueves, 19 de enero de 2017

Número de descomposiciones en diferencia de cuadrados


Después de jugar bastante con los números naturales, he observado la disparidad existente entre ellos respecto al número de sus descomposiciones en diferencias de cuadrados. Nos referimos al número de soluciones de

N=a^2-b^2 con a y b enteros y a>0 y b>=0 para un N dado.

Hay números que no admiten ninguna descomposición de este tipo, como el 22, y otros que admiten muchas. Un ejemplo es el 1680, que admite doce:

1680=421^2-419^2=212^2-208^2=143^2-137^2=109^2-101^2=89^2-79^2=
=76^2-64^2=67^2-53^2=52^2-32^2=47^2-23^2=44^2-16^2=43^2-13^2=41^2-1^2

En la tabla lo puedes comprobar:



¿De qué depende esto? Lo iremos viendo más adelante.

Obtención del número de descomposiciones

Al igual que hemos procedido con otros temas, comenzaremos encontrando las soluciones sin apoyo de la teoría, para más adelante fundamentarlo y después extraer propiedades y curiosidades.

Cualquier algoritmo para resolver esta cuestión se puede basar en el hecho de que una descomposición de este tipo equivale a expresar N mediante un producto de factores con la misma paridad, ambos pares o ambos impares. En efecto, si N=b^2-a^2 resultaría N=(a+b)(a-b), y ambos factores tienen la misma paridad, como se comprueba estudiando los cuatro casos posibles para a y b: par-par, par-impar, impar-par e impar-impar. Es fácil. A la inversa: si N=m*n, ambos de la misma paridad resultaría que (m+n)/2 y (m-n)/2 serían enteros, cumpliéndose que N=((m+n)/2)^2-((m-n)/2)^2

El número de descomposiciones de N en diferencia de cuadrados coincide con el de productos con factores de la misma paridad y resultado N.

Para hacer más fáciles los trabajos, admitiremos que b pueda ser nulo, o de forma equivalente, que los dos factores sean iguales.

Establecida esta propiedad, la búsqueda se efectuaría encontrando divisores d de N tales que d y N/d tengan la misma paridad. Esto lo podemos lograr fácilmente, usando divisiones enteras y aritmética modular. En el siguiente ejemplo lo implementamos como función Basic para hojas de cálculo:

Public Function numdifcuad(n)
Dim i, m
m = 0 ‘Contador de casos
For i = 1 To Sqr(n) ‘Se llega a la raíz cuadrada para evitar repeticiones de divisores
If n / i = n \ i Then ‘Si el cociente es igual al cociente entero, es que es divisor
If Abs(i - n / i) Mod 2 = 0 Then m = m + 1 ‘el divisor i y el cociente n/i tienen la misma paridad, porque su diferencia da resto 0 módulo 2, luego incrementamos el contador.
End If
Next i
numdifcuad = m
End Function

Así de sencillo. Con esta función podemos contar las descomposiciones posibles para cada número. En la tabla puedes observar las correspondientes a los primeros números:



Verás que se dan muchos casos, desde el 22 o el 14, que no admiten descomposiciones, hasta el 16 o el 15, que admiten 2. Si siguiéramos leyendo hacia abajo descubriríamos que 45 es el primero que admite 3 casos: 45=23^2-22^2=9^2-6^2=7^2-2^2.que se corresponden con las descomposiciones en producto 45=45*1=15*3=9*5, todas con factores de igual paridad. El primero con cinco casos es el 96, y así podríamos seguir hasta 1680 que vimos presenta doce.



Estos resultados están ya publicados en https://oeis.org/A034178



Podemos comprobar que coinciden con nuestros resultados

Análisis teórico de la situación

Es importante distinguir en principio si N es par o impar.

Caso 1: N es impar

Si N es impar, todos los posibles pares de factores N=m*n tendrán la misma paridad, luego sólo tenemos que contarlos. Recordamos que la función TAU, que cuenta los divisores, viene dada por la fórmula


En ella a, b c,…representan los factores primos y p, q, r,…sus exponentes. Como los divisores han de formar pares, deberemos encontrar la mitad de la función (si esta es par), por tanto, si expresamos el número de descomposiciones mediante la función ND, tendremos:


Lo comprobamos. Según la tabla el 21 admite dos descomposiciones. Como 21=3*7, ?(21)=(1+1)(1+1)=4, y su mitad entera es dos, como indica la tabla.

Si TAU es impar, será porque todos los exponentes serán pares, con lo que N será un cuadrado, apareciendo entonces un par nuevo al multiplicar la raíz cuadrada por sí misma. Podemos unificar ambas situaciones:


El segundo paréntesis valdrá cero en el caso par y uno en el impar.

Es el caso del número 3969:



Los factores son todos impares, los exponentes, 4 y 2, pares. TAU valdrá en este caso (1+4)(1+2)=15. Su mitad entera es 7, y añadimos 1 por su raíz cuadrada. Coincide entonces con el valor 8 que nos da NUMDIFCUAD.

Caso 2: N es par

En este caso aparece el factor 2, lo que obliga a que los dos factores que buscamos sean ambos pares y deban contener el 2 como factor. Esto no es posible si el factor 2 es único, con exponente 1, y esa es la causa de que 6, 10, 14 o 22 no presenten descomposiciones.

No admiten descomposiciones en diferencias de cuadrados los múltiplos de 2 que no lo sean de 4, es decir, los del tipo 4k+2

N es potencia de 2

Siguiendo un razonamiento similar al caso impar, deberemos encontrar la función TAU. Como tratamos con un solo exponente, sea p, el número de divisores será 1+p, pero al ser el caso par, el divisor 1 no nos interesa, y nos quedarían tan solo p divisores. Así, para p=5 dispondríamos de 2, 4, 8, 16 y 32. La potencia completa, en este caso 32, no nos valdría, porque tendría que formar par con el 1, que es impar, luego sólo nos quedan p-1 divisores disponibles. Esto vuelve a confirmar que el caso p=1 no produce pares de la misma paridad.

Los pares engendrados por el 2 serán pues, (p-1)/2 si p es impar y Int((p-1)/2)+1 si es par, por el par que se gana por su raíz cuadrada. Unificando:

N no es potencia de 2

Si el número es potencia de 2, sin factores impares, esta expresión vale, pero, en caso contrario, estas posibilidades del factor 2 se deberán multiplicar por las correspondientes al factor impar. La complicación surge del hecho de cada par puede producir productos idénticos, que se contarían dos veces, y hay que tener en cuenta los pares de factores repetidos en el caso de que p sea par. Por ello, la única forma de encajar todo es:


Multiplicamos el número de pares con factores desiguales de la potencia de 2 contenida en N por todos los factores de la parte impar de N, y después, en otro producto, los factores con repetición se multiplican sólo por los pares que se forman en la parte impar. Algo complicado, pero funciona.

Hemos plasmado los tres casos en una única función, a la que llamaremos NUMDIFCUAD2:

Public Function numdifcuad2(n)
Dim p, q, r, s, t, m

m = n: p = 0
While m Mod 2 = 0: m = m / 2: p = p + 1: Wend ‘Extraemos la potencia de 2
If p = 1 Then numdifcuad2 = 0: Exit Function ‘Caso imposible
q = n / 2 ^ p ‘q es la parte impar
If q = 1 And p > 1 Then numdifcuad2 = Int((p - 1) / 2) + (p - 1) Mod 2
‘Es potencia de 2 pura
If p = 0 And q > 1 Then t = fsigma(q, 0): numdifcuad2 = Int(t / 2) + t Mod 2
‘Es un número impar
If p > 1 And q > 1 Then t = fsigma(q, 0): numdifcuad2 = t * Int((p - 1) / 2) + ((p - 1) Mod 2) * (Int(t / 2) + t Mod 2)
‘Tiene parte par y parte impar
End Function

La complejidad del cálculo nos ha aconsejado comprobar mediante tablas si las dos versiones de NUMDIFCUAD coinciden. Lo hemos probado desde 1 hasta 100000, sin que aparecieran discrepancias. Como ejemplo, adjuntamos los valores de algunos números de seis cifras entre los que hay impares, pares y una potencia de 2 pura:



Otra interpretación

Como todo cuadrado es suma de números impares consecutivos, como por ejemplo 16=1+3+5+7, al restar dos cuadrados se eliminarán sumandos impares, quedando tan sólo aquellos que no sean comunes. Es, por ejemplo, el caso de 44=12^2-10^2=21+23 o el de 72=9^2-3^2=7+9+11+13+15+17

Así que el número de descomposiciones que estamos estudiando coincide con el de formas de expresar el número como suma de impares.





lunes, 9 de enero de 2017

Sandwich de semiprimos


Unos comentarios de James Tanton (@jamestanton) en Twiter me han hecho interesarme por aquellos números tales que tanto su anterior como su posterior son semiprimos. No los recorreremos todos, porque son muchos, sino que los clasificaremos por tipos.

Todos los números que poseen esta propiedad han de ser pares, salvo el 5, porque si fueran impares, los semiprimos adyacentes deberían ser dobles de primos, que serían números consecutivos, lo que salvo el caso de 2 y 3 es imposible.

Consecuencia inmediata es que un número primo mayor que 5 no puede estar encerrado entre dos semiprimos.

Los comentarios citados arriba se referían a los cuadrados, y estos ya están publicados en http://oeis.org/A108278

Cuadrados “sandwicheados”

En realidad, en esa página figuran las bases, pero elevando al cuadrado nos resultarán los cuadrados pedidos:

144, 900, 1764, 3600, 10404, 11664, 39204, 97344, 213444, 272484, 360000, 656100, 685584, 1040400, 1102500, 1127844, 1633284, 2108304, 2214144,…

En todos ellos el anterior y el posterior son semiprimos. Tomemos el cuadrado 213444=462^2. Su anterior 213443=461*463  y el posterior 213445=5*42689, ambos semiprimos. El primero nos lleva a una situación interesante, y es que si el cuadrado central es n^2, el anterior será n^2-1=(n-1)(n+1), y al ser semiprimo el producto ambos factores serán primos y más aún, primos gemelos. Es lo que ha ocurrido con el 144.

Si un cuadrado está rodeado por dos semiprimos, el anterior es producto de dos primos gemelos.

Respecto a los factores de n^2+1, han de ser del tipo 4k+1, según estudiamos hace tiempo. Puedes seguir el razonamiento en el apartado dedicado a “Un cuadrado y una unidad” en el documento

http://www.hojamat.es/publicaciones/hojanum1.pdf

Al deber ser pares, estos cuadrados serán todos múltiplos de 4.

Si deseas reproducirlos con PARI, este puede ser el código:

for(i=2,2000,n=i*i;if(bigomega(n-1)==2&&bigomega(n+1)==2,print1(n,"; ")))

Triangulares entre semiprimos

También los triangulares pueden estar comprendidos entre dos semiprimos. Los primeros están publicados en http://oeis.org/A121898

120, 300, 528, 780, 2628, 3240, 3828, 5460, 13530, 18528, 19110, 22578, 25878, 31878, 32640, 37128, 49770, 56280, 64980, 72390, …

En PARI:

for(i=2,2000,n=i*(i+1)/2;if(bigomega(n-1)==2&&bigomega(n+1)==2,print1(n,"; ")))

Además de ser pares, son múltiplos de 3, y por tanto de 6. En efecto, si un triangular no es múltiplo de 3 sólo puede ser porque su orden sea del tipo 3k+1, ya que entonces el triangular tendría la expresión (3k+1)(3k+2)/2, que no contiene ningún factor 3 (en los otros casos sí). Pero en este caso, al restarle 1 no obtendríamos un semiprimo. Lo desarrollamos:


Al tener el factor 9 no puede ser semiprimo. Además hemos descubierto que es nueve veces el triangular de orden k.

Por ejemplo, el número triangular de orden 13, 91=13*14/2, no es múltiplo de 3, y su anterior, 90, no puede ser semiprimo, y es igual a 9*10=9*T(4)

El desarrollo anterior se puede invertir, es decir, que si multiplicamos por 9 un triangular y sumamos 1, obtenemos otro triangular no múltiplo de 3 o 6.
Sólo los números triangulares N múltiplos de 6 pueden tener semiprimos N-1 anteriores a ellos.

Oblongos entre semiprimos

¿Ocurrirá algo similar con los oblongos? La respuesta es negativa.

Recordemos que los números oblongos son los dobles de los triangulares, es decir, los que se pueden expresar como N=k(k+1). Así, 56=7*8 es oblongo, y su anterior 55=5*11 y el posterior 57=3*19 son semiprimos. Cumple la condición de estar entre semiprimos, pero no es múltiplo de 3 (par sí tiene que ser).

Los primeros oblongos con la propiedad requerida son:

56, 552, 870, 1056, 1190, 1640, 1892, 2652, 4032, 5256, 5402, 6806, 8372, 9120, 9506, 9702, 10920, 11772, 12656, 12882, 15006, 15252, 15500, 16256, 16770, 17556, 18632, 23256, 24492, 27722, 29070, 30800, 33306, 33672, 34410, 36290, 40200, 40602, 44310, 45582, 46872, 49506,…

Esta sucesión estaba inédita y la hemos publicado en https://oeis.org/A276565

Tomamos uno de ellos, por ejemplo 16256=127*128, y por tanto, oblongo. Su anterior es semiprimo, ya que 16255=5*3251, y 16257=3*5419, también lo es.

El posterior no puede ser múltiplo de 5, porque los oblongos terminan todos en 0, 2 o 6, y al sumar no obtendremos ni 5 ni 0 como última cifra.

El anterior no puede ser múltiplo de 3. Si lo es el oblongo, es claro que al restar 1 deja de serlo. Si no lo es, sería del tipo

(3k+1)(3k+2)-1 = 9k^2+9k+1 y tampoco.

Si deseas obtener más términos, puedes adaptar este código en PARI:

for(i=2,2000,n=i*(i+1);if(bigomega(n-1)==2&&bigomega(n+1)==2,print1(n,"; ")))

Números de Fibonacci

Están publicados los números de la sucesión de Fibonacci comprendidos entre semiprimos. Sólo hay cuatro con pocas cifras: 5, 34, 144, 46368.  Se conjetura que no hay infinitos.

Puedes estudiarlos en http://oeis.org/A167023

Cubos perfectos

Los cubos rodeados de semiprimos son muy escasos. El primero es 216=6^3, con 215=5*43 y 217=7*31.

Los siguientes llegan a ser casi inabordables: 216, 1302170688, 7211429568, 20346417000, 71887512312, 281268868608, 1394417360448, 17571944311992, 28350304855488, 170400029184000, 450335804625000, 504966851923968, 616121259098688, 1064394808685208, 3442267015299000, 3517494650695368, 3540860163178632, …

Es preferible tratar con sus bases. Las tienes publicadas en http://oeis.org/A268043

6, 1092, 1932, 2730, 4158, 6552, 11172, 25998, 30492, 55440, 76650, 79632, 85092, 102102, 150990, 152082, 152418, 166782, 211218,…

Estos números tienen una propiedad importante, y es que su anterior y posterior han de formar un par de primos gemelos. La idea es sencilla: si n^3-1 ha de ser semiprimo, al ser múltiplo de n-1, este ha de ser primo, pues en caso contrario el otro no sería semiprimo. Igual ocurre con n+1. Por ejemplo, 166782 está rodeado por los primos gemelos 166781 y 166783.

Los puedes encontrar con PARI:

for(i=2,2000,n=i^3;if(bigomega(n-1)==2&&bigomega(n+1)==2,print1(i,"; ")))

Potencias enteras

Hemos estudiado los cuadrados y cubos entre semiprimos, pero podríamos generalizar a todas las potencias de base y exponente enteros mayores que 1.
No es muy difícil encontrarlos si se dispone de una función ESPOTENCIA o similar. En nuestro equipo disponemos de ella, y hemos podido emprender la búsqueda, consiguiendo esta lista de los primeros:

144, 216, 900, 1764, 2048, 3600, 10404, 11664, 39204, 97344, 213444, 248832, 272484, 360000, 656100, 685584, 1040400, 1102500, 1127844, 1633284, 2108304, 2214144, 3504384, 3802500, 4112784, 4536900, 4588164, 5475600, 7784100, 7851204, 8388608, 8820900, 9000000, 9734400…

Los hemos publicado en https://oeis.org/A276564

Puedes reproducirla con PARI:

for(i=2,10^7,if(ispower(i)&&bigomega(i-)==2&&bigomega(i+1)==2,print1(i,", ")))

Con base prima hay muy pocos. Los primeros son 2048 y 8388608

Otros casos

Podríamos seguir el estudio con pentagonales (los primeros serían 5, 92, 590, 1080, 1820, 8400,…) o hexagonales (120, 780, 3828, 19110,…), pero por hoy ya está bien. Lo dejamos como propuesta.

miércoles, 4 de enero de 2017

Bienvenida al 2017



Para aquellos de nuestros lectores que no nos sigan en las redes sociales, resumimos los cálculos sobre el 2017 publicados en Twitter (o de próxima aparición)

Poco se puede decir del próximo año 2017, salvo que es primo, pero los siguientes años serán semiprimos: 2018=2×1009 y 2019=3×673.
Entre todos los desarrollos del 2017 que he estudiado, me quedo con este, por su sencillez y doble simetría:
2017=(47-4)×47-4

No es fácil engendrar un año con sólo dos capicúas:
2017=686+1331

Me gusta la expresión del 2017 como suma de cubos:
2017=7^3+11^3+7^3
que da lugar a esta triple simetría:
2017=343+1331+343

El nuevo año también admite una suma simétrica de tres cuadrados:
2017 =18^2+37^2+18^2

Todo número es suma de tres números triangulares. El nuevo año también:
2017=28×29/2+29×30/2+48×49/2
2017=31×32/2+38×39/2+39×40/2

Curiosa expresión del 2017 como suma de potencias con bases crecientes y exponentes decrecientes:
2017=1^5+2^4+4^3+44^2

El año 2017 “pitagorea” con cuadrado y sin él:
2017^2=792^2+1855^2
2017=44^2+9^2

Las cifras de los números importantes dan la bienvenida al 2017. Hoy con π:
2017=314×1×5+(92+6)×5-35-8

Hoy reciben al nuevo año las cifras del número “e”:
2017=2718-((2+81)×(8+2)-84-5×9)

El número áureo también aporta sus cifras al nuevo año:
2017=1618+0+339+88-7-4-9-8

Las cifras de la raíz de dos se unen a la fiesta del 2017:
2017=1414+213+562-(3+7)×(3+0)×9+50+48

También contamos en el 2017 con las cifras de la raíz cuadrada de tres:
2017=1732+0+50×8+0-7-5-6!+8+87×7

Para el nuevo año, cuadrados de cuatro en cuatro:
2017=14^2+19^2+26^2+28^2=16^2+20^2+20^2+31^2=17^2+24^2+24^2+24^2=18^2+21^2+24^2+26^2

Les toca el turno a los desarrollos monocifra:
2017=(1111-111)×(1+1)+((1+1)^(1+1))^(1+1)+1    
2017=2222-222+(2+2)^2+2/2     
2017= (333+3)×(3+3)+3/3            
2017= (444+44+4)×4+44+4+4/4 
2017= (555+-55+5)×(5-5/5)-(5+5+5)/5     
2017=666×(6+6+6)/6+6+6+6+6/6             
2017= (77+77-7-7)×(7+7)+7×7+7+7/7      
2017=88×(8+8+8-8/8)-8+8/8      
2017=999+999+9+9+9/9              

El año 17 engendrado por el 17:
2017=17×17×1×7+1-7

No podía faltar un pandigital para el nuevo año:
2017=(5+80)×4×6-9-1-7-3×2       

El 2017 subiendo del 1 al 11                       
2017=1×(2+3)×(4+5)×67-8-9×10×11         

Con tres potencias de dos se engendra todo un año:        
2017=2^11-2^5+2^0      

Hacemos bajar al 2017 tres anchos escalones                    
2017=3×(3+332+2)×2-2-1×(1+1+1)           

Las cifras de los ocho primeros primos despiden estos desarrollos:
2017=(2+3)×5×(7+1)×11-(3+171+9)          
Hasta el 2018.

lunes, 19 de diciembre de 2016

Numeradores de los números armónicos


Un número armónico Hn es un número racional formado mediante la suma de los inversos de los números naturales hasta n, es decir

La sucesión de los números armónicos tiene límite infinito, por ser divergente la serie armónica a la que pertenecen los sumandos. Es fácil ver que los denominadores de los números armónicos son, salvo simplificaciones, los primeros factoriales, pero aquí nos van a interesar los numeradores.

Con hoja de cálculo no es difícil encontrarlos, pues situamos en columna los distintos valores hasta n y en otra columna vamos dividiendo n! entre 1, 2, 3,… para después sumar la columna de cocientes (que serán todos enteros). Aquí tienes el desarrollo para n=6



Se ha dividido el factorial de 6, 720, entre 1,2,…6, sumando los cocientes hasta conseguir 1764, luego el sexto número armónico tendrá la expresión 1764/720.

Nos basamos en este cálculo en la operación de sacar el denominador común.

Si nos atenemos a la definición primera, vemos que H(n+1)=H(n)+1/(n+1). Si  imaginamos mentalmente el resultado de esa suma entenderemos que el numerador de H(n) quedará multiplicado por n+1 con la suma posterior del primitivo denominador, es decir, n!. Esto nos lleva a una fórmula de recurrencia fácil para los numeradores: a(n+1)=a(n)*(n+1)+n!

Si aplicamos la fórmula de recurrencia al 1, obtenemos la sucesión 1, 3, 11, 50, 274, 1764, 13068, 109584,…que coincide con los números de Stirling de primera clase, publicados en http://oeis.org/A000254 y que puedes estudiar allí. Nosotros simplificaremos las fracciones de los números armónicos, con lo que los valores anteriores se verán sustituidos por

 1, 3, 11, 25, 137, 49, 363, 761, 7129, 7381, 83711, 86021, 1145993,...

Estos son los numeradores que consideraremos en principio. Están también publicados (http://oeis.org/A001008)

Según las ideas contenidas en los párrafos anteriores, la primera suma que consideramos y la operación de simplificar, el algoritmo adecuado para generar estos numeradores simplificados puede ser el siguiente, que lo hemos orientado a una función de VBA:

Public Function numearmonico(n)
Dim i, v, f
v = 0
f = 1
For i = 1 To n: f = f * i: Next i ‘Se construye el factorial
For i = 1 To n
v = v + f / i ‘Se suman los cocientes del factorial con 1..n
Next i
v = v / mcd(v, f) ‘Simplificación de la fracción
numearmonico = v
End Function

Con esta función es fácil reproducir de nuevo los valores de los numeradores:



Todo lo anterior tenía por objeto que practicáramos un poco con cálculos y algoritmos, pero si hubiéramos acudido al lenguaje PARI, la generación de los numeradores es tan sencilla que se limita a su definición:

for(n=1,30,print1(numerator(sum(i=1, n, 1/i)),", "))

Obtenemos así rápidamente los primeros 30 términos:



Propiedades

Estos numeradores son también interesantes por sus propiedades. Algunas sorprenden, como que el primo p al cuadrado divide al término a(p-1) para p>3 (teorema de Wolstenholme)

Lo vemos en la tabla:


En ella se observa cómo los cuadrados de 5, 7, 11, 13 y 17 dividen a los términos de orden 4, 6, 10, 12 y 16 respectivamente. Sorprendente.

En http://oeis.org/A001008 puedes leer otras propiedades similares.

Desarrollo por Taylor

La función generatriz de los números armónicos es log(1-x)/(x-1). Podemos aprovecharla para generar de nuevo estos numeradores mediante PARI. Basta desarrollar la F.G usando el desarrollo de Taylor:

print(taylor(log(1-x)/(x-1),x,20))

Obtenemos:



Y basta leer los denominadores.

Números armónicos generalizados

Todo lo que hemos estudiado se puede extender al caso en el que la suma de fracciones se efectúe con potencias en los denominadores:

Si repetimos el algoritmo que hemos usado, pero con potencias, resultarían los numeradores generalizados, que se obtienen cambiando los sumandos 1/n por 1/n^m en el algoritmo presentado más arriba.

En Hoja de Cálculo:

Public Function numearmonico(n, m) ‘Valores de n y el exponente m
Dim i, v, f
v = 0
f = 1
For i = 1 To n: f = f * i ^ m: Next i
For i = 1 To n
v = v + f / i ^ m
Next i
v = v / mcd(v, f)
numearmonico = v
End Function

Y en PARI:

for(n=1,30,print1(numerator(sum(i=1, n, 1/i^m)),", "))

Resultan una serie de sucesiones, a las que OEIS da el nombre de “Wolstenholme numbers”, y que tienen propiedades similares a los numeradores estudiados. Las primeras, según su exponente son:

Con el 2 

http://oeis.org/A007406

1, 5, 49, 205, 5269, 5369, 266681,… y en ellas el primo p divide a a(p-1). Por ejemplo, el 5 divide a a(4)

Con el 3

http://oeis.org/A007408

1, 9, 251, 2035, 256103, 28567, 9822481,…

Con 4

http://oeis.org/A007410

1, 17, 1393, 22369, 14001361, 14011361,…

Con 5

http://oeis.org/A099828

1, 33, 8051, 257875, 806108207, 268736069,…

Todos tienen propiedades similares, según puedes consultar en los enlaces. Podemos resumir todos en una tabla de hoja de cálculo usando la función antes definida. Aquí tienes los primeros:




lunes, 5 de diciembre de 2016

Trocear y desplazar


En esta entrada estudiaremos algunos números enteros que mantienen una propiedad (ser primos, cuadrados, triangulares,…) si rotamos algunas de sus cifras, es decir, que algunas situadas a la derecha las desplazamos a los primeros lugares de la izquierda, o la operación inversa, rotar de izquierda a derecha. Por brevedad, sólo consideraremos el desplazamiento de derecha a izquierda.

Por ejemplo, 622521 es cuadrado perfecto, porque 622521=789^2.  Si ahora rotamos las cifras de la derecha 21 y las situamos al principio, resulta otro cuadrado: 216225=465^2. Recorreremos algunos casos concretos, diferenciándolos por la propiedad que conservan y por el número de cifras rotadas

Cuadrados con rotación de una cifra

El cuadrado de 13516 es 182682256. Si desplazamos el 6 al primer lugar queda otro cuadrado: 618268225=24865^2. Todos los números que poseen una propiedad similar están recogidos en http://oeis.org/A045877 y los primeros son: 1, 2, 3, 16, 21, 31, 129, 221, 247, 258, 1062, 1593, 1964, 2221, 13516, 17287, 18516,… En esta lista figuran las bases de los cuadrados que presentan esa propiedad, no los cuadrados.

Este ejemplo ya publicado nos servirá de guía para presentar la función rotate(m,n), que rota n cifras de la derecha en el número m. La idea de la función es calcular en primer lugar el módulo de m respecto a 10^n, porque ese será el número formado por las cifras que rotan. Una vez obtenido, lo multiplicamos por 10^(número de cifras –n), lo que sitúa esas cifras en la parte izquierda. Para rellenar la parte derecha, tomaremos las otras cifras, que son el cociente entero del número entre 10^n. Esas no las multiplicamos porque quedan a la derecha. La expresión completa sería a = (m Mod 10 ^ n) * 10 ^ (k - n) + m \ 10 ^ n.
Adjuntamos la función completa en Basic VBA:

Public Function rotate(m, n)
Dim a, k

If m < 10 ^ n Then rotate = m: Exit Function
k = Int(Log(m) / Log(10) + 1)
a = (m Mod 10 ^ n) * 10 ^ (k - n) + m \ 10 ^ n
rotate = a
End Function

Hemos usado módulo y división entera. Otra técnica sería convertir el número en una cadena, trocearla, rotarla después y volverla a convertir en número. Esa técnica no nos viene bien en hoja de cálculo, porque introduce unos molestos espacios en blanco. Sí es adecuada en PARI. El listado siguiente lo podrás entender con facilidad:

rotate(a,p)=if(a<10^p, a, eval(Str(a%10^p, Str(a\10^p))))

La función Str convierte en cadena tanto el módulo a%10^p como el cociente entero a\10^p. El paréntesis que los une sirve para concatenar las cadenas en orden inverso, y eval las convierte de nuevo en número.

Con el siguiente código podemos reproducir la sucesión presentada arriba, de números cuyos cuadrados siguen siéndolo si rotamos una cifra:

rotate(a,p)=if(a<10^p, a, eval(Str(a%10^p, Str(a\10^p))))
for(i=1,10^5,s=i*i;r=rotate(s,1);if(issquare(i)&&issquare(r),print1(i,", ")))

Aquí tienes el resultado, que coincide con el de arriba.


Cuadrados con rotación de dos cifras

Si la función rotate la usamos con dos cifras podemos descubrir qué cuadrados siguen siéndolo si trasladamos sus dos últimas cifras de derecha a izquierda. Para simplificar el problema eliminaremos de la lista los que terminen en dos ceros, que introducirían demasiados elementos en el listado. Basta concretar el parámetro de cifras en dos y exigir que el módulo de cada candidato respecto a 100 no sea cero.

Con hoja de cálculo descubrimos sólo tres casos de tres cifras con esas condiciones: 144, 196 y 625, que se convierten en los cuadrados 441, 961 y 256 respectivamente. Con cuatro cifras no aparecen, y con cinco 21609, 38416, 60025 y 86436.

Para mayor número de cifras usaremos PARI. El código es similar al que usamos más arriba, con dos ligeros retoques:

rotate(a,p)=if(a<10^p, a, eval(Str(a%10^p, Str(a\10^p))))
for(i=10,10^5,if(i%10>>0,n=i*i;r=rotate(n,2);if(issquare(r),print1(n,", "))))

Con él aparecen todos los casos hasta 10^10:

144, 196, 625, 21609, 38416, 60025, 86436, 120409, 161604, 236196, 272484, 363609, 436921, 481636, 622521, 638401, 646416, 870489, 904401, 1560001, 6240004, 14010049, 20602521, 29822521, 77422401, 82410084, 1632321604, 3672723609, 6023622544, 6529286416, 9089524921,…

Es atractivo 60025=245^2, que se convierte en 25600=160^2
Puedes adaptar el código para encontrar cuadrados que con la rotación de tres cifras sigan siendo cuadrados. Los primeros son estos:

16384, 25600, 36864, 61009, 2896804, 7049025,…

Y aquí tienes los primeros que rotan cuatro cifras. Los incluimos como un reto para reproducirlos:

166464, 210681, 214369, 216225, 364816, 690561, 842724, 898704, 962361, 1560001, 2010724, 6240004, 8042896, 14010049, 20412324, 23242041, 32410249, 56040196, 81649296, 82410084, 92968164,…

Triangulares con rotación de cifras

Se puede repetir el estudio con los triangulares. Con rotación de una cifra ya están publicados en http://oeis.org/A068071:

1, 3, 6, 10, 55, 66, 210, 253, 666, 780, 4851, 8911, 15400, 17766, 125751, 303810, 407253, 906531, 1125750, 2648451, 5925403, 6706953, 15772536, 22207780, 50085036, 78406503, 438094800, 1623331710, 1764803755,…

Por ejemplo, 407253 es el triangular número 902, porque  407253=902*903/2, y al desplazar la cifra 3 queda 340725, triangular número 825.

Para que podamos reproducirlos necesitamos recordar que un número N es triangular si 8*N+1 es un cuadrado. Todo lo que hemos presentado en PARI hasta ahora nos vale, sustituyendo issquare(i) por issquare(8*i+1)

Así, la sucesión anterior la podemos reproducir así (en sus primeros términos):

rotate(a,p)=if(a<10^p, a, eval(Str(a%10^p, Str(a\10^p))))
for(i=1,10^4,n=i*(i+1)/2;r=rotate(n,1);if(issquare(8*r+1),print1(n,", ")))



Con rotación de dos cifras, a partir de 100 resultan:

300, 325, 666, 5050, 5151, 5565, 6105, 6555, 13530, 14196, 36046, 78606, 187578, 395605, 508536, 837865, 975106, 1400301, 7340196, 9135675, 14426506, 24580566, 36265386, 38434528, 209582101, 338715378, 495070311, 2936270028, 4201373611, 5664257830, 5735794065, 8997105153, 9652787040, …

Con el código en PARI

rotate(a,p)=if(a<10^p, a, eval(Str(a%10^p, Str(a\10^p))))
for(i=24,10^5,n=i*(i+1)/2;r=rotate(n,2);if(issquare(8*r+1),print1(n,", ")))

Prueba a reproducir esta sucesión, y comprueba algún caso, como el de 187578, que es el triangular número 612, mientras 781875 es el número 1250.

Y con tres a partir de 1000:


1035, 1485, 1891, 30135, 46360, 60031, 75078, 96141

Rotación de primos

Con una cifra están publicados en A234901:

2, 3, 5, 7, 11, 13, 17, 31, 37, 71, 73, 79, 97, 113, 131, 173, 197, 199, 271, 277, 311, 313, 337, 373, 379, 397, 419, 479, 491, 571, 577, 593, 617, 631, 673, 719, 733, 811, 839, 877, 911, 919, 971, 977…

Para no cansar, reproducimos los primos que rotan dos cifras:

Con dos: 101, 103, 107, 113, 127, 131, 149, 157, 163, 181, 191, 197, 199, 307, 311, 317, 331, 337, 359, 367, 373, 701, 709, 719, 727, 733, 739, 757, 761, 787, 797, 907, 919, 937, 941, 947, 971, 983, 991

Podemos usar la siguiente rutina en PARI para ir descubriendo rotaciones de primos con más cifras. Basta ir cambiando el valor del segundo parámetro de rotate. La siguiente sirve para obtener las rotaciones de dos cifras:

rotate(a,p)=if(a<10^p, a, eval(Str(a%10^p, Str(a\10^p))))
forprime(i=101,10^5,r=rotate(i,2);if(isprime(r),print1(i,", ")))

Comprueba estos listados:

Primos con rotación de tres cifras:

1103, 1109, 1123, 1163, 1181, 1193, 1301, 1303, 1319, 1321, 1327, 1361, 1777, 1783, 1907, 1913, 1931, 1933, 1949, 1951, 1979, 1987, 1993, 1997, 2113, 2131, 2161, 2311, 2333, 2339, 2347, 2377, 2381, 2389, 2393, 2399, 2707, 2713, 2729, 2741, 2777, …

Primos con rotación de cuatro cifras:

10007, 10069, 10091, 10099, 10103, 10151, 10181, 10193, 10211, 10253, 10259, 10267, 10271, 10273, 10301, 10333, 10337, 10357, 10369, 10391, 10427, 10459, 10487, 10501, 10559, 10601, 10613, 10627, 10631, 10657, 10687, 10691, 10733,…

Y así podemos seguir.


miércoles, 23 de noviembre de 2016

Máximo producto en la partición de un número (2)


En la entrada anterior estudiamos la función que asigna a cada número entero positivo el máximo producto formado entre los sumandos de todas su particiones. Mediante recurrencia calculamos los valores de la función correspondiente, a la que llamamos MPP, pudiéndose formar tablas similares a la siguiente:



El objetivo de esta segunda entrada sobre el tema es demostrar que todos los valores de la función se calculan con una de estas tres expresiones: 3^n, 2*3^n, 4*3^n, de forma cíclica salvo el valor inicial 1.

Factores del producto máximo

Razonaremos a continuación que los valores de mpp(n) (n>4) solo poseen los factores 2 y 3. Todo se basa en que si un número N se descompone en dos sumandos a y b, el máximo producto a*b se produce en el interior del intervalo (1,N), debido a que la función x(N-x) presenta el máximo en  N/2. Consecuencia de esto es que si una partición recibe un sumando nuevo para crear otra partición superior, ese sumando se puede ir descomponiendo en sumandos 1, 2 y 3, de forma que el producto aumente. Por ejemplo, si el nuevo sumando es 6, se puede sustituir por 3+3, cuyo producto es 9, superior al 6. Como el 6 es factor del producto, también lo será el 9, y el producto aumentará.

Con sumandos mayores se puede proceder de igual forma, Por ejemplo, el 10 se puede sustituir por 3+3+3+1, con lo que el producto final, en lugar de ser multiplicado por 10, lo hará en 27 (obsérvese que el 1 no interviene en el proceso).

Todos los valores del producto máximo presentarán la descomposición del tipo 2^p*3^q.

Concretemos algo más:

Caso 1: números del tipo 3k (múltiplos de 3)

Los primeros valores, según la tabla anterior, son: mpp(3)=3, mpp(6)=9, mpp(9)=27, es decir, las primeras potencias de 3. Si pasamos al valor 12, según el razonamiento de los párrafos anteriores, el producto máximo recibirá el factor 3, luego se dará que mpp(12)=3^4=81 y mpp(15)=3^5=243.

Si N=3*k, su función mpp(N)=3^k

Caso 2: Números del tipo 3k+2

Para ellos es inútil acudir al anterior 3k+1, pues el producto sólo perdería el factor 1 que no incrementa su valor, pero si acudimos al máximo de 3k=3^k, el nuevo factor máximo es el 2, y queda mpp(3k+2)=2*3^k.

Por ejemplo, mpp(8)=mpp(3*2+2)=2*3^2=18, como ya sabemos por la entrada anterior. Otro: mpp(14)=mpp(3*4+2)=2*3^4=162.

Si N=3*k+2, su función mpp(N)=2*3^k

Caso 3: Números del tipo 3k+1

Si acudimos al valor 3(k-1)+2, basta usar como nuevo sumando el 2 para ver que el máximo que buscamos es 4*3^(k-1).

Vemos algún ejemplo: mpp(13)=mpp(3*4+1)=4*3^(4-1)=4*27=108, como ya sabemos por casos anteriores. El valor de mpp(31) = mpp(3*10+1) = 4*3^9 = 4*19683 = 78732.

Lo hemos comprobado con la definición de función de la entrada anterior:


Si N=3*k+1, su función mpp(N)=4*3^(k-1)

Nueva definición de la función

Estas propiedades que acabamos de estudiar nos permiten simplificar mucho la definición de la función MPP. La presentamos en varias versiones para que las analices. La primera está construida con las funciones de Excel o Calc. Es esta, que la hemos copiado actuando sobre la celda D11:

=3^(ENTERO((D11)/3)-(RESIDUO(D11;3)=1))*2^RESIDUO(-D11;3)

Parece difícil. Intenta comprenderla. La primera parte decide a qué exponente elevaremos el 3, y la segunda el mismo problema para el 2. Para el primero hallamos el valor de k que acompaña al 3 y le restamos un 1 en el caso 3k+1. Para el segundo elegimos 1, 2 o 4 según el residuo respecto al 3, pero con signo menos. Ya, es complicado.

Con ella hemos construido esta tabla, para los números 40 a 45:



Código en Basic VBA:

Es quizás más fácil de entender la definición (segunda ya) de mpp en VBA. Como ya tenemos una versión, le llamaremos a esta MPP_2. Su desarrollo puede ser este:

Public Function mpp_2(n)
Dim k, m, p

k = Int(n / 3) ‘se determina el valor de k en n=3k+b
If n Mod 3 = 1 Then m = 1 Else m = 0 ‘si es múltiplo de 3, el exponente de 3 disminuye en 1
p = 3 ^ (k - m) ‘parte correspondiente al factor 3
m = 3 - n Mod 3: If m = 3 Then m = 0 ‘se prepara el exponent del 2
mpp_2 = p * 2 ^ m ‘se ensambla la función
End Function

Este código es más simple y rápido que el anterior. Hemos comprobado su equivalencia, y se demuestra aquí el poder simplificador del razonamiento matemático.

Para quienes entiendan el lenguaje PARI, copio aquí la solución de M. Somos:

 mpp(n) = floor( 3^(n - 4 - (n - 4) \ 3 * 2) * 2^( -n%3))

Es muy parecida a la que proponemos.

Otra recurrencia

Si repasas la página http://oeis.org/A000792 que nos viene ayudando en este estudio, podrás leer más propiedades y fórmulas sobre estos productos máximos.

Una muy sencilla y curiosa es la de Ivan Neretin:

a(n) = a(n-1) + largest proper divisor of a(n-1), n > 2

En efecto, estudiamos los tres casos:

Si a(n-1)=3*k, sabemos que mpp(a(n-1))=3^k, y su mayor divisor propio md=3^(k-1). Sumamos y obtenemos mpp(a(n))=mpp(3*k+1)=3^k+3^(k-1)=4*3^(k-1), como ya sabemos.

Si a(n-1)=3*k+1, se tendrá: mpp(a(n-1))=4*3^(k-1), su mayor divisor propio md=2*3^(k-1). Sumamos: mpp(a(n))=mpp(3*k+2)=4*3^(k-1)+2*3^(k-1)=2*3^k, expresión ya conocida.

Si a(n-1)=3*k+2, mpp(a(n-1))=2*3^k, su mayor divisor propio md=3^k y al sumar obtenemos la esperada 3^(k+1).

Si dispones de la función en VBA de “mayor divisor”, en una columna de hoja de cálculo puedes engendrar la sucesión completa. En este blog disponemos de la función MAYORDIV (no la desarrollamos porque acude a otras más complicadas), lo que nos permite desarrollar la recurrencia:

Se ha adosado la recurrencia a la derecha de una tabla de MPP, para ver la equivalencia.

Si no dispones de la función MAYORDIV, puedes copiar la que sigue:

Public Function mayordiv(n)
Dim f, m
Dim es As Boolean

f = 2
m = 1
If n / 2 = n \ 2 Then es = True: m = 2 Else es = False
While f * f <= n And Not es
If n / f = n \ f Then es = True: m = f
f = f + 1
Wend
If m = 1 Then mayordiv = 1 Else mayordiv = n / m
End Function

lunes, 14 de noviembre de 2016

Máximo producto en la partición de un número (1)


Ya sabemos que una partición de un número entero positivo N es una suma de números también enteros positivos cuyo resultado es ese número. Ya hemos tratado en este blog el tema de las particiones, y lo volveremos a desarrollar próximamente. Si no tienes claro el concepto puedes acudir a las direcciones

http://hojaynumeros.blogspot.com.es/2011/02/particiones-de-un-numero.html

https://es.wikipedia.org/wiki/Partici%C3%B3n_(teor%C3%ADa_de_n%C3%BAmeros)

Está muy estudiado el tema del desarrollo de las particiones y el del cálculo de su número. Aquí nos interesará el máximo producto que se puede lograr multiplicando los sumandos de cada partición. Lo introducimos con ejemplos:

Máximo producto logrado con particiones

Tomemos el número 6. Mentalmente se pueden escribir sus particiones (no se tiene en cuenta el orden): 6=5+1=4+2=3+3=4+1+1=… En cada partición calculamos el producto entre sumandos: 6, 5, 8, 9, 4,… y nos quedamos con el máximo. En el esquema lo verás mejor:



Figuran las once particiones del 6 y en la columna de la derecha los productos de sumandos. En la partición de sumando único lo elegimos como producto. Se observa que el máximo producto es 9. Como el resultado es único, constituye una función del número elegido, que podríamos escribir como MPP (Máximo Producto en Particiones) y se tendría que MPP(6)=9.

Otro ejemplo: En PARI las particiones se obtienen con la función partitions. Si pedimos las particiones del número 8 las obtenemos como vectores independientes:



Si multiplicas los sumandos dentro de cada vector descubrirás que el máximo producto es 18=2*3*3, luego MPP(8)=18

Estos resultados figuran en la página de OEIS http://oeis.org/A000792 con una definición recursiva que ya trataremos:

1, 1, 2, 3, 4, 6, 9, 12, 18, 27, 36, 54, 81, 108, 162, 243, 324, 486, 729, 972, 1458, 2187, 2916, 4374, 6561, 8748, 13122, 19683, 26244, 39366, 59049, 78732, 118098, 177147, 236196, 354294, 531441, 708588, 1062882, 1594323, 2125764, 3188646, 4782969,…

Como era de esperar, los valores son crecientes (cada uno es un máximo que se apoya en los anteriores) y pronto adquieren una buena tasa de crecimiento. La página citada contiene una fórmula recursiva que es fácil de entender.

a(n) = max{ (n-i)*a(i) : i<n}; a(0) = 1

Se define a(0) como 1, y también es fácil entender las siguientes: a(1)=1, a(2)=2, a(3)=3,… La recursividad tampoco es difícil de captar. Se trata de multiplicar cada valor menor que n por el máximo correspondiente a su diferencia con n. En efecto, las particiones de n se forman eligiendo esos valores i:1…n-1 para luego unirlos con las particiones de n-i. Por ejemplo, en las particiones del 7, si elegimos el valor 4, se deberá combinar con las particiones del 3 para formar las particiones del 7 que contengan un 4. Igual ocurre con todos los valores: 5 se añadirá a las particiones del 2 y 3 se unirá a las particiones de 4.

Si conservamos los valores máximos de cada partición de n-i, al multiplicarlos por i resultarán productos de la partición superior, candidatos a ser máximos. Al recorrer todos los valores menores que n dispondremos de n-1 posibles máximos, y uno de ellos será el MPP.

Algoritmo de construcción de la función MPP

Las ideas anteriores nos permitirán construir la función  MPP. En VBA de Excel se puede usar esta definición de función:

Public Function mpp(n)
Dim mx, i, m, j, mm
Dim a(50) ‘Está preparado para n<=50. Se puede ampliar a otro número

If n = 0 Or n = 1 Then mpp = 1: Exit Function
If n = 2 Then mpp = 2: Exit Function ‘Casos particulares
a(0) = 1: a(1) = 1: a(2) = 2: mx = 2
If n > 2 Then
For i = 1 To n ‘Se recorren los valores anteriores para la recursión
m = 1 ‘Valor provisional del máximo
For j = 1 To i
mm = j * a(i - j)
If m < mm Then m = mm ‘Se busca un máximo nuevo mediante los productos con los anteriores
Next j
mx = m: a(i) = m ‘Se incorpora el máximo a la lista
Next i
End If
mpp = mx ‘Máximo final
End Function

Con esta función podemos encontrar el máximo producto entre particiones de cualquier número. Si es mayor que 50 bastará cambiar la dimensión del vector de máximos. En la tabla siguiente hemos recogidos los valores de MPP para los números comprendidos entre 30 y 40:


Para quienes conozcan el lenguaje PARI (gratuito y muy recomendable para estos temas) se inserta un código para esta función, que también devuelve los valores entre 30 y 40:

mpp(n)=my(a=vector(50), m, mm, mx=2,mp=1);a[1]=1;a[2]=2;if(n<2,mp=1,if(n==2,mp=2,for(i=3,n,m=1;for(j=1,i, d=i-j;if(d>0,mm = j * a[d],mm=j);if(m<mm,m=mm));mx=m;a[i]=m));mp=mx;mp)
for(k=30,40,print1(mpp(k),", "))

Aquí tienes el resultado, que coincide con el de Excel:



Interpretación algebraica

Estos valores coinciden con los cardinales máximos de los subgrupos del grupo simétrico S(n). Usando la descomposición en ciclos se les puede dar un significado (ver http://hojaynumeros.blogspot.com.es/2013/10/ciclos-2-descomposicion-en-ciclos.html)

Por ejemplo, en el caso del 8 visto más arriba, mpp(8)=18, y puede dársele el significado de que es el cardinal máximo de un subgrupo propio del grupo de permutaciones de 8. En concreto, usando la descomposición en ciclos, podría ser G=(1,2)(3,4,5)(6,7,8) o cualquiera de sus isomorfos. En el caso del 6 es fácil ver que el subgrupo maximal es el GM=(1,2,3)(4,5,6), de cardinal 3*3=9, que es el valor de mpp(6).

Existe una forma directa y simple para calcular mpp(n), sin recurrencias ni algoritmos. Como es un cambio importante en el desarrollo que hemos llevado hasta ahora, lo dejamos para la próxima entrada.