Sistemas Dinámicos Discretos

Importante: Abrir y evaluar el documento definiciones.nb que contiene las funciones necesarias

Determinación de ciclos

f[x_] := x^2 - 2

Plot[{x, f[x]}, {x, -2, 2}, 
        PlotRange -> {{-2.1, 2.1}, {-2.1, 2.1}}, 
        AspectRatio -> 1, 
        PlotLabel -> "f(x)=x^2-2"];

Plot[{x, f[f[x]]}, {x, -2, 2}, 
        PlotRange -> {{-2.1, 2.1}, {-2.1, 2.1}}, 
        AspectRatio -> 1, PlotLabel -> "f^2(x)"];

Plot[{x, f[f[f[x]]]}, {x, -2, 2}, 
        PlotRange -> {{-2.1, 2.1}, {-2.1, 2.1}}, 
        AspectRatio -> 1, PlotLabel -> "f^3(x)"];

Iterate[f, x, 3] // Expand

Series temporales e histogramas

h[x_] := x^2 - 1.3

Orbit[h, 0, 5]

orb1 = Orbit[h, 0, 49];

TimeSeries[orb1];

orb2 = Orbit[logis[3.9], 0.5, 49];

TimeSeries[orb2];

Histogramas

orb3 = Orbit[logis[3.83], 0.123, 300];

Histogram[orb3, HistogramCategories -> 100];

orb4 = Orbit[logis[3.9], 0.123, 2000];

Histogram[orb4, HistogramCategories -> 100];

Iteración gráfica

Table[WebDiagram[logis[3.5], 0.1, i];, {i, 1, 9}];

Orbit[logis[3.5], 0.1, 9]

WebDiagram[logis[3.5], 0.38282, 4];

Orbit[logis[3.5], 0.38282, 8]

web2 = WebDiagram[logis[3.5], 0.1, 9];

Puntos fijos atractores y repulsores

Solve[logis[2.8][x] == x, x]

WebDiagram[logis[2.8], 0.07, 20];

Solve[logis[3.2][x] == x, x]

web4 = WebDiagram[logis[3.2], 0.032, 40];

Solve[Iterate[logis[3.2], x, 2] == x, x]

La función logística y el caos

Un ejemplo sin caos:  f(x)=2 x (1-x)

pfijos1 = Solve[logis[2][x] == x, x]

D[logis[2][x], x] /. pfijos1[[1]]

D[logis[2][x], x] /. pfijos1[[2]]

Orbit[logis[2], 0.01, 10]

WebDiagram[logis[2], 0.01, 10];

Orbit[logis[2], 0.99, 10]

WebDiagram[logis[2], 0.99, 10];

Un ejemplo intermedio: f(x)=3.2 x (1-x)

pfijos2 = Solve[logis[3.2][x] == x, x]

D[logis[3.2][x], x] /. pfijos2[[1]]

D[logis[3.2][x], x] /. pfijos2[[2]]

pciclos2 = Solve[Iterate[logis[3.2], x, 2] == x, x]

(D[logis[3.2][x], x] /. 
          pciclos2[[2]])*(D[logis[3.2][x], x] /. 
          pciclos2[[4]])

WebDiagram[logis[3.2], 0.3, 30];

WebDiagram[logis[3.2], 0.011232580472832009, 4];

Orbit[logis[3.2], 0.011232580472832009, 5]

Cambiamos la séptima cifra decimal de 5 a 9

Orbit[logis[3.2], 0.011232980472832009, 10]

Orbit[logis[3.2], 0.011232980472832009, 30]

Orbit[logis[3.2], 0.011232980472832009, 90]

WebDiagram[logis[3.2], 0.011232980472832009, 100];

Cálculo órbita hacia atrás

atras1 = 
      Solve[logis[3.2][x] == Evaluate[x /. pfijos2[[2]]], 
        x]

atras2 = 
      Solve[logis[3.2][x] == Evaluate[x /. atras1[[1]]], 
        x]

atras3 = 
      Solve[logis[3.2][x] == Evaluate[x /. atras2[[1]]], 
        x]

atras4 = 
      Solve[logis[3.2][x] == Evaluate[x /. atras3[[1]]], 
        x]

Un ejemplo con caos: f(x)=4 x (1-x)

pfijos3 = Solve[logis[4][x] == x, x]

D[logis[4][x], x] /. pfijos3[[1]]

D[logis[4][x], x] /. pfijos3[[2]]

Orbit[logis[4], 0.009607359798384785, 7]

Orbit[logis[4], 0.0010705383806982505, 7]

WebDiagram[logis[4], 0.123, 300];

orbcaos = Orbit[logis[4], 0.123, 1000];

Histogram[orbcaos, HistogramCategories -> 100];

orbcaos2 = Orbit[logis[4], 0.123, 100];

TimeSeries[orbcaos2];

El método de Euler, la ecuación logística y el caos

g[r_][x_, y_] = logis[r][y]

reglas = 
      DSolve[{y'[x] == r y[x] (1 - y[x]), y[0] == y00}, 
        y[x], x]

exacta[r_][y0_][x_] := 
      Evaluate[y[x] /. reglas[[1]] /. {y00 -> y0}]

y'=2.5 y (1-y),   y(0)=2, x en [0,1],  h=0.1,  n=10

puntos1 = EulerValues[g[2.5], 0, 2, 0.1, 10]

ej1 = 
        Show[
Plot[exacta[2.5][2][x], {x, 0, 1}, 
            PlotRange -> All, PlotStyle -> RGBColor[0, 0, 1], 
            DisplayFunction -> Identity], 
          ListPlot[puntos1, PlotJoined -> True, PlotRange -> All, 
            DisplayFunction -> 
              Identity], 
DisplayFunction -> $DisplayFunction];

y'=2.5 y (1-y),   y(0)=2, x en [0,1],  h=0.05,  n=20

puntos2 = EulerValues[g[2.5], 0, 2, 0.05, 20];

ej2 = 
        Show[
Plot[exacta[2.5][2][x], {x, 0, 1}, 
            PlotRange -> All, PlotStyle -> RGBColor[0, 0, 1], 
            DisplayFunction -> Identity], 
          ListPlot[puntos2, PlotJoined -> True, PlotRange -> All, 
            DisplayFunction -> 
              Identity], 
DisplayFunction -> $DisplayFunction];

y'=53 y (1-y),   y(0)=1.15, x en [0,1], h=0.05,  n=20

puntos3 = EulerValues[g[53], 0, 1.15, 0.05, 20];

ej3 = 
        Show[
Plot[exacta[53][1.15][x], {x, 0, 1}, 
            PlotRange -> All, PlotStyle -> RGBColor[0, 0, 1], 
            DisplayFunction -> Identity], 
          ListPlot[puntos3, PlotJoined -> True, PlotRange -> All, 
            DisplayFunction -> 
              Identity], 
DisplayFunction -> $DisplayFunction];

y'=26.5 y (1-y),   y(0)=1.2,  h=0.1,  n=2000  (próximo al caos   1+h r = 26.5 * 0.1 = 3.65)

puntoscaos = EulerValues[g[26.5], 0, 1.2, 0.1, 2000];

grafpuntoscaos = 
        ListPlot[puntoscaos, PlotRange -> All];

grafpolcaos = 
        ListPlot[puntoscaos, PlotRange -> All, 
          PlotJoined -> True];

Eliminamos esta situación  si 1+h r no es próximo a cuatro  (1+h r = 1+26.5 * 0.01 = 1.265)

puntos5 = EulerValues[g[26.5], 0, 1.2, 0.01, 30];

ej5 = 
        Show[
Plot[exacta[26.5][1.2][x], {x, 0, 
              0.3}, PlotRange -> All, 
            PlotStyle -> RGBColor[0, 0, 1], 
            DisplayFunction -> Identity], 
          ListPlot[puntos5, PlotJoined -> True, PlotRange -> All, 
            DisplayFunction -> 
              Identity], 
DisplayFunction -> $DisplayFunction];