Informática geek, matemáticas, pensamiento crítico y alguna otra cosa de vez en cuando.

2009-12-26

Fraccionar archivos usando Reed-Solomon

Hace tiempo que necesitaba una utilidad capaz de dividir un archivo en varias partes. Pero no que hiciera sólo eso, que para eso ya existen muchas (Hacha, por ejemplo), sino una tal que si se pierden algunas partes (con un valor máximo de «algunas» decidido de antemano) todavía se pueda recuperar el archivo original.

El uso que le quería dar es una copia de seguridad distribuida. Se divide el archivo que se quiere guardar en múltiples partes y se suben a diferentes lugares: MediaFire, MegaUpload, FileFactory, 4shared, ADrive... debidamente cifradas de antemano si queremos preservar nuestra intimidad, claro. Si a la hora de recuperar las partes hay alguno de los servicios que ha cerrado, o que ha borrado el archivo, o cuya clave no nos funciona por alguna razón, siempre es posible reconstruir el original usando los demás, suponiendo que el número de bajas no sea excesivo.

La solución obvia es subir el archivo completo a cada uno de los alojamientos remotos. El problema está en el tamaño máximo soportado por cada alojamiento y en el ancho de banda total que hay que usar para subir tantas veces el mismo archivo.

Esto sugiere una pregunta: ¿Es posible obtener, a partir de un archivo dado, n+m archivos, de manera que si de ellas se pierden m partes cualesquiera, se pueda reconstruir el archivo original?

Sí, se puede, y además de forma óptima, es decir, que si el archivo tiene B bytes, el tamaño de cada bloque sería exactamente B/n si B es múltiplo de n. Si no lo es, no importa: nos guardamos la longitud original del archivo y añadimos bytes (por ejemplo, ceros) hasta llegar al siguiente múltiplo de n. Después, tras reconstruirlo, truncamos el archivo por su longitud original.

De hecho, el principio para hacerlo es el que se emplea en el RAID. No me refiero a la marca de insecticida, sino a un sistema genérico para agrupar varios discos (Redundant Array of Inexpensive —o Independent— Disks).

Según la forma de usar los discos, existen varios tipos de RAID, a los que se les asigna números. Algunos de esos tipos contienen información redundante para prevenir fallos de las unidades de disco. Los más básicos son RAID 1, RAID 3 y RAID 6. RAID 1 es el sistema más obvio: todos los discos tienen el mismo contenido. Es como la idea de subir el mismo archivo a todas partes. RAID 3 lo hace mejor: si hay N discos, la capacidad efectiva es N-1, y si se estropea uno, se puede reconstruir a partir de los demás. El sistema que emplea es muy sencillo: cada bit del disco extra es el XOR de todos los demás. Así, si R, S y T son los discos que contienen los datos normales, el disco U contiene R⊕S⊕T. Si se estropea R, se puede recuperar su contenido calculando S⊕T⊕U, ya que puesto que U=R⊕S⊕T, S⊕T⊕U = S⊕T⊕(R⊕S⊕T) = (quitando paréntesis y reordenando) S⊕S⊕T⊕T⊕R = 0⊕0⊕R = R. Lo mismo se aplica al resto de discos.

Sin embargo, no se puede extender el mismo truco a otro disco más y conseguir que si fallan dos discos cualesquiera, ambos se puedan recomponer. El truco del XOR está limitado a un disco fallido como máximo. Esto quiere decir que si falla un disco durante la recomposición de los datos originales (o hasta que se reemplaza el disco fallido), se pierden los datos. Para conseguir tolerancia a dos fallos o más, hace falta recurrir a los códigos Reed-Solomon. Esto es lo que hace el RAID 6.

No es trivial operar con estos códigos; hacen falta algunos conceptos de álgebra lineal y de estructuras algebraicas. Puesto que no encontré ningún programa ya hecho para dividir y recomponer datos, tuve que hacerme uno propio. Lo que sigue es una síntesis de lo que he averiguado durante su elaboración.

Hay un tutorial sobre códigos Reed-Solomon aplicados a sistemas RAID, aquí: [1] (PDF). La descripción de la matriz necesaria para generar los códigos (y su inversa para decodificarlos) es errónea; hay un artículo que la enmienda aquí: [2] (PDF).

Así que tenemos un archivo conceptualmente dividido en n bloques de datos, y queremos generar a partir de ellos un total de n+m bloques de datos, tales que si se pierden cualesquiera m de los n+m bloques, podamos reconstruir los que se han perdido. Para conseguirlo, dividiremos los bloques en palabras, es decir, conjuntos de bits de una longitud fija predeterminada w. El número de bits por palabra w es un valor escogido de tal forma que 2w sea mayor que el máximo valor de n+m que vamos a emplear. Por ejemplo, si escogemos w=8 (por tanto las palabras son bytes), podremos dividir el archivo original en un máximo de 255 partes contando las que se pueden perder; por ejemplo, n=200, m=55.

Pero para empezar, supongamos que trabajamos en el dominio de los números reales; luego ya veremos cómo hacer lo mismo sólo con enteros de w bits. El problema consiste entonces en encontrar, a partir de una serie de datos d1, d2, ..., dn, una serie de elementos c1, c2, ..., cn+m que tengan la propiedad deseada. Veamos cómo hallarlos.

Dispondremos los datos en forma de una matriz columna, D. Para resolver el problema hemos de encontrar dos matrices A y B. La matriz A debe cumplir esta condición: que el producto AD sea una matriz columna C que contiene los elementos c1, c2, ..., cn+m. La matriz B, por su parte, tiene que cumplir esta condición: si cogemos n elementos distintos cualesquiera de los n+m que tiene C y formamos con ellos una matriz columna C’, entonces BC’=D.

Si tenemos esas dos matrices A y B, el problema estará resuelto: si podemos usar A para calcular C, obtendremos c1, c2, ..., cn+m y por tanto generar los bloques; si podemos usar B para calcular D a partir de cualquier subconjunto de n elementos distintos de C, podemos reconstruir los datos originales.

El problema, pues, radica en cómo construir las matrices A y B. Veamos la forma.

Para poder multiplicarla por D, la matriz A tiene que tener n columnas; para dar como resultado una matriz columna con n+m elementos, tiene que tener n+m filas. La característica que deseamos que tenga A, ya veremos después por qué, es que al tachar m filas cualesquiera, la matriz resultante sea siempre invertible (lo que equivale a requerir que su determinante no sea nulo). ¿Cómo lo conseguimos?

La respuesta viene de la mano de una familia de matrices llamadas matrices de Vandermonde. Estas matrices se construyen como sigue: sean a1, a2, ..., ak escalares arbitrarios; entonces la fila i está compuesta por una sucesión geométrica que comienza en ai0=1:

Definición de la matriz de Vandermonde

(Se define 00=1 a este efecto). El determinante de una matriz cuadrada de Vandermonde nunca es cero si ninguno de los escalares está repetido (concretamente, su valor es ∏(aj - ai) para todos los posibles pares ordenados ai, aj de escalares que intervienen, que no será cero a menos que haya escalares repetidos), así que, si formamos una matriz rectangular de Vandermonde más alta que ancha con todos los escalares distintos, entonces al tachar cualesquiera filas que la hagan cuadrada, el determinante de la matriz resultante no será cero.

De modo que ese es el truco para construir nuestra matriz A: como una matriz de Vandermonde. Usaremos como escalares los números de fila empezando en cero, y se elevarán a los números de columna empezando en cero. La primera fila será pues: 00=1, 01=0, 02=0, ..., 0n-1=0; la segunda será 10=1, 11=1, 12=1, ..., 1n-1=1; la tercera será 20=1, 21=2, 22=4, 23=8, ..., 2n-1, etcétera. Todas sus submatrices cuadradas de orden n serán entonces invertibles, que es la propiedad que buscábamos.

Ahora veamos cómo construir B para reconstruir los datos originales a partir de C’, que como hemos dicho, es un subconjunto de n elementos de C. Está claro que B tiene que ser una matriz cuadrada de orden n, ya que al multiplicarla por C’, que tiene n elementos, tiene que darnos D, que también tiene n elementos.

Fijemos nuestra atención por un momento en cómo se construye la matriz C. Cada elemento ci es el producto de la fila i de A por la única columna de D. Construyamos ahora una matriz A’ cuyas filas sean únicamente aquellas filas de A cuyo índice i coincide con el correspondiente de cada elemento ci en C’. Es decir, construyamos A’ de tal manera que se cumpla la relación A’D=C’.

Entonces ahora es obvio que para hallar D, basta con invertir la matriz A’ puesto que A’-1A’D=A’-1C’, luego D=A’-1C’, por tanto la matriz B buscada es simplemente A’-1. Sabemos además que A’ es invertible porque la hemos construido así a propósito.

Recapitulando lo que llevamos hasta ahora, se trata de construir un par de matrices A y B tales que al tachar cualesquiera m filas de A, la submatriz resultante siga siendo invertible; para reconstruir los datos a partir de un subconjunto C’ de n bloques de los n+m que hemos generado, construimos A’ como la submatriz de A resultante de dejar únicamente las filas que se corresponden con las presentes en C’, y al invertirla obtenemos B.

Veamos un ejemplo. Tomemos n=3 y m=2. Tenemos tres datos d1, d2 y d3 a partir de los cuales queremos obtener cinco elementos c1, c2, c3, c4 y c5. Construimos primero la matriz A y la usamos para obtener los bloques:

Ejemplo del esquema de codificación Reed-Solomon con n=3, m=2

Ahora supongamos que sólo tenemos los bloques 1, 3 y 5 y nos faltan los bloques 2 y 4, y queremos reconstruir los datos d1, d2 y d3. Construimos la matriz C’ que consta simplemente de c1, c3 y c5 y construimos también A’, usando las filas 1, 3 y 5 de A, la invertimos y la multiplicamos por C’ para obtener D:

Ejemplo del esquema de decodificación Reed-Solomon basada en la codificación anterior, con n=3, m=2 y teniendo los bloques 1, 3 y 5.

Como cortesía adicional, podemos refinar la matriz A para que cada ci para i entre 1 y n sea igual a di. Para ello, se pueden convertir las n primeras filas de A en una matriz identidad mediante transformaciones elementales, que son unas transformaciones que no afectan a la invertibilidad de la matriz. El algoritmo a aplicar, que también nos sirve para invertir la matriz, se llama algoritmo de Gauss-Jordan. Está explicado en la referencia [2].

Muy bien, muy bonito, pero estamos trabajando con reales de tamaño arbitrario, no con enteros de w bits. ¿Cómo hacemos para que no nos aparezcan números racionales o demasiado grandes?

Los cuerpos de Galois

Un cuerpo finito o de Galois (abreviado GF) es un cuerpo con un número de elementos finito. Hasta aquí, obvio. El más simple es GF(2) que tiene dos elementos, llamémosles 0 y 1. En ese cuerpo, la suma se define igual que en la aritmética en módulo dos (es decir, equivale a la operación XOR) y el producto es el aritmético habitual, que al ser realizado sólo con los elementos 0 y 1, equivale a la operación AND.

Como curiosidad, existen cuerpos de Galois de orden primo y de orden la potencia de un primo, pero no de ningún otro. Por ejemplo, no existen cuerpos de Galois de orden 6, porque 6 no es primo ni es la potencia de un primo, pero sí de orden 9, porque 9 = 3² y 3 es primo. Para construir un cuerpo de Galois de orden primo p (como el caso p = 2 recién nombrado), se puede usar la aritmética en módulo p. Sin embargo, la adición y la multiplicación módulo pi para i > 1 no forman un cuerpo, por lo que hay que emplear otro sistema para construirlo.

Para crear cuerpos de Galois de orden 2w, se usan como elementos polinomios de grado w con coeficientes en GF(2), es decir, los coeficientes pueden ser sólo 0 y 1, y al operar entre ellos se aplican las reglas de GF(2). La suma es entonces suma de polinomios, pero al ser los coeficientes de GF(2), al sumar o restar dos coeficientes se hace en realidad un XOR entre ellos. El producto es... bueno, el producto es un dolor de muelas. Para calcularlo hay que hallar el producto de los polinomios, módulo un polinomio especial llamado generador, característico del cuerpo con el que estamos trabajando.

Para traducir los polinomios a números enteros, se guarda el coeficiente del término de grado i en el bit i del número. Así, por ejemplo, el polinomio x2+1 en GF(24) se almacena en un entero con el valor binario 0101, es decir, el valor decimal 5. Para sumar (o restar) dos elementos del cuerpo, basta entonces con hacer el XOR bit a bit. En la práctica, para multiplicar de forma eficiente se usan una tabla de logaritmos y una de antilogaritmos, es decir, a*b=antilog(log(a)+log(b)), considerando a=0 y b=0 como casos especiales. La división se traduce igualmente como el antilogaritmo de la resta de logaritmos. (Esa suma y esa resta concretas se hacen en aritmética convencional módulo 2w-1).

Lo más interesante es que las propiedades de las matrices con las que aquí trabajamos se mantienen cuando operamos en GF(2w), con tal de que 2w sea mayor que el número de filas y de columnas de la matriz. Esto nos permite trabajar todo el tiempo con enteros mientras las operaciones las realicemos en el cuerpo. Obsérvese que el RAID 3 equivale a trabajar con w=1 y construir la matriz A de la siguiente forma: las n primeras filas forman una submatriz identidad, y la fila n+1 consiste en todo unos. Aunque en ese caso no se cumple que 2w > n+m, sí que se cumple que todas las submatrices de orden n de A son invertibles, que es el requisito principal.

El programa

Con todos esos datos ya se puede realizar un programa como el que yo requería. He escogido w=8 porque, además de trabajar con bytes, que es rápido, existe el problema de que el algoritmo de Gauss-Jordan es de orden O(n³) y el almacenamiento de las matrices requiere O(n²). Si hubiera escogido w=16, la matriz requeriría un mínimo de 8 Gb de memoria para almacenarla completa y el algoritmo de Gauss-Jordan prefiero no pensar cuántas operaciones requeriría. Con w=8, la matriz necesita a lo sumo 64 Kb y en el peor caso Gauss-Jordan requiere del orden de 8 millones de operaciones. El inconveniente es que no se puede dividir el archivo en más de 255 partes en total, pero era una solución de compromiso necesaria. De todas formas, para mis propósitos es suficiente. El polinomio generador que he usado es el que viene en el tutorial para w=8, o sea, x8 + x4 + x3 + x2 + 1.

He probado a dividir un archivo de 4,4 Gb en 255 partes de las cuales 8 se podían perder. Lo he reconstruido tras quitar ocho partes elegidas a dedo. El resultado coincidía con el original tal y como se esperaba.

Aún no he decidido la licencia (probablemente será GPL2 o GPL3), y está pendiente de una limpieza para hacerlo presentable. Si hay demanda, tal vez lo haga público. ¿Alguien más quiere el programa?

Referencias

[1] http://www.cs.utk.edu/~plank/plank/papers/CS-96-332.pdf
[2] http://www.cs.utk.edu/~plank/plank/papers/CS-03-504.pdf

2 comments:

Julio said...

El algoritmo de Gauss-Jordan diagonaliza la matriz entera. Hay otros métodos, como la descomposición LU, aunque no estoy seguro de si hay alguna ganancia de tiempo en el total de operaciones.

No se si lo habrás probado. El paquete GSL trae algoritmos para tratar con matrices, y hacer estas descomposiciones.

pgimeno said...

Estuve buscando un algoritmo adecuado para invertir la matriz. Vi el LU, pero preferí implementar Gauss-Jordan, que era muy sencillote, y comprobar si era lo suficientemente rápido. Resulta que prácticamente no se nota, de modo que lo dejé así.