Lo injusto de los sorteos por letras

Todos hemos visto alguna vez algún sorteo en el que se escoge una letra aleatoriamente y los ganadores se seleccionan por orden alfabético, empezando por la primera persona cuyo apellido tenga esa letra como inicial. Todos nos hemos dado cuenta de que es un sistema injusto; por ejemplo, es un chollo llamarse Abad, porque probablemente conseguirás aquello que se sortee tanto si sale la W como si sale la X, Y, Z o la A.

La mayoría de nosotros, incluido yo, hemos pensado que la diferencia de probabilidades sería pequeña... bueno, pensar, lo que se dice pensar, no lo hemos pensado, simplemente lo hemos dado por supuesto.

Bueno, pues no, la diferencia no es pequeña. Yo me he dado cuenta después de que mi amigo Raúl Corvillo me llamase la atención sobre dos blogs, "Un dato vale más que mil palabras" y "La ciencia para todos". También está este artículo de una mujer cabreantemente apellidada Grima.

El problema es que las iniciales de los apellidos en España están distribuidas de una forma más caprichosa de lo que se podría esperar; estos datos están sacados del primer blog:

	frecuencia	% sobre el total
A	2.884.390	6,7%
B	2.263.664	5,2%
C	3.969.992	9,2%
D	1.747.696	4,0%
E	781.910	1,8%
F	1.877.528	4,3%
G	4.857.351	11,2%
H	992.297	2,3%
I	424.730	1,0%
J	722.854	1,7%
K	55.885	0,1%
L	2.250.441	5,2%
M	5.291.515	12,2%
N	699.534	1,6%
O	803.973	1,9%
P	3.042.595	7,0%
Q	185.195	0,4%
R	3.565.620	8,2%
S	3.201.882	7,4%
T	1.425.424	3,3%
U	171.705	0,4%
V	1.631.083	3,8%
W	48.578	0,1%
X	14.690	0,0%
Y	92.553	0,2%
Z	269.539	0,6%
TOTAL	43.272.624	99,8%

Por si algún día le pasase algo al blog de Eduardo, también hago una copia de su histograma, que ilustra maravillosamente lo irregular que es la distribucion de las iniciales:

Para echar las cuentas, he escrito este programita:

#include "stdafx.h" #include <iostream> using namespace std; int _tmain(int argc, _TCHAR* argv[]){ int i, winner, letter, imax, imin; double freq[26] = {6.7, 5.2, 9.2, 4.0, 1.8, 4.3, 11.2, 2.3, 1.0, 1.7, 0.1, 5.2, 12.2, 1.6, 1.9, 7.0, 0.4, 8.2, 7.4, 3.3, 0.4, 3.8, 0.1, 0.0, 0.2, 0.6 }; double pfirst[26]; double plast [26]; double sum = 0.0, p, noreach=0.0, percent, max, min, nochance = 0.0; for (i=0; i<26; i++) { sum += freq[i]; pfirst[i] = 0; plast [i] = 0; } cout << "sum = " << sum << "\n"; cout << "escriba la probabilidad de ganar como porcentaje:\n"; cin >> p; cout << "\n"; for (winner=0; winner<26; winner++) { percent = p; letter = winner; if (freq[letter]>p) nochance += freq[letter]-p; while (percent > freq[letter]) { pfirst[letter] += 100.0/26.0; plast [letter] += 100.0/26.0; percent -= freq[letter]; letter = (letter+1)%26; } pfirst[letter] += 100.0/26.0; } max = 0; min = 10000; for (i=0; i<26; i++) { if (pfirst[i]>max) { max = pfirst[i]; imax = i; } if (plast[i]<min) { min = plast[i]; imin = i; } } cout << "% de poblacion que no puede ganar: " << nochance << "\n"; cout << "la probabilida de que gane la primera " << (char)(65+imax) << " es " << pfirst[imax] << "\n"; cout << "la probabilida de que gane la ultima " << (char)(65+imin) << " es " << plast [imin] << "\n"; }

Empecemos imaginando que se sortea algo que va a ganar el 1% de los participantes; por ejemplo, se sortean 100 plazas de campamentos para niños y se presentan 10.000. Resulta que el 79% de los participantes no tiene ninguna posibilidad de conseguir una plaza, mientras que los 100 primeros solicitantes cuyo apellido empiece por la A tienen nada más y nada menos que el 19,23% de probabilidades de conseguirla.

Hasta aquí no hay ninguna sorpresa, era obvio que el sistema era especialmente injusto si la probabilidad de ganar era pequeña.

La auténtica sorpresa empieza a aparecer cuando descubrimos que, al aumentar el número de premios, el sistema no se hace justo rápidamente. Por ejemplo, si la probabilidad "global" de ganar es un razonable 15%, hay gente que tiene un 31% de posibilidades de ganar (las primeras A) mientras que otros tienen sólo un 4% (las últimas G).

¿Y si el porcentaje de premiados fuese el 50%? Ahora las primeras C ganan con probabilidad 62%, mientras que las últimas M sólo tienen un 38% de probabilidad de ganar. Es decir, Cabaretero tiene un 60% más de probabilidades de ganar que Mutante.

Si vamos al extremo de dar premios al 90% de los solicitantes, resulta que las primeras A, G y M tienen garantizado el premio, mientras que las últimas M sólo tienen un 69,23% de probabilidades de conseguirlo.

Supongamos que estuviésemos dispuestos a aceptar que un sorteo es "tolerablemente injusto" (he aquí un oxímoron) si al más beneficiado le da sólo un 25% más de probabilidades de ganar que al más perjudicado. ¿Cuándo es tolerable un sorteo basado en las iniciales de los apellidos? Sorpresa, sólo si al menos el 97% de los participantes van a conseguir el premio. No era esto lo que yo me esperaba.

Es previsible que este sistema funcione incluso peor en pueblos pequeños donde algunos apellidos se repiten mucho.

Francamente, me parece un sistema inaceptable, especialmente cuando hay una solución obvia y al alcance de cualquier notario del siglo XXI: se numeran las solicitudes (ya sea alfabéticamente, usando otro criterio, o sin usar ningún criterio), se escoge un número al azar, y se otorgan los premios a partir de ahí.

Simulación de olas

Le estamos poniendo olas al simulador de aerogeneradores de Gamesa para ver qué pasa cuando el oleaje aporrea la base de la torre de un aerogenerador en el mar. En esta simulación, las olas llegan a tener hasta 7 metros de altura y el agua tiene 20 metros de profundidad. Las flechas indican la velocidad del agua. La mayor parte del trabajo es de David Campillo.

He intentado publicar aquí un video en formato .avi, pero no he conseguido que funcione, así pongo un enlace en Facebook.

Sintaxis bizarra

En el trabajo, al hacer un interfaz entre MATLAB y una librería nuestra en C++, hemos seguido unas instrucciones quizás un poco anticuadas, y nos hemos encontrado con unas líneas bizarras de C++. Tras una investigación sintáctica, las rarezas se pueden simplificar a esto:

double a[2][2][2][2];
double (*data)[2][2][2][2];
data = new double[110][2][2][2][2];
data = (double (*)[2][2][2][2]) new double[32];

La única línea que resulta familiar es la primera; declara que a es una variable de un tipo peculiar, está formada por 16 doubles en un array unidimensional (sin los cuatro niveles usuales de indirección) al que se accede a través de no un índice, sino 4; es decir, el decimocuarto double no es a[13] sino a[1][1][0][1] (o quizás a[1][0][1][1], no estoy seguro ahora).

La segunda línea, que ya empieza a ser rara, simplemente declara que data es un puntero a un objeto del tipo descrito antes.

La tercera línea construye un objeto del tipo anterior pero de dimensiones [110][2][2][2][2]. Obsérvese que esto es posible porque todas las dimensiones son estáticas, y por tanto el compilador puede saber qué tamaño tendrá el objeto construido. La gracia está en que, aunque este objeto no tiene las dimensiones de los objetos a los que apunta data, en realidad se puede ver como un array de 110 objetos del tipo de data, luego el puntero del new se puede convertir al tipo de puntero de data. O algo parecido.

La cuarta línea construye un array de 32 doubles de los de toda la vida, y entonces convierte el puntero al tipo de puntero adecuado para poder asignarlo a data, de forma que los 32 doubles se ven como dos objetos del tipo double[2][2][2][2]. Una de las cosas llamativas de esta sintaxis es que, si le quitas los paréntesis al (*), el compilador no te entiende.

Para quien se haya quedado intrigado; Visual Studio es capaz de compilar esta línea, pero ¿cómo la interpreta?

double** b = (double**)(double* (***)[2][2][10]) new (double* (**)[2][50][2]) ();

Cómo ordenar una baraja cortándola

Un problema planteado por un compañero de trabajo. Hacemos un mazo con una baraja y lo desordenamos; sabiendo dónde están las cartas, ¿es posible cortar el mazo varias veces de alguna forma para reordenarla?

Los cortes "normales" consisten en separar el mazo en dos submazos y poner encima el que estaba debajo. Este tipo de cortes no cambia el orden relativo de las cartas, simplemente las rota, cambiando la que está encima, y por tanto no sirve para ordenar mazos salvo en raros casos.

Ahora bien, si se permiten cortes en los que el mazo se divida en varios submazos, ¿se puede conseguir ordenarlo?

Un ejemplo. Tenemos un mazo de 8 cartas, de la A a la H. Originalmente están en el orden BCEAFGHD, siendo la primera carta, B, la que está arriba. Podemos cortar el mazo por dos sitios, digamos (BC)(EA)(FGHD), dividiéndolo en tres bloques de cartas: BC, EA, FGHD. Primero dejaremos sobre la mesa el bloque FGHD, luego el EA, y luego el BC. Entonces recogeremos primero el primer bloque, FGHD, luego el segundo y finalmente el tercero, y tendremos las cartas del mazo en este orden: (FGHD)(EA)(BC). Ahora podría volver a cortar este mazo en tres bloques, (FGH)(DE)(ABC), y al recogerlos tendría (ABC)(DE)(FGH), que ya está ordenado.

La cuestión es si siempre podemos hacer esto, y cuántos cortes de qué tipo son necesarios.

Pensando en mazos pequeños, enseguida se encuentra un ejemplo curioso. Si tenemos un mazo de tres cartas, y sólo permitimos hacer cortes en tres bloques, no podremos ordenarlo, porque lo único que haremos será intercambiar las cartas de arriba y abajo. Es fácil ver que, si tenemos un mazo con n cartas que siempre se corta en n bloques, lo único que haremos será invertir el orden de las cartas.

Sin embargo, tenemos el siguiente teorema:

Teorema: si se permiten cortes de dos y tres bloques, siempre es posible ordenar un mazo de n cartas con 3n-3 cortes.

Demostración: supongamos que tenemos un grupo de cartas ya ordenado, D, y una carta B que queremos poner al principio o al final de D. Empezamos con un mazo ABCDE, donde alguno de los grupos A, C o E podría estar vacío. Si queremos colocar B delante de D, cortamos ABCDE en (A)(BC)(DE) para obtener DEBCA, que partimos en (DE)(B)(CA) y obtenemos CABDE, con B justo delante de D como queríamos. Si quisiéramos colocar B detrás de D, lo que podríamos hacer es partir en mazo en (AB)(CD)(E), luego partimos ECDAB en (E)(CDA)(B), y finalmente partimos BCDAE en (BC)(D)(AE) para llegar a AEDBC. Podría parecer que estos cortes usan tres bloques, pero como alguno de los grupos podría estar vacio, en realidad a veces usaremos cortes con sólo dos bloques. Por tanto, en dos o tres cortes podemos colocar cualquier carta delante o detrás de un bloque ya ordenado, de mode que en n-1 movimientos podemos ordenar el mazo, así que como mucho son necesarios 3n-3 cortes.

Por ejemplo, usando cortes de dos o tres bloques, un baraja española de 40 cartas siempre se puede reordenar en 117 cortes o menos. No es ésta una cota demasiado astuta. Es fácil ver que una baraja de 3 cartas siempre se puede reordenar en dos cortes de dos o tres bloques, y una de cuatro cartas en tres. Por cierto, si permitimos cortes de dos, tres o cuatro bloques, las barajas de cuatro cartas siempre se pueden reordenar en dos cortes.

¿Cuántos cortes de dos o tres bloques son necesarios para una baraja con n cartas en el peor caso? Me temo que es una pregunta demasiado difícil, pero se puede encontrar una cota trivial. Por ejemplo, un mazo de 40 cartas se puede cortar en dos bloques de 39 formas, y en tres bloques de comb(39,2) = 39*38/2 = 741 formas; incluso si todas las combinaciones de estos 780 movimientos diesen siempre resultados diferentes, no podríamos conseguir las 40! permutaciones diferentes usando 16 cortes, ya que 780¹⁶<40! ; por tanto, el peor caso necesitará al menos 17 cortes.

En cambio, si permitimos cortes en dos, tres, ... , cuarenta bloques, tenemos que el número total de cortes que podemos usar es 2³⁹-1 (porque podemos cortar o no entre cada par de cartas, pero no podemos no cortar en ningún sitio). Entonces, incluso si todas las combinaciones de estos cortes fuesen diferentes, necesitaremos 5 cortes en el peor caso, ya que (2³⁹-1)⁴<40!.

El problema de determinar los cortes óptimos para un mazo parece difícil, pero estuve pensando en este problema y en una tarde de hackeo compulsivo escribí este programa, que usa cortes de hasta n bloques y encuentra una solución con n-1 cortes o menos. En promedio consigue ordenar un mazo de 40 cartas ordenadas aleatoriamente en unos 14 cortes; AVISO, es un programa feo de narices:

//#include <string> //#include <fstream> #include "stdafx.h" #include <iostream> #include "math.h" #include <time.h> using namespace std; // Un tramo es una secuencia de cartas ya ordenadas, como 12 13 14. Quizás // contiene una única carta. Nunca se cortará en mitad de un tramo. // // Un bloque es una serie de cartas que formará uno de los cortes. Un bloque está // formado por uno o más tramos. // // Un superbloque está formado exactamente por dos bloques que al intercambiar // su orden relativo juntarán dos cartas. // // Por ejemplo: un mazo de 8 cartas como 2 8 6 7 1 3 4 5 se puede partir en // cinco tramos [2] [8] [6 7] [1] [3 4 5], cuatro bloques ([2]) ([8]) ([6 7] [1]) ([3 4 5]) // y un superbloque ([2]) ([8]) <([6 7] [1]) ([3 4 5])>. Al romper el mazo en los // cuatro bloques y cambiar su orden quedará ([3 4 5]) ([6 7] [1]) ([8]) ([2]), y se juntan // las cartas 5 y 6 que antes eran las cartas en los extremos del superbloque. // // También se juntan las cartas 1 y 2; a esto le he llamado optimización, porque en // realidad los tramos [2] y [8] podrían haber formado un bloque ([2] [8]) en vez de // dos bloques ([2]) ([8]), pero entonces no se habrían juntado. // // El algoritmo eurístico es el siguiente: // 1. formar los tramos (existe una forma única de hacerlo) // 2. formar los superbloques decrecientes. Cosas como (x+1 y) (z x) que al // cortar resultarán en (z x x+1 y). Uso una strategia codiciosa, hago una lista // con todos los posibles superbloques (definidos por las parejas de tramos de // sus extremos), elimino los superbloques que contienen completamente a otros // superbloques, y cojo todos los superbloques que puedo empezando por la // izquierda. // 3. entre los tramos que me quedan sin asignar a superbloques, busco superbloques // crecientes; es decir, cosas como (x y) (z x+1) que al reordenar se convertirán // en (z x+1 x y) y que aunque ahora no mejorarán nada, en la siguiente iteración // se podrán usar para hacer un superbloque decreciente. // 4. buscar optimizaciones - por dónde partir los superbloques para conseguir alguna // mejora. // 5. si quedan superbloques sin dividir en bloques, los rompo por la izquierda. int debug = 0; // escribe todo const int NCARTAS = 40; int carta[NCARTAS]; int usada[NCARTAS]; int buffer[NCARTAS]; // tramos int nTramos; int tp[NCARTAS]; // tramo principio int tf[NCARTAS]; // tramo final int tb[NCARTAS]; // bloque al que pertenece este tramo int sbat[NCARTAS]; // superbloque al que ha sido asignado este tramo // parejas para hacer superbloques int nParejas; int tpp[NCARTAS]; // tramo principio pareja int tfp[NCARTAS]; // tramo final pareja // super bloques int nSB; int sbp[NCARTAS]; // principio int sbf[NCARTAS]; // final int sbd[NCARTAS]; // direccion int sbl[NCARTAS]; // libre para optimizar; -1 no, positivo partir ahí // bloques int nBloques; // random values #include <time.h> struct RANDOM{ unsigned long super_z0; unsigned long super_w0; }; void InitRandoms(struct RANDOM *random) { srand( (unsigned)time( NULL ) ); random->super_z0=rand();//2345; random->super_w0=rand();//6542; } unsigned long SuperRand(struct RANDOM *random, unsigned long max) { unsigned long super_znew0 = ((random->super_z0=36969*(random->super_z0&65535)+(random->super_z0>>16))<<16); unsigned long super_wnew0 = ((random->super_w0=18000*(random->super_w0&65535)+(random->super_w0>>16))&65535); unsigned long SUPER_RAND0 = (super_znew0+super_wnew0); return ((SUPER_RAND0)%max); } RANDOM randomSeed; int myRand(int max) { return SuperRand(&randomSeed, 1000000) %max; } double myDoubleRand() { return (double)SuperRand(&randomSeed, 1048576)/1048576.0; } void getRandomPermutation() { int i, j, iCarta, contador; for (i=0; i<NCARTAS; i++) usada[i] = 0; for (j=0; j<NCARTAS; j++) { iCarta = (int) (myDoubleRand()*((double)(NCARTAS-j))); contador = 0; for (i=0; i<NCARTAS; i++) { if (usada[i]==0) { contador++; if (contador>iCarta) { carta[i] = j+1; usada[i] = 1; break; } } } } } void parteCartasEnTramos() { int i; // añadir al primer tramo la primera carta tp[0] = 0; tf[0] = 0; nTramos = 1; // decidir si cada carta se añade al último tramo o si empieza tramo nuevo for (i=1; i<NCARTAS; i++) { if (carta[i]==carta[tf[nTramos-1]]+1) { tf[nTramos-1]++; }else{ tp[nTramos] = i; tf[nTramos] = i; nTramos++; } } } void printTramos() { if (debug==1) { int i; for (i=0; i<nTramos; i++) { if (tp[i]==tf[i]) { cout << carta[tp[i]] << " "; }else{ cout << carta[tp[i]] << "-" << carta[tf[i]] << " "; } } cout << "\n"; } } const int DECRECIENTE = 0; const int CRECIENTE = 1; void buscaSuperBloques(int direccion, int primerTramo, int ultimoTramo) { int i, p1, p2, minimo, primeraPareja; bool cambios; // lista de todas las parejas nParejas = 0; for (p1=primerTramo; p1<=ultimoTramo-1; p1++) { for (p2=p1+1; p2<=ultimoTramo; p2++) { if (direccion==DECRECIENTE && carta[tp[p1]]==carta[tf[p2]]+1) { if (debug==1) cout << "pareja " << p1 << "-" << p2 << "\n"; tpp[nParejas] = p1; tfp[nParejas] = p2; nParejas++; } if (direccion==CRECIENTE && carta[tp[p1]]+1==carta[tf[p2]]) { if (debug==1) cout << "pareja " << p1 << "-" << p2 << "\n"; tpp[nParejas] = p1; tfp[nParejas] = p2; nParejas++; } } } // Eliminar aquellas parejas que contengan completamente a otras parejas cambios = true; while (cambios) { cambios = false; for (p1=0; p1<nParejas; p1++) { for (p2=0; p2<nParejas; p2++) { if (p1!=p2) { // si p1 incluye a p2, eliminar p1 if (tp[tpp[p1]] <= tp[tpp[p2]] && tf[tfp[p1]] >= tf[tfp[p2]]) { if (debug==1) cout << "eliminando pareja " << tpp[p1] << "-" << tfp[p1] << " porque incluye pareja " << tpp[p2] << "-" << tfp[p2] << "\n"; tpp[p1] = tpp[nParejas-1]; tfp[p1] = tfp[nParejas-1]; nParejas--; cambios = true; break; } } } if (cambios) break; } } if (debug==1) { cout << "Lista limpia de parejas:\n"; for (i=0; i<nParejas; i++) { cout << i << " " << tpp[i] << "-" << tfp[i] << "\n"; } } // mientras haya parejas, convertir la primera en superbloque y borrar las que solapen while (nParejas>0) { minimo = 10000; for (i=0; i<nParejas; i++) { if (tpp[i]<minimo) { minimo = tpp[i]; primeraPareja = i; } } if (debug==1) cout << "conviertiendo pareja " << tpp[primeraPareja] << "-" << tfp[primeraPareja] << " en superbloque\n"; sbp[nSB] = tpp[primeraPareja]; sbf[nSB] = tfp[primeraPareja]; sbd[nSB] = direccion; sbl[nSB] = -1; // marcar libre para optimizar for (i=tpp[primeraPareja]; i<=tfp[primeraPareja]; i++) sbat[i] = nSB; nSB++; i = 0; while (i<nParejas) { if (tpp[i]<=sbf[nSB-1] && tfp[i]>=sbp[nSB-1]) { if (debug==1) cout << "borrando pareja " << tpp[i] << "-" << tfp[i] << " porque solapa\n"; tpp[i] = tpp[nParejas-1]; tfp[i] = tfp[nParejas-1]; nParejas--; }else{ i++; } } } } void marcaBloques() { int i, j; // dividir superbloques en bloques; for (i=0; i<nSB; i++) { if (sbl[i]==-1) sbl[i]=sbp[i]; // si todavía están libres para optimizar, partir en el primer sitio posible for (j=sbp[i]; j<=sbl[i]; j++) { tb[j] = nBloques; } nBloques++; for (j=sbl[i]+1; j<=sbf[i]; j++) { tb[j] = nBloques; } nBloques++; } // agrupar los tramos que queden en bloques for (i=0; i<nTramos; i++) { if (tb[i]==-1) { tb[i] = nBloques; j = i+1; while (j<nTramos && tb[j]==-1) { tb[j] = nBloques; j++; } nBloques++; } } if (debug==1) { cout << "partición de tramos en bloques:\n"; for (i=0; i<nTramos; i++) cout << tb[i] << " "; cout << "\n"; } } void printBloques() { int bloque, i; bloque = -1; for (i=0; i<nTramos; i++) { if (tb[i]!=bloque) { // nuevo bloque. Cerrar el anterior si bloque != -1. if (bloque!=-1) cout << ") "; cout << "("; if (tp[i]==tf[i]) { cout << carta[tp[i]]; }else{ cout << carta[tp[i]] << "-" << carta[tf[i]]; } bloque = tb[i]; }else{ // mismo bloque. cout << " "; if (tp[i]==tf[i]) { cout << carta[tp[i]]; }else{ cout << carta[tp[i]] << "-" << carta[tf[i]]; } } } cout << ")\n"; } // Como printBloques, pero leo al revés, escribo cartas en vez de tramos en buffer, y luego copio a cartas. void mueveCartas() { int i, j, k, bloque, finalBloque, write; bloque = -1; write = 0; for (i=nTramos-1; i>=0; i--) { if (tb[i]!=bloque) { if (bloque != -1) { // escribir en buffer las cartas de los tramos i+1 a finalBloque for (j=i+1; j<=finalBloque; j++) { for (k=tp[j]; k<=tf[j]; k++) { buffer[write++] = carta[k]; } } } bloque = tb[i]; finalBloque = i; } } // escribir en buffer las cartas de los tramos 0 a finalBloque for (j=0; j<=finalBloque; j++) { for (k=tp[j]; k<=tf[j]; k++) { buffer[write++] = carta[k]; } } for (i=0; i<NCARTAS; i++) { carta[i] = buffer[i]; } } void printEstructura() { int iTramo, iCarta; if (debug==1) { cout << "pos carta tramo bloque superbloque\n"; for (iCarta=0; iCarta<NCARTAS; iCarta++) { // buscar tramo al que pertenece esta carta for (iTramo=0; iTramo<nTramos; iTramo++) { if (tp[iTramo]<=iCarta && tf[iTramo]>=iCarta) { printf("%2d %2d %2d %2d %2d\n", iCarta, carta[iCarta], iTramo, tb[iTramo], sbat[iTramo]); } } } } } void buscaSuperBloquesCrecientes() { int i, primero, ultimo; // buscar tramos sin asignar a superbloques for (i=0; i<nTramos; i++) { if (sbat[i]==-1) { primero = i; ultimo = i; while (i+1 < nTramos && sbat[i+1]==-1) { i++; ultimo = i; } buscaSuperBloques(CRECIENTE, primero, ultimo); } } } // poner en primero y ultimo los tramos delante de tramo que no están asignados // a ningún bloque ni a ningún superbloque void anterior(int tramo, int & primero, int & ultimo) { ultimo = tramo - 1; primero = tramo; while (primero-1 >= 0 && tb[primero-1]==-1 && sbat[primero-1]==-1) { primero--; } } // poner en primero y ultimo los tramos detrás de tramo que no están asignados // a ningún bloque ni a ningún superbloque void posterior(int tramo, int & primero, int & ultimo) { primero = tramo + 1; ultimo = tramo; while (ultimo < nTramos-1 && tb[ultimo+1]==-1 && sbat[ultimo+1]==-1) { ultimo++; } } // marcar los tramos de primero a ultimo como pertenecientes a un nuevo bloque void nuevoBloque(int primero, int ultimo) { int i; for (i=primero; i<=ultimo; i++) { tb[i] = nBloques; } if (ultimo>=primero) nBloques++; } void optimizaSuperbloques() { int i, primero, ultimo, j, k; // buscar superbloques sin optimizar y mirar si los puedo partir de forma que // formen una pareja creciente con los tramos anteriores que no están en un superbloque // y que no están asignados a ningún tramo. loop1: for (i=0; i<nSB; i++) { if (sbl[i]==-1 && sbp[i]+1<sbf[i]) { anterior(sbp[i], primero, ultimo); for (j=primero; j<=ultimo; j++) { for (k=sbp[i]+1; k<sbf[i]; k++) { if (carta[tf[k]]==carta[tp[j]]-1) { cout << "opt: unire " << carta[tf[k]] << " con " << carta[tp[j]] << "\n"; sbl[i] = k; // partir este superbloque desde sbp[i] hasta tramo k; nuevoBloque(primero, j-1); nuevoBloque(j, ultimo); goto loop1; } } } } } // Idem pareja creciente con tramos posteriores loop2: for (i=0; i<nSB; i++) { if (sbl[i]==-1 && sbp[i]+1<sbf[i]) { posterior(sbf[i], primero, ultimo); for (j=primero; j<=ultimo; j++) { for (k=sbp[i]+1; k<sbf[i]; k++) { if (carta[tp[k]]==carta[tf[j]]+1) { cout << "opt: unire " << carta[tp[k]] << " con " << carta[tf[j]] << "\n"; sbl[i] = k-1; // partir este superbloque desde sbp[i] hasta tramo k-1; nuevoBloque(primero, j); nuevoBloque(j+1, ultimo); goto loop2; } } } } } // Idem pareja decreciente con tramos anteriores loop3: for (i=0; i<nSB; i++) { if (sbl[i]==-1 && sbp[i]+1<sbf[i]) { anterior(sbp[i], primero, ultimo); for (j=primero; j<=ultimo; j++) { for (k=sbp[i]+1; k<sbf[i]; k++) { if (carta[tf[k]]==carta[tp[j]]+1) { cout << "opt: unire " << carta[tf[k]] << " con " << carta[tp[j]] << "\n"; sbl[i] = k; // partir este superbloque desde sbp[i] hasta tramo k; nuevoBloque(primero, j-1); nuevoBloque(j, ultimo); goto loop3; } } } } } // Idem pareja decreciente con tramos posteriores loop4: for (i=0; i<nSB; i++) { if (sbl[i]==-1 && sbp[i]+1<sbf[i]) { posterior(sbf[i], primero, ultimo); for (j=primero; j<=ultimo; j++) { for (k=sbp[i]+1; k<sbf[i]; k++) { if (carta[tp[k]]==carta[tf[j]]-1) { cout << "opt: unire " << carta[tp[k]] << " con " << carta[tf[j]] << "\n"; sbl[i] = k-1; // partir este superbloque desde sbp[i] hasta tramo k-1; nuevoBloque(primero, j); nuevoBloque(j+1, ultimo); goto loop4; } } } } } } int _tmain(int argc, _TCHAR* argv[]){ clock_t time = clock(); int i, nCortes; InitRandoms(&randomSeed); getRandomPermutation(); nCortes = 0; while (true) { for (i=0; i<NCARTAS; i++) cout << carta[i] << " "; cout << "\n"; parteCartasEnTramos(); if (nTramos==1) { cout << "FIN; usados " << nCortes << " cortes.\n"; break; } printTramos(); // partir en superbloques. nSB = 0; // ningún superbloque for (i=0; i<NCARTAS; i++) sbat[i]=-1; // ningún tramo asignado a ningún superbloque buscaSuperBloques(DECRECIENTE, 0, nTramos-1); buscaSuperBloquesCrecientes(); for (i=0; i<nTramos; i++) tb[i] = -1; nBloques = 0; cout << "\n"; optimizaSuperbloques(); // parte superbloques no optimizados en cualquier sitio for (i=0; i<nSB; i++) { if (sbl[i]==-1) { sbl[i] = sbp[i]; } } marcaBloques(); nCortes++; cout << "corte numero " << nCortes << ":\n"; printBloques(); mueveCartas(); } time = clock() - time; cout << "time spent = " << time << " ms\n"; cout << "press any key to end\n"; char cinget; cin.get(cinget); return 0; }

La idea es la siguiente.

Definamos un tramo como una serie de cartas (posiblemente sólo una) que ya están ordenadas crecientemente. Al hacer cortes, nunca cortaremos por la mitad un tramo (no tendría sentido deshacer el trabajo ya hecho), así que los bloques (grupos de cartas por donde corto) están formados por tramos. Obviamente, si un tramo acaba en la carta x, otro tramo empieza con la x+1, a menos que la carta x sea la última del mazo. El mazo está ordenado si y sólo si contiene un único tramo.

Si el mazo no está ordenado, entonces existe un tramo que empieza por x+1 y que aparece antes que un tramo que acaba en x. Porque, si no, todos los tramos ya estarían ordenados entre sí, y por tanto todo el mazo estaría ordenado.

Vale, pues entonces cortamos por delante de x+1, por detrás de x, y en algún sitio intermedio (como x y x+1 están en tramos diferentes sé que hay al menos una separación de tramos entre ellos por donde puedo cortar). Al reordenar los bloques habremos juntado dos tramos en un solo tramo. Así que cada vez que hacemos un corte podemos reducir en uno el número de tramos. Como empezamos con como mucho n tramos, en n-1 cortes o menos llegaremos a tener un único tramo, y entonces el mazo estará ordenado.

Un ejemplo: tengo el mazo 2 8 6 7 1 3 4 5; sus tramos son [2] [8] [6 7] [1] [3 4 5]; como el tramo [6 7] empieza por 6 y el [3 4 5] acaba en 5, parto el mazo en tres bloques ([2] [8]) ([6 7]) ([1] [3 4 5]) y al reordenarlo me quedará 1 3 4 5 6 7 2 8, que ahora tiene sólo 4 tramos: [1] [3 4 5 6 7] [2] [8].

Obviamente, hay varias formas de agrupar los tramos en bloques; si ese mazo lo hubiese partido en ([2] [8]) ([6 7] [1]) ([3 4 5]), al reordenar no sólo empalmaría [3 4 5] con [6 7], sino también [1] con [2], y me quedaría con un mazo con sólo tres tramos: [3 4 5 6 7] [1 2] [8]. Y además, hay otros pares de tramos que podría haber usado para partir el mazo; [2] y [1], y [8] y [6 7], también me habrían servido porque cumplen que la primera carta del primer tramo es la siguiente (x+1) a la última carta del segundo tramo (x).

El programa es una caca escrita en una tarde, y la estrategia heurística que sigue no está demasiado pensada, pero a pesar de todo, parece que siempre consigue reordenar un mazo de 40 cartas en 17 cortes o menos - el promedio está en unos 14 cortes. He hecho un par de pruebas con mazos de 1000 cartas, y los reordena en unos 200 cortes. Parece que en una iteración con n tramos este programa consigue juntar unos log(n) tramos, así que para reordenar un mazo de n cartas necesita algo más de n/log(n) cortes. Estoy convencido de que esto se debe de poder mejorar bastante.

Por qué el descubridor de la distribución t de Student no se llamaba Student

Resulta que, a principios del siglo XX, Claude Guinness se dedicaba a contratar a los mejores licenciados de Oxford y Cambridge para que trabajasen mejorando los procesos industriales de su cervecería.

Uno de sus fichajes fue William Sealy Gosset, un joven químico y matemático que trabajó en el control de calidad de stout, un tipo de cerveza negra. Gosset atacó el problema de minimizar el número de pruebas necesario para decidir si el stout era bueno. Los aficionados a la cerveza pensarán que este problema es estúpido, así que para justificarlo diré que cuanto menos stout se gaste en la fábrica haciendo pruebas, más barato se podrá vender el resto.

Por no entrar en detalles, cuando se hace un número grande de cierto tipo de pruebas, el resultado sigue aproximadamente una distribución normal. Pero Gosset quería hacer precisamente un número pequeño de pruebas, y la aproximación habitual no le servía. Así que echó sus cuentas, descubrió su distribución, y la encontró bastante útil, así que quiso publicarla.

Pero entonces se topó con la política de publicaciones del señor Claude Guinness, al que no le hacía nada de gracia que sus trabajadores revelasen a su competencia lo que él consideraba información confidencial sobre sus métodos de elaboración. Guinnes ya había tenido unos cuantos problemas con estos temas, y no estaba muy dispuesto a regalar al mundo un método que le ahorraba bastante dinero. Francamente, esto era un poco ponerle puertas al campo; no puedes reunir a los mejores científicos del mundo para que descubran la mejor forma de hacer cerveza y esperar que no publiquen sus resultados. Por otro lado, lo que Gosset quería publicar no eran las temperaturas óptimas para fermentar variedades de lúpulo, sino una distribución estadística...

Así que Gosset y Guinness llegaron a un acuerdo: Gosset podría publicar sus descubrimientos, pero usando un pesudónimo científico, para no dar pistas a los otros fabricantes de cerveza sobre cómo abaratar sus procesos de calidad. Este tipo de apaños era relativamente frecuente allá por 1908.

Esta es la razón por la que Gosset publicó sus artículos bajo el nombre de Student, y así hoy en día todos conocemos la distribución t de Student pero casi nadie sabe quién fue Gosset.

Incidentalmente, es posible que Guinness se preocupase demasiado. Es verdad que en aquellos tiempos ya estaba de moda aplicar metodologías científicas a todo tipo de industrias, pero la industria cervecera en concreto resultó ser muy conservadora, y la mayoría de sus competidores continuó con sus métodos tradicionales hasta que el éxito de Guinness les apabulló y les obligó a modernizarse.

Gosset siguió publicando y alcanzó renombre, a pesar de su nombre falso. Tuvo unas cuantas oportunidades de trabajar en la universidad, pero las rechazó porque tenía una familia grande y en Guinness le pagaban muy bien. Debió de ser un tío muy majo, porque consiguió la proeza de hacerse amigo a la vez de Pearson y de Fisher. Estos dos matemáticos se profesaban un intensísimo desprecio por sus desacuerdos sobre la filosofía de las pruebas de significancia, y sus agrias discusiones abarcaron dos generaciones; como esta historia es larga, los interesados podrán alucinar leyendo, por ejemplo, http://www.econ.uba.ar/www/institutos/epistemologia/marco_archivos/ponencias/actas%20xiii/trabajos%20episte/urbisaia-brufman_trabajo.pdf o http://stats.org.uk/statistical-inference/Lenhard2006.pdf (la verdad es que no sé si hay mejores artículos donde cuenten bien esta historia, no he buscado mucho).