Inicio

fractal.c

/*
Compilar con la opcion c99:
gcc -std=c99 -c filename.c -o filename.o

* Enlazar con la opcion de librerías matematicas y con allegro:
gcc -lm -o filename filename.o `allegro-config --libs`

* En Anjuta, ejecutar en terminal allegro-config --libs, y su salida son las librerías que
tenemos que pasar al gcc. Así que vamos a Opciones-> Comandos -> Construir Archivo y le ponemos
ahí todas las opciones.
*/

/*
fractal.c Este fractal es la portada del libro: Análisis Numérico Las matemáticas del Cálculo científico, de Kincaid, y que por desgracia está
descatalogado. En él se explica cómo programarlo y aquí también claro, ahí va:
Dibuja los fractales de las cuencas de atracción de las raíces cúbicas de la unidad. El
polinomio es p(z)=z**3-1. Tiene una raíz real (1) y dos complejas: (-1+sqrt(3)*I)/2 y
(-1-sqrt(3)*I)/2 donde I es sqrt(-1).
Mediante el método de Newton realizamos 40 iteraciones y si queda a menos de 0.25 de una raíz,
consideramos que converge a ella, y le asignamos un color.
Hay puntos que no convergen a ninguna raíz y se llaman conjunto de Julia (1918). Estos puntos
forman un fractal. Así que el límite entre cuencas es donde se muestra la geometría fractal.
*/
#include "allegro.h"
#include <math.h>
#include <complex.h>

int main(void)
{
PALETTE palette;
int i,j,n;
int amarillo,verde, rojo, azul;
double x, y ,a,id,jd;
_Complex double z,r2, r3;
/*Iniciar modo grafico*/ if (allegro_init() != 0)
return 1;
install_keyboard();
set_color_depth(16);
if (set_gfx_mode(GFX_AUTODETECT, 1024, 768, 0, 0) !=0) {
if (set_gfx_mode(GFX_AUTODETECT,640, 480, 0, 0) != 0) { allegro_message("Error setting graphics mode\n%s\n",allegro_error);
return 1;
}
}
set_palette(black_palette);
set_palette(palette);

amarillo=makecol(255,255,50);
verde=makecol(0,255,0);
azul=makecol(0,0,255);
rojo=makecol(255,0,0);
/* Código del fractal*/ r2=(-1+sqrt(3)*I)/2;
r3=(-1-sqrt(3)*I)/2;
/* Para una monitor de 1024 por 768 píxeles*/ for (i=0; i<=1024;i++) {
for (j=0; j<=768; j++) {
/* Le restamos 512 y 384 para que el punto (0,0) del plano complejo quede en la mitad de la
pantalla */
id=i-512;
jd=j-384;
z=(id+jd*I); for (n=1;n<=40;n++) z=z-(z/3)+(1/(3*z*z));
x=cabs(z-r2);

/* Aunque haga los cálculos para todas las raíces no merece la pena poner else en los if
porque tarda lo mismo.*/
if (x<=0.25) putpixel(screen,i,j, rojo);
x=cabs(z-r3);
if (x<=0.25) putpixel(screen,i,j, amarillo);
x=cabs(z-1);
if (x<=0.25) putpixel(screen,i,j, azul);
}
}
do {} while (!keypressed());
return 0;
}
END_OF_MAIN();

Inicio

Valid XHTML 1.0 Strict