Function Repository Resource:

IteratedMap2D

Source Notebook

Obtain the orbit of a 2D iterated map from given initial conditions

Contributed by: E. Chan-López, Erika Yuridia Flores Chan, Héctor Argote Morales & Víctor Castellanos

ResourceFunction["IteratedMap2D"][map,icv,depvar,trange,stepsize]

gives the orbit corresponding to the recurrence equations represented by map with initial conditions icv, dependent variables depvars, an iterator specification trange and a stepsize.

Details

The map must be a list of recurrence equations.
The icv must be a list of initial conditions.
The depvars must be a list containing the dependent variables.
The trange must be a list containing the independent variable and a matrix.
The stepsize corresponds to the steps taken to reach the maximum number of iterations.
ResourceFunction["IteratedMap2D"][map,icv,depvars,trange,stepsize,"TimeSeriesData"] can be used to obtain the time series data from orbit data.

Examples

Basic Examples (6) 

Define the Gingerbreadman map:

In[1]:=
Gingerbreadman[{x_[n_], y_[n_]}] := {x[n + 1] == 1 - y[n] + Abs[x[n]],
    y[n + 1] == x[n]};
X = {x[n], y[n]}
Out[2]=

Initial condition:

In[3]:=
X0 = {{-0.1, 0.}}
Out[3]=

Interval of the maximum number of iterations:

In[4]:=
ninterval = {{1, 10000}}
Out[4]=

Step size of iterations for initial condition X0:

In[5]:=
nstep = ConstantArray[1, Length[X0]]
Out[5]=

Generate the orbit of the Gingerbreadman map:

In[6]:=
ResourceFunction[
  "IteratedMap2D", ResourceSystemBase -> "https://www.wolframcloud.com/objects/resourcesystem/api/1.0"][Gingerbreadman[X], X0, X, {n, ninterval}, nstep] // Short
Out[6]=

Visualize the orbit:

In[7]:=
ListPlot[%, AspectRatio -> Automatic, PlotStyle -> Black]
Out[7]=

Scope (2) 

Define the following dynamical system:

In[8]:=
eqns[{x_[n_], y_[n_]}] := {x[n + 1] == (4/5)*x[n] + y[n], y[n + 1] == 41/50 - x[n]^2};
X = {x[n], y[n]}
Out[9]=

Initial condition:

In[10]:=
X0 = {{0.2, 0.2}}
Out[10]=

Interval of the maximum number of iterations:

In[11]:=
ninterval = {{1, 100000}}
Out[11]=

Step size of iterations for initial condition X0:

In[12]:=
nstep = ConstantArray[1, Length[X0]]
Out[12]=

Orbit data:

In[13]:=
ResourceFunction[
  "IteratedMap2D", ResourceSystemBase -> "https://www.wolframcloud.com/objects/resourcesystem/api/1.0"][eqns[X], X0, X, {n, ninterval}, nstep] // Short
Out[13]=

Visualize the orbit:

In[14]:=
ListPlot[%, PlotStyle -> {{Black, Opacity[0.5]}}]
Out[14]=

Use IteratedMap2D for several initial conditions with the map:

In[15]:=
eqns[{x_[k_], y_[k_]}] := With[{a = 1/13, b = 1/17, c = 1/23}, {x[k + 1] == y[k] - Sign[a*x[k]]*Sqrt[Abs[a b c x[k]*y[k] + c x[k]]], y[k + 1] == a/13 - x[k]}];
X = {x[k], y[k]}
Out[16]=

Initial conditions:

In[17]:=
X0 = {{0.17, 0.17}, {0.099, 0.044}, {0.07, 0.14}, {-0.14, 0.067}, {-0.117, 0.039}, {-0.134, -0.029}, {-0.0705, -0.0298}, {-0.036, 0.00968}, {-0.0186, -0.0135}, {0.053, 0.0329}, {0.0237, -0.024}, {0.088, 0.091}, {0.1058, 0.193}, {0.0802, 0.143}, {-0.0294, -0.0182}, {0.15969, 0.1061}, {0.0856, 0.1671}, {-0.0799, -0.0288}, {-0.00859, -0.0052}, {0.1112, -0.035}, {0.214889, 0.0757}}
Out[17]=

Interval of the maximum number of iterations for each initial condition:

In[18]:=
nintervals = ConstantArray[{1, 2500}, Length[X0]]
Out[18]=

Step size of iterations for each initial condition:

In[19]:=
nsteps = ConstantArray[1, Length[X0]]
Out[19]=

Orbits data:

In[20]:=
ResourceFunction[
  "IteratedMap2D", ResourceSystemBase -> "https://www.wolframcloud.com/objects/resourcesystem/api/1.0"][eqns[X], X0, X, {k, nintervals}, nsteps] // Short
Out[20]=

Visualize the orbits:

In[21]:=
Graphics[
 MapIndexed[{PointSize[0.004], ColorData["Rainbow"][#2[[1]]/Length[X0]], Point[#1]} &, %], Axes -> True, AspectRatio -> 1/GoldenRatio]
Out[21]=

Applications (28) 

Neimark-Sacker bifurcation (13) 

Use IteratedMap2D considering three initial conditions with the following two-dimensional discrete-time dynamical system :

In[22]:=
vfield = {r x (1 - y), x};
vars = {x, y}
Out[23]=

Jacobian matrix:

In[24]:=
J[{x_, y_}][r_] := Evaluate[D[vfield, {vars}]]
MatrixForm@J[vars][r]
Out[25]=

The non-trivial fixed point:

In[26]:=
fp = SolveValues[vfield - vars == 0 && x > 0 && y > 0 && r > 0, vars][[1]];
MatrixForm@List@fp
Out[27]=

The Neimark-Sacker critical bifurcation value:

In[28]:=
r0 = SolveValues[Det[J[Normal[fp]][r]] == 1, r][[1]]
Out[28]=

The transversality condition:

In[29]:=
ReplaceAll[r -> r0][D[Sqrt[Det[J[Normal[fp]][r]]], r]]
Out[29]=

The non-trivial fixed point is unstable for 𝓇>𝓇0 and locally stable for 𝓇<r0, with the stable closed invariant curve bifurcating from for 𝓇>𝓇0. Define the system for 𝓇<r0 and 𝓇>𝓇0:

In[30]:=
seqns[{x_[n_], y_[n_]}] := With[{r = 1.98}, {x[n + 1] == r x[n] (1 - y[n]), y[n + 1] == x[n]}];
ueqns[{x_[n_], y_[n_]}] := With[{r = 2.012}, {x[n + 1] == r x[n] (1 - y[n]), y[n + 1] == x[n]}];
X = {x[n], y[n]}
Out[32]=

Initial conditions for 𝓇<r0 and 𝓇>𝓇0:

In[33]:=
X0s = {{0.30, 0.25}, {0.44, 0.22}}
Out[33]=
In[34]:=
X0u = {{0.495, 0.495}, {0.485, 0.485}, {0.44, 0.22}}
Out[34]=

Intervals of the maximum number of iterations for 𝓇<r0 and 𝓇>𝓇0:

In[35]:=
snintervals = {{1, 200}, {1, 500}}
Out[35]=
In[36]:=
unintervals = {{1000, 1500}, {1, 500}, {1, 200}}
Out[36]=

Step size of iterations for each initial condition when 𝓇<r0 and 𝓇>𝓇0:

In[37]:=
snsteps = ConstantArray[1, Length[X0s]]
Out[37]=
In[38]:=
unsteps = ConstantArray[1, Length[X0u]]
Out[38]=

Orbits data for 𝓇<r0:

In[39]:=
Short[sdata = ResourceFunction[
   "IteratedMap2D", ResourceSystemBase -> "https://www.wolframcloud.com/objects/resourcesystem/api/1.0"][seqns[X], X0s, X, {n, snintervals}, snsteps]]
Out[39]=

Orbits data for 𝓇>r0:

In[40]:=
Short[udata = ResourceFunction[
   "IteratedMap2D", ResourceSystemBase -> "https://www.wolframcloud.com/objects/resourcesystem/api/1.0"][ueqns[X], X0u, X, {n, unintervals}, unsteps]]
Out[40]=

Visualize the phase portrait graphics for the two conditions:

In[41]:=
{Graphics[
   MapIndexed[{PointSize[0.006], ColorData["Rainbow"][#2[[1]]/Length[X0s]], Point[#1]} &, sdata], {Frame -> True, FrameTicksStyle -> GrayLevel[0], FrameLabel -> {
TraditionalForm[r < Subscript[r, 0]], None}, Axes -> False, AspectRatio -> GoldenRatio^(-1)}], Graphics[
   MapIndexed[{PointSize[0.006], ColorData["Rainbow"][#2[[1]]/Length[X0u]], Point[#1]} &, udata], {Frame -> True, FrameTicksStyle -> GrayLevel[0], FrameLabel -> {
TraditionalForm[r > Subscript[r, 0]], None}, Axes -> False, AspectRatio -> GoldenRatio^(-1)}]} // GraphicsRow
Out[41]=

Generate time series data for the discrete-time dynamical system:

In[42]:=
ResourceFunction[
  "IteratedMap2D", ResourceSystemBase -> "https://www.wolframcloud.com/objects/resourcesystem/api/1.0"][ueqns[X], X0u, X, {n, unintervals}, unsteps, "TimeSeriesData"] // Short
Out[42]=
In[43]:=
GraphicsGrid[
 Partition[
  Flatten[MapIndexed[{ListPlot[#1, PlotStyle -> {PointSize[0.01], ColorData["Rainbow"][#2[[1]]/Length[X0u]]}]} &, %]], 3]]
Out[43]=

Multiple Attractors in a Predator-Prey Model (5) 

Define a predator-prey model using the Crowley-Martin functional response and initial conditions:

In[44]:=
eqns[{x_[n_], y_[n_]}] := With[{a = 24647595988461150985968273113/665860725091715916800000000,
     b = 247579691107/12800000000, c = 6400000000/80393230369, d = 323537/40000, k = 2828296/1090611, r = 1414148/323537, \[Gamma] = 323537/60000}, {x[n + 1] == r x[n] (1 - x[n]/k) - (
      a x[n]*y[n])/((1 + b x[n]) (1 + c y[n])), y[n + 1] == (\[Gamma] a x[n]*y[n])/((1 + b x[n]) (1 + c y[n])) - d y[n]}];
X0 = {{0.995, 0.995}, {0.995, 0.995}, {0.77, 0.82}, {0.9, 0.77}}
Out[45]=

Interval of the maximum number of iterations for each initial condition:

In[46]:=
nintervals = {{8000, 10000}, {1, 5000}, {1, 800}, {1, 800}}
Out[46]=

Step size of iterations for each initial condition:

In[47]:=
nsteps = {1, 8, 1, 1}
Out[47]=

Orbits data:

In[48]:=
ResourceFunction[
  "IteratedMap2D", ResourceSystemBase -> "https://www.wolframcloud.com/objects/resourcesystem/api/1.0"][eqns[X], X0, X, {n, nintervals}, nsteps] // Short
Out[48]=

Visualize the orbits:

In[49]:=
Graphics[
 MapIndexed[{PointSize[0.004], ColorData["Rainbow"][#2[[1]]/Length[X0]], Point[#1]} &, %], Axes -> True, AspectRatio -> 1/GoldenRatio]
Out[49]=

A chaotic attractor:

In[50]:=
eqns[{x_[n_], y_[n_]}] := With[{a = 17275862227/622054400, b = 6883/512, c = 256/2209, d = 55/8, k = 488/189, r = 244/55, \[Gamma] = 55/12}, {x[n + 1] == r x[n] (1 - x[n]/k) - (
      a x[n]*y[n])/((1 + b x[n]) (1 + c y[n])), y[n + 1] == (\[Gamma] a x[n]*y[n])/((1 + b x[n]) (1 + c y[n])) - d y[n]}];
X = {x[n], y[n]}
Out[51]=
In[52]:=
X0 = {{0.995, 0.995}}
Out[52]=

Interval of the maximum number of iterations:

In[53]:=
ninterval = {{1, 200000}}
Out[53]=

Step size of iterations for the initial condition X0:

In[54]:=
nstep = ConstantArray[1, Length[X0]]
Out[54]=

Orbit data:

In[55]:=
ResourceFunction[
  "IteratedMap2D", ResourceSystemBase -> "https://www.wolframcloud.com/objects/resourcesystem/api/1.0"][eqns[X], X0, X, {n, ninterval}, nstep] // Short
Out[55]=

Visualize the orbit:

In[56]:=
ListPlot[%, PlotStyle -> {{Black, Opacity[0.3]}}, AspectRatio -> 1]
Out[56]=

First strange attractor:

In[57]:=
eqns[{x_[n_], y_[n_]}] := With[{a = 4592791259697/181297792000, b = 152923/12800, c = 6400/48841, d = 261/40, k = 776/301, r = 388/87, \[Gamma] = 87/20}, {x[n + 1] == r x[n] (1 - x[n]/k) - (a x[n]*y[n])/((1 + b x[n]) (1 + c y[n])), y[n + 1] == (\[Gamma] a x[n]*y[n])/((1 + b x[n]) (1 + c y[n])) - d y[n]}]

Orbit data:

In[58]:=
ResourceFunction[
  "IteratedMap2D", ResourceSystemBase -> "https://www.wolframcloud.com/objects/resourcesystem/api/1.0"][
  eqns[X], {{0.995, 0.995}}, {x[n], y[n]}, {n, ninterval}, {1}] // Short
Out[58]=

Visualize the orbit:

In[59]:=
ListPlot[%, PlotStyle -> {{Black, Opacity[0.15]}}, AspectRatio -> 1]
Out[59]=

Second strange attractor:

In[60]:=
eqns[{x_[n_], y_[n_]}] := With[{a = 8592991153063/334815411200, b = 155587/12800, c = 6400/49729, d = 263/40, k = 2344/909, r = 1172/263, \[Gamma] = 263/60}, {x[n + 1] == r x[n] (1 - x[n]/k) - (a x[n]*y[n])/((1 + b x[n]) (1 + c y[n])), y[n + 1] == (\[Gamma] a x[n]*y[n])/((1 + b x[n]) (1 + c y[n])) - d y[n]}]

Orbit data:

In[61]:=
ResourceFunction[
  "IteratedMap2D", ResourceSystemBase -> "https://www.wolframcloud.com/objects/resourcesystem/api/1.0"][
  eqns[X], {{0.995, 0.995}}, {x[n], y[n]}, {n, ninterval}, {1}] // Short
Out[61]=

Visualize the orbit:

In[62]:=
ListPlot[%, PlotStyle -> {{Black, Opacity[0.3]}}, AspectRatio -> 1]
Out[62]=

Third strange attractor with multiple periodic orbits:

In[63]:=
eqns[{x_[n_], y_[n_]}] := With[{a = 24686881771241562225059473113/665860725091715916800000000,
     b = 247579691107/12800000000, c = 6400000000/80393230369, d = 323537/40000, k = 2828296/1090611, r = 1414148/323537, \[Gamma] = 323537/60000}, {x[n + 1] == r x[n] (1 - x[n]/k) - (
      a x[n]*y[n])/((1 + b x[n]) (1 + c y[n])), y[n + 1] == (\[Gamma] a x[n]*y[n])/((1 + b x[n]) (1 + c y[n])) - d y[n]}];
X = {x[n], y[n]}
Out[64]=

Initial condition:

In[65]:=
X0 = {{0.77, 0.82}}
Out[65]=

Orbit data:

In[66]:=
ResourceFunction[
  "IteratedMap2D", ResourceSystemBase -> "https://www.wolframcloud.com/objects/resourcesystem/api/1.0"][eqns[X], X0, X, {n, ninterval}, nstep] // Short
Out[66]=

Visualize the orbit:

In[67]:=
ListPlot[%, PlotStyle -> {{Black, Opacity[0.3]}}, AspectRatio -> 1]
Out[67]=

Zoom into an "interesting" region:

In[68]:=
Show[%, PlotStyle -> {Black, Opacity[0.45], PointSize[0.0025]}, PlotRange -> {{1.19, 1.29}, {0.7, 0.77}}, Frame -> True, FrameStyle -> Directive[Black]]
Out[68]=

Standard Map (5) 

The standard map is defined by the equation 𝓍n+1=𝓍n+ε sin(𝓍n) , 𝓎n+1=𝓎n+𝓍n+1 where 𝓍n and 𝓎n are taken modulo 2π:

In[69]:=
StandardMap[{x_[n_], y_[n_]}][\[CurlyEpsilon]_] := {x[n + 1] == Mod[x[n] + y[n] + \[CurlyEpsilon] Sin[x[n]], 2 \[Pi]], y[n + 1] == Mod[y[n] + \[CurlyEpsilon] Sin[x[n]], 2 \[Pi]]};
X0 = {{0.0952, 0.2801}, {2.9067, 4.5882}, {0.6671, 5.8718}, {0.8733, 4.5413}, {0.9015, 1.03502}, {5.8875, 5.6839}, {5.3908, 0.2054}, {1.48266, 0.1584}, {1.6982, 0.1427}, {1.9887, 0.0645}, {2.3261, 0.0801}, {2.6635, 0.0958}, {3.0009, 0.0801}, {5.9812, 5.6842}, {0.1986, 4.0092}, {0.1612, 2.0112}, {1.4451, 3.1748}, {0.1237, 1.646}, {0.4704, 1.6512}, {0.0581, 2.5017}, {1.8856, 2.6791}, {2.148, 2.7313}, {2.551, 2.893}, {0.0581, 4.448}, {0.3673, 4.3436}, {5.3346, 4.6776}, {1.642, 4.8237}, {1.164, 5.8307}, {1.2296, 5.4759}, {0.208, 1.7242}, {0.3673, 2.2147}, {2.9166, 3.0809}, {0.6766, 1.6303}};

Interval of the maximum number of iterations for each initial condition:

In[70]:=
nintervals = ConstantArray[{1, 5000}, Length[X0]]
Out[70]=

Step size of iterations for each initial condition:

In[71]:=
nsteps = ConstantArray[1, Length[X0]]
Out[71]=

Orbits data:

In[72]:=
ResourceFunction[
  "IteratedMap2D", ResourceSystemBase -> "https://www.wolframcloud.com/objects/resourcesystem/api/1.0"][StandardMap[X][1.023], X0, X, {n, nintervals}, nsteps] // Short
Out[72]=

Visualize the orbits:

In[73]:=
Graphics[
 MapIndexed[{PointSize[0.003], ColorData["DarkBands"][#2[[1]]/Length[X0]], Point[#1]} &, %], Axes -> True, AspectRatio -> 1, Background -> Darker[Gray, 0.9], AxesStyle -> White]
Out[73]=

Gumowski-Mira Attractor (5) 

The Gumowski-Mira map is defined by 𝓍n+1=𝒷 𝓎n+𝒻(𝓍n) , 𝓎n+1=𝒻(𝓍n+1)-𝓍n, where :

In[74]:=
f[x_[n_]][a_] := a x[n] + (2 (1 - a) x[n]^2)/(1 + x[n]^2)^2;
GMAttractor[{x_[n_], y_[n_]}][{a_, b_}] := {x[n + 1] == b y[n] + a x[n] + (2 (1 - a) x[n]^2)/(1 + x[n]^2)^2, y[n + 1] == ReplaceAll[
       x[n] -> b y[n] + a x[n] + (2 (1 - a) x[n]^2)/(1 + x[n]^2)^2][
      f[x[n]][a]] - x[n]};
X0 = CompressedData["
1:eJwBwQE+/iFib1JlAgAAABsAAAACAAAAmpmZmZmZ2T+amZmZmZnpP3im89iz
m+Q/pJU/eXFZ8D8A0nv1WufkPzbXUirsFPM/aMFoEnKH5j/z+xHTiVTxP6Bd
al/YS9G/hXUPYVefxj9gSs8tzSLJvxRij9pU8rc/uO105ZUw1z/qfuPiMW68
P8CI2ReBAsI/1DYRTok4sT/AQggq7N6ZP7wI/P6Mv6I/gO5KipLBqD/eWnqY
XdO2PyAOF/XqkuM/Gmfl3GhhxD9wPESZBozUv6A6cfTdIek/BE+hRKfV8z9t
fPbEm2npP97QD+C9EvU/VFGZd+Ph6D+2NLWPaS/3P4gIkRSQbug/8InNFKbF
5z/U57zzY4H5P0R5mB7VtgBAonow32bU+D90IQtGnBfyvxjyAfT56PI/WP8E
YwqD4z9k3ou/CqsBQI5FyTBZpgRAfUZjNnez+7+K2HL/SE8FQKih69vvk/q/
EFXx6uZj7z9Mg/MnlfkCwEJI7EwxPvI/8pVS2/rTA8BesHjJA27zPw2ppQJV
CwXAQqcGdOBe9T+CS379GG4GwKrpYsPETv0/gjM36v2cBcDCSx2XAvf/vzRg
qa74owJAZtDkwQ==
"]
Out[76]=

Interval of the maximum number of iterations for each initial condition:

In[77]:=
nintervals = ConstantArray[{1, 5000}, Length[X0]]
Out[77]=

Step size of iterations for each initial condition:

In[78]:=
nsteps = ConstantArray[1, Length[X0]]
Out[78]=

Orbits data:

In[79]:=
ResourceFunction[
  "IteratedMap2D", ResourceSystemBase -> "https://www.wolframcloud.com/objects/resourcesystem/api/1.0"][GMAttractor[X][{0.266, 1}], X0, X, {n, nintervals}, nsteps] // Short
Out[79]=

Visualize the orbits:

In[80]:=
Graphics[
 MapIndexed[{PointSize[0.003], ColorData["DarkBands"][#2[[1]]/Length[X0]], Point[#1]} &, %], Axes -> True, AspectRatio -> 1, Background -> Darker[Gray, 0.9], AxesStyle -> White]
Out[80]=

Properties and Relations (4) 

The map exhibit a Neimark-Sacker bifurcation for 𝓇>𝓇0 at . IteratedMap2D can be used with the resource function PlotGrid to exhibit in a clear fashion a simple bifurcation diagram for the above map:

Orbits data for 𝓇<r0:

In[81]:=
X1s = {{0.30, 0.25}, {0.44, 0.22}}
Out[81]=
In[82]:=
Short[sdata1 = ResourceFunction[
   "IteratedMap2D", ResourceSystemBase -> "https://www.wolframcloud.com/objects/resourcesystem/api/1.0"][{x[n + 1] == 1.98 x[n] (1 - y[n]), y[n + 1] == x[n]}, X1s, {x[n], y[n]}, {n, {{1, 200}, {1, 500}}}, ConstantArray[1, Length[X1s]]]]
Out[82]=

Orbits data for 𝓇>r0:

In[83]:=
X1u = {{0.495, 0.495}, {0.485, 0.485}, {0.44, 0.22}}
Out[83]=
In[84]:=
Short[udata1 = ResourceFunction[
   "IteratedMap2D", ResourceSystemBase -> "https://www.wolframcloud.com/objects/resourcesystem/api/1.0"][{x[n + 1] == 2.012 x[n] (1 - y[n]), y[n + 1] == x[n]}, X1u, {x[n], y[n]}, {n, {{1000, 1500}, {1, 500}, {1, 200}}}, ConstantArray[1, Length[X1u]]]]
Out[84]=

Visualize the bifurcation diagram:

In[85]:=
ResourceFunction[
  "PlotGrid"][{{Graphics[
    MapIndexed[{PointSize[0.006], ColorData["Rainbow"][#2[[1]]/Length[X1s]], Point[#1]} &, sdata1], {Frame -> True, FrameTicksStyle -> GrayLevel[0], FrameLabel -> {
TraditionalForm[r < Subscript[r, 0]], None}, Axes -> False, AspectRatio -> GoldenRatio^(-1)}], Graphics[
    MapIndexed[{PointSize[0.006], ColorData["Rainbow"][#2[[1]]/Length[X1u]], Point[#1]} &, udata1], {Frame -> True, FrameTicksStyle -> GrayLevel[0], FrameLabel -> {
TraditionalForm[r > Subscript[r, 0]], None}, Axes -> False, AspectRatio -> GoldenRatio^(-1)}]}}, ImageSize -> 700]
Out[85]=

Publisher

Ramón Eduardo Chan López

Version History

  • 1.0.0 – 16 September 2022

Source Metadata

Related Resources

License Information