TEMAS DE PROGRAMACIÓN

Los programas pueden ser implementados en cualquiera de los lenguajes más o menos standard: Matlab, Maple, C++, Fortran, etc. Importante: todo programa tiene que estar acompañado de una "documentación", soporte escrito donde se explique claramente qué se hizo, como funciona cada una de las componentes, gráficos o tablas de los experimentos realizados, conclusiones obtenidas, dificultades encontradas, etc..

El trabajo puede ser hecho tanto individualmente o como de a dos. Hay flexibilidad en los temas: ustedes también pueden proponer o discutir conmigo otros temas de programación que les interese hacer. Ante cualquier duda o inquietud, no duden en contactarme por e-mail, o mejor aún, personalmente hasta principios de Feb/2006.

1. Algoritmo de Strassen para la inversión de matrices densas

Este algoritmo fue explicado en el curso, muestra que la complejidad teórica de la inversión es igual a la complejidad de la multiplicación de matrices, y por lo tanto la inversión puede hacerse en mejor complejidad que el N^3 del algoritmo de eliminación. El objetivo es estudiar este algoritmo desde un punto de vista práctico.

  1. Implementar un algoritmo de multiplicación de matrices. En una primera etapa, se puede usar la multiplicación standard, pero luego hay que implementar alguna multiplicación rapida, por ejemplo el algoritmo clásico de Strassen o bien el del articulo de Kaporin, que tiene complejidad N^2.8  en lugar del N^3 del algoritmo standard. 

  2. Implementar el algoritmo de tipo "divide-and-conquer" de Strassen, explicado en el libro de Pan y en las notas del curso.

  3. Determinar el rango de aplicación práctica, es decir la máxima dimensión N que puede resolver en un tiempo razonable (e.g. 10').

  4. Estudiar su estabilidad numérica (= pérdida de precisión por errores de redondeo). 

  5. Comparar los resultados (velocidad y estabilidad) con los del algoritmo de eliminación (implementarlo basándose en el pseudo-código en las notas del curso (capítulo 4) o bien usar una versión ya implementada).

Para hacer más: la versión presentada no es estable ya que no tiene pivotage. Testear formas de estabilizarlo, por ejemplo usando precondicionadores estructurados, ver [Pan, S 5.6].

Referencias:

  1. Aho, Hopcroft, Ullman, The design and analysis of computer algorithms. Addison-Wesley 1975.

  2. Pan, Structured matrices and polynomials: unified superfast algorithms (Chapter 4). Birkhauser, 2001.

  3. Kaporin, A practical algorithm for matrix multiplication. Num. Linear Alg. Appl. 1999.

  4. Notas del curso (Capítulo 5)

2. Algoritmo rápido unificado de descomposición PLU de matrices estructuradas

El algoritmo rápido presentado en el curso, suponía que los operadores de desplazamiento S,T eran respectivamente diagonal y triangular superior. Esto restringe la aplicabilidad del algoritmo, por ejemplo la estructura tipo Toeplitz no puede ser tratada directamente, sino que hay que reducirla al tipo Vandermonde usando un precondicionador circulante.

En el libro de Pan (y parcialmente en las notas del curso) se explica que esta restricción no es necesaria, y cómo el algoritmo se adapta a operadores quasi-diagonales, lo cual cubre todos las estructuras clásicas en forma unificada. Las ideas son muy similares a las presentadas en el curso, pero cambian los detalles. El objetivo es implementar la versión rápida (complejidad O(\alpha*N^2)), estudiar su velocidad y estabilidad numérica.

  1. Implementar el algoritmo. Usar el algoritmo rápido presentado en el curso, pero usando la actualización de generadores para el caso quasi-diagonal siguiendo [Pan, Ch. 4].

  2. Incorporar pivotage parcial.

  3. Con esta base, programar su aplicación a las cuatro estructuras clásicas (Toeplitz, Hankel, Vandermonde, Cauchy). En cada caso hay que implementar fórmulas de reconstrucción de la primera línea y columna de una matriz a partir de sus generadores, ver [Pan, S 4.4].

  4. Determinar el rango de aplicación práctica, es decir la máxima dimensión N que puede resolver en un tiempo razonable (e.g. 10').

  5. Estudiar su estabilidad numérica (= pérdida de precisión por errores de redondeo). 

  6. Comparar los resultados (velocidad y estabilidad) con el algoritmo de eliminación (implementarlo basándose en el pseudo-código en las notas del curso (capítulo 4) o bien usar una versión ya implementada).

Referencias:

  1. Pan, Structured matrices and polynomials: unified superfast algorithms (Chapter 4). Birkhauser, 2001.

  2. Notas del curso (Capítulo 5)

3. Algoritmo super-rápido de decodificación de códigos de Reed-Solomon

El algoritmo super-rápido para matrices estructuradas explicado en el curso puede aplicarse a la decodificación de códigos de Reed-Solomon de longitud N en O(N*log^3(N)) ops de F_q. Esta implementación requiere un cierto trabajo a nivel de detalles (que hasta donde yo sé, no están escritos en ningún lado) pero es un trabajo accesible y sumamente interesante.

  1. Implementar el algoritmo super-rapido adaptándolo específicamente para la estructura tipo Vandermonde. Lo que se pide es la parte de  álgebra lineal, es decir el cálculo del núcleo de la matriz de Berlekamp-Massey. Se puede obviar la división final de polinomios, para la cual sería necesaria implementar aparte (o conseguir) la división rápida basada en la TFR.

  2. Implementar la aritmética de F_q para q=2^r, o mejor aún, usar alguna implementación ya disponible.

  3. Determinar el rango de aplicación práctica, es decir la máxima dimensión N que puede resolver en un tiempo razonable (e.g. 10').

  4. Déterminar experimentalemente la probabilidad de fracaso del algoritmo (es decir que alguno de los pivots encontrados sea 0).

Referencias:

  1. Olshevsky, Shokrollahi, A displacement approach to decoding algebraic codes. Contemporary Math. 1999.

  2. Pan, Structured matrices and polynomials: unified superfast algorithms (Chapter 5). Birkhauser, 2001.

  3. Notas del curso (Capítulo 5).

4. Multigrilla para la ecuación de Hemholtz

Se trata de resolver numéricamente la ecuación de Hemholtz en el cuadrado:

-u_xx-u_yy+a*u_x=f        para (x,y) en el cuadrado [0,1]^2;
                          u=0        para (x,y) en el borde;

via el multigrilla. Se tratará de determinar el comportamiento del problema con respecto a los distintos parámetros. Se recomienda programar previamente la ecuación de Hemholtz en dimensión 1

  1. Discretizar la ecuación con diferencias finitas en una malla uniforme.

  2. Implementar las distintas componentes del ciclo multigrilla y del multigrilla completo.

  3. Experimente cambiando el parámetro de regularización, la cantidad de pasos de regularización por ciclo, ciclo V y W. Considere en particular la ecuación -e(u_xx+u_yy)+a*u_x=x(1-x)sin(lx) que tiene solución exacta u=(Ax^2+Bx+C)sin(lx) con A=-1/s, B=(1-2aA)/s, (2eA-aB)/s, y s=1/(el^2). Compare los algorithmos para e=0.01, 0.1, 1; a=0.1, 1, 10.

Referencias:

  1. Briggs, Henson, McCormick: A multigrid tutotial. SIAM 2000.


5. Distribución del condicionamiento de matrices aleatórias con estructuras particulares

Generar matrices aleatórias con estructura dada: simétrica, banda, sparse con soporte fijado de antemano, Toeplitz,
Toeplitz-block-Toeplitz, Hankel, Vandermonde, Cauchy.  Se busca obtener experimentalmente la distribución del  condicionamiento, y compararla con los resultados de  Edelman  para matrices densas y coeficientes siguiendo una distribución gaussiana.

  1. Aprender a generar coefficientes siguiendo una distribución dada, aprender a calcular el condicionamiento de una matriz (se puede y se recomienda usar las rutinas ya implementadas en Matlab o en librerías específicas).

  2. Manejar el programa eficientemente para calcular numerosos ejemplos en la mayor dimensión posible. Sería ideal disponer de máquinas para hacer grandes cálculos, no sé (pero puedo averiguarlo) si hay en la UBA ni son disponibles para los estudiantes.

  3. Comparar las distribuciones encontradas con las obtenidas teóricamente por Edelman. En primer lugar verificar experimentalmente los resultados teóricos: el caso denso (real o complejo) con el condicionamiento de Demmel k_D; para luego pasar a las distintas estructuras.

Referencias:

  1. Edelman, Ph.D. Thesis. MIT 1989. 

  2. Schatzman, Les bonnes matrices sont rares, aparecerá en Matapli.