The Henon Map

The classical Henon map:

a=14/10
b=3/10
F[{x_,y_}]:={1-a (x^2)+y,b x}

The map F is invertible:

Finv[{x_,y_}]:={(y/b),-1+x+((a (y^2))/(b^2))}
Finv[F[{x,y}]]=={x,y}
F[Finv[{x,y}]]=={x,y}

Here is a plot of an orbit (forgetting the initial part):

We first set q=(F^1000)(0,0). Then we plot the next 10 000 points in the orbit.

q=Nest[F,{0.,0.},1000]
Orbit=NestList[F,q,20000];
plot1=ListPlot[Orbit,PlotRange->All,PlotStyle->AbsolutePointSize[1]]

A failed attempt at finding a Trapping Region.

I plot a polygon and its inverse image under the map. If the inverse image contains the polygon in its interior, then the inverse image is a trapping region.

polygon=Polygon[{{1.1,-0.15},{1.31,-0.05},{1.31,0.1},{0.5,0.25},{-1.32,0.42},{-1.05,-0.5}}]
Vertex[n_,polygon_]:=polygon[[1,Mod[n-1,Length[polygon[[1]]]]+1]]
edges=Table[ParametricPlot[Vertex[i-1,polygon]+t*(Vertex[i,polygon]-Vertex[i-1,polygon]),{t,0,1},PlotStyle->Red,PlotRange->{{-2,1.5},{-1,0.5}}],{i,1,Length[polygon[[1]]]}];
imageEdges=Table[ParametricPlot[Finv[Vertex[i-1,polygon]+t*(Vertex[i,polygon]-Vertex[i-1,polygon])],{t,0,1},PlotStyle->Green],{i,1,Length[polygon[[1]]]}];

Show[Join[{plot1},edges,imageEdges]]

So, the green region above is “almost a trapping region.”

A fixed point

The point p below is a fixed point.

p={(1/28) (-7+(Sqrt[609])),(3/280) (-7+(Sqrt[609]))}
p==F[p]

The following is the Jacobian Matrix at the fixed point.

M=N[D[F[{x,y}],{{x,y}}]/.{x->p[[1]],y->p[[2]]}];
M//MatrixForm
Eigenvalues[M]

So this is a saddle point!

We change coordinates to make the Jacobian diagonal, and the fixed point equal to zero.

{PM,DM}=JordanDecomposition[M]

So, M is approximately equal to the following expression, where D is a diagonal matrix:

PM.DM.Inverse[PM]

ToDiagonalCoordinates[{x_,y_}]:=Inverse[PM].({x,y}-p)
ToOriginalCoordinates[{x_,y_}]:=PM.{x,y}+p
ChangedMap[{x_,y_}]:=ToDiagonalCoordinates[F[ToOriginalCoordinates[{x,y}]]]
ChangedMapInverse[{x_,y_}]:=ToDiagonalCoordinates[Finv[ToOriginalCoordinates[{x,y}]]]

In the new map ChangedMap, the stable manifold is nearly vertical and the unstable manifold is nearly horizontal.

The following shows orbits of points near the stable manifold.

ListPlot[NestList[ChangedMap,{0,0.1},8],PlotRange->All,Joined->True]
ListPlot[NestList[ChangedMap,{0,-0.1},8],PlotRange->All,Joined->True]

A nice box:

xr=0.1;
yr=0.1;
box=Polygon[{{-xr,-yr},{xr,-yr},{xr,yr},{-xr,yr}}]
edges=Table[ParametricPlot[Vertex[i-1,box]+t*(Vertex[i,box]-Vertex[i-1,box]),{t,0,1},PlotStyle->Red,PlotRange->{{-2.2*xr,2.2*xr},{-1.1*yr,1.1*yr}}],{i,1,4}];
imageEdges=Table[ParametricPlot[ChangedMap[Vertex[i-1,box]+t*(Vertex[i,box]-Vertex[i-1,box])],{t,0,1},PlotStyle->Green],{i,1,4}];

Show[Join[edges,imageEdges]]

The following routine decides if a point in the box lies to the right or the left of the stable manifold:

Side[{x_,y_}]:=Block[{out=NestWhile[{ChangedMap[#[[1]]],-#[[2]]}&,{{x,y},1},Abs[#[[1,1]]]<xr&]},Sign[out[[1,1]]*out[[2]]]]

This plots the “left or right” function. The boundary between the colors is the stable manifold.

ContourPlot[Side[{x,y}],{x,-0.101,0.1},{y,-0.101,0.1}]

We can use this to find points on the stable manifold with any given y-coordinate, to a given number of bits of accuracy:

StablePoint[y_,bits_]:=Block[{xleft=-0.1,xright=0.1},
For[i=0,i<bits,i++,Block[{xmiddle=(xleft+xright)/2},
If[Side[{xmiddle,y}]==1,xright=xmiddle,xleft=xmiddle]]];{(xleft+xright)/2,y}]

We compile a list of points on the stable manifold.

NumPoints=1000;
stableManifold=Table[StablePoint[i*0.1/NumPoints,30],{i,-NumPoints+0.5,NumPoints}];

The local stable manifold:

ListPlot[stableManifold,Joined->True,PlotRange->{{-xr,xr},{-yr,yr}},AspectRatio->Automatic]

We move it into our original coordinates, and plot it twice. The second time we draw it with our orbit approximating the attractor.

stableManifold2=Table[ToOriginalCoordinates[stableManifold[[i]]],{i,1,Length[stableManifold]}];
stableManifoldPlot=ListPlot[stableManifold2,Joined->True,PlotStyle->Red]
Show[plot1,stableManifoldPlot]

Now we can iterate it to obtain more of the stable manifold:

stableManifold2=Table[Finv[stableManifold2[[i]]],{i,1,Length[stableManifold]}];
stableManifoldPlot=ListPlot[stableManifold2,Joined->True,PlotStyle->Red,PlotRange->{{-1.5,1.5},{-0.5,0.5}}];
Show[plot1,stableManifoldPlot]

stableManifold2=Table[Finv[stableManifold2[[i]]],{i,1,Length[stableManifold]}];
stableManifoldPlot=ListPlot[stableManifold2,Joined->True,PlotStyle->Red,PlotRange->{{-1.5,1.5},{-0.5,0.5}}];
Show[plot1,stableManifoldPlot]

stableManifold2=Table[Finv[stableManifold2[[i]]],{i,1,Length[stableManifold]}];
stableManifoldPlot=ListPlot[stableManifold2,Joined->True,PlotStyle->Red,PlotRange->{{-1.5,1.5},{-0.5,0.5}}];
Show[plot1,stableManifoldPlot]

stableManifold2=Table[Finv[stableManifold2[[i]]],{i,1,Length[stableManifold]}];
stableManifoldPlot=ListPlot[stableManifold2,Joined->True,PlotStyle->Red,PlotRange->{{-1.5,1.5},{-0.5,0.5}}];
Show[plot1,stableManifoldPlot]

The Unstable manifold

We do the analogous thing to find the unstable manifold.

The following decides if a point lies above or below the unstable manifold.

Side2[{x_,y_}]:=Block[{out=NestWhile[ChangedMapInverse,{x,y},Abs[#[[2]]]<yr&]},Sign[out[[2]]]]

We plot this function below:

ContourPlot[Side2[{x,y}],{x,-0.101,0.1},{y,-0.101,0.1}]

We can use this to find points with a given x-coordinate on the unstable manifold.

UnstablePoint[x_,bits_]:=Block[{yleft=-0.1,yright=0.1},
For[i=0,i<bits,i++,Block[{ymiddle=(yleft+yright)/2},
If[Side2[{x,ymiddle}]==1,yright=ymiddle,yleft=ymiddle]]];{x,(yleft+yright)/2}]

We find some points on the unstable manifold:

NumPoints=5000;
unstableManifold=Table[UnstablePoint[i*0.1/NumPoints,30],{i,-NumPoints+0.5,NumPoints}];

Now we plot the unstable manifold. It hugs the x-axis very closely.

ListPlot[unstableManifold,Joined->True,PlotRange->{{-xr,xr},{-yr,yr}}]

We move the points back into our original coordinates:

unstableManifold2=Table[ToOriginalCoordinates[unstableManifold[[i]]],{i,1,Length[unstableManifold]}];
unstableManifoldPlot=ListPlot[unstableManifold2,Joined->True,PlotStyle->Green]
Show[plot1,unstableManifoldPlot]

Then, we can iterate forward to obtain more of the unstable manifold.

unstableManifold2=Table[F[unstableManifold2[[i]]],{i,1,Length[unstableManifold]}];
unstableManifoldPlot=ListPlot[unstableManifold2,Joined->True,PlotStyle->Green,PlotRange->{{-1.5,1.5},{-0.5,0.5}}];
Show[unstableManifoldPlot]

unstableManifold2=Table[F[unstableManifold2[[i]]],{i,1,Length[unstableManifold]}];
unstableManifoldPlot=ListPlot[unstableManifold2,Joined->True,PlotStyle->Green,PlotRange->{{-1.5,1.5},{-0.5,0.5}}];
Show[unstableManifoldPlot]

unstableManifold2=Table[F[unstableManifold2[[i]]],{i,1,Length[unstableManifold]}];
unstableManifoldPlot=ListPlot[unstableManifold2,Joined->True,PlotStyle->Green,PlotRange->{{-1.5,1.5},{-0.5,0.5}}];
Show[unstableManifoldPlot]

unstableManifold2=Table[F[unstableManifold2[[i]]],{i,1,Length[unstableManifold]}];
unstableManifoldPlot=ListPlot[unstableManifold2,Joined->True,PlotStyle->Green,PlotRange->{{-1.5,1.5},{-0.5,0.5}}];
Show[unstableManifoldPlot]

unstableManifold2=Table[F[unstableManifold2[[i]]],{i,1,Length[unstableManifold]}];
unstableManifoldPlot=ListPlot[unstableManifold2,Joined->True,PlotStyle->Green,PlotRange->{{-1.5,1.5},{-0.5,0.5}}];
Show[unstableManifoldPlot]

unstableManifold2=Table[F[unstableManifold2[[i]]],{i,1,Length[unstableManifold]}];
unstableManifoldPlot=ListPlot[unstableManifold2,Joined->True,PlotStyle->Green,PlotRange->{{-1.5,1.5},{-0.5,0.5}}];
Show[unstableManifoldPlot]

unstableManifold2=Table[F[unstableManifold2[[i]]],{i,1,Length[unstableManifold]}];
unstableManifoldPlot=ListPlot[unstableManifold2,Joined->True,PlotStyle->Green,PlotRange->{{-1.5,1.5},{-0.5,0.5}}];
Show[unstableManifoldPlot]

unstableManifold2=Table[F[unstableManifold2[[i]]],{i,1,Length[unstableManifold]}];
unstableManifoldPlot=ListPlot[unstableManifold2,Joined->True,PlotStyle->Green,PlotRange->{{-1.5,1.5},{-0.5,0.5}}];
Show[unstableManifoldPlot]

unstableManifold2=Table[F[unstableManifold2[[i]]],{i,1,Length[unstableManifold]}];
unstableManifoldPlot=ListPlot[unstableManifold2,Joined->True,PlotStyle->Green,PlotRange->{{-1.5,1.5},{-0.5,0.5}}];
Show[unstableManifoldPlot]

unstableManifold2=Table[F[unstableManifold2[[i]]],{i,1,Length[unstableManifold]}];
unstableManifoldPlot=ListPlot[unstableManifold2,Joined->True,PlotStyle->Green,PlotRange->{{-1.5,1.5},{-0.5,0.5}}];
Show[unstableManifoldPlot]

unstableManifold2=Table[F[unstableManifold2[[i]]],{i,1,Length[unstableManifold]}];
unstableManifoldPlot=ListPlot[unstableManifold2,Joined->True,PlotStyle->Green,PlotRange->{{-1.5,1.5},{-0.5,0.5}}];
Show[unstableManifoldPlot]