Basics of Non Linear Dynamics

Renan Cabrera Lafuente
Academia Nacional de Ciencias de Bolivia

renanbo@hotmail.com
rencabla@ceibo.entelnet.bo
La Paz Bolivia -2000

Abstract 
The purpose for this document is to introduce to some topic of Non Linear Dynamics by a practical aproach through a series of computer programs in that can be used both to satisfy personal curiosity and also to face real problems in a direct way.

Needs["Transmat`" ,"f:\\renan\\fiscom\\book\\Transmat.m"]

Needs["nltools`" ,"f:\\renan\\fiscom\\book\\nltools.m"]

<< Graphics`PlotField`

Needs["ProgrammingInMathematica`RungeKutta`","f:\\renan\\fiscom\\book\\RUNGEKUTTA.m"]

Today Non Linear Dynamics is regarded as one of the most prolific fields of study because of the broad area of application and the facilities that arises from the increace of  computer's speed  as well  as the development of versatile softwares that are  able to make not only numerical but also symbolic calculus.

Non Linear Finite Difference Equations

A finite Diference Equation is a relation that must be obeyed between two or more consecutive discrete elements of a series.  A Non linear finite Difference Equation  is such that there must exist at least one non linear term in the equation.
A particular interesting one is:

(Subscript[x, t + 1]) = R (Subscript[x, t]) (1 - (Subscript[x, t]))

Which can be expresed as follows

f[R_][x_] = R x (1 - x);

And its graph for a particular R value is:

Plot[ f[2][x] ,  {x, 0, 1}, Frame -> True, 
      FrameLabel -> {(Subscript[x, t]), (Subscript[x, "t+1"])}, RotateLabel -> False]

Now we can solve this equation by nesting f[R] for a  initial x value.

x[1]  = 0.1 with R = 1.5 and for n=25 steps.

\(\(\ solf\  = NestList[\ f[1.5], 0.1, 25]; 


    ListPlot[solf\ , Frame -> True, 
      FrameLabel -> {"t", \*""\!\(x\_t\)""}, 
	
      RotateLabel -> False, 
	PlotStyle -> {PointSize[0.02], Hue[0]}, 
      PlotLabel -> "Steady Aproach Solution"]\)\)

The solution is greatly changed just by changing the parameter R.  
  x[1]  = 0.1 with R = 3.3

\(\(\ solf\  = NestList[\ f[3.3], 0.1, 25]; 

    ListPlot[solf\ , Frame -> True, 
      FrameLabel -> {"t", \*""\!\(x\_t\)""}, 
	
      RotateLabel -> False, PlotStyle -> {Hue[0], PointSize[0.02]}
	, 
      PlotLabel -> "Periodic Solution"]\)\)

Another interesting solution for R =4.

\(\(\ solf\  = NestList[\ f[4], 0.1, 40]; 

    ListPlot[solf\ , Frame -> True, PlotJoined -> True, 
      FrameLabel -> {"t", \*""\!\(x\_t\)""}, 
	
      RotateLabel -> False, PlotLabel -> "Chaotic Solution"]\)\)

Stability

The question is whether or not there exist a steady solution for a particular parameter value. The answer can be found by finding fixed points (Subscript[x, t]) such as

(Subscript[x, t]) = f ((Subscript[x, t]))

Which can be understood as the point where a straight line with slope =1 cuts the f ((Subscript[x, t])) function.

\(Show[

      Plot[\ \(f[2.9]\)[x]\ , \ {x, 0, 1}, Frame -> True, 
        FrameLabel -> {x\_t, x\_"t+1"}, RotateLabel -> False, 
        DisplayFunction -> Identity, 
		PlotStyle -> {Hue[0]}]
, 
	
      Plot[x, {x, 0, 1}, DisplayFunction -> Identity]
, 
      DisplayFunction -> $DisplayFunction, GridLines -> Automatic]\)

Put the cursor over the graphic and push  Crtl to see in the down bar the values pointed out by the cursor.

In order to find an analytical fixed point we have to solve equation (2)

FP  =  Drop[Solve[x == f[2.9][x], x] , 1]

But we drop the obvious  x = 0 solution.

Despite of the existence of fixed points we dont know yet if they are  able to allow a steady solution. The called local stability depends of the slope in the particular fixed point.  The local stability requires a slope "m" be  |m| < 1. 
For instance, the later fixed point is stable as we can see below

D[ f[2.9][x] , x ]  /.  FP

\(solf\  = NestList[\ f[2.9], 0.5, 25]; 

    ListPlot[solf\ , Frame -> True, 
      FrameLabel -> {"t", \*""\!\(x\_t\)""}, 
	
      RotateLabel -> False, PlotRange -> {0, 1}, 
      PlotStyle -> {Hue[0], PointSize[0.02]}
	, 
      PlotLabel -> "Steady Solution around the fixed point"]\)

Bifurcation

As we saw, solutions may vary suddenly as we change the parameter. The value of the parameter at theses  points are called critical points. A importan type of critial point is the bifurcation point. In our case a bifurcation point is where the solution duplicate its period. To see this we can make a Feigenbaum graph.

\(\(Rmin\  = \ 1; 
Rmax\  = \ 4; 
dR\  = \ 0.01; 

Nit\  = \ 200; 

    SampleNit = 50; 


    Feig\  = \ 
      Flatten[Table[
	Drop[NestList[\ f[R], 0.1\ , Nit], Nit - SampleNit], 
          
	{R, Rmin, Rmax, dR}]
	]; 

\)\)

\(\(FeigPair\  = \ 
      Transpose[{

          Flatten[
	
            Table[\ R\ , {R, Rmin, Rmax, dR}, {SampleNit + 1}]\ \ \ \ \ 

            		\ \ \ \ \ \ ]
		, Feig}]; \)\)

\(ListPlot[\ 
FeigPair\ 	, PlotStyle -> {PointSize[0.002]}, 
	
      Frame -> True, FrameLabel -> {"R", " "}\ ]\)

\(Rmin\  = \ 2.9; 
Rmax\  = \ 4; 
dR\  = \ 0.002; 
Nit\  = \ 200; 

    SampleNit = 70; 


    Feig\  = \ 
      Flatten[Table[
	Drop[NestList[\ f[R], 0.1\ , Nit], Nit - SampleNit], 
          
	{R, Rmin, Rmax, dR}]
	]; \)

\(\(FeigPair\  = \ 
      Transpose[{

          Flatten[
	
            Table[\ R\ , {R, Rmin, Rmax, dR}, {SampleNit + 1}]\ \ \ \ \ 

            		\ \ \ \ \ \ ]
		, Feig}]; \)\)

\(ListPlot[\ 
FeigPair\ 	, PlotStyle -> {PointSize[0.002]}\ , 
      Frame -> True, FrameLabel -> {"R", " "}]\)

One Dimensional Differential Equations

The study of a differential equation can be summarized to the knoledge of its behavior around its stability points and the dynamics as the equation parameter goes to the infinity  t -> ∞.
Lets take the following case.

(d/(dt)) x  =  f(x)  =  x  +  α   (x^3)

f[x_]  =  x + α  (x^3);

The fixed points x are those which make f(x) = 0

FixedP = Solve[f[x] == 0, x]

Then α must be negative (α < 0) in order to have fixed points others that the origin

Next, in order to find out the stability of the fixed points we just examine ((df(x))/(dx))

D[f[x], x]  /.  FixedP

Then origin is the only unstable fixed point and the others are stable

Now lets find a solution for a particular value of α = -1

\(\(	{\ D[x[t], t]\  == x[t] + f[\ x[t]\ ] /. \ α -> \(-1\)\ , 
      x[0] == 1}\)\)

\(\(
\(Solx\  = 
      NDSolve[\ 

        	{\ D[x[t], t]\  == f[\ x[t]\ ] /. \ α -> \(-1\)\ , 
            x[0] == 1/2}\  // Evaluate, x, {t, 0, 6}\ \ ]; \)\)\)

\(Plot[\ \ Evaluate[x[t] /. Solx]\ \ \ \ , {t, 0, 5}, Frame -> True, 
      PlotRange -> {0, 1.5}, GridLines -> Automatic, 
	
      FrameLabel -> {"t", "x[t]"}, RotateLabel -> False, 
      PlotStyle -> {Thickness[0.008]}]\)

But for this particular differential equation  we certainly have a exact solution

\(DSolve[\ 

      	{\ D[x[t], t]\  == f[\ x[t]\ ] /. \ α -> \(-1\)\ , 
        x[0] == 1/2}, x[t], t]\ 
\)Plot[x[t] /. %, {t, 0, 5}, Frame -> True, PlotRange -> {0, 1.5}, 
      GridLines -> Automatic, FrameLabel -> {"t", "x[t]"}, 
      RotateLabel -> False]

Two and Three Dimensional Differential Equations

Solving Systems of Linear Differential Equations

Introduction

Solving Systems of Linear Differential Equations is a very important task because nonlinear differential equations ca be approximed by Linear ones in small neighboorhoods where we would be interested on. A two dimensional Linear Differential Equation ca be expresed as:

((d )/(dt))(({{(Subscript[X, 1])}, {(Subscript[X, 2])}}))= (({{(Subscript[a, 1, 1]), (Subscript[a, 1, 2])}, {(Subscript[a, 2, 1]), (Subscript[a, 2, 2])}}))(({{(Subscript[X, 1])}, {(Subscript[X, 2])}}))

This type of differential equations is able to be solved exactly with Mathematica software. The direct way to face this equation is by using the builted in DSolve Mathematica function as it shows below.

\(\(DSolve[\ 

      	{Dt[x[t], t]\ \  == \ x[t]\ \  - 2\ \ y[t]\ , 
        Dt[y[t], t]\ \  == \ x[t]\ \  + 3\ \ y[t]\ }, 
	{x[t], y[t]}, t]\ 
    \)\)

Although we can solve this type of equation as shown above, it could be useful another method by solving it step by step. 
The whole system is  characterized by the matrix:

A =(({{(Subscript[a, 1, 1]), (Subscript[a, 1, 2])}, {(Subscript[a, 2, 1]), (Subscript[a, 2, 2])}}))

One elegant method  says that it could be solved by finding its eigen-equation solution.

A (Subscript[(Overscript[ξ, ⇀]), n])  =  
      (Subscript[λ, n]) (Subscript[(Overscript[ξ, ⇀]), n])

Then forming a transformation Matrix.

(Subscript[U, i, n])  =  (Subscript[(((Subscript[ξ, i]))), n])

And the new system equation (in the eigenvector coordinate system ) to be solved shall be.

((d )/(dt))(({{(Subscript[Y, 1])}, {(Subscript[Y, 2])}}))= (({{(Subscript[λ, 1]), 0}, {0, (Subscript[λ, 2])}}))(({{(Subscript[Y, 1])}, {(Subscript[Y, 2])}}))

Which is easy to solve (more easy when all eigenvalues are different ) and the old variables are:

(Overscript[X, ⇀])  =  U (Overscript[Y, ⇀])

A little different method (it is easy to see that it works out for more dimensions than two) can be achieved by defining the Transmisión Matrix Φ which obeys the following  relations:

((d )/(dt))(Overscript[X, ⇀]) =  A (Overscript[X, ⇀])
                  (Overscript[X, ⇀]) =  Φ  (Subscript[(Overscript[X, ⇀]), 0]) 
                  Φ =  (e^( A t))

Which (Subscript[(Overscript[X, ⇀]), 0]) stands for the  initial conditions.

Then     Φ = U \(Φ\&^\) (U^(-1))

Which  \(Φ\&^\) is the Transmision Matrix in the eigenvectors system coordinate, of course  only when there are n different eigenvectors with n equal to the matrix dimension.
If the later were not the case then Φ is the Transmision Matrix in the Jordan decompostion matrix system and U is the Jordan similarity decompostion matrix.

All the described procedure can be made all by the LinearTransmissionMatrix[A,t] function implemented in the package Transmat.m loaded at the begining of this document.

Needs["Transmat`" ,"f:\\renan\\fiscom\\book\\Transmat.m"]

Example 1

A  =  {{-1, 0, 0}, {0, -4, 4}, {0,  -1, 0}}

LinearTransmissionMatrix[A, t]

Example 2

B  = {{1, -2}, {1, 3}}

FB  = LinearTransmissionMatrix[B, t]

We can now graph ((d )/(dt))(({{(Subscript[X, 1])}, {(Subscript[X, 2])}}))

PlotVectorField[ B.{x, y} ,   {x, -1, 1}, {y, -1, 1} , 
      ScaleFactor -> 0.5 ]

Example 3

\(Clear[\ B\ , FB]; 
B\  = {{\(-1\), 3}, {\(-1\), 1}}; 

    FB\  = LinearTransmissionMatrix[B, t]\)

PlotVectorField[ B.{x, y} ,   {x, -1, 1}, {y, -1, 1} , 
      ScaleFactor -> 0.5 ]

Example 4

\(Clear[\ B\ , FB]; 
B\  = {{\(-2\), \(-1\)}, {3, \(-1\)}}; 

    FB\  = LinearTransmissionMatrix[B, t]\)

PlotVectorField[ B.{x, y} ,   {x, -1, 1}, {y, -1, 1}  , 
      ScaleFactor -> 0.5]

Example 5

\(Clear[\ B\ , FB]; 
B\  = {{\(-2\), 5/7}, {7, \(-3\)}}; 

    FB\  = LinearTransmissionMatrix[B, t]\)

PlotVectorField[ B.{x, y} ,   {x, -1, 1}, {y, -1, 1} ,  
      ScaleFactor  ->  0.5]

Example 6

\(Clear[\ B\ , FB]; 
B\  = {{\(-1\), 2}, {1, 1}}; 

    FB\  = LinearTransmissionMatrix[B, t]\)

PlotVectorField[ B.{x, y} ,   {x, -1, 1}, {y, -1, 1} ,   
      ScaleFunction  ->  (1  # &), ScaleFactor -> 0.5]

The Lotka-Volterra equations

Let x be the population of a prey species, and  y the populaion of a predator species. Both of them are in close and entagled relation by their interrelation.   Lotka and Volterra had proposed a system of differential equations.

(d/(dt))x = α x  - β x y 

(d/(dt))y = γ x y  -  δ y

Where  α , β , γ and  δ are positive constants.

The system is characterized by the vector

LV  = 
       { α x  - β x y ,     
        γ x y  -  δ y  };

The fixed points are

LVFixedP  = Solve[ LV  == {0, 0} , {x, y}]

Now we have to expand the Lotka Volterra equation around the fixed Points by taking the linear part of its Taylor series.

Taylor[LV, {x, y}, {0, 0}, 1]

\(\(LVM\  = {{α, 0}, {0, \(-δ\)}}; 

    LinearTransmissionMatrix[LVM, t]
\)\)

Then the origin is a saddle point

\(\(	\ \ 
    Taylor[LV, {x, y}, {0, 0}, 1]\  /. {α -> 1, \ δ -> 1}\)\)

\(PlotVectorField[

      Evaluate[\ 
        Taylor[LV, {x, y}, {0, 0}, 1]\  /. {α -> 1, \ δ -> 1}], 
      \ \ {x, 0, 1}, {y, 0, 1}\ , \ \ ScaleFunction\  -> \ \((1  #\ &)\), 
      ScaleFactor -> 0.2, PlotLabel -> {"x=Prey", "y=Predator"}]\)

Taylor[LV, {x, y}, {(δ/γ), (α/β)}, 1]

\(PlotVectorField[

      Evaluate[\ 
        Taylor[LV, {x, y}, {δ\/γ, α\/β}, 1] /. {
            α -> 1, \ δ -> 1, β -> 1, γ -> 1}], 
      \ \ {x, 0, 2}, {y, 0, 2}\ , \ \ ScaleFunction\  -> \ \((1  #\ &)\), 
      ScaleFactor -> 0.3, PlotLabel -> {"x=Prey", "y=Predator"}]\)

Then the point is a focus center

The Lorentz Equations

The most famous nonlinear system of ordinary differential equations is the Lorentz equations system. As in the later case the system can be characterized by a vector.

Lorentz  =  {a (y - x),  x (b - z) - y, x*y - c*z}

The fixed Points are then

Solve[Lorentz  == 0 , {x, y, z}]

This 3 dimensional equation is best understood by finding a specific numerical solution

\(Lorentz0\  = \ Lorentz /. {a -> 10, b -> 28, c -> 8/3}; 


    LorentzList = 
	
      RKSolve[
		Lorentz0, {x, y, z}, {0.5, 0.5, 0.}, 

        		{t, 0, 60, 0.02}]; \)

\(Show[Graphics3D[Line[LorentzList]], 
	
      PlotRange -> {{\(-30\), 40}, {\(-30\), 30}, {0, 60}}]\)

The 3 dimensional equations may show chaos as in this case.

The Rossler Equations

The most famous nonlinear system of ordinary differential equations is the Lorentz equations system. As in the later case the system can be characterized by a vector.

Rossler  = 
       {-a*(y + z),  x  +  b*y ,  c +  z*(x - d) }

The fixed Points are then

Solve[Lorentz  == 0 , {x, y, z}]

This 3 dimension equation is best understood by finding a specific numerical solution

\(Rossler0\  = \ Rossler /. {a -> 1, b -> 0.3, c -> 0.4, d -> 5.7}; 


    RosslerList = 
	
      RKSolve[
		Rossler0, {x, y, z}, {0.5, 0.5, 0.5}, 

        		{t, 0, 60, 0.03}]; \)

\(Show[Graphics3D[Line[RosslerList]], 
	
      PlotRange -> {{\(-15\), 20}, {\(-15\), 15}, {0, 40}}]\)

Bibliography

Daniel Kaplan and Leon Glass, (1995), Understanding Non Linear Dynamics ,Springer Verlag.

Kiseliov, Krasnov, Makarenko (1979),  Problems of Differential Equations, Mir Publishing