Recocido Simulado

Simulated Annealing - Threshold algorithms (threshold = umbral)

El algoritmo de recocido simulado (Simulated Annealing Algorithm - SAA) pertenece a una clase de Algoritmos de busqueda local (Local Search Algorithms – LSA) comunmente llamada Algoritmos de Umbral (Threshold Algorithm - TA). Hay dos razones por las cuales los TA resultan interesantes dentro de los LSA:

1) parecen andar bien en una amplia gama de problemas reales (practicos)

2) algunos TA, como el SAA, tienen caracteristicas que permiten hacer un analisis de la convergencia.

Esqueleto de un TA

Sea (S,c) una instancia de un problema de optimizacion combinatoria, donde:

• S es el conjunto de soluciones factibles

• c es la funcion costo (a valores reales positivos)

El problema es hallar un i en S que minimize c.

Para implementar un TA son necesarios ademas:

• Una funcion entorno N de S en partes de S.

• Una sucesion tk (los llamados threshold)

La manera de elegir los tk y el criterio de aceptacion de una nueva solucion definen 3 tipos de TA:

Dado i en S en la iteracion k

Genero j en N(i)

Utilizo los valores c(j) – c(i) y tk para decidir aceptar o no la solucion j

Local Search Improvement (mejora continua): tk = 0 (para todo k)

Si c(j) – c(i) < tk = 0 entonces acepto j

Threshold Accepting (umbral de aceptacion): se fija la sucesion tk tal que tk  tk+1, tk > 0, y tk tiende a 0 cuando k tiende a infinito.

Si c(j) – c(i) < tk entonces acepto j

En este caso, todas las soluciones que disminuyen el costo son aceptadas, y las que incrementan el costo son aceptadas en forma limitada. A medida que aumenta k (progresa el algoritmo) solo se aceptan incrementos pequeños, hasta que eventualmente solo se aceptan mejoras.

Simulated Annealing (recocido simulado): los tk se toman como en el threshold accepting pero el criterio de aceptacion es probabilistico

Si c(j) – c(i)  0 entonces acepto j

Si c(j) – c(i) > 0 entonces acepto j con probabilidad exp [(c(i) – c(j)) / tk]

(en la iteracion k se genera un numero al azar r y se acepta j si r < exp [(c(i) – c(j)) / tk])

En este caso, cada vecino de una solucion tiene una probabilidad positiva de reemplazar a la solucion actual. Los tk se eligen de forma tal que a medida que avanzan las iteraciones, aceptar soluciones con grandes incrementos en el costo es menos probable (pero sigue existiendo una probabilidad positiva de aceptarlos).

Analogia Fisica

El metodo del recocido se utiliza en la industria para obtener materiales mas resistentes, o mas cristalinos, en general, para mejorar las cualidades de un material.

El proceso consiste en “derretir” el material (calentarlo a muy alta temperatura). En esa situacion, los atomos adquieren una distribucion “azarosa” dentro de la estructura del material y la energia del sistema es maxima. Luego se hace descender la temperatura muy lentamente por etapas, dejando que en cada una de esas etapas los atomos queden en equilibrio (es decir, que los atomos alcancen una configuracion optima para esa temperatura). Al final del proceso, los atomos forman una estructura cristalina altamente regular, el material alcanza asi una maxima resistencia y la energia del sistema es minima.

Experimentalmente se comprueba que si la temperatura se hace descender bruscamente o no se espera suficiente tiempo en cada etapa, al final la estructura del material no es la optima.

La rama de la Fisica llamada Mecanica Estadistica se encargo de desarrollar una serie de metodos para estudiar el comportamiento de grandes cantidades de atomos de un sistema. Debido a que en promedio, en un sistema hay 1023 atomos por cm3, solamente puede estudiarse el comportamiento mas probable del sistema en equilibrio a una dada temperatura. La experimentacion mostro que los atomos de un sistema en un proceso de recocido se comportan según el factor de probabilidad de Boltzman. En 1953 Metropolis modelo el proceso de recocido: en cada paso del algoritmo se le da al atomo un desplazamiento azaroso y se mide el cambio de energia E. Si E  0 se acepta el desplazamiento. Si E > 0, se acepta el desplazamiento con probabilidad exp (-E / T.K), donde T es la temperatura del sistema y K es la constante de Boltzman.

The simulated annealing method (El metodo del Recocido Simulado)

Sea S el conjunto de soluciones posibles del sistema (a las que identificamos con los diferentes “estados del sistema”) y tenemos dada una funcion costo sobre los elementos de S (a la que identificamos con la “energia del sistema”). Se quiere encontrar un elemento en S que minimize la funcion costo (analogamente, se trata de encontrar un estado en el cual la energia del sistema sea minima). Asumimos que los estados del sistema tienen la funcion de distribucion de probabilidad de Boltzman, i.e., la probabilidad de que el sistema se encuentre en el estado j es:

P (j) = (1/Zt) exp [- c(j) / t]

Donde Zt =  exp [- c(i) / t] (suma sobre todos los elementos i de S)

t es la temperatura del sistema y c(i) es el costo de la solucion i.

Sea S* el subconjunto de S de las soluciones que minimizan c (globalmente, es decir soluciones optimas del problema). Para t suficientemente chico:

exp [- c(j*) / t] >> exp [- c(j) / t]

para todo j != j*

Entonces:

P (j) = exp [- c(j) / t] / {S* exp [- c(j*) / t]} = 0 si j != j*

1/S* si j = j*

(Esto se obtiene de tomar t tendiendo a 0)

Por lo tanto:  P (j*) = 1 (suma sobre todos los j* en S*).

Simulated Annealing Algorithm (Algoritmo de Recocido Simulado)

Version monotona

El algoritmo se divide en etapas. A cada etapa le corresponde una temperatura menor que la que tenia la etapa anterior (a esto hace referencia la monotonia: despues de cada etapa la temperatura baja, se enfria el sistema). Por lo tanto hace falta un criterio de cambio de la temperatura (“cuanto tiempo” se espera en cada etapa para dar lugar a que el sistema alcance su “equilibrio termico”).

Datos iniciales y parametros a ser definidos para poder inicializar el algoritmo:

Temperatural inicial (T0)

La temperatura inicial T0 debe ser una temperatura que permita casi (o todo) movimiento, es decir que la probabilidad de pasar del estado i al j (en N(i)) sea muy alta, sin importar la diferencia c(j) – c(i). Esto es que el sistema tenga un alto grado de libertad. En problemas como TSP, donde el input son los nodos de un grafo y las soluciones posibles son distintas formas de recorrer estos nodos, pude tomarse T0 proporcional a la raiz cuadrada de la cantidad de nodos. En general se toma un valor T0 que se cree suficientemente alto y se observa la primera etapa para verificar que el sistema tenga un grado de libertad y en funcion de esta observacion se ajusta T0.

Solucion inicial (i0)

En todas las versiones, el sistema debe ser “derretido” antes de implementar el algoritmo. Esto es que la solucion factible incial que llamamos i0 deberia ser una solucion tomada al azar del conjunto de soluciones factibles. En algunos problemas esto puede hacerse utilizando pseudo-random numbers provistos por una maquina (ejemplo de esto es el bandwidth problem). Pero en muchos casos ya es problemático encontrar una solucion, por lo que es imposible tomar una al azar. En estos casos se implementa un algoritmo “greedy” tipo local search para buscar una solucion factible y se toma esta como i0 (ejemplo de esto es el TSP).

Funcion entorno (N)

Factor de enfriamiento

Tnext = T (factor de enfriamiento geometrico,  < 1, muy cercano a 1)

Tnext = 1 / (1 + T) (donde  es un real positivo cercano a cero)

Criterio de cambio de la temperatura

Se usan dos parametros: K = cantidad de iteraciones que estamos dispuestos a hacer en cada etapa (equivalente a la cantidad de tiempo que vamos a esperar a que el sistema alcance su equilibrio termico para una dada temperatura T); A = cantidad de aceptaciones que se permiten hacer en cada etapa.

A medida que T disminuye se supone que al sistema le resulta mas dificil alcanzar un equilibrio porque es mas dificultoso el movimiento, entonces hay que esperar mas tiempo, esto se traduce en aumentar K.

Parametro de aumento de K (, se usan valores alrededor de 1,05)

Criterio de STOP

a) Lundy and Mees: si el algoritmo se detiene cuando T <  / [ln (#S – 1)/]

Donde #S es el cardinal del conjunto de soluciones (debe tenerse un metodo de estimar este valor).

Entonces, si i es la solucion que da el algoritmo e i* en un optimo global,

P(|c(i) – c(i*)| < ) = 

b) En general se utiliza un parametro de congelamiento (frozen: FRZN). Como a medida que disminuye la temperatura, aumenta el parametro K y A permanece constante, la proporcion A/K se hace pequeña. Asumimos que si A/K < FRZN el sistema esta congelado (la cantidad de aceptaciones respecto de la cantidad de iteraciones es muy chica, esto da la idea de que cambiar de configuracion es muy dificil).

El algoritmo:

1. i = i0

2. T = T0

3. K = K0

...