Numérico búsqueda de órbitas periódicas

Question

Supongamos que $\phi$ es un mapa definido en todos los de $R^n$, con la imagen en $R^n$. Hay un caso conocido de un buen algoritmo para encontrar el periódico puntos de $\phi$ orden $m$, (o puntos fijos de $\phi^m$) además de:

1. Aplicando el método de Newton en algunos (al azar) inicial $m$ $x_0,x_1, ..., x_{m-1}$ y buscando las soluciones del sistema $\phi(x_k)=x_{k+1}$;
2. O aplicar el método de Newton en una sola (al azar) punto inicial $x$, y buscando la solución de $\phi^m(x)=x$ (y por tanto para calcular el producto de la Jacobians de $\phi$ en su primer $m$ poderes).

Yo también agradecería cualquier referencia acerca de esta cuestión. Gracias!

Answer 1

Si estás interesado en el atractivo de las órbitas luego, por supuesto, sólo se puede iterar la función. Una de las principales teorema en una dimensión dinámica de los estados del atractivo de las órbitas siempre atraen al menos un punto crítico. Esto significa que sólo hay un par de semillas iniciales que usted necesita para ejecutar la prueba. No estoy seguro de si hay un análogo a esto en las dimensiones superiores de la dinámica, sin embargo. Además, no te ayuda para no atractiva órbitas de todos modos.

En el caso más general de (potencialmente) no atractiva de las órbitas, la primera cosa a tener en cuenta es que probablemente no es una buena idea para tratar de encontrar puntos fijos de $f^{m}$, como parece que sabes. La razón es que el método de Newton, al igual que muchos algoritmos numéricos, es bastante sensible a las condiciones iniciales. Para evitar esto, supongamos que estamos buscando un punto de plazo, $m$ y definir un mapa de $F:{\mathbb R}^{n m} \rightarrow {\mathbb R}^{n m}$ por

$$F(x_1,x_2,\ldots,x_m) = (f(x_n), f(x_1), \ldots, f(x_{n-1})).$$

Note that each $x_i$ here is a vector in $\mathbb R^n$ and that an orbit of period $m$ for $f$ forms a vector that is a fixed point of $F$. Shifting to a higher dimensional map like this lightens the problem of sensitivity.

Next, if you want all periodic points, there are algorithms that can simultaneously find all solutions to many types of equations and systems of equations. One simple technique, that works for a univariate polynomial, is to realize the polynomial as the characteristic polynomial of a matrix and then use techniques from numerical linear algebra to find the eigenvalues. A system of polynomials or even algebraic equations can be reduced using Groebner basis techniques, which is somewhat analogous to Gaussian elimination; the eigen-approach can then be used on this system.

Of course these techniques and others are implemented in many modern numerical tools. My personal favorite is Mathematica, which can also solve many systems of transcendental equations, if you're will to restrict the domain of search.

Example

As an example, let's explore the chaotic Henon map defined by

$$f(x,y) = (y+1-1.4x^2, 0.3x).$$

This is a fairly well known 2D chaotic map with a strange attractor that looks like so:

f[{x_, y_}] = {y + 1 - 1.4 x^2, 0.3 x};
hpic = ListPlot[NestList[f, {0, 0}, 1000]]

Note that I'll work this example out in Mathematica but will try to keep the code to a minimum. I think the code for the important steps is fairly clear. For example, we clearly see the definition of $f$ in the code above. And, again, there is certainly other numerical software out there that could do this.

Next, we define the function $F$:

F[xys:{{_, _} ..}] := f /@ RotateRight[xys];

The {{_,_}..} business defines $F$ for lists of pairs; given such a list, the list is rotated to the right and then $f$ se transforma en el resultado. Ahora, supongamos que queremos que todas las órbitas de período 7.

vars = Table[{x[k], y[k]}, {k, 1, 7}];
orbits = Select[Chop[
  NSolve[F[vars] == vars, Flatten[vars]]], 
  FreeQ[#, _Complex] &];
firstOrbit = First[orbits]

(* Out: {x[1] -> -0.285392, y[1] -> -0.326186, x[2] -> 0.559786, 
    y[2] -> -0.0856176, x[3] -> 0.475678, y[3] -> 0.167936, 
    x[4] -> 0.851159, y[4] -> 0.142703, x[5] -> 0.128443, 
    y[5] -> 0.255348, x[6] -> 1.23225, y[6] -> 0.038533, 
    x[7] -> -1.08729, y[7] -> 0.369675}
*)

La parte clave es simplemente la NSolve paso, que encontrar las órbitas. Vamos a comprobar la órbita de la muestra.

{x0, y0} = {-0.2853920648214598, -0.3261858137547133}; 
seventhIterate = Nest[f, {x0, y0}, 7]
seventhIterate - {x0, y0}

(* Out: 
  {-0.285392, -0.326186} *)
  {1.04361*10^-14, 6.10623*10^-16}
*)

No está mal! He aquí lo que la órbita se ve como en el atractor de sí mismo:

Hay mucho que hacer aquí, sin embargo. En particular, de orden superior, se itera va a ser más difícil. A veces un homotopy enfoque puede ayudar allí. Creo que necesita un poco más de conocimiento y de la inteligencia en relación con el sistema, sin embargo. Hice un homotopy enfoque a una pregunta diferente en mathematica.stackexchange.com: http://mathematica.stackexchange.com/questions/19480/

Answer 2

Una caja de herramientas de software como GAIO podría ayudar.