Particiones y números de Stirling de segunda clase

Un número de Stirling de segunda clase es el número de formas de particionar un conjunto de n elementos en k subconjuntos no vacíos, denotándose como S2(n, k) o {nk}.

Se denota S2(n, k) para diferenciarlo de los números de Stirling de primera clase que, sin signo, lo denotábamos como S(n, k) o [nk], mientras que con signo eran s(n, k). En lo que sigue vemos las relaciones de las números de Stirling de primera y segunda clase con el factorial ascendente y descendente, también llamado símbolo de Pochhammer o rising factorial, denotado como x⁽ⁿ⁾ o xⁿ para el ascendente y falling factorial, denotado como (x)_n o xⁿ, para el descendente:

Los números de Stirling de segunda clase obedecen a la siguiente recurrencia:

que usando la denotación S2(n, k) y agregando los casos bases queda así:

En la aplicación Combine Tools, en la pestaña "Numéricos", implementamos este algoritmo:

// S2(n, k) = S2(n-1, k-1) + k * S2(n-1, k)
function stirlingS2(n, k) {
    if (n===k){
        return 1;
    } else if (k>n || k===0) {
        return 0;
    } else {
        return stirlingS2(n-1, k-1) + k * stirlingS2(n-1, k);
    }
}

obteniéndose esta tabla de primeros valores:

Se definen recurrencias para obtener los números de Stirling de segunda clase:

Y también una fórmula cerrada, es decir, que no es recursiva:

Dado que un número de Stirling S2(n, k) es una partición de orden k-ésima, podemos contar todas las particiones de un conjunto de n elementos así:

En la aplicación lo implementamos con este simple algoritmo bell2(n) que usa el algoritmo anterior stirlingS2(n, k):

function bell2(n) {
    let sum = 0;
    for (let j=0; j<=n; j++) {
        sum += stirlingS2(n, j);
    }
    return sum;
}

Obtenemos estos valores para el número de particiones con los primeros valores de n:

Algoritmo para generar particiones usando números de Stirling de segunda clase

Usando la expresión [1] que vimos antes B(n) = ∑_j=0..n S2(n, j), vamos a implementar un algoritmo que genere las particiones, observando que el código implementa la recurrencia S2(n, k) = k × S2(n-1, k) + S2(n-1, k-1), donde la variable k está implícita en la longitud del subresultado res:

/* B(n) = sum_(j=0)^(n) S2(n, j)
   S2(n, j) = S2(n-1, j-1) + j × S2(n-1, j) */
function partition4(list=[], n=list.length, res=[], result=[]) {
    if (n===0) {
        //iter += 1;
        //if (withCopyCost) iter += list.length;
        result.push(copy2(res));
    } else {
        //iter += 1;
        let pivot = list[list.length-n];
        for (let j=0; j<res.length; j++) {
            //iter += 1;
            res[j].push(pivot);
            partition4(list, n-1, res, result);
            res[j].pop();
        }
        res.push([pivot]);
        partition4(list, n-1, res, result);
        res.pop();
   }
   return result;
}

Sobre el conjunto con tres elementos {a, b, c} obtenemos estas particiones:

que de forma compacta, tal como explicamos en el tema sobre particiones, también puede representarse así:

En cualquier caso la aplicación nos los devuelve como Arrays, aunque deben entenderse que son conjuntos. La vista de la salida en JSON sería esta:

[[["a","b","c"]],
[["a","b"],["c"]],
[["a","c"],["b"]],
[["a"],["b","c"]],
[["a"],["b"],["c"]]]

Observamos la simplicidad del algoritmo con respecto a los que vimos en los dos temas anteriores. El siguiente esquema expone la evolución del algoritmo para la ejecución de partition4(["a", "b", "c"]),

list=[a,b,c] res=[]
// FIRST CALL
partition(n=3, res=[])
n=3  pivot=a j=0..-1
// EMPTY LOOP
res.push([a]) => [[a]]
partition(n=2, res=[[a]])
    n=2  pivot=b  j=0..0
    j=0
        res[0].push(b) => [[a,b]]
        partition(n=1, res=[[a,b]])
            n=1 pivot=c  j=0..0
            j=0
                res[0].push(c) => [[a,b,c]]
                partition(n=0, res=[[a,b,c]]) ==>> [[a,b,c]]
                res[0].pop() = [[a,b]]
            res.push([c]) => [[a,b],[c]]
            partition(n=0, res=[[a,b],[c]])   ==>> [[a,b],[c]]
            res.pop() => [[a,b]]
        res.pop() => []
    res.push([b]) => [[a],[b]]
    partition(n=1, [[a],[b]])
        n=1 pivot=c j=0..1
        j=0
            res[0].push(c) => [[a,c],[b]]
            partition(n=0, res=[[a,c],[b]])   ==>> [[a,c],[b]]
            res[0].pop() => [[a],[b]]
        j=1
            res[1].push(c) => [[a],[b,c]]
            partition(n=0, res=[[a],[b,c]])   ==>> [[a],[b,c]]
            res[1].pop() => [[a],[b]]]
        res.push([c]) => [[a],[b],[c]]
        partition(n=0, [[a],[b],[c]])         ==>> [[a],[b],[c]]
        res.pop() => [[a],[b]]
    res.pop() => [[a]]
res.pop() => []

Con => indicamos el estado del subresultado res en cada momento, observando que cuando hacemos una llamada con n=0 llegamos al caso base, devolviendo una partición que ponemos en una columna a la derecha precedida por ==>>.

Para calcular el coste observamos el algoritmo, donde el bucle for (let j=0; j<res.length; j++) empieza con el array res=[] vacío con longitud 0 y puede llegar a tener la longitud máxima n, pues para una lista inicial [a, b, c] tendremos 1 partición con longitud 1: [a,b,c], 3 particiones con longitud 2: [[a,b],[c]], [[a,c],[b]], [[a],[b,c]] y 1 partición de longitud 3: [[a],[b],[c]]. Por lo tanto la longitud de res varía entre k=0 inicialmente y un máximo de k=n=3 en este ejemplo. Por lo tanto el primer planteamiento del coste sería el siguiente:

Esto es una recurrencia con dos variables, pues aparte de n tenemos k. Dada la dificultad saber cuánto vale k en cada caso, intentaremos obtener el coste de un forma indirecta. Por un lado si ejecutamos el algoritmo contando iteraciones obtendremos los valores Iter que figuran en la siguiente tabla, usando también los valores B(n) ya vistos en la tabla al final del apartado anterior, obteniendo las diferencias Iter - B(n):

Si prescindimos del primer valor para n=0, la secuencia a(n) = 1, 3, 8, 23, 75, 278, 1155, 5295 se encuentra en OEIS A024716 y obedece a la siguiente fórmula para n≥1

También obedece a la siguiente fórmula alternativa según OEIS, en este caso para n≥0:

Así que tomando [3] tenemos que

de donde obtenemos el coste del algoritmo partition4(list)

Si usamos [4] ajustándolo al mismo offset n≥1, para lo cual usamos C(n, i+1) llegamos a una segunda fórmula alternativa para el coste de partition4(list)

Donde B(n) son los números de Bell que se calcula con esta recurrencia

y S2(n, k) son los números de Stirling de segunda clase que tiene la siguiente recurrencia

El coste calculado puede compararse con el recuento de iteraciones en la aplicación obteniendo una tabla de valores, donde se separan con "#" comprobándose que son iguales:

En caso de coste de copia hay que agregar n B(n) pues hay B(n) subresultados y cada uno tiene una longitud n:

Particiones parciales y algoritmo para generarlas

En el apartado anterior vimos que no fue fácil desarrollar la recurrencia del coste debido a que el algoritmo partition4(list) hacía uso de la recurrencia de los números de Stirling de segunda clase S2(n, k) = S2(n-1, k-1) + k × S2(n-1, k). Si esta recurrencia depende de dos variables n, k, cualquier algoritmo que la use también tendrá implícito esas dos variables.

Intentaremos implementar directamente S2(n, k) = S2(n-1, k-1) + k × S2(n-1, k) en algo que podríamos denominar particiones parciales de un conjunto de n elementos tomados en k particiones. Para los primeros valores de S2(n, k) tenemos:

Desarrollando los resultados de esos primeros valores:

Este algoritmo lo denominaremos kPartition(list, k), en este caso la primera versión kPartition1:

// S2(n, k) = S2(n-1, k-1) + k × S2(n-1, k)
function kPartition1(list=[], k=0) {
    let result = [];
    let n = list.length;
    if (n<k || k===0) {
        //iter += 1;
        result.push([]);
    } else if (n===k) {
        //iter += n;
        result.push([...list.map(v => [v])]);
    } else if (k===1) {
        //iter += n;
        result.push([[...list]]);
    } else {
        //iter += 1;
        let pivot = list[n-1];
        //iter += (n-1);
        let sublist = list.slice(0, n-1);
        let res = kPartition1(sublist, k-1);
        /* res.length = S2(n-1, k-1) */
        for (let r of res) {
            //iter += 1;
            r.push([pivot]);
            /* n+k = iter copy2(r) =  iter copyArray(r) */
            //if (withCopyCost) iter += n+k;
            result.push(copy2(r));
            r.pop();
        }
        res = kPartition1(sublist, k);
        /* res.length = S2(n-1, k) */
        for (let r of res) {
            //iter += 1;
            /* res.length = k */
            for (let rr of r) {
                //iter += 1;
                rr.push(pivot);
                /* n+k = iter copy2(r) =  iter copyArray(r) */
                //if (withCopyCost) iter += n+k;
                result.push(copy2(r));
                rr.pop();
            }
        }
    }
    return result;
}

En cuanto a las dos sentencias if (withCopyCost) iter += n+k, es necesario copiar los subresultados para pasarlos al array final result, pues en otro caso estaríamos pasasando la misma referencia. Cada subresultado parcial es como, por ejemplo, el subresultado [[a, b], [c]], por lo que hay que copiar con una profundidad de 2. Podemos usar una función copyArray(arr=[]):

function copyArray(arr=[]) {
    let arrCopy = [];
    for (let i=0; i<arr.length; i++) {
        iter += 1;
        arrCopy[i] = [];
        for (let j=0; j<arr[i].length; j++) {
            iter += 1;
            arrCopy[i].push(arr[i][j]);
        }
    }
    return arrCopy;
}

O bien usar copy2(arr=[[]]) que en la aplicación se utiliza para copiar particiones seleccionando algún método. Por defecto se usa slice. El método loop equivale a la anterior función copyArray():

function copy2(arr=[[]]) {
    if (copyMethod==="slice") {
        return arr.slice(0).map(v => v.slice(0));
    } else if (copyMethod==="map") {
        return arr.map(v => v.map(w => w));
    } else if (copyMethod==="concat") {
        return [].concat(arr).map(v => [].concat(v));
    } else if (copyMethod==="spread"){
        return [...arr].map(v => [...v]);
    } else if (copyMethod==="json") {
        return JSON.parse(JSON.stringify(arr));
    } else if (copyMethod==="loop") {
        let arrCopy = [];
        for (let i=0; i<arr.length; i++){
            arrCopy[i] = [];
            for (let j=0; j<arr[i].length; j++) {
                arrCopy[i].push(arr[i][j]);
            }
        }
        return arrCopy;
    } else {
        return arr;
    }
}

El planteamiento de la recurrencia del coste T'(n, k), considerando el coste de copia, es el siguiente:

Si hacemos n+k = 0 que equivale a ignorar el coste de copiar los subresultados, la recurrencia sin coste de copia queda así:

donde se observa claramente que S2(n-1, k-1) + S2(n-1, k) k = S2(n, k) con lo que nos queda por ahora así el coste de kPartition1(n, k) en la aplicación, pues no he podido encontrar algo mejor:

Podemos usar la fórmula cerrada (no recursiva) [1] para calcular los números de Stirling de segunda clase:

Considerando el coste de copiar los subresultados, que según vemos en el código del algoritmo se produce cuando insertamos los subresultados en el array result que contiene los resultados finales, hemos de suma el coste n+k que es del copiar un subresultado:

...
//if (withCopyCost) iter += n+k;
result.push(copy2(r));
...

La fórmula con coste de copia para kPartition1(list, k) es, entonces, lo mismo que [5] que vimos antes al plantear la recurrencia del coste del algortimo, pues como ya dije, no pude obtener una mejor fórmula:

Algoritmo para generar particiones usando particiones parciales

Usando las particiones parciales vistas antes y sabiendo que el total de particiones es el número de Bell, entonces B(n) = ∑_j=0..n S2(n, j), implementamos un quinto algoritmo no recursivo para generar particiones:

// B(n) = sum_(j=0)^(n) S2(n, j)
function partition5(list=[]) {
    let result = [];
    let n = list.length;
    //iter += 1;
    for (let k=1; k<=n; k++){
        //iter += 1;
        let res = kPartition1(list, k);
        /* res.length = S2(n, k) */
        for (let j=0; j<res.length; j++) {
            //iter += 1;
            result.push(res[j]);
        }
    }
    return result;
}

El coste de este algoritmo iterativo se deduce fácilmente:

siendo T_k el coste de las particiones parciales kPartition1(list, k) que vimos en el apartado anterior. Con copia es lo mismo, usando el coste de copia de las particiones parciales:

Comparando algoritmos para generar particiones

Comparamos coste de los algoritmos vistos para generar particiones: partition1(list), partition2(list), partition3(list) en los temas anteriores y partition4(list) y partition5(list) en este tema:

Con partition5(list) se consiguen los peores resultados, pues hace uso de las particiones parciales. El mejor resultado es para partition4(list).