Introducción a dynamics
=======================

El paquete adicional dynamics incluye varias funciones para crear diversas representaciones gráficas de sistemas dinámicos y fractales, y además una implementación del método numérico de Runge-Kutta de cuarto orden, para resolver sistemas de ecuaciones diferenciales.

Para usar las funciones en este paquete será necesario primero que todo cargarlo con load("dynamics"), y las funciones que crean gráficas necesitan que Xmaxima esté instalado.

Funciones y variables para dynamics
===================================

 -- Función: chaosgame ([[X1, Y1]...[XM, YM]], [X0, Y0], B, N, ...opciones...);

     Usa el método llamado juego del caos, para producir fractales: se dibuja un punto inicial (X0, Y0) y luego se elije aleatoriamente uno de los M puntos [X1, Y1]...[XM, YM]. Después se dibuja un nuevo punto que estará en el segmento entre el último punto dibujado y el punto que se acabó de elegir aleatoriamente, a una distancia del punto elegido que será B veces la longitud del segmento. El proceso se repite N veces.

 -- Función: evolution (F, Y0, N,...opciones...);

     Dibuja N+1 puntos en una gráfica bidimensional (serie de tiempo), en que las coordenadas horizontales de los puntos son los números enteros 0, 1, 2, ..., N, y las coordenadas verticales son los valores Y(N) correspondientes, obtenidos a partir de la relación de recurrencia

                  y(n+1) = F(y(n))

     Con valor inicial Y(0) igual a Y0. F deberá ser una expresión que dependa únicamente de la variable Y (y no de N), Y0 deberá ser un número real y N un número entero positivo.

 -- Función: evolution2d ([F, G], [X0, Y0], N, ...opciones...);

     Muestra, en una gráfica bidimensional, los primeros N+1 puntos de la sucesión definida a partir del sistema dinámico discreto con relaciones de recurrencia:

                  x(n+1) = F(x(n), y(n))    y(n+1) = G(x(n), y(n))

     Con valores iniciales X0 y Y0. F y G deben ser dos expresiones que dependan únicamente de X y Y.

 -- Función: ifs ([R1,...,RM],[A1,...,AM], [[X1,Y1]...[XM, YM]], [X0,Y0], N, ...opciones...);

     Usa el método del Sistema de Funciones Iteradas (IFS, en inglés Iterated Function System). Ese método es semejante al método descrito en la función chaosgame, pero en vez de aproximar el último punto al punto elegido aleatoriamente, las dos coordenadas del último punto se multiplican por una matriz 2 por 2 AI correspondiente al punto que fue elegido aleatoriamente.

     La selección aleatoria de uno de los M puntos atractivos puede ser realizada con una función de probabilidad no uniforme, definida con los pesos R1,...,RM. Esos pesos deben ser dados en forma acumulada; por ejemplo, si se quieren 3 puntos con probabilidades 0.2, 0.5 y 0.3, los pesos R1, R2 y R3 podrían ser 2, 7 y 10, o cualquier otro grupo de números que tengan la misma proporción.

 -- Función: julia (X, Y, ...OPCIONES...)

     Crea un fichero gráfico con la representación del conjunto de Julia del número complejo (X + i Y). Los parámetros X y Y deben ser reales. El fichero se crea en el directório actual o en el directório del usuario, usando el formato gráfico XPM. El programa puede demorar varios segundos a ser ejecutado y cuando termina imprime un mensaje con el nombre del fichero creado.

     Se asignan diferentes colores a los puntos que no pertenecen al conjunto de Julia, de acuerdo con el número de iteraciones que demore la secuencia, comenzando en ese punto, a salir fuera del círculo de convergencia con radio 2. El número máximo de iteraciones se define con la opción LEVELS; después de ejecutadas ese número de iteraciones, si la secuencia aun está dentro del círculo de convergencia, el punto será coloreado con el color definido por la opción COLOR.

     Todos los colores usados para los puntos que no pertenecen al conjunto de Julia tendrán los mismos valores de saturación (SATURATION) y valor (VALUE), pero con diferentes ángulos de tonalidad, distribuidos uniformemente en el intervalo entre HUE y (HUE + HUERANGE).

     Se puede dar a la función una secuencia de opciones. La lista de posibles opciones aparece en una sección más al frente.

 -- Función: mandelbrot (OPTIONS)

     Crea un fichero gráfico con la representación del conjunto de Mandelbrot. El fichero se crea en el directório actual o en el directório del usuario, usando el formato gráfico XPM. El programa puede demorar varios segundos a ser ejecutado y cuando termina imprime un mensaje con el nombre del fichero creado.

     Se asignan diferentes colores a los puntos que no pertenecen al conjunto de Mandelbrot, de acuerdo con el número de iteraciones que demore la secuencia generada por ese punto a salir fuera del círculo de convergencia con radio 2. El número máximo de iteraciones se define con la opción LEVELS; después de ejecutadas ese número de iteraciones, si la secuencia aun está dentro del círculo de convergencia, el punto será coloreado con el color definido por la opción COLOR.

     Todos los colores usados para los puntos que no pertenecen al conjunto de Mandelbrot tendrán los mismos valores de saturación (SATURATION) y valor (VALUE), pero con diferentes ángulos de tonalidad, distribuidos uniformemente en el intervalo entre HUE y (HUE + HUERANGE).

     Se puede dar a la función una secuencia de opciones. La lista de posibles opciones aparece en una sección más al frente.

 -- Función: orbits (F, Y0, N1, N2, [X, X0, XF, XSTEP], ...opciones...);

     Dibuja el diagrama de órbitas de una familia de sistemas dinámicos discretos unidimensionales, con un parámetro X; ese tipo de diagrama se usa para mostrar las bifurcaciones de un sistema discreto unidimensional.

     La función F(Y) define una secuencia que comienza con un valor inicial Y0, igual que en el caso de la función evolution, pero en este caso la función también dependerá del parámetro X, el cual tomará valores comprendidos en el intervalo de X0 a XF, con incrementos XSTEP. Cada valor usado para el parámetro X se muestra en el eje horizontal. En el eje vertical se mostrarán N2 valores de la sucesión Y(N1+1),..., Y(N1+N2+1), obtenidos después de dejarla evolucionar durante N1 iteraciones iniciales.

 -- Función: rk (EDO, var, inicial, dominio)
 -- Función: rk ([EDO1,...,EDOm], [v1,...,vm], [inic1,...,inicm], dominio)

     La primera forma se usa para resolver numéricamente una ecuación diferencial ordinaria de primer orden (EDO), y la segunda forma resuelve numéricamente un sistema de M de esas ecuaciones, usando el método de Runge-Kutta de cuarto orden. VAR representa la variable dependiente. EDO debe ser una expresión que dependa únicamente de las variables independiente y dependente, y define la derivada de la variable dependiente en función de la variable independiente.

     La variable independiente se representa con DOMINIO, que debe ser una lista con cuatro elementos, como por ejemplo:

          [t, 0, 10, 0.1]

     el primer elemento de la lista identifica la variable independiente, el segundo y tercer elementos son los valores inicial y final para esa variable, y el último elemento da el valor de los incrementos que deberán ser usados dentro de ese intervalo.

     Si se van a resolver M ecuaciones, deberá haber M variables dependientes V1, V2, ..., VM. Los valores iniciales para esas variables serán INIC1, INIC2, ..., INICM. Continuará existiendo apenas una variable independiente definida por la lista DOMAIN, como en el caso anterior. EDO1, ..., EDOM son las expresiones que definen las derivadas de cada una de las variables dependientes en función de la variable independiente. Las únicas variables que pueden aparecer en cada una de esas expresiones son la variable independiente y cualquiera de las variables dependientes. Es importante que las derivadas EDO1, ..., EDOM sean colocadas en la lista en el mismo orden en que fueron agrupadas las variables dependientes; por ejemplo, el tercer elemento de la lista será interpretado como la derivada de la tercera variable dependiente.

     El programa intenta integrar las ecuaciones desde el valor inicial de la variable independiente, hasta el valor final, usando incrementos fijos. Si en algún paso una de las variables dependientes toma un valor absoluto muy grande, la integración será suspendida en ese punto.  El resultado será una lista con un número de elementos igual al número de iteraciones realizadas. Cada elemento en la lista de resultados es también una lista con M+1 elementos: el valor de la variable independiente, seguido de los valores de las variables dependientes correspondientes a ese punto.

 -- Función: staircase (F, Y0, N, ...opciones...);

     Dibuja un diagrama de escalera (o diagrama de red) para la sucesión definida por la ecuación de recurrencia

                  y(n+1) = F(y(n))

     La interpretación y valores permitidos de los parámetros de entrada es la misma que para la función evolution. Un diagrama de escalera consiste en una gráfica de la función F(Y), junto con la recta G(Y) = Y. Se comienza por dibujar un segmento vertical desde el punto (Y0, Y0) en la recta, hasta el punto de intersección con la función F. En seguida, desde ese punto se dibuja un segmento horizontal hasta el punto de intersección con la recta, (Y1, Y1); el procedimiento se repite N veces hasta alcanzar el punto (YN, YN).

Opciones

Cada opción es una lista con dos o más elementos. El primer elemento en la lista es el nombre de la opción y el resto consiste en los argumentos para esa opción.

Las opciones aceptadas por las funciones evolution, evolution2, staircase, orbits, ifs y chaosgame son las siguientes:

  * domain especifica los valores mínimo y máximo de la variable independiente para la gráfica de la función F representada por staircase.

  * pointsize define el radio de cada punto dibujado, en unidades de puntos.  El valor por omisión es 1.

  * xaxislabel es la etiqueta que se escribirá en el eje horizontal.

  * xcenter es la coordenada x del punto que deberá aparecer en el centro de la gráfica. Esta opción no es usada por la función orbits.

  * xradius es mitad de la longitud del intervalo de valores de x que será representado. Esta opción no es usada por la función orbits.

  * yaxislabel es la etiqueta que se escribirá en el eje vertical.

  * ycenter es la coordenada y del punto que deberá aparecer en el centro de la gráfica.

  * yradius es mitad de la longitud del intervalo de valores de y que será representado.

Las opciones aceptadas por los programas julia y mandelbrot son las siguientes:

  * size acepta uno o dos argumentos. Si se da solo un argumento, el ancho y la altura del fichero gráfico creado serán iguales a ese valor en pixels. Si se dan dos argumentos, esos dos valores serán usados para el ancho y la altura. El valor por omisión es 400 pixels tanto para el ancho como para la altura. Si los dos valores no son iguales, el conjunto aparecerá distorcionado.

  * levels define el número máximo de iteraciones, que es también el número de colores usado para los puntos que no pertenecem al conjunto. El valor por omisión es 12; valores mayores implican tiempos de procesamiento más elevados.

  * huerange define el intervalo de ángulos usados para la tonalidad de los puntos que no pertenecen al conjunto. El valor por omisión es 360, que hace que los colores usados abarcarán todo el rango de tonalidades. Valores mayores que 360 implican repetición de algunos valores de la tonalidad, y pueden usarse valores negativos para que el ángulo de tonalidad sea menor a medida que el número de iteraciones aumente.

  * hue define la tonalidad, en grados, del primer color usado para representar los puntos que no pertenecen al conjunto.  Su valor por omisión es 300 grados, que corresponde al color magenta; el color correspondiente para otros valores estandar es rojo 0, naranja 45, amarillo 60, verde 120, cian 180 y azul 240. Consulte también la opción HUERANGE.

  * saturation define el nivel de saturación que será usado para los colores de los puntos que no pertenecen al conjunto. Debe ser un valor entre 0 y 1. El valor por omisión es 0.46.

  * value define el valor de los colores usados para puntos que no pertenezcan al conjunto. Debe estar comprendido entre 0 y 1; cuanto mayor sea, mas brilhantes serán los colores. Su valor por omisión es 0.96.

  * color debe ir seguido de tres parámetros que definen la tonalidad, saturación y valor del color que será usado para los puntos del conjunto. El valor por omisión es 0 para los tres parámetros, que corresponde al negro. Consulte las explicación sobre el rango de valores aceptados en la explicación de las opciones HUE, SATURATION y VALUE.

  * center deberá tener dos parámetros reales que dan las coordenadas, en el plano complejo, del punto en el centro de la región representada. El valor por omisión es 0 para las dos coordenadas (el origen).

  * radius es el radio de el mayor círculo que cabe dentro de la región quadrada que será representada. El valor por omisión es 2.

  * filename da el nombre del fichero donde se guardará la gráfica producida. A ese nombre se le acrecentará la terminación .xpm. Si el fichero ya existe, será substituido por el fichero producido por la función. El valor por omisión es julia para el conjunto de Julia y mandelbrot para el conjunto de Mandelbrot.

Ejemplos

Representación gráfica y diagrama de escalera de la secuencia: 2, cos(2), cos(cos(2)),...

     (%i1) load("dynamics")$
     (%i2) evolution(cos(y), 2, 11, [yaxislabel, "y"], [xaxislabel,"n"]);
     (%i3) staircase(cos(y), 1, 11, [domain, 0, 1.2]);

Si su procesador es lento, tendrá que reducir el número de iteraciones usado en los ejemplos siguientes. Y el valor de POINTSIZE que da mejores resultados depende del monitor y de la resolución que use. Tendrá que experimentar con diferentes valores.

Diagrama de órbitas para el mapa cuadrático

             y(n+1) = x + y(n)^2

     (%i4) orbits(y^2+x, 0, 50, 200, [x, -2, 0.25, 0.01], [pointsize, 0.9]);

Para ampliar la región alrededor de la bifurcación en la parte de abajo, cerca de x = -1.25, use el comando:

     (%i5) orbits(x+y^2, 0, 100, 400, [x,-1,-1.53,-0.001], [pointsize,0.9],
                  [ycenter,-1.2], [yradius,0.4]);

Evolución de un sistema en dos dimensiones, que conduce a un fractal:

     (%i6) f: 0.6*x*(1+2*x)+0.8*y*(x-1)-y^2-0.9$
     (%i7) g: 0.1*x*(1-6*x+4*y)+0.1*y*(1+9*y)-0.4$
     (%i8) evolution2d([f,g],[-0.5,0],50000,[pointsize,0.7]);

Y una ampliación de una pequeña región en el fractal:

     (%i9) evolution2d([f,g],[-0.5,0],300000,[pointsize,0.7], [xcenter,-0.7],
                       [ycenter,-0.3],[xradius,0.1],[yradius,0.1]);

Una gráfica del triangulo de Sierpinsky, obtenida con el juego del caos:

     (%i9) chaosgame([[0, 0], [1, 0], [0.5, sqrt(3)/2]], [0.1, 0.1], 1/2,
                      30000, [pointsize,0.7]);

El helecho de Barnsley, obtenido con el Sistema de Funciones Iteradas:

     (%i10) a1: matrix([0.85,0.04],[-0.04,0.85])$
     (%i11) a2: matrix([0.2,-0.26],[0.23,0.22])$
     (%i12) a3: matrix([-0.15,0.28],[0.26,0.24])$
     (%i13) a4: matrix([0,0],[0,0.16])$
     (%i14) p1: [0,1.6]$
     (%i15) p2: [0,1.6]$
     (%i16) p3: [0,0.44]$
     (%i17) p4: [0,0]$
     (%i18) w: [85,92,99,100]$
     (%i19) ifs(w,[a1,a2,a3,a4],[p1,p2,p3,p4],[5,0],50000,[pointsize,0.9]);

Para crear un fichero llamado dinamica9.xpm con la representación gráfica del conjunto de Mandelbrot, con 12 colores, use el comando:

     mandelbrot([filename,"dinamica9"])$

y el conjunto de Julia del número (-0.55 + i 0.6) puede ser obtenido con:

     julia(-0.55, 0.6, [levels, 36], [center, 0, 0.6], [radius, 0.3],
           [hue, 240], [huerange, -180], [filename, "dinamica10"])$

la gráfica se guardará en el fichero dinamica10.xpm y mostrará la región desde -0.3 hasta 0.3 en la dirección x, y desde 0.3 hasta 0.9 en la dirección y. Serán usados 36 colores, comenzando con azul e terminando con amarillo.

Para resolver numéricamente la ecuación diferencial

               dx/dt = t - x^2

Con valor inicial x(t=0) = 1, en el intervalo de t desde 0 hasta 8, y con incrementos de 0.1, se usa:

     (%i20) resultados: rk(t-x^2,x,1,[t,0,8,0.1])$

los resultados quedarán guardados en la lista resultados.

Para resolver numéricamente el sistema:

             dx/dt = 4-x^2-4*y^2     dy/dt = y^2-x^2+1

para t entre 0 y 4, con valores iniciales -1.25 y 0.75 para (x, y) en t=0:

     (%i21) sol: rk([4-x^2-4*y^2,y^2-x^2+1],[x,y],[-1.25,0.75],[t,0,4,0.02])$