PERCOLACIÓN
PERCOLACIÓN
Julian David Toro Medina Departamento de sistemas e industrial Facultad de ingeniería Universidad Nacional de Colombia
El problema hace referencia al efecto en el que el agua se filtra desde la superficie hasta las capas subterráneas de la tierra, en donde conforma ríos que fluyen hacia el mar. Haremos, pues, cuenta de que tenemos una malla que representa la “roca sólida”, a medida que vamos abriendo celdas o “retirando rocas” el agua avanza hasta llegar a las capas inferiores. Solo cuando haya una conexión
entre la primera y la última fila, diremos que está percolado. Para desarrollar un algoritmo que represente este caso, se diseñó una “malla”, que en realidad es una matriz llena de ceros, en la que se van abriendo espacios, o llenando de unos, para formar un camino que conecte la primera y la última fila, pero estos espacios se abren generando valores aleatorios que determinan la posición en la matriz en la que se va a cambiar el cero por un uno. Al final, después de que se consigue que el sistema esté percolado, se procede a determinar el umbral de percolación, que consiste en la razón o proporción de celdas “abiertas”, o unos, que hay en la malla, o matriz. Al final se realizan un número determinado de experimentos, que varían entre cien y quinientos, para determinar el valor promedio del umbral de percolación y, desde ahí, hallar el intervalo de confianza para este umbral, con un 95% de confianza. EL CÓDI GO
Universidad Nacional de Colombia
Página 1
PERCOLACIÓN
Para desarrollar el ejercicio es necesario tener en cuenta la malla con que se va a trabajar, cómo se asignan los lugares para “abrir” espacios y las herramientas necesarias para desarrollarlo. Es bueno
revisar una matriz 8x8 que facilite la asimilación de estos conceptos. Inicialmente se crea una matriz A con las dimensiones deseadas, que en un principio, está llena de ceros. Tenemos:
[
]
Ahora generaremos dos valores, uno i y el otro j, que serán los que se usarán para identificar la posición en que se cambiará el valor de un cero por un uno, es decir, abrir un espacio en la malla. Estos valores i y j son generados aleatoriamente. i=aleatorio(entre 1 y 8);
i=3
j=aleatorio(entre 1 y 8);
j=6
[
]
Hay que tener en cuenta evitar que el programa genere aleatorios que hagan referencia a un espacio en el que previamente ya hay un uno, ya que la maquina realiza una actividad mayor cambiando celdas que ya están abiertas, además, si después de esto se verifica si está percolado, se pierde tiempo comprobando una y otra vez la misma matriz, ya que no se le ha hecho ningún cambio. Después de “abrir y abrir” espacios, asignar unos, se procede a verificar si la malla está percolada,
si hay un camino de unos que me conecte la primera y la última fila de la matriz. Tenemos:
[
]
Esta matriz está percolada, podemos verificarlo siguiendo el camino señalado. Universidad Nacional de Colombia
Página 2
PERCOLACIÓN
Ahora se incluirán unos conceptos que funcionan bastante bien para entender cómo solucionarlo. Antes de continuar, esbozaremos la idea del algoritmo a aplicarse. Primero que todo hay más de una vía “ posible”, vías que en un momento dado se cierran y no con tinúan. Puede ser necesario, además de tener en cuenta las celdas abiertas, revisar las “conexiones” existentes, ya que si se guía la
solución a través de estas conexiones se facilita el trabajo, debido a que no solo revisa unos, sino que comprueba las conexiones que hay aumentando la posibilidad, de que al seguirlas, se consiga el camino correcto. Si llegado el caso las conexiones no nos dan una pista de hacia dónde continuar, se puede revisar el camino de unos cercano para decidir qué rumbo tomar dentro de nuestro recorrido. Por último, es importante tener en cuenta cuándo hay caminos ciegos, que al ser tentadores, lleven a nuestro recorrido hacia lugares sin salida, además hay que impedir que este viaje sobre sus propios pasos y genere bucles en el recorrido. La primera herramienta a utilizarse será una matriz “auxiliar de conexión”, la cual indicará las
conexiones existentes dentro de nuestra matriz original. Esta funciona, o se genera, recorriendo la matriz A y poniendo unos en la posición en donde encuentre que hay dos unos en la misma columna y en filas sucesivas. A esta matriz se le denominará como B. Entonces:
[
]
[
]
Si se observan las posiciones de los unos en la matriz A que han sido resaltados, estos son sucesivos y están en la misma columna, por lo tanto en la matriz B son representados mediante un uno, porque están conectados. Esta matriz es más pequeña, tiene una fila menos, ya que la fila superior no está conectada más arriba con otra, es decir, se recorre la celda colocando los unos iniciando desde abajo, por lo tanto, la primera, al no tener filas encima, no tendrá valores en la matriz auxiliar. Ya partiendo desde aquí se procede a avanzar y recorrer la malla para encontrar si está percolada, en el camino se irá entendiendo la utilidad de esta matriz auxiliar B. 1. Primero se recorre la matriz B buscando unos, ya que si hay alguna conexión se podrá avanzar y encontrar un posible camino, en lugar de recorrer la matriz A buscando unos y
Universidad Nacional de Colombia
Página 3
PERCOLACIÓN
luego verificando si se puede continuar. En este caso se prefirió recorrer la matriz de abajo hacia arriba, por lo tanto la fila a verificarse será la última.
a.
[
]
En este caso hay tres caminos posibles, en las columnas 1, 3 y 6. 2. Después de encontrar estos valores, se procede a iniciar una cadena. Desde aquí se busca verificar si la fila superior en la matriz B, es decir, la penúltima, también tiene conexión. En este caso, ninguna de las tres columnas cumple este requisito. Cabe aclarar que primero nos centramos en la columna 1, y a medida que tengamos intentos fallidos de encontrar un camino, pasamos a la siguiente columna haciendo la verificación inicial. 3. Si llegado el caso no se consigue que haya una conexión inmediata, se busca si hay un valor de uno en la fila superior en la matriz A, permitiendo por ahí continuar abriéndonos camino en la matriz. Si este caso se da, se procede a verificar si hay formas de irse hacia la derecha o hacia la izquierda, en la fila superior. El problema es decidir cuál camino tomar si tenemos ambas posibilidades, he aquí la magia de la matriz auxiliar, ya que se puede calcular la densidad de conexiones que haya en cualquiera de las direcciones haciendo una sumatoria de los valores que haya en la matriz B de acuerdo a un parámetro de qué tan extensa debe ser esta malla de densidad, y por medio de este valor decidir qué camino tomar, ya que el camino con mayor densidad de conexiones probablemente será el que permita que haya una cadena y sea el que determine si está percolado.
Universidad Nacional de Colombia
[
]
[
] Página 4
PERCOLACIÓN
En este caso las columnas 1 y 3 tienen unos en la fila superior, pero al fin y al cabo, no hay otro camino que seguir, por lo tanto no continúa con la cadena. En cambio, la columna 6 si tiene posibilidades de continuar, ya que tiene un uno en la fila superior y hay dos posibilidades de camino por tomar, ahora queda decidir cuál de los dos hay que tomar. Para un valor de dos para el cálculo de densidad, es decir, dos filas más arriba y dos columnas a la derecha o izquierda según corresponda, se obtienen valores de densidad iguales a uno y dos para izquierda y derecha, respectivamente. Es de aquí que se toma la decisión de continuar hacia la derecha.
[
]
[
]
Ahora cabe revisar qué pasaría si el valor de este cálculo de densidad es tres. Se obtienen valores iguales a tres para ambas direcciones, ¿qué se debe hacer? En este caso se calcula un nuevo valor de densidad, pero en la matriz A. Se obtiene valores iguales a cuatro y cinco para izquierda y derecha, respectivamente, por lo que se escoge continuar hacia la derecha.
Se pueden tomar otros criterios, como lo es escoger irse hacia el lado contrario de donde esté más cerca por ejemplo, para así decidir qué camino tomar cuando tenemos estas dos posibilidades.
Universidad Nacional de Colombia
Página 5
PERCOLACIÓN
4. Si se sigue aplicando este algoritmo se llega a la primera fila encontrando, entonces, que la matriz está percolada. Pero hay que tener en cuenta muchos detalles al momento de ponerlo en funcionamiento. Primero que todo, cuando se asignen los valores que se van a utilizar para recorrer la matriz hay que cuidarse de que estos nunca se salgan de sus límites para evitar que genere errores. Por ejemplo, en nuestra matriz de ocho, al generar el valor de densidad llegamos a tomar en cuenta, como este puede ser de tres, posiciones como B(4,9) por lo que se generará un error. También hay que tener en cuenta que este se podría llegar a devolver a posiciones que ya previamente había recorrido, por lo que hay que evitar que lo haga y que, más bien, siga en la dirección que llevaba, es decir, si iba hacia la derecha que continúe hacia la derecha. 5. Ahora un caso particular del que hay que prestar mucha atención.
[
]
Si fijamos nuestra atención en los unos azules podemos ver que estos representan un “hueco”, es un camino ciego al que el algoritmo llevará y que no representa la solución.
Para esto se verifica, cuando se revisan las conexiones, si más arriba hay algo con la forma
( )
y se le advierte que debe devolverse y continuar en la dirección que llevaba.
CÁLCUL O DEL UM BRAL DE PERCOLA CIÓN
Conociendo cómo se determina si está percolado, se procede a calcular el umbral de percolación. Para esto se halla la prop orción de celdas “abiertas” que hay en la malla, o más fácil aún, sumar los valores que hay en la matriz y dividirla entre sus dimensiones al cuadrado. En nuestro ejemplo:
[ Universidad Nacional de Colombia
] Página 6
PERCOLACIÓN
Esta suma sería igual a tres en la primera fila, a seis en la segunda, a cinco en la tercera y así sucesivamente hasta tener un total de treinta y ocho. Como es una matriz 8x8, esta tiene un total de sesenta y cuatro posiciones. Por lo tanto:
Esta matriz tiene un umbral de percolación de 0,59375. Hay que recordar que este se calcula siempre después de encontrar que, a medida que vamos creando aleatorios y “abriendo” celdas, la matriz o malla quedan percoladas. Al realizar este experimento muchas veces se puede determinar la media de estos umbrales y su desviación estándar para así poder establecer un intervalo de confianza para estos umbrales. APLICANDO ESTO A MATRICES DE 20x20, 40x40, 50x50 Y 100x100 Los respectivos valores para estas mallas, teniendo en cuenta que el intervalo de confianza es del 95%, son:
Matriz 20x20 Número de experimentos: 500 Umbral de percolación: 0,6154 Intervalo de confianza: 0,6109 – 0,6199
Matriz 40x40 Número de experimentos: 500 Umbral de percolación: 0,6503 Intervalo de confianza: 0,6470 – 0,6537
Matriz 50x50 Número de experimentos: 100 Umbral de percolación: 0,6559 Intervalo de confianza: 0,6497 – 0,6622
Matriz 100x100 Número de experimentos: 10 Umbral de percolación: 0,6860 Intervalo de confianza: 0,6711 – 0,7009
Universidad Nacional de Colombia
Página 7
PERCOLACIÓN
CÓDIGO EN MATLAB® clear all clc close all %Determinar el tamaño de la matriz s=input('Ingrese el tamaño de la red: '); %Determinar el tamaño de la matriz de densidad rr=input('Ingrese el tamaño de la red de densidad: '); g=1; %contador de experimento actual exp=input('Ingrese el numero de experimentos: '); h=1; %contador que controla el numero de experimentos datos=zeros(1,exp); for h=1:exp %% Generando nuevo problema %Creando las matrices A=zeros(s,s); B=zeros(s-1,s); NotPer=1; %Variable que controla que no este percolado while NotPer==1 nuevo=1; while nuevo==1 i=randi([1,s],1); j=randi([1,s],1); if A(i,j)~=1 nuevo=0; end end for i=1:s-1 for j=1:s if A(s+1-i,j)==1 && A(s-i,j)==1 B(s-i,j)=1; else B(s-i,j)=0; end end end %%Analiza si esta percolado suma=0; for j=1:s suma=suma+B(s-1,j); end
Universidad Nacional de Colombia
Página 8
PERCOLACIÓN
r=rr; while NotPer==1 && r>=2 m=1; while m<=s && suma>0 if B(s-1,m)==1 %Hay posibilidad de una cadena i=s; j=m; cadena=1; while i>=1 n=2; af=zeros(2,n); bf=zeros(2,n); hueco=0; while cadena==1 af=zeros(2,n); for z=1:n-1 for zz=1:2 af(zz,z)=bf(zz,z); end end af(1,n)=i; af(2,n)=j; if B(i-1,j)==1 if i-2>=1 && B(i-2,j)==1 %Verifica si hay una columna conectada if hueco==1; %Comprueba que hay un hueco, busca salida lateral if i>=1 && j>1 && j=1 if j-l>=1 IB=IB+B(i-k,j-l); IA=IA+A(i-k+1,j-l); end if j+l<=s DB=DB+B(i-k,j+l); DA=DA+A(i-k+1,j+l); end end end end if DB>IB && j1 %Mayor densidad a la izquierda j=j-1;
Universidad Nacional de Colombia
Página 9
PERCOLACIÓN
if af(1,n-1)==i && af(2,n-1)==j j=j+2; end else if DA>IA && j1 %Mayor densidad 1 a la izquierda j=j-1; if af(1,n-1)==i && af(2,n-1)==j j=j+2; end else if s-j>=j %Cercano a la izquierda j=j-1; if af(1,n-1)==i && af(2,n-1)==j j=j+2; end elseif s-j=1 && j==1%Pegado a la izquierda if A(i,j+1)==1 %1 a la derecha j=j+1; else %no hay cadena m=m+1; cadena=0; i=0; end elseif i>=1 && j==s %Pegado a la derecha if A(i,j-1)==1 %1 a la izquierda j=j-1; else %no hay cadena
Universidad Nacional de Colombia
Página 10
PERCOLACIÓN
m=m+1; cadena=0; i=0; end else %no hay cadena if i~=0 m=m+1; cadena=0; i=0; end end end else i=i-1; end elseif A(i-1,j)==1 if i>2 && i1 && A(i-2,j)==0 && A(i-1,j+1)==0 && A(i-1,j-1)==0 devolver=1; hueco=1; while devolver==1 if A(i,j)==1 && A(i,j+1)==0 && A(i,j-1)==0 i=i+1; end if A(i,j)==1 && A(i,j-1)==1 && A(i,j+1)==1 DB=0; IB=0; DA=0; IA=0; for k=1:r for l=1:r if i-1-k>=1 if j-l>=1 IB=IB+B(i-k,j-l); IA=IA+A(i-k+1,j-l); end if j+l<=s DB=DB+B(i-k,j+l); DA=DA+A(i-k+1,j+l); end end end end if DB>IB && j1 %Mayor densidad coneccion a la izquierda j=j-1; devolver=0; hueco=0; else if DA>IA && j
Universidad Nacional de Colombia
Página 11
PERCOLACIÓN
devolver=0; hueco=0; elseif DA1 %Mayor densidad 1 a la izquierda j=j-1; devolver=0; hueco=0; else if s-j>=j %Cercano a la izquierda j=j-1; devolver=0; hueco=0; elseif s-j=s i=0; m=m+1; cadena=0; else i=i+1; j=j+1; end devolver=0; hueco=0; end end end else i=i-1; end else %Busca 1 laterales if i>=1 && j>1 && j
Universidad Nacional de Colombia
Página 12
PERCOLACIÓN
IA=0; for k=1:r for l=1:r if i-1-k>=1 if j-l>=1 IB=IB+B(i-1-k,j-l); IA=IA+A(i-k,j-l); end if j+l<=s DB=DB+B(i-1-k,j+l); DA=DA+A(i-k,j+l); end end end end if DB>IB && j1 %Mayor densidad coneccion a la izquierda i=i-1; j=j-1; else if DA>IA && j1 %Mayor densidad 1 a la izquierda i=i-1; j=j-1; else if s-j>=j %Cercano a la izquierda i=i-1; j=j-1; elseif s-j
Universidad Nacional de Colombia
Página 13
PERCOLACIÓN
if i>=1 && j==1%Pegado a la izquierda if A(i-1,j+1)==1 %1 a la derecha i=i-1; j=j+1; else %no hay cadena m=m+1; cadena=0; i=0; end elseif i>=1 && j==s %Pegado a la derecha if A(i-1,j-1)==1 %1 a la izquierda i=i-1; j=j-1; else %no hay cadena m=m+1; cadena=0; i=0; end else %no hay cadena m=m+1; cadena=0; i=0; end end end else %Busca 1 laterales for a=2:4 if n>a && af(1,n-a)==af(1,n) && af(2,n-a)==af(2,n) cadena=0; i=0; m=m+1; end end if i>=1 && j>1 && j=1 if j-l>=1 IB=IB+B(i-k,j-l); IA=IA+A(i-k+1,j-l); end if j+l<=s DB=DB+B(i-k,j+l); DA=DA+A(i-k+1,j+l); end end end end
Universidad Nacional de Colombia
Página 14
PERCOLACIÓN
if DB>IB && j1 %Mayor densidad a la izquierda j=j-1; if af(1,n-1)==i && af(2,n-1)==j j=j+2; end else if DA>IA && j1 %Mayor densidad 1 a la izquierda j=j-1; if af(1,n-1)==i && af(2,n-1)==j j=j+2; end else if s-j>=j %Cercano a la izquierda j=j-1; if af(1,n-1)==i && af(2,n-1)==j j=j+2; end elseif s-j=1 && j==1%Pegado a la izquierda if A(i,j+1)==1 %1 a la derecha j=j+1; else %no hay cadena m=m+1;
Universidad Nacional de Colombia
Página 15
PERCOLACIÓN
cadena=0; i=0; end elseif i>=1 && j==s %Pegado a la derecha if A(i,j-1)==1 %1 a la izquierda j=j-1; else %no hay cadena m=m+1; cadena=0; i=0; end else %no hay cadena if i~=0 m=m+1; cadena=0; i=0; end end end end if i==1 cadena=0; i=0; NotPer=0; m=s+1; end bf=zeros(2,n); for z=1:n for zz=1:2 bf(zz,z)=af(zz,z); end end n=n+1; end end else m=m+1; end end r=r-1; end end disp 'Percolado'; a=0; for i=1:s for j=1:s a=a+A(i,j); end end disp 'Abiertas:'; disp(a); disp 'Porcentaje:'; disp(a/(s*s));
Universidad Nacional de Colombia
Página 16
PERCOLACIÓN
datos(1,h)=a/(s*s); %almaceno umbral g=g+1; disp 'Esperimento numero:'; disp(g); end disp 'Numero de experimentos:'; disp(exp); disp 'Umbral:'; disp(mean(datos)); Min=(mean(datos))-1.96*(std(datos)/(sqrt(h))); %limite inferior intervalo Max=(mean(datos))+1.96*(std(datos)/(sqrt(h))); %limite superior intervalo disp 'Intervalo de confianza:'; disp(Min); disp(Max);
Universidad Nacional de Colombia
Página 17
