Fractales de Mandelbrot y Julia

Fractales de Mandelbrot y Julia Por Ariel S.

Los fractales más conocidos son los conjuntos de Julia y de Mandelbrot. Es muy fácil generar dichos conjuntos usando números complejos. Para cada número complejo a definamos la función como $f_{c}(z)=z^{2}+c$ a para todo $z\in \U{2102}$ .

El conjunto de Julia asociado al número complejo se define como el conjunto a de los números complejos a tales que la sucesión definida por

está acotada (algunos autores llaman conjunto de Julia a la frontera de ). Es decir

El conjunto de Mandelbrot es el conjunto de los números complejos tales que la sucesión definida por

está acotada. Es decir,

En otros términos, el conjunto de Mandelbrot es el conjunto de los números complejos tales que $0\in J(z)$ .

Para representar la sucesión suele usarse la notación , con el convenio de que . Esta sucesión se llama la órbita de por $f_{c}$ y sus elementos se obtienen iterando sucesivamente la función $f_{c}$ partiendo del valor inicial . La notación $f_{c}^{(n)}$ representa la composición de la función a consigo misma n veces.

Implementación de los algoritmos

Lo que se explicará a continuación es una manera de implementar un algoritmo generador de imagenes fractales a partir de las definiciones del conjnto de Mandelbrot y de Julia antes mencionados. Los códigos fuentes exhibidos están en lenguaje C y fueron testeados en el compilador Borland C++ Builder 6.

Lo primero que podemos hacer es crear algunas funciones que serán utiles a la hora de realizar cálculos con números complejos, para las cuales debemos definir un nuevo tipo de datos, el de números complejos, de la siguiente manera

typedef struct
{double x; double y;}complejo;

Ahora una función para calcular el conjugado de un complejo

complejo Cconj(complejo a)
{    complejo c;
     c.x = a.x;
     c.y = -a.y;
     return(c);}

Para calcular la norma o módulo del complejo

double Cmod(complejo a)
{ double m = a.x*a.x + a.y*a.y;
return(m);
}

Para calcular la suma de dos números complejos

complejo Csuma(complejo a, complejo b)
{    complejo c;
     c.x = a.x + b.x;
     c.y = a.y + b.y;
     return(c);
}

Para calcular la suma de un número complejo y un número real

complejo CsumaR(complejo a, double r)
{    complejo c;
     c.x = a.x + r;
     c.y = a.y;
     return(c);
}

Para calcularel producto de dos números complejos

complejo Cprod(complejo a, complejo b)
{    complejo c;
     c.x = a.x*b.x - a.y*b.y;
     c.y = a.x*b.y + a.y*b.x;
     return(c);
}

Para calcular el producto de un número complejo y un número real

complejo CprodR(complejo a, double r)
{    complejo c;
     c.x = r*a.x;
     c.y = r*a.y;
     return(c);
}

Para calcular el cociente de dos números complejos

complejo Cdiv(complejo a, complejo b)
{    complejo c;
     c = CprodR(Cprod(a, Cconj(b)), 1/Cmod(b));
     return(c);
}

Para calcular el cuadrado de un número complejo

complejo Csqr(complejo a)
{    complejo c;
     c.x = a.x*a.x - a.y*a.y;
     c.y = 2*a.x*a.y;
     return(c);
}

Para calcular la n-ésima potencia de un número complejo, donde n es un número entero

complejo Cnpot(complejo a, int n)
{    complejo c, u;
     bool nNegativo;
     c = a;
     if (n < 0) {
     nNegativo = true; n = -n;
     };
     if (n == 0) {
     c.x = 1; c.y = 0;
     }
     else if (n == 1) {c = a;}
     else {
     for(int k = 1; k < n; k++)
     {c = Cprod(a, c);}
     }
     if (nNegativo == true){
     u.x = 1;
     u.y = 0;
     Cdiv(u, c);
     }
     return(c);
}

Ahora si se puede definir una función que llamamos Julia la cual realiza lo siguiente: dado un complejo c fijo, le damos como parámetro de entrada también otro complejo z, entonces calcula el número de iteraciones que se necesitan para que la sucesión se salga del disco de centro 0 y radio 2. En caso de que se necesiten muchas iteraciones para que ocurra esto se le impone también otra condición de salida del bucle, que se llegue al número maximo de iteraciones el cual en este caso llame itera_max, por ejemplo con un valor de 250. Como se observa en la definición, los valores de z tales que su n-ésima iteracion esta dentro del disco de radio 2, o que hayan superado ese número de iteraciones los consideramos como pertenecientes al conjunto de Julia y se los representará con un punto en su correspondiente posición en el plano complejo iluminado con un color dependiente de su rapidez de salir del bucle (parámetro de retorno k).

int Julia(complejo z, complejo c, int itera_max)
{   complejo w;
    int k = 0;
    w = z;
    float m = Cmod(w);
    while ((m ≤ 2) && (k ≤ itera_max))
    {
        w = Csqr(w);
        w = Csuma(w, c);
        m = Cmod(w);
        k++;
    }
    return(k);
}

Para representar los valores del conjunto de Mandelbrot hacemos algo similar utilizando la correspondiente función recursiva para calcular el tiempo de escape del bucle de cada uno de los valores de plano complejo. Nuevamente las condiciones serán, que el valor de la iteración esté dentro del disco de centro 0 y radio 2 y que no se haya superado el número máximo de iteraciones. Análogamente al caso anterior, le asignamos un color dependiente de la rapidez con que se salió del bucle, el valor k. Y de esa manera repreentamos el conjunto de Mandelbrot en el plano complejo.

int Mandelbrot(complejo z,int itera_max)
{    complejo w;
     int k=0;
     w=z;
     float m=Cmod(w);
     while ((m<=2) && (k<=itera_max))
     {    w=Csqr(w);
          w=Csuma(w,z);
          m=Cmod(w);
          k++;
     }
     return(k);
}

Una vez que tenemos implementados los dos algoritmos que nos devuelven los colores de cada punto correspodiente a los conjuntos que querramos representar, programamos un rutina que realice esto. En este caso he utilizado un componente TPaintBox con nombre g que representará al plano complejo. Sobre él se pasan a calcular los valores de k (color asociado al punto del plano) y se lo representa como un pixel de color número k sobre g, seleccionado de una paleta de colores previamente definida a la cual llamé pal, concretamente pal es un vector de itera_max números enteros los cuales representan el valor hexadecimal de los colores de la paleta.

Los valores de xmin, ymin, xmax e ymax son las coordenadas de los vértices del rectángulo que representará la región del plano complejo para la cual se calcularán los valores de k.

Esta rutina lo que hace es mostrar una imagen del conjunto de Julia para un cirerto valor asociado c y cieras coordenadas de visualización que se les de como parámetros de entrada.

void graf(float xmin, float ymin, float xmax, float ymax, complejo c)
{    int col;
     complejo z;
     float px, py, dx, dy;
     px = xmin;
     py = ymin;
     dx = (xmax - xmin)/tamX;
     dy = (ymax - ymin)/tamY;
     for(int cx = 0; cx ≤ tamX; cx++)
     {    for(int cy = 0; cy ≤ tamY; cy++)
          {
             z.x = px;
             z.y = py;
             col = Julia(z, c, itmax);
             Form1 -> g -> Canvas -> Pixels[cx][cy] = pal[col];
             py = py + dy;
          }
          px = px + dx;
          py = ymin;
     }
}

Si en vez de querer graficar el conjunto de Julia queremos que sea el de Mandelbrot simplemente debemos reemplazar la linea del código correspondiente a la función del color,

col = Mandelbrot(z, itmax);

Al ir moviendo el cursor del mouse sobre el área del PaintBox g, podemos ir obteniendo las coordenadas correspodientes al punto en el plano complejo que representa realmente en la pantalla, esto dependerá del nivel de acercamiento que tenga el conjunto, es decir, de los valores de xmin, ymin, xmax e ymax. Sobre el evento OnClick de PaintBox g podemos obtener las coordenadas de x e y correspondientes al sistema de coordenadas relativo a g, cuyo origen es el vértice superior izquierdo. Haciendo un par de obvios cálculos de conversión podemos obtener las coordenadas reales del punto en el plano complej que representa dicho pixel en la pantalla:

float dx = xmax - xmin, dy = ymax - ymin;
float xx = X/tamX, yy = Y/tamY;
g_posx = xmin + xx*dx;
g_posy = ymin + yy*dy;

donde xmin, ymin, xmax e ymax son las coordenadas de la región del plano complejo que se está visualizando en pantalla, tamX y tamY son las dimensiones (medidas en pixeles) del PaintBox g, muestra area de visualización en la pantalla.

De esta manera podemos programar que cada vez que se haga click sobre una zona de la imagen fractal, se haga un acercamiento o zoom con un cierto nivel de magnificacion, al cual llame nivel_zoom, en global_zoom se guarda el acercamiento total que se realizó (global porque es una variable global) tomando como acercamiento de nivel 1 al disco de centro 0 y radio 2.

global_zoom = global_zoom + nivel_zoom;
xmin = g_posx - 2/global_zoom;
ymin = g_posy - 2/global_zoom;
xmax = g_posx + 2/global_zoom;
ymax = g_posy + 2/global_zoom;

Podemos obtener de esta forma imagenes como las siguientes:

Dejá un comentario, sugerencia u opinión