Algoritmos de representación gráfica para el conjunto de Mandelbrot

En ambos algoritmos de tiempo de escape, optimizados y no optimizados, las coordenadas x e y de cada punto se utilizan como valores iniciales en un cálculo iterativo (que se describe en detalle más adelante). El resultado de cada iteración se utiliza como valor inicial para la siguiente. Durante cada iteración, se comprueba si los valores han alcanzado una condición crítica de "escape" o "salida". Si se alcanza dicha condición, el cálculo se detiene, se dibuja el píxel y se examina el siguiente punto x , y . Para algunos valores iniciales, el escape se produce rápidamente, tras un número reducido de iteraciones. Para valores iniciales muy cercanos al conjunto, pero que no pertenecen a él, puede requerir cientos o miles de iteraciones. Para valores dentro del conjunto de Mandelbrot, el escape nunca se producirá. El programador o usuario debe elegir cuántas iteraciones —o cuánta "profundidad"— desea analizar. Cuanto mayor sea el número máximo de iteraciones, mayor será el detalle y la sutileza que se aprecien en la imagen final, pero mayor será el tiempo necesario para calcular la imagen fractal.

Las condiciones de escape pueden ser simples o complejas. Dado que ningún número complejo con una parte real o imaginaria mayor que 2 puede formar parte del conjunto, una solución común es escapar cuando cualquiera de los coeficientes supera 2. Un método computacionalmente más complejo que detecta los escapes antes consiste en calcular la distancia desde el origen utilizando el teorema de Pitágoras , es decir, determinar el valor absoluto , o módulo , del número complejo. Si este valor supera 2, o equivalentemente, cuando la suma de los cuadrados de las partes real e imaginaria supera 4, el punto ha alcanzado el escape. Variantes de representación computacionalmente más intensivas incluyen el método Buddhabrot , que encuentra los puntos de escape y grafica sus coordenadas iteradas.

El color de cada punto representa la rapidez con la que los valores alcanzaron el punto de escape. Generalmente, se usa el negro para mostrar los valores que no logran escapar antes del límite de iteraciones, y se utilizan colores cada vez más brillantes para los puntos que sí lo logran. Esto proporciona una representación visual de cuántos ciclos fueron necesarios para alcanzar la condición de escape.

Para generar dicha imagen, la región del plano complejo que estamos considerando se subdivide en un cierto número de píxeles . Para colorear cualquiera de esos píxeles, hagamos lo siguiente:${\displaystyle c}$sea ​​el punto medio de ese píxel. Ahora iteramos el punto crítico 0 bajo${\displaystyle P_{c}}$, comprobando en cada paso si el punto de la órbita tiene un módulo mayor que 2. Cuando esto ocurre, sabemos que${\displaystyle c}$Si el parámetro no pertenece al conjunto de Mandelbrot, lo coloreamos según el número de iteraciones utilizadas para determinarlo. En caso contrario, continuamos iterando hasta un número fijo de pasos, tras lo cual decidimos que nuestro parámetro "probablemente" pertenece al conjunto de Mandelbrot, o al menos está muy cerca de él, y coloreamos el píxel de negro.

En pseudocódigo , este algoritmo se vería así. El algoritmo no utiliza números complejos y simula manualmente operaciones con números complejos usando dos números reales, para aquellos que no disponen de un tipo de dato complejo . El programa puede simplificarse si el lenguaje de programación incluye operaciones con tipos de datos complejos.

para cada píxel (Px, Py) en la pantalla hacer x0 := coordenada x escalada del píxel (escalada para que se encuentre en la escala X de Mandelbrot (-2.00, 0.47)) y0 := coordenada y escalada del píxel (escalada para que se encuentre en la escala Y de Mandelbrot (-1,12, 1,12)) x := 0.0 y := 0.0 iteración := 0 max_iteration := 1000 mientras (x*x + y*y ≤ 2*2 Y iteración < iteración_máxima) hacer xtemp := x*x - y*y + x0 y := 2*x*y + y0 x := xtemp iteración := iteración + 1 color := paleta[iteración] trazar(Px, Py, color)

Aquí, relacionando el pseudocódigo con${\displaystyle c}$,${\displaystyle z}$y${\displaystyle P_{c}}$:

• ${\displaystyle z=x+iy\ }$
• ${\displaystyle z^{2}=x^{2}+2ixy}$-${\displaystyle y^{2}\ }$
• ${\displaystyle c=x_{0}+iy_{0}\ }$

y así, como puede verse en el pseudocódigo en el cálculo de x e y :

• ${\displaystyle x=\matop {\mathrm {Re} } (z^{2}+c)=x^{2}-y^{2}+x_{0}}$y${\displaystyle y=\mathop {\mathrm {Im} } (z^{2}+c)=2xy+y_{0}.\ }$

Para obtener imágenes a color del conjunto, se puede asignar un color a cada valor del número de iteraciones ejecutadas mediante diversas funciones (lineales, exponenciales, etc.). Una forma práctica, sin ralentizar los cálculos, es usar el número de iteraciones ejecutadas como entrada para una paleta inicializada al inicio. Si la tabla de colores tiene, por ejemplo, 500 entradas, la selección de color es n  mod  500, donde n es el número de iteraciones.

El código de la sección anterior utiliza un bucle while interno no optimizado para mayor claridad. En la versión no optimizada, se deben realizar cinco multiplicaciones por iteración. Para reducir el número de multiplicaciones, se puede utilizar el siguiente código para el bucle while interno:

x2:= 0 y2:= 0 w:= 0 mientras (x2 + y2 ≤ 4 y iteración < iteración_máxima) hacer x:= x2 - y2 + x0 y:= w - x2 - y2 + y0 x2:= x * x y2:= y * y w:= (x + y) * (x + y) iteración:= iteración + 1

El código anterior funciona mediante una simplificación algebraica de la multiplicación compleja:

${\displaystyle {\begin{aligned}(iy+x)^{2}&=-y^{2}+2iyx+x^{2}\\&=x^{2}-y^{2}+2iyx\end{aligned}}}$

Utilizando la identidad anterior, el número de multiplicaciones se puede reducir a tres en lugar de cinco.

El bucle while interno anterior se puede optimizar aún más expandiendo w a

${\displaystyle w=x^{2}+2xy+y^{2}}$

Sustituyendo w en${\displaystyle y=wx^{2}-y^{2}+y_{0}}$rendimientos ${\displaystyle y=2xy+y_{0}}$y por lo tanto ya no es necesario calcular w .

El pseudocódigo optimizado para lo anterior es:

x:= 0 y:= 0 x2:= 0 y2:= 0 mientras (x2 + y2 ≤ 4 y iteración < iteración_máxima) hacer x2:= x * x y2:= y * y y:= 2 * x * y + y0 x:= x2 - y2 + x0 iteración:= iteración + 1

Tenga en cuenta que en el pseudocódigo anterior,${\displaystyle 2xy}$parece aumentar el número de multiplicaciones en 1, pero dado que 2 es el multiplicador, el código se puede optimizar mediante${\displaystyle (x+x)y}$.

Un método de coloración más complejo implica el uso de un histograma que asocia cada píxel con el número máximo de iteraciones de dicho píxel antes de la salida/salida. Este método distribuirá los colores de manera uniforme en la misma área general y, lo que es importante, es independiente del número máximo de iteraciones elegido. ^{[ 1 ]}

Este algoritmo tiene cuatro pasadas. La primera pasada consiste en calcular el número de iteraciones asociadas a cada píxel (pero sin que se muestre ningún píxel). Estos se almacenan en un array IterationCounts[x][y], donde xy yson las coordenadas x e y de dicho píxel en la pantalla, respectivamente.

El primer paso de la segunda pasada es crear una matriz NumIterationsPerPixel[n], donde el tamaño de la matriz nes el número máximo de iteraciones. A continuación, se debe iterar sobre la matriz de pares de píxeles y recuento de iteraciones IterationCounts[x][y], y recuperar el recuento de iteraciones guardado de cada píxel, i , por ejemplo . Después de recuperar el recuento de iteraciones ii = IterationCounts[x][y] de cada píxel , es necesario indexar la matriz en i e incrementar el valor indexado (que inicialmente es cero), por ejemplo .NumIterationsPerPixelNumIterationsPerPixel[i] = NumIterationsPerPixel[i] + 1

para (x = 0; x < ancho; x++) hacer para (y = 0; y < alto; y++) hacer i:= IterationCounts[x][y] NumIteracionesPorPíxel[i]++

La tercera pasada recorre el array NumIterationsPerPixel y suma todos los valores almacenados, guardándolos en total . El índice del array representa el número de píxeles que alcanzaron ese número de iteraciones antes de la cancelación.

total: = 0 para (i = 0; i < max_iterations; i++) hacer total += NumIterationsPerPixel[i]

Después de esto, comienza la cuarta pasada y se indexan todos los valores en el array IterationCounts. Para cada recuento de iteración i , asociado a cada píxel, se suma el recuento a la suma global de todos los recuentos de iteración desde 1 hasta i en el array NumIterationsPerPixel . Este valor se normaliza dividiendo la suma por el valor total calculado previamente.

hue[][]:= 0.0 para (x = 0; x < ancho; x++) hacer para (y = 0; y < alto; y++) hacer iteración:= RecuentoDeIteraciones[x][y] para (i = 0; i <= iteración; i++) hacer hue[x][y] += NumIterationsPerPixel[i] / total /* Debe ser una división de punto flotante. */ ... color = paleta[tono[m, n]] ...

Finalmente, el valor calculado se utiliza, por ejemplo, como índice para una paleta de colores.

Este método se puede combinar con el método de coloración suave que se describe a continuación para obtener imágenes más estéticas.

Esta imagen se generó con el algoritmo de tiempo de escape. Hay "bandas" de color muy evidentes.

Esta imagen se generó con el algoritmo de conteo de iteraciones normalizado. Las bandas de color se reemplazaron por un degradado suave. Además, los colores presentan el mismo patrón que se observaría si se utilizara el algoritmo de tiempo de escape.

El algoritmo de tiempo de escape es popular por su simplicidad. Sin embargo, crea bandas de color que, como un tipo de aliasing , pueden restar valor estético a una imagen. Esto se puede mejorar utilizando un algoritmo conocido como "conteo de iteraciones normalizado" ^{[ 2 ] [ 3 ]} , que proporciona una transición suave de colores entre iteraciones. El algoritmo asocia un número real${\displaystyle \nu }$con cada valor de z mediante la conexión del número de iteración con la función potencial . Esta función viene dada por

${\displaystyle \phi (z)=\lim _{n\to \infty }(\log |z_{n}|/P^{n}),}$

donde z _n es el valor después de n iteraciones y P es la potencia a la que se eleva z en la ecuación del conjunto de Mandelbrot ( z _{n +1}  = z _n ^P + c , P es generalmente 2).   

Si elegimos un radio de rescate grande N (por ejemplo, 10 ¹⁰⁰ ), tenemos que

${\displaystyle \log |z_{n}|/P^{n}=\log(N)/P^{\nu (z)}}$

para algún número real${\displaystyle \nu (z)}$y esto es

${\displaystyle \nu (z)=n-\log _{P}(\log |z_{n}|/\log(N)),}$

y como n es el primer número de iteración tal que | z _n | > N , el número que restamos de n está en el intervalo [0,  1) .

Para la coloración debemos tener una escala cíclica de colores (construida matemáticamente, por ejemplo) que contenga H colores numerados del 0 al H  −  1 ( H =  500, por ejemplo). Multiplicamos el número real${\displaystyle \nu (z)}$Mediante un número real fijo que determina la densidad de los colores en la imagen, tome la parte entera de este número módulo H y úsela para buscar el color correspondiente en la tabla de colores. 

Por ejemplo, modificando el pseudocódigo anterior y utilizando también el concepto de interpolación lineal se obtendría:

para cada píxel (Px, Py) en la pantalla hacer x0:= coordenada x escalada del píxel (escalada para que se encuentre en la escala X de Mandelbrot (-2,5, 1)) y0:= coordenada y escalada del píxel (escalada para que se encuentre en la escala Y de Mandelbrot (-1, 1)) x:= 0.0 y:= 0.0 iteración:= 0 max_iteration:= 1000 // Aquí se elige N = 2^8 como un radio de rescate razonable.mientras x*x + y*y ≤ (1 << 16) y iteración < max_iteration hacer xtemp:= x*x - y*y + x0 y:= 2*x*y + y0 x:= xtemp iteración:= iteración + 1 // Se utiliza para evitar problemas de punto flotante con puntos dentro del conjunto. Si iteración < max_iteration entonces // se elimina la raíz cuadrada del término interno utilizando reglas de simplificación de logaritmos. log_zn:= log(x*x + y*y) / 2 nu:= iniciar sesión(log_zn / iniciar sesión(2)) / iniciar sesión(2) // Reorganizando la función potencial. // Dividiendo log_zn por log(2) en lugar de log(N = 1<<8) // porque queremos que toda la paleta abarque desde el // centro hasta el radio 2, NO nuestro radio de rescate. iteración:= iteración + 1 - nu color1:= paleta[piso(iteración)] color2:= paleta[floor(iteración) + 1] // iteración % 1 = parte fraccionaria de la iteración. color:= linear_interpolate(color1, color2, iteración % 1) trazar(Px, Py, color)

Normalmente, al renderizar un fractal, la distribución de los colores de una paleta determinada a lo largo del fractal es estática. Si deseamos desplazar la ubicación con respecto al borde del fractal o ajustar la paleta para que se reproduzca de forma cíclica, podemos realizar algunos cambios sencillos al calcular el número final de iteraciones antes de utilizarlo para seleccionar un elemento de nuestra paleta.

Una vez que hayamos obtenido el número de iteraciones, podemos hacer que el rango de colores sea no lineal. Elevar un valor normalizado al rango [0,1] a una potencia n transforma un rango lineal en un rango exponencial, lo que en nuestro caso puede modificar la apariencia de los colores en el exterior del fractal y permitirnos resaltar otros colores o acercar toda la paleta al borde.

${\displaystyle v=((\mathbf {i} /max_{i})^{\mathbf {S} }\mathbf {N} )^{1.5}{\bmod {\mathbf {N} }}}$

donde i es nuestro contador de iteraciones después del rescate, max_i es nuestro límite de iteraciones, S es el exponente al que elevamos las iteraciones y N es el número de elementos en nuestra paleta. Esto escala el contador de iteraciones de forma no lineal y escala la paleta para que se repita aproximadamente en proporción al zoom.

Luego podemos introducir v en cualquier algoritmo que deseemos para generar un color.

Una opción que podríamos considerar es evitar por completo el uso de paletas o la mezcla de colores. De hecho, existen varios métodos que podemos aprovechar para generar colores suaves y uniformes, creando el color directamente en la superficie.

f(v) se refiere a la función de transferencia sRGB.

Un método ingenuo para generar un color de esta manera es escalar directamente v a 255 y pasarlo a RGB como tal.

rgb = [v * 255, v * 255, v * 255]

Un inconveniente de esto es que RGB no es lineal debido a la gamma; en su lugar, considere sRGB lineal. La conversión de RGB a sRGB utiliza una función de compresión inversa en los canales. Esto hace que la gamma sea lineal y nos permite sumar correctamente los colores para el muestreo.

srgb = [v * 255, v * 255, v * 255]

La coloración HSV se puede lograr mapeando el número de iteraciones de [0,max_iter) a [0,360), elevándolo a la potencia de 1.5 y luego módulo 360. $${\displaystyle h=((i/max_{i})360)^{1.5}{\bmod {3}}60}$$ Entonces podemos simplemente tomar el recuento de iteraciones mapeado exponencialmente en el valor y devolver

hsv = [powf((i / max) * 360, 1.5) % 360, 100, (i / max) * 100]

Este método también se aplica a HSL, excepto que en este caso aplicamos una saturación del 50%.

hsl = [powf((i / max) * 360, 1.5) % 360, 50, (i / max) * 100]

Uno de los métodos de coloración más uniformes desde el punto de vista perceptual consiste en pasar el número de iteraciones procesadas al canal LCH. Si utilizamos el método cíclico y de mapeo exponencial descrito anteriormente, podemos aplicar el resultado a los canales Luma y Chroma. También podemos aplicar un mapeo exponencial al número de iteraciones, escalarlo a 360 y pasar este valor módulo 360 al canal de tono.

${\textstyle {\begin{array}{lcl}x&\in &\mathbb {Q+} \\s_{i}&=&(i/max_{i})^{\mathbf {x} }\\v&=&1.0-cos^{2}(\pi s_{i})\\L&=&75-(75v)\\C&=&28+(75-75v)\\H&=&(360s_{i})^{1.5}{\bmod {3}}60\end{array}}}$

Un problema que queremos evitar son los colores fuera de la gama cromática. Esto se puede lograr con un pequeño truco basado en el cambio de los colores dentro de la gama cromática en relación con la luminancia y la crominancia. A medida que aumentamos la luminancia, debemos disminuir la crominancia para mantenernos dentro de la gama cromática.

s = iteraciones/max_i; v = 1,0 - powf(cos(pi * s), 2,0); LCH = [75 - (75 * v), 28 + (75 - (75 * v)), powf(360 * s, 1.5) % 360];

Se puede calcular la distancia desde el punto c (en el exterior o en el interior ) hasta el punto más cercano en el límite del conjunto de Mandelbrot. ^{[ 4 ] [ 5 ]}

Estimación de distancia exterior

La prueba de la conexidad del conjunto de Mandelbrot proporciona, de hecho, una fórmula para el mapa uniformizador del complemento de${\displaystyle M}$(y la derivada de este mapa). Mediante el teorema del cuarto de Koebe , se puede estimar la distancia entre el punto medio de nuestro píxel y el conjunto de Mandelbrot con una precisión de hasta un factor de 4.

En otras palabras, siempre que el número máximo de iteraciones sea suficientemente alto, se obtiene una imagen del conjunto de Mandelbrot con las siguientes propiedades:

1. Cada píxel que contiene un punto del conjunto de Mandelbrot está coloreado de negro.
2. Cada píxel de color negro está cerca del conjunto de Mandelbrot.

El límite superior b para la estimación de distancia de un píxel c (un número complejo) del conjunto de Mandelbrot viene dado por ^{[ 6 ] [ 7 ] [ 8 ]}

${\displaystyle b=\lim _{n\to \infty }{\frac {2\cdot |{P_{c}^{n}(c)|\cdot \ln |{P_{c}^{n}(c)}}|}{|{\frac {\partial }{\partial {c}}}P_{c}^{n}(c)|}},}$

dónde

• ${\displaystyle P_{c}(z)\,}$representa polinomio cuadrático complejo
• ${\displaystyle P_{c}^{n}(c)}$representa n iteraciones de${\displaystyle P_{c}(z)\to z}$o${\displaystyle z^{2}+c\to z}$, comenzando con${\displaystyle z=c}$: ${\displaystyle P_{c}^{0}(c)=c}$, ${\displaystyle P_{c}^{n+1}(c)=P_{c}^{n}(c)^{2}+c}$;
• ${\displaystyle {\frac {\partial }{\partial {c}}}P_{c}^{n}(c)}$es el derivado de ${\displaystyle P_{c}^{n}(c)}$con respecto a c. Esta derivada se puede encontrar comenzando con ${\displaystyle {\frac {\partial }{\partial {c}}}P_{c}^{0}(c)=1}$y luego${\displaystyle {\frac {\partial }{\partial {c}}}P_{c}^{n+1}(c)=2\cdot {}P_{c}^{n}(c)\cdot {\frac {\partial }{\partial {c}}}P_{c}^{n}(c)+1}$Esto se puede verificar fácilmente utilizando la regla de la cadena para la derivada .

La idea detrás de esta fórmula es simple: Cuando las líneas equipotenciales para la función potencial${\displaystyle \phi (z)}$yacen cerca, el número${\displaystyle |\phi '(z)|}$es grande y, a la inversa, por lo tanto, las líneas equipotenciales para la función${\displaystyle \phi (z)/|\phi '(z)|}$debería acostarse aproximadamente a la misma hora.

Desde el punto de vista de un matemático, esta fórmula solo funciona en el límite donde n tiende a infinito, pero se pueden encontrar estimaciones muy razonables con solo unas pocas iteraciones adicionales después de que finaliza el bucle principal.

Una vez que se encuentra b , por el teorema de Koebe 1/4, sabemos que no hay ningún punto del conjunto de Mandelbrot cuya distancia a c sea menor que b/4 .

La estimación de distancia puede utilizarse para dibujar el límite del conjunto de Mandelbrot; véase el artículo « Conjunto de Julia ». En este método, los píxeles suficientemente cercanos a M se dibujan con un color diferente. Esto crea dibujos donde los finos «filamentos» del conjunto de Mandelbrot se aprecian con claridad. Esta técnica se utiliza con éxito en las imágenes en blanco y negro de conjuntos de Mandelbrot que aparecen en los libros «La belleza de los fractales » ^{[ 9 ]} y «La ciencia de las imágenes fractales» ^{[ 10 ] .}

Aquí hay una imagen de ejemplo en blanco y negro generada utilizando estimaciones de distancia:

La estimación de distancias también se puede utilizar para generar imágenes 3D de conjuntos de Mandelbrot y Julia.

Estimación de la distancia interior

También es posible estimar la distancia de un punto límite periódico (es decir, hiperbólico ) al límite del conjunto de Mandelbrot. El límite superior b para la estimación de la distancia viene dado por ^{[ 4 ].}

${\displaystyle b={\frac {1-\left|{{\frac {\partial }{\partial {z}}}P_{c}^{p}(z_{0})}\right|^{2}}{\left|{{\frac {\partial }{\partial {c}}}{\frac {\partial }{\partial {z}}}P_{c}^{p}(z_{0})+{\frac {\partial }{\partial {z}}}{\frac {\partial }{\partial {z}}}P_{c}^{p}(z_{0}){\frac {{\frac {\partial }{\partial {c}}}P_{c}^{p}(z_{0})}{1-{\frac {\partial }{\partial {z}}}P_{c}^{p}(z_{0})}}}\right|}},}$

dónde

• ${\displaystyle p}$es el período,
• ${\displaystyle c}$es el punto a estimar,
• ${\displaystyle P_{c}(z)}$es el polinomio cuadrático complejo${\displaystyle P_{c}(z)=z^{2}+c}$
• ${\displaystyle P_{c}^{p}(z_{0})}$es el${\displaystyle p}$iteración de pliegues de${\displaystyle P_{c}(z)\to z}$, comenzando con${\displaystyle P_{c}^{0}(z)=z_{0}}$
• ${\displaystyle z_{0}}$es alguno de los${\displaystyle p}$puntos que hacen el atractor de las iteraciones de${\displaystyle P_{c}(z)\to z}$comenzando con ${\displaystyle P_{c}^{0}(z)=c}$;${\displaystyle z_{0}}$Satisface${\displaystyle z_{0}=P_{c}^{p}(z_{0})}$,
• ${\displaystyle {\frac {\partial }{\partial {c}}}{\frac {\partial }{\partial {z}}}P_{c}^{p}(z_{0})}$, ${\displaystyle {\frac {\partial }{\partial {z}}}{\frac {\partial }{\partial {z}}}P_{c}^{p}(z_{0})}$, ${\displaystyle {\frac {\partial }{\partial {c}}}P_{c}^{p}(z_{0})}$y ${\displaystyle {\frac {\partial }{\partial {z}}}P_{c}^{p}(z_{0})}$son varios derivados de ${\displaystyle P_{c}^{p}(z)}$, evaluado en${\displaystyle z_{0}}$.

De forma análoga al caso exterior, una vez que se encuentra b , sabemos que todos los puntos que se encuentran a una distancia de b /4 de c están dentro del conjunto de Mandelbrot.

Hay dos problemas prácticos con la estimación de la distancia interior: primero, necesitamos encontrar${\displaystyle z_{0}}$precisamente, y segundo, necesitamos encontrar${\displaystyle p}$precisamente. El problema con${\displaystyle z_{0}}$es que la convergencia a${\displaystyle z_{0}}$mediante iteración${\displaystyle P_{c}(z)}$requiere, teóricamente, un número infinito de operaciones. El problema con cualquier dado${\displaystyle p}$En ocasiones, debido a errores de redondeo, un período se identifica erróneamente como un múltiplo entero del período real (por ejemplo, se detecta un período de 86, cuando el período real es solo 43 = 86/2). En tal caso, la distancia se sobreestima, es decir, el radio reportado podría contener puntos fuera del conjunto de Mandelbrot.

Una forma de mejorar los cálculos es determinar de antemano si el punto dado se encuentra dentro de la cardioide o en el bulbo de periodo 2. Antes de pasar el valor complejo por el algoritmo de tiempo de escape, compruebe primero que:

${\displaystyle p={\sqrt {\left(x-{\frac {1}{4}}\right)^{2}+y^{2}}}}$,

${\displaystyle x\leq p-2p^{2}+{\frac {1}{4}}}$,

${\displaystyle (x+1)^{2}+y^{2}\leq {\frac {1}{16}}}$,

donde x representa el valor real del punto e y el valor imaginario. Las dos primeras ecuaciones determinan que el punto se encuentra dentro de la cardioide, y la última, dentro del bulbo de periodo 2.

La prueba cardioide se puede realizar de forma equivalente sin la raíz cuadrada:

${\displaystyle q=\left(x-{\frac {1}{4}}\right)^{2}+y^{2},}$

${\displaystyle q\left(q+\left(x-{\frac {1}{4}}\right)\right)\leq {\frac {1}{4}}y^{2}.}$

Los brotes de tercer orden y superiores no tienen pruebas equivalentes, porque no son perfectamente circulares. ^{[ 11 ]} Sin embargo, es posible determinar si los puntos están dentro de círculos inscritos dentro de estos bulbos de orden superior, lo que impide que muchos, aunque no todos, los puntos del bulbo sean iterados.

Para evitar tener que realizar un gran número de iteraciones para los puntos dentro del conjunto, se puede realizar una comprobación de periodicidad, que verifica si un punto alcanzado al iterar un píxel ya se había alcanzado antes. Si es así, el píxel no puede divergir y debe estar en el conjunto. La comprobación de periodicidad implica una compensación, ya que la necesidad de recordar puntos consume instrucciones de gestión de datos y memoria, pero ahorra instrucciones de cálculo. Sin embargo, al comparar con una sola iteración anterior se pueden detectar muchos periodos con una mínima sobrecarga de rendimiento. Por ejemplo, dentro del bucle while del pseudocódigo anterior, realice las siguientes modificaciones:

xold := 0 yold := 0 período := 0 mientras (x*x + y*y ≤ 2*2 y iteración < max_iteration) hacer xtemp := x*x - y*y + x0 y := 2*x*y + y0 x := xtemp iteración := iteración + 1 si x ≈ xold e y ≈ yold entonces iteración := iteración_máxima /* Establecer al máximo para la representación gráfica del color */ break /* Estamos dentro del conjunto de Mandelbrot, salimos del bucle while */ período:= período + 1 si el período > 20 entonces período := 0 xold := x yold := y

El código anterior almacena un nuevo valor de x e y en cada vigésima iteración, por lo que puede detectar periodos de hasta 20 puntos de longitud.

Dado que el conjunto de Mandelbrot es completo , ^{[ 12 ]} cualquier punto encerrado por una figura cerrada cuyos bordes se encuentren completamente dentro del conjunto de Mandelbrot debe pertenecer a este. El trazado de bordes funciona siguiendo las lemniscatas de los distintos niveles de iteración (bandas de color) alrededor del conjunto y luego rellenando toda la banda a la vez. Esto también proporciona un aumento de velocidad, ya que ahora se puede omitir un gran número de puntos. ^{[ 13 ]}

En la animación mostrada, los puntos fuera del conjunto se colorean con un algoritmo de tiempo de escape de 1000 iteraciones. Trazar el borde del conjunto y rellenarlo, en lugar de iterar los puntos interiores, reduce el número total de iteraciones en un 93,16 %. Con un límite de iteraciones mayor, el beneficio sería aún mayor.

La comprobación de rectángulos es un método más antiguo y sencillo para representar el conjunto de Mandelbrot. La idea básica es que si cada píxel del borde de un rectángulo tiene la misma cantidad de iteraciones, entonces el rectángulo se puede rellenar de forma segura con esa cantidad de iteraciones. Existen varias variantes de este método; sin embargo, todas son más lentas que el método de trazado de bordes, ya que terminan calculando más píxeles. Una variante calcula solo los píxeles de las esquinas de cada rectángulo; sin embargo, esto provoca imágenes dañadas con más frecuencia que calcular el borde completo, por lo que solo funciona razonablemente bien si se utilizan cajas pequeñas de aproximadamente 6x6 píxeles y no se recurre a partir de cajas más grandes ( método Fractint ).

El método más sencillo para comprobar rectángulos consiste en verificar los bordes de rectángulos de igual tamaño, que se asemejan a un patrón de cuadrícula. (Algoritmo de Mariani.) ^{[ 14 ]}

Una variante más rápida y ligeramente más avanzada consiste en calcular primero un recuadro más grande, por ejemplo, de 25x25 píxeles. Si todo el borde del recuadro tiene el mismo color, simplemente se rellena con ese color. Si no, se divide el recuadro en cuatro recuadros de 13x13 píxeles, reutilizando los píxeles ya calculados como borde exterior y compartiendo los píxeles interiores en forma de cruz entre los recuadros interiores. De nuevo, se rellenan los recuadros que tienen un solo color de borde. Y los que no lo tienen, se dividen ahora en cuatro recuadros de 7x7 píxeles. Y luego los que "fallan" en recuadros de 4x4. (Algoritmo de Mariani-Silver).

Aún más rápido es dividir las cajas por la mitad en lugar de en cuatro. En ese caso, lo óptimo sería usar cajas con una relación de aspecto de 1,4:1 , de modo que se puedan dividir como se doblan las hojas A3 en hojas A4 y A5 (el método DIN).

Al igual que con el trazado de bordes, la comprobación de rectángulos solo funciona en áreas con un color discreto. Pero incluso si el área exterior utiliza una coloración suave/continua, la comprobación de rectángulos seguirá acelerando el costoso área interior del conjunto de Mandelbrot. A menos que el área interior también utilice algún método de coloración suave, por ejemplo, la estimación de distancia interior .

La simetría horizontal del conjunto de Mandelbrot permite omitir partes del proceso de renderizado cuando el eje real está presente en la imagen final. Sin embargo, independientemente de la parte que se refleje, se renderizará la misma cantidad de puntos.

Los conjuntos de Julia presentan simetría respecto al origen. Esto significa que los cuadrantes 1 y 3 son simétricos, al igual que los cuadrantes 2 y 4. Para que los conjuntos de Mandelbrot y de Julia presenten simetría, es necesario tratarla de forma diferente para cada tipo de gráfico.

La representación de conjuntos de Mandelbrot y Julia mediante el método de tiempo de escape se presta extraordinariamente bien al procesamiento paralelo. En máquinas multinúcleo, el área a representar se puede dividir en una serie de áreas rectangulares que luego se pueden proporcionar como un conjunto de tareas para ser renderizadas por un grupo de hilos de renderizado. Este es un problema de computación fácilmente paralelizable ^{[ 15 ]} . (Cabe destacar que se obtiene la mayor aceleración excluyendo primero las áreas simétricas del gráfico y luego dividiendo las regiones únicas restantes en áreas rectangulares). ^{[ 16 ]}

Aquí hay un breve video que muestra el conjunto de Mandelbrot siendo renderizado usando multihilo y simetría, pero sin seguimiento de límites:

Por último, aquí tenéis un vídeo que muestra la misma imagen del conjunto de Mandelbrot renderizada utilizando multihilo, simetría y seguimiento de límites:

Las imágenes muy ampliadas requieren más de los 64-128 bits de precisión estándar que proporcionan la mayoría de las unidades de punto flotante de hardware , lo que obliga a los renderizadores a utilizar bibliotecas matemáticas lentas de "BigNum" o de " precisión arbitraria " para realizar los cálculos. Sin embargo, esto se puede acelerar mediante la explotación de la teoría de perturbaciones . Dado

${\displaystyle z_{n+1}=z_{n}^{2}+c}$

como la iteración, y un pequeño épsilon y delta, es el caso que

${\displaystyle (z_{n}+\epsilon )^{2}+(c+\delta )=z_{n}^{2}+2z_{n}\epsilon +\epsilon ^{2}+c+\delta ,}$

o

${\displaystyle =z_{n+1}+2z_{n}\epsilon +\epsilon ^{2}+\delta ,}$

entonces si uno define

${\displaystyle \epsilon _{n+1}=2z_{n}\epsilon _{n}+\epsilon _{n}^{2}+\delta ,}$

Se puede calcular un único punto (por ejemplo, el centro de una imagen) usando aritmética de alta precisión ( z ), dando una órbita de referencia , y luego calcular muchos puntos alrededor de él en términos de varios desplazamientos iniciales delta más la iteración anterior para épsilon, donde épsilon-cero se establece en 0. Para la mayoría de las iteraciones, épsilon no necesita más de 16 cifras significativas, y en consecuencia se puede usar punto flotante de hardware para obtener una imagen bastante precisa. ^{[ 17 ]} A menudo habrá algunas áreas donde las órbitas de los puntos divergen lo suficiente de la órbita de referencia que se necesita precisión adicional en esos puntos, o bien se necesitan órbitas de referencia locales adicionales calculadas con alta precisión. Al medir la distancia de la órbita entre el punto de referencia y el punto calculado con baja precisión, se puede detectar que no es posible calcular el punto correctamente, y se puede detener el cálculo. Estos puntos incorrectos se pueden recalcular más tarde, por ejemplo, a partir de otro punto de referencia más cercano.

Además, es posible aproximar los valores iniciales para los puntos de baja precisión con una serie de Taylor truncada , lo que a menudo permite omitir una cantidad significativa de iteraciones. ^{[ 18 ]} Existen renderizadores que implementan estas técnicas y ofrecen mejoras de velocidad para imágenes altamente ampliadas de aproximadamente dos órdenes de magnitud. ^{[ 19 ] [ 20 ]}

Una explicación alternativa de lo anterior:

Para el punto central del disco${\displaystyle c}$y sus iteraciones${\displaystyle z_{n}}$y un punto arbitrario en el disco${\displaystyle c+\delta }$y sus iteraciones${\displaystyle z'_{n}}$Es posible definir la siguiente relación iterativa:

${\displaystyle z'_{n}=z_{n}+\epsilon _{n}}$

Con${\displaystyle \epsilon _{1}=\delta }$. Iteraciones sucesivas de${\displaystyle \epsilon _{n}}$Se puede encontrar utilizando lo siguiente:

${\displaystyle z'_{n+1}={z'_{n}}^{2}+(c+\delta )}$

${\displaystyle z'_{n+1}=(z_{n}+\epsilon _{n})^{2}+c+\delta }$

${\displaystyle z'_{n+1}={z_{n}}^{2}+c+2z_{n}\epsilon _{n}+{\epsilon _{n}}^{2}+\delta }$

${\displaystyle z'_{n+1}=z_{n+1}+2z_{n}\epsilon _{n}+{\epsilon _{n}}^{2}+\delta }$

Ahora, partiendo de la definición original:

${\displaystyle z'_{n+1}=z_{n+1}+\epsilon _{n+1}}$,

Resulta que:

${\displaystyle \epsilon _{n+1}=2z_{n}\epsilon _{n}+{\epsilon _{n}}^{2}+\delta }$

Como la relación iterativa relaciona un punto arbitrario con el punto central mediante un cambio muy pequeño${\displaystyle \delta }$, entonces la mayoría de las iteraciones de${\displaystyle \epsilon _{n}}$También son pequeños y se pueden calcular utilizando hardware de punto flotante.

Sin embargo, para cada punto arbitrario del disco es posible calcular un valor para un valor dado.${\displaystyle \epsilon _{n}}$sin tener que iterar a través de la secuencia desde${\displaystyle \epsilon _{0}}$, expresando${\displaystyle \epsilon _{n}}$como una serie de poder de${\displaystyle \delta }$.

${\displaystyle \epsilon _{n}=A_{n}\delta +B_{n}\delta ^{2}+C_{n}\delta ^{3}+\dotsc }$

Con${\displaystyle A_{1}=1,B_{1}=0,C_{1}=0,\dotsc }$.

Ahora, dada la ecuación de iteración de${\displaystyle \epsilon }$, es posible calcular los coeficientes de la serie de potencias para cada${\displaystyle \epsilon _{n}}$:

${\displaystyle \epsilon _{n+1}=2z_{n}\epsilon _{n}+{\epsilon _{n}}^{2}+\delta }$

${\displaystyle \epsilon _{n+1}=2z_{n}(A_{n}\delta +B_{n}\delta ^{2}+C_{n}\delta ^{3}+\dotsc )+(A_{n}\delta +B_{n}\delta ^{2}+C_{n}\delta ^{3}+\dotsc )^{2}+\delta }$

${\displaystyle \epsilon _{n+1}=(2z_{n}A_{n}+1)\delta +(2z_{n}B_{n}+{A_{n}}^{2})\delta ^{2}+(2z_{n}C_{n}+2A_{n}B_{n})\delta ^{3}+\dotsc }$

Por lo tanto, se deduce que:

${\displaystyle A_{n+1}=2z_{n}A_{n}+1}$

${\displaystyle B_{n+1}=2z_{n}B_{n}+{A_{n}}^{2}}$

${\displaystyle C_{n+1}=2z_{n}C_{n}+2A_{n}B_{n}}$

${\displaystyle \vdots }$

Los coeficientes de la serie de potencias se pueden calcular como una serie iterativa utilizando únicamente los valores de las iteraciones del punto central.${\displaystyle z}$y no cambian para ningún punto arbitrario en el disco. Si${\displaystyle \delta }$es muy pequeño,${\displaystyle \epsilon _{n}}$debería poder calcularse con suficiente precisión utilizando solo unos pocos términos de la serie de potencias. Como los contornos de escape de Mandelbrot son "continuos" sobre el plano complejo, si se ha calculado el tiempo de escape de un punto, entonces el tiempo de escape de los vecinos de ese punto debería ser similar. La interpolación de los puntos vecinos debería proporcionar una buena estimación de dónde empezar en el${\displaystyle \epsilon _{n}}$serie.

Además, la interpolación separada de los puntos del eje real y del eje imaginario debería proporcionar un límite superior e inferior para el punto que se está calculando. Si ambos resultados son iguales (es decir, ambos escapan o no escapan), entonces la diferencia${\displaystyle \Delta n}$se puede utilizar para recusarse hasta que se pueda establecer un límite superior e inferior. Si se puede utilizar hardware de punto flotante para iterar el${\displaystyle \epsilon }$serie, entonces existe una relación entre cuántas iteraciones se pueden lograr en el tiempo que lleva usar el software BigNum para calcular un valor dado${\displaystyle \epsilon _{n}}$Si la diferencia entre los límites es mayor que el número de iteraciones, es posible realizar una búsqueda binaria utilizando el software BigNum, reduciendo sucesivamente a la mitad la diferencia hasta que resulte más eficiente en términos de tiempo encontrar el valor de escape utilizando hardware de punto flotante.

• 
  Teoría de perturbaciones con nivel de zoom 2.4e-57
• 
  Teoría de perturbaciones con nivel de zoom 1,53e-62
• 
  Teoría de perturbaciones con nivel de zoom 2,84e-82
• 
  Teoría de perturbaciones con nivel de zoom 7.17e-36
• 
  Teoría de perturbaciones con nivel de zoom 4.31e-51
• 
  Teoría de perturbaciones con nivel de zoom 1.01e-55