Physics 213 Wednesday 10-10-03 PHASE SPACE AND DISSIPATION - We have now seen numerous examples of nonlinear behavior, different kinds of complex dynamics (periodic, quasiperiodic, nonperiodic); transitions or bifurcations between these dynamics; and some ways to characterize the nontransient dynamics with Poincare sections and power spectra. I would like to start giving you the mathematical insight to understand where much of this remarkable behavior comes from. As you will see, the mathematics is not so much hard as new and different. Along the way, we will discover specific items that we can compute analytically or, more typically, numerically to get much insight about many nonlinear systems. Much of this material is what active researchers in nonlinear dynamics regard as the "classical" stuff needed to follow most current research in nonlinear dynamics. - We will use a modern GEOMETRIC way of thinking about nonlinear systems, an approach pioneered by the French mathematician and theoretical physicist Henri Poincare at the turn of the century. The idea is to give up trying to find specific analytical solutions or approximations to nonlinear systems, but instead to consider families of nonlinear solutions. As we will see, one can often make powerful statements about families of solutions without knowing any particular solution. For example, we will prove for some flows that all solutions are bounded for all time without knowing any particular solution. We will be able to show that some flows can never have periodic or quasiperiodic behavior without knowing a any explicit solutions. - This geometric approach will lead us to several concepts that are useful in many different contexts in the sciences and in modeling, and that are quite exciting as novel ways of thinking about objects it the world aroundu usP - phase space, - dissipation, - fractals, - fractal dimensions. - embeddings of time series, This geometric viewpoint is essential for understanding how deterministic dynamics (flows or maps) can produce random time series with broad-band power spectra. - Here is a quick outline of what is to come. A mathematical theorem for existence and uniqueness of solutions for flows suggests the concept of a phase space or state space. Solutions to a given flow (orbits) are then geometric objects in this phase space, e.g., a constant solution is a point, a periodic solution is a closed loop, and a quasiperiodic dynamics with two incommensurate frequencies is a topological donut. We will see that DISSIPATION IN PHYSICAL MODELS CORRESPONDS TO AVERAGE SHRINKING OF PHASE SPACE VOLUMES ALONG AN ORBIT. Dissipation leads to the concept of attractors and basins of attraction, that initial conditions in phase space always ends up in a special region of phase space where non-transient behavior is observed. Because phase space contracts, these attractors have zero volume (in phase space). Because phase space volumes do not contract uniformly, the attractors can be FRACTAL sets of points and this is what specifically leads to chaotic behavior (average exponential magnification of small perturbations). The degree of complexity of an attractor can be determined analytically and empirically by calculating one of many kinds of fractal dimensions for the attractor. - The key mathematical elements that will enter into our arguments are relatively simple, e.g., - consider phase space and orbits in phase space - continuity or smoothness of the flow or map; - boundedness of solutions; - nonuniform contractions of phase space volumes (dissipation) It is rather astounding how much these simple elements restrict the possible behaviors of general equations, leading to interesting dynamics. In particular, continuity, boundedness, dissipation and local instability are all that is needed to have chaos. AUTONOMOUS FIRST-ORDER FLOWS: - We will begin by discussing the "easiest" case of flows (continuous state, continuous time), then show later how maps are related to such flows. The discussion for maps will also provide a good review of the concepts we will be working out for flows and the repetition should help you to appreciate the math and the ideas more. - As we discussed earlier in the semester, all flows can be written in a FIRST-ORDER AUTONOMOUS form dX/dt = F(X) where X is an N-dimensional vector of time: X = ( x1(t) , ..., xN(t) ) , and where F(X) is an N-dimensional vector function of the vector X: F(X) = ( f1(X) , ..., fN(X) ) . In this form, mathematicians have proven a basic theorem of existence and uniqueness of solutions depending on the properties of the component functions of F, and most software writers assume this first-order autonomous form for integrating arbitrary sets of ordinary differential equations (odes). - Recall that the system dX/dt=F(X) is a FIRST-ORDER set of equations because at most first-order time derivatives appear. The odes dX/dt=F(X) are also called AUTONOMOUS because the right side vector function, F(X), does not depend on the integration variable t (time) explicitly. PHASE SPACE, PHASE SPACE DIMENSION, AND ORBITS - Given a N-dimensional autonomous flow dX/dt = F(X), the vector X(t) can be considered as living in a space, called PHASE SPACE, with orthogonal coordinate axes that are labeled by the individual variables x1, x2, ..., xN. A better phrase than "phase space" (which is an historic term from several centuries ago) would be STATE SPACE, since the basic theorem we describe below says that a point in phase space gives a complete description of the system, i.e., constitutes a dynamical state. - Phase space generally has different physical units for each axis, e.g., for an harmonic oscillator, one variable might be position with units of meters, the other variable might be speed with units of m/s. value in time t defines a vector X(t) which identifies a point in the phase space. - The number of components N of the vector X appearing in the first-order autonomous flow is called the PHASE-SPACE DIMENSION of the flow and is one of many dimensions that we will discuss throughout this course, e.g., capacity, topological dimension, Hausdorff dimension, correlation dimension, and Lyapunov dimension. - In a N-dimensional phase space, a solution X(t) to the flow traces out a N-dimensional curve called an ORBIT. As a set of points, an orbit does not have time or a direction associated with it. Simple examples of orbits are points (time-independent solutions), closed loops (periodic solutions), and distorted donuts (quasiperiodic solutions with two incommensurate frequencies). A mathematical theorem states that a donut can not be embedded into a plane. This suggests that QUASIPERIODIC EVOLUTION REQUIRES AN AUTONOMOUS SYSTEM WITH AT LEAST THREE DIMENSIONS, This statement is a simple example of how strong statements can be made without knowing explicit solutions. This statement is, in fact correct: you can not obtain quasiperiodic dynamics from ANY flow with just two variables, no matter how complicated the right side functions. - Some people instead plot what are called INTEGRAL CURVES instead of orbits, which is the solution X(t;X0) against time t, in a N+1 dimensional space. This is not too useful except in one- and two-dimensional phase spaces, but does emphasize the time needed to trace out an orbit. - The above different kinds of orbits hint at what a mathematical theory of bifurcations must look like. As some parameter is varied, a particular kind of orbit in phase space "disappears" and is replaced by a new kind of orbit. E.g., the S -> P transition of Gollub and Benson might correspond to a fixed point disappearing and a limit cycle (closed loop or periodic orbit) appearing. As a simple example, consider the one-dimensional ode: x' = r - x^2 , where r is some real-valued parameter. For r < 0, there are no real-valued fixed points in the one-dimensional phase space. For r=0, there is a single fixed point x=0 (is it stable?), and for r=0 there are TWO fixed points x=Sqrt[r] and x=-Sqrt[r] of which only one is stable (PROVE THIS!). So as the parameter r is varied, the phase space has no fixed points, then one fixed point (best thought of as degenerate), then two fixed points. EXAMPLES OF PHASE SPACE - The logistic flow x' = c x (1 - x) is autonomous and first order, and so its solutions x(t) live in a one-dimensional phase space, the real line. Orbits in this phase space look like line segments or like points (x=0, x=1). - The harmonic oscillator x'' + x = 0, can be put in the autonomous form: x1' = x2 , x2' = - x1 , so solutions live in a two-dimensional phase space whose axes x1 and x2 have different physical units, position and velocity. Orbits look like a point (x1=x2=0) or like circles centered on the origin. - If we add some friction and so obtain a damped harmonic oscillator, x'' + d x' + x = 0 , d > 0, we still have a two-dimensional phase space but now most orbits look like spirals or lines that slowly wind in towards the origin. As we saw earlier it the course, you get spirals if d < 2. - The nonlinear oscillator x'' = - Sin(x) also lives in a two-dimensional phase space but of a different topology. Since Sin(x) is periodic with period 2Pi, the flow can not tell the difference between x values separated by values of 2Pi. The phase space is thus equivalent to a CYLINDER, where height along the cylinder is x2 and angle around the cylinder is x1. - The Lorenz equations are already autonomous and live in a three-dimensional phase space. Lorenz orbits look like points, closed loops, and weird fractal objects such as the famous Lorenz attractor. We will show soon that, for the typical parameters of interest for the Lorenz equations, there can not be any quasiperiodic solutions. THE EXISTENCE-UNIQUENESS THEOREM FOR FIRST-ORDER AUTONOMOUS FLOWS: - The importance of first-order autonomous flows and the significance of writing arbitrary flows in this form derives from a fundamental mathematical theorem giving existence and uniqueness for solutions. Most of you have presumably seen this theorem mentioned in your introductory course on odes, although perhaps not a full proof. - In informal terms, this theorem says the following: - Assume a first-order autonomous flow of the form X' = F(X), with N components. - Assume that a point in the flow's N-dimensional phase space is specified, i.e., you give N numbers as an initial condition, X0 = (x0[1], ..., x0[N]) at time t0 (this is some point in phase space); - Assume that the vector field F(X) is smooth in the vicinity of this initial condition X0. Smooth means F and the N^2 partial derivatives Df_i/Dx_j defining the Jacobian matrix are continuous at X0 and close to X0. - Then: there is an open interval of time containing the initial time t0, that is the inequality | t - t0 | < c , holds for some positive constant c, such that there EXISTS a UNIQUE vector function X(t) with the property that it passes through the initial condition at time t0: X(t0) = X0 , and such that X(t) is a solution of the flow on this open time interval, dX/dt = F(X) for t in (t0-c,t0+c) - A more precise formulation can be found in many textbooks, e.g., "Introduction to Numerical Analysis" by J. Stoer and R. Bulirsch, Springer-Verlag, 1980. The precise theorem allows a weaker condition than smoothness (continuous zero and first derivatives), namely a so-called LIPSCHITZ CONDITION. For most physical and modeling problems, F(X) is smooth and hence Lipschitz. Nonsmooth functions can arise in mechanical engineering through slip-stick friction forces. I don't know of practical ways to test the Lipschitz continuity of a N-dimensional function so the definition, while mathematically satisfying since it is a weaker and so broader requirement than differentiability, is too hard to apply in practice. - The theorem is essential in that it tells us what constitutes a dynamical system: a N-dimensional vector X is a state for the dynamical rule X'=F(X) PROVIDED that F(X) is smooth. The vector X is sufficient information to use all the rules, and the rules can be used to generate a new state. - Please appreciate that this theorem only gives LOCAL existence, that a solution exists for some time just before and just after the given initial time. It says nothing about whether solutions exist long before or long after the initial time, which is a GLOBAL statement. Further, the theorem says nothing about how big the constant c is, i.e., how large is the time interval over which the solution exists. The theorem can not be made stronger in that there are many flows for which solutions exist for short times (locally) but not for long times (globally). E.g., the 1d flow dx/dt = 1 + x^2 , has a right side f(x) = 1 + x^2 that is smooth for any initial condition x0 on the real axis. Yet its solution through the initial value x=0 at time t=0 is easily seen to be the function x(t) = tan(t) , which blows up after a finite time, at t=Pi/2. So this is an example of a solution that can not be extended indefinitely beyond the initial time t=0. Many odes have such so-called REAL-TIME SINGULARITIES and such singularities are generally difficult to predict before hand, for the same reason that most nonlinear equations are hard to analyze. - As we discussed earlier in the course, the class of linear dynamical systems, for which F(X) = M X for some matrix M, is one of the few classes of problems for which GLOBAL solutions can be established. Given any initial condition X0, there exists a unique solution X(t) for all time that passes through X0 at time t=t0. We have proved this earlier by simply writing down an explicit solution in terms of eigenvalues and eigenvectors of the matrix M. We conclude that linear flows are better behaved than nonlinear flows, but of course they are also more boring mathematically and unphysical scientifically since physical systems can't grow exponetially for all time. - The assumption of some kind of smoothness or differentiability is also essential for the uniqueness part of the theorem to hold. For example, consider the autonomous equation dx/dt = x^(1/3) with x(0) = 0 . The function f(x) is continuous but not differentiable (nor Lipschitz continuous) close to x0=0, and in fact there are TWO different solutions, x = 0 and x = (2t/3)^(3/2) , that pass through the same initial condition x=0 and that are solutions of the ode. Thus for this flow, a single number x0 is NOT a state of the system, additional information is needed to pick out a solution. Fortunately, most physical models are sufficiently smooth for the fundamental theorem to apply. IMPLICATIONS OF THE BASIC THEOREM FOR NONLINEAR DYNAMICS: - The first important point of this theorem is that by putting a flow in first-order autonomous form, you know exactly how many numbers must be specified to obtain a unique solution, namely the phase-space dimension, N, of the flow. - Since one obtains a unique solution for any given initial condition X0, it becomes important to specify the dependence on X0 explicitly. A common notation is: X(t) = X(t;X0) which indicates that X(t) is the unique solution satisfying dX/dt=F(X) such that X(t0) = X0. The statement X(t0;X0) = X0 is an identity. - If you understand this notation, it should also be clear to you that the following statement is a tautology: X(t2;X(t1;X0)) = X(t2+t1;X0) , for all times t2 > t1 > t0. That is, the solution that goes through the initial vector X(t1;X0), which, in turn, is the value of the solution X(t;X0) at time t1 which goes through the initial state X0 at time t0, has to be the same as the solution that goes through the initial condition X0 at time t0. This mathematical statement is an obscure but precise way of stating that any point on the solution of a flow determines the entire flow. - As an explicit example of this notation, the harmonic oscillator x'' = - x , can be written in autonomous form as (x1,x2)' = (x2,-x1). An initial condition involves specifying a 2-vector X0 = ( x10, x20 ). Now a general solution to the harmonic oscillator is: x(t) = A Cos[t] + B Sin[t] , since Cos[t], Sin[t] are two independent solutions to a second-order ode. We can solve the equations: x10 = x1(t0) = A Cos[t0] + B Sin[t0] , x20 = x1'(t0) = - A Sin[t0] + B Cos[t0] , for the constants A and B in terms of the initial data x10, x20. This gives us a general solution explicitly in terms of the initial data: x1(t;X0) = ( x10 Cos[t0] - x20 Sin[t0] ) Cos[t] + ( x20 Cos[t0] + x10 Sin[t0] ) Sin[t] and similarly for x2(t;X0), the other component of X. You can see how x1(t0;X0) becomes the correct initial value x10, through the trigonometric relation, Cos^2 + Sin^2 = 1. - In general, such an explicit solution is impossible for most nonlinear flows. But the fundamental theorem says that this kind of dependence exists, whether we can find it or not. DEPENDENCE OF SOLUTIONS ON PARAMETER SPACE - The notation X(t;X0) suggested by the basic theorem is incomplete. Most dynamical systems depend on some k-dimensional vector of parameters P = (p1,p2,...,pk), so that the vector field F(X) should really be written as F(X;P). Then the solution X(t;X0) is also a function of parameters and we should write: X(t) = X(t;X0;P) . We have seen several examples earlier of dynamical systems that depend on parameters. We have also discussed how it is best to transform variables and parameters to dimensionless form so as to reduce the number of effective parameters, i.e., the dimension of the parameter space P. - A technical aside that you can skip on a first reading: If the vector field of the flow F(X;P) is continuously differentiable r times w.r.t. to the various parameters (mathematicians would say that F(X;P) is C^r w.r.t. the parameter vector P), there is another theorem that the solution vector X(t;X0;P) is also C^r w.r.t the components of the vector P. So typically, the solution is a smooth function of parameters, even if we don't know this function explicitly. You can verify this theorem for some of the elementary examples that we have studied, e.g., a damped linear oscillator. - We see that a full description of dynamical systems requires THREE different spaces: - The space of initial conditions, X0, of dimension the phase-space dimension N; - Phase space itself, in which the evolution takes place; - The parameter space, P, which can have an arbitrary finite dimension. Engineers often call parameter space CONTROL SPACE, since it contains the mathematical knobs that "control" the model. - The dependence of nonlinear models on these three spaces leads to a more specific goal for the subject of dissipative nonlinear dynamics: TO UNDERSTAND HOW SOLUTIONS VARY OVER LONG TIMES (STATISTICALLY STEADY, NONTRANSIENT CASE) AS ONE VARIES PARAMETERS AND INITIAL CONDITIONS. We might also include the experimental or computational goals of testing any mathematical theory. UNIQUENESS IMPLIES NON-CROSSING OF ORBITS - I would now like to discuss how the existence and uniqueness theorem for odes imposes strong geometric constraints on orbits. There are two separate implications that we will discuss: 1. Orbits can not cross each other in finite time. 2. The vector field F(X) is always tangent to the orbit X(t). - A first important consequence of the basic theorem is that DIFFERENT ORBITS X(T;X0) CAN NOT CROSS IN PHASE SPACE IN FINITE TIME. This follows from the uniqueness part, since an intersection of two orbits at a point X* in phase space would violate the uniqueness theorem for the initial condition chosen to be X0=X*. - The non-crossing does not hold for maps, X[i+1]=G(X[i]), which explains why maps can have more much more complex dynamics than flows in 1d and 2d phase spaces. Remember again the difference between the logistic flow x' = x (1 - x), and the logistic map, x[n+1] = r x[n] (1 - x[n]). The former has only transients to the solution x=1, the latter has periodic cycles of all order, chaos, and an extremely complex bifurcation diagram as a function of the single parameter "r". - WARNING: Non-crossing of orbits does not hold for PROJECTIONS of orbits on lower-dimensional planes. Thus an orbit (x1(t),x2(t),x3(t)) in a 3d phase space may look like some non-intersecting knot (e.g., a trefoil knot), but the projection on the x1-x2 plane will have self-intersections. Don't get confused by projections. NON-CROSSING COROLLARY: NO PERIODIC ORBITS IN 1D FLOWS - For 1d autonomous flows, dx/dt=f(x), the phase space is simply the real line (all possible values of the variable x(t)). The non-crossing theorem then implies that all bounded dynamics are boring: the solution can only go monotonically up or down in value since a reversal in direction (e.g., x(t) going to the right and then going to the left) would constitute a forbidden crossing in phase space. So no matter how complex the right side function f(x) for the autonomous ode, if f(x) is smooth so that unique solutions exist, then there can never be periodic or quasiperiodic or chaotic behavior. The only non-transient bounded orbits are fixed points, and all bounded orbits look like line segments, the orbits move from the initial point x0 towards some fixed point and take infinitely long to do so. Unbounded orbits look like infinite half-lines. - You should be able to see that it will always take an infinite amount of time for any orbit x(t) to reach a fixed point. This follows from a linearity argument: as the orbit x(t) gets sufficiently close to a fixed point x*, one can replace the dynamics x'=f(x) by the linearized dynamics around x*, y' = f'(x*) y. The solution of this is exponential decay Exp[-t f'(x*)] and it takes an infinite amount of time for an exponential to decay to zero. - The absence of periodic orbits (or any nontrivial nontransient orbit) in 1d can also be understood by another argument. Any flow of the form dx/dt = f(x) can be written in the form: dx/dt = - dE/dx , where an energy-like quantity E(x) is defined by integrating the right side function: E[x] = - Integral[ f(x') dx' , {x', x0, x} ] . Here the integral goes from the initial data x0 to some arbitrary value x. How does E[x] change along an orbit? Let's differentiate E[x] w.r.t. time with the assumption that dx/dt=f(x) is given by following the given dynamical system: dE/dt = dE/dx dx/dt = f(x(t)) dE/dx . But the derivative d/dx of an integral from x0 to x is simply the integrand evaluated at x (PROVE THIS!) and so: dE/dt = f(x(t)) (- f(x)) = - f(x)^2 <= 0 , the time derivative is negative except at a fixed point at which point it vanishes. Thus the quantity E[x] must decrease along all orbits of the dynamics and will continue to decrease until f(x)=0, i.e., until a fixed point is reached. The energy-like quantity is called a LYAPUNOV FUNCTIONAL, named after a Russian engineer Lyapunov at the turn of the century, who showed how energy-like quantities could be defined in many cases and used to establish stability of some dynamical system. We will come back to Lyapunov functionals again when we discuss mathematical neural nets and the Hopfield model. NON-CROSSING COROLLARY FOR ECOLOGICAL MODELS - The 3-species problem that you investigated in an earlier homework assignment had a flow of a common form for population models: d(n1)/dt = n1 g1(n1,n2,n3) , d(n2)/dt = n2 g2(n1,n2,n3) , d(n3)/dt = n3 g3(n1,n2,n3) . where the functions g1, g2, g3 are polynomials and so smooth. An set of odes of this form have to have the nice property that if at initial condition N0 = (n10, n20, n30) starts in the positive octant n1 > 0, n2 > 0, n3 > 0, then it stays it that octant for all time. Perhaps the easiest way to see this is that the planes n1=0, n2=0, and n3=0 are all filled with solutions that stay within those planes, and so it is impossible by the uniqueness theorem for an orbit N(t) = (n1(t), n2(t), n3(t)) to cross into negative values since it would mean colliding with a preexisting solution. POINCARE-BENDIXSON THEOREM: NO QP2 OR CHAOS IN THE PLANE - Let's now consider the implication of bounded unique orbits in the plane, i.e., for two-dimensional phase spaces. You can think of a bounded region as some large irregular circle and you can think of an orbit as being a line that is drawn longer and longer in this fixed bounded region and that can not cross itself. If the orbit is not a fixed point or periodic loop or does not converge to one of these objects, it becomes hard to see how to make the curve run around forever in a bounded region while avoiding itself all the time. Now you can draw such curves with a fine pencil but how does a fixed right side vector F(x,y) = (f1(x,y),f2(x,y)) know to avoid intersecting itself, especially after going away for some long time and coming back? A local Taylor series suggests that it is unlikely that there is enough information in the 2d vector field F(x) to guide an orbit away from itself for all time. - For 2d flows, there is a powerful and famous theorem due to Poincare and Bendixson that effectively says that the only BOUNDED non-transient orbits are fixed points or limit cycles. This, too, is a consequence of non-crossing since if any orbit in the plane crossed itself in finite time, it must be a periodic orbit. You can try to draw nonperiodic bounded orbits and convince yourself that it is difficult to do. See S. Wiggins, "Introduction to Applied Nonlinear Dynamics Systems and Chaos", Springer-Verlag, for a fairly readable statement and proof of this theorem. - The exact statement of the Poincare-Bendixson theorem is a bit more obscure (see Guckenheimer and Holmes, Section 1.8): Any orbit confined to a bounded closed region without fixed points is either periodic or converges to a periodic orbit inside that region. EXAMPLE 7.3.1 FROM STROGATZ Consider 2d dynamical system in polar coordinates: r' = r (1 - r^2) + a r Cos[t] , t' = 1 , with parameter a satisfying 0 <= a < 1. Then we can ask if there is some region such that all orbits stay in that region and there are no fixed points in the region. E.g., let's ask if we can find a radius rmin such that r' > 0 for r > rmin: 1 - r^2 + a Cos[t] > 0 . Since Cos[t] >= -1, we can satisfy this equality for all values of angle t provided: rmin < Sqrt[1 - a] , e.g., we can choose rmin = 0.999 Sqrt[1-a]. Similarly, we can ask if there is a maximum radius rmax such that r' < 0 for all r > rmax and for all angles t: 1 - r^2 + a Cos[t] < 0 . Now we can use the fact that Cos[t] <= 1 which leads to the inequality: rmax > Sqrt[1 + a] , and we could choose a specific value rmax = 1.001 Sqrt[1+a]. So we now have a trapping region: 0.999 Sqrt[1 - a] < r < 1.001 Sqrt[1 + a] , for all angles t. Further, there are no fixed points in this region since t'=1 and t is always increasing linearly with time. So there must be at least one periodic orbit in this region. See Strogatz Section 7.3 for further discussion. - In three- or higher-dimensional phase spaces, it is common to find orbits that are confined to a closed bounded region that does not converge to a fixed point or periodic loop. Examples include tori (quasiperiodic motion) or strange attractors (sets traced out by chaotic orbits, e.g., the butterfly wings traced out by the Lorenz orbit). Useful extensions of the Poincare-Bendixson theorem to more than two phase-space dimensions don't exist. Further, there does not yet exist a classification theorem of what possible nontransient behavior can occur in the higher dimensional phase spaces, the possibilities rapidly become too rich and difficult. REVIEW OF KEY POINTS SO FAR: - Uniqueness theorem for first-order autonomous flows X'=F(X) implies that orbit X(t;X0) can never cross itself. This has strong implications for possible dynamics in 1d and 2d phase spaces: - 1d: only bounded nontransient behavior is a fixed point - 2d: only bounded nontransient behavior is a fixed point or limit cycle (periodic closed loop). This is one way of stating the Poincare-Bendixson theorem. - Any orbit that approaches a fixed point will take an infinite amount of time to do so. This is implied by the uniqueness theorem for odes (orbit can't pass through a fixed point, that would give two distinct orbits going through the same initial condition) but also by linearization of any flow near the fixed point, in which case various coefficients of eigenvectors are decaying exponentially rapidly, which takes infinitely long to reach zero. - Poincare-Bendixson can be used to establish existence of a limit cycle in some bounded closed region IF you can prove that that there are no fixed points in that region AND IF you can show that some given orbit stays within the bounded closed region forever. F(X) OF X'=F(X) IS TANGENT TO ORBITS: - A second important geometric insight about orbits in phase space is this: The vector field F(X) defining a flow dX/dt = F(X) is tangent to the orbit at the point X in phase space. By tangent I mean that the orbit always moves for an infinitesimal distance in the direction of the vector F(X). This is almost obvious if you consider integrating the flow X'=F(X) over a short time step h with the forward Euler algorithm: X(t+h) = X(t) + h F(X) . So to get to a point just a bit further in time (t+h where the time step h is small), you need to move in the direction F(X). More carefully, if a path X(t) is differentiable w.r.t time t, by definition dX/dt points in the tangent direction. But if the path X(t) is given by some flow, then dX/dt = F(X) and so F(X) is the tangent vector to the flow at point X. CHALLENGE: If X(t) is a differentiable path then T(t)=dX/dt lies in the tangent direction. But T(t) is itself a path in phase space since it is a N-dimensional vector evolving in time. So what direction does N(t)=dT/dt (tangent to the tangent) lie in w.r.t the original path X(t)? What direction does dN/dt lie in compared to T(t) and X(t)? How many different directions can you produce by successive differentiations in a N-dimensional phase space? TANGENCY TO ORBITS AND QUALITATIVE SHAPES OF ORBITS - Tangency to orbits will help us to understand stability of fixed points and of periodic orbits: all arrows defined by F(X) point inwards towards a stable object, while at least one arrow points outwards for unstable objects. - One use of tangency is the ability to sketch approximate orbits without knowing any explicit solutions. One simply follows the local arrow heads F(X) in small increments, X goes to X + F(X)dt. This illustrates one way to get qualitative information without a analytical solution. - A simple example: the flow for the harmonic oscillator, F(X) = (x2, -x1) , is normal to (perpendicular to) the vector X=(x1,x2), which immediately gives circular motion in phase space. Similarly, the vector field for a damped harmonic oscillator F(X) = (x2, -a x2 - x1) , with positive damping constant a > 0 points in towards the origin which explains the spiraling motion of the resulting orbit for any nonzero initial condition. - Another simple example is d(x1,x2)/dt = (-x1,-x2). Here the vector is always pointing in towards the origin and all orbits move in to the origin. Since x1=x2=0 is a fixed point of this flow, any orbit takes an infinite amount of time to approach this the origin. This is consistent with the known exponential decay: it takes an infinite amount of time for Exp[-t] to decay to zero. - Mathematica provides a function that allows you to visualize two- and three-dimensional vector fields. For example, you can type: << Graphics` PlotVectorField[ (* vector field for oscillator *) { x2, -x1 }, { x1, -2, 2 }, { x2, -2, 2 } ] PlotVectorField[ (* vector field for strongly damped oscillator *) { x2, - 4. x2 - x1 }, { x1, -2, 2 }, { x2, -2, 2 } ] PlotVectorField[ { x2, - x1 - (2. - x1^2) x2 } , (* Van der Pol field *) { x1, -2, 2 }, { x2, -2 ,2 }, ColorFunction -> Hue (* add some color *) ] Mathematica determines a grid of points. At each point, Mathematica then determines a direction and appropriate length for a small vector that will be drawn at that point. It is much harder to visualize three-dimensional vector fields and is rarely insightful to do so. NULLCLINES FOR TWO-DIMENSIONAL FLOWS - A variation on this theme of following tangent vectors that is useful only for two-dimensional phase spaces is to sketch the NULLCLINES of a two-dimensional vector field F(X). These are the two curves for which x'=0 and y'=0 respectively if the 2d flow is written as x'=f(x,y) and y'=g(x,y). Three observations: - Along the nullcline x'=0, the tangent vector F(X)=(0,g) is always pointing in the y-direction, i.e., the flow is vertical (up or down). - Similarly, along the nullcline y'=0, the tangent vector F(X)=(f,0) is always pointing in the x-direction, i.e., the flow is horizontal (right or left). - Nullclines intersect in fixed points since x'=y'=0 imply F(x,y) = 0. - Here is an illustration of how to use nullclines. Consider the two-dimensional flow: u' = v - g(u) , v' = - u , where g(u) is a CUBIC nonlinearity: g(u) = - u + (1/3) u^3 . What kind of dynamics does this system have? Let's look at the two nullclines: u' = 0 iff v = g(u) , y' = 0 iff u = 0 . So one nullcline is the cubic polynomial v=g(u) and the vector field is strictly vertical on this curve. The vector field has the value: F(u,v) = ( 0, -u ) , along nullcline v=g(u) on this cubic nullcline and so points up for u negative and points down for u positive. The second nullcline is just the vertical axis u=0 along which the flow is strictly horizontal: F(u,v) = ( v, 0 ) , along null cline u=0 The flow points horizontally to the right above the u axis (v > 0) and points horizontally to the left below the u axis (v < 0). So these nullclines suggest some kind of swirling motion going clockwise around the origin, perhaps some kind of periodic behavior. This is indeed the case and it is not too hard to apply the Poincare-Bendixson theorem here and show rigorously that there must be a limit cycle surrounding the origin. This two-dimensional model with a cubic nullcline is a quite important class of models in chemistry and in physiology since it leads to OSCILLATORY media (local oscillations) or to EXCITABLE media in which small perturbations decay but larger perturbations grow to a large extent before decaying. Muscle and brain tissue are excitable media and so would be modeled locally with equations of this kind. - CHALLENGE: Why are nullclines not a useful concept for a three-dimensional phase space? - Readable discussion about nullclines and applications to the Fitzhugh-Nagumo model of electrical conduction in a neuron is given in Leah Edelstein-Keshet's book "Mathematical Models in Biology". TANGENCY OF ORBITS AND BOUNDEDNESS - Nearly all useful models of nature assume that solutions are bounded for all time. (An exception would be hyperbolic wave equations, say partial{U(x,t)}/partial{t} = c partial{U}/partial{x} , in which information like a sound wave can propagate infinitely far away at speed c). It would be nice to have some ways to prove boundedness of solutions for a given mathematical model, then we would know for sure that a certain behavior can't occur, e.g., in numerical or experimental studies. Furthermore, boundedness is an implicit assumption in analyzing many time series, one assumes that the series does not diverge for any observation time. - The fact that F(X) is tangent to orbits provides an opportunity to prove boundedness of ALL solutions for many equations. The idea is to find some N-1 dimensional closed surface that encloses some N-dimensional region in phase space, and then show that all initial conditions on this surface have tangent vectors F(X) that point into the surface. Then any orbit starting on the surface will move inside the surface and can never get out, proving boundedness of all orbits starting on or inside the surface. - To say that a vector F(X) "points inside a surface" is the same as the mathematical condition that F(X) . N < 0 , i.e., the dot product between the vector F(X) and the unit vector N that is normal to the surface at the phase point X, has a negative value. Then F(X) points somewhat in the opposite direction of the normal and any orbit going through X must enter the surface. - BOUNDEDNESS OF THE LORENZ EQUATIONS: In his famous paper, Lorenz was able to prove boundedness of solutions associated with his chaotic orbits for his equations: x' = s ( y - x ) , y' = - y + r x - x z , z' = - b z + x y . He did this by finding an ellipsoidal surface f(x,y,z): f(x,y,z) = (1/(2s)) x^2 + (1/2) y^2 + (1/2) z^2 - (r + 1) z - mu = 0 , where mu >= 0 is some large constant. The normal vector to the surface f at any point is given by Grad[f] (do you remember why?): N propto Grad[f] and Grad[f] = ( f_x, f_y, f_z ) = ( x/s , y , z - (r+1) ) , and so F . Grad[f] = - x^2 - y^2 - b z^2 + (r + 1) b Z . But for large enough mu (size of the ellipsoid), the quadratic terms will dominate the linear z term and this dot product will always be negative. This means that every orbit starting on the surface must enter the ellipsoid and can't get out, which proves boundedness for all ics on or in the ellipsoid. You may like to examine this argument carefully to study for which parameter values of r, sigma, and b the boundedness works. CHALLENGE: Can you convince yourself whether ALL starting points sufficiently far from the origin x=y=z=0 will enter into the ellipsoid? - This trick of boundedness is not always easy to apply, the necessary surface may be not exist or may have a complex shape. As an example, see if you can prove boundedness of all orbits within some region for the famous van der Pol oscillator x'' + r (x^2 - 1) x' + x = 0 , or in first-order autonomous form: x1' = x2 , x2' = r (1 - x1^2) x2 - x1 , where r is some nonnegative constant bifurcation parameter. You should try plotting the vector field F(X) and tracking some orbits along vectors to see if there is a suitable surface in phase space to establish boundedness. - Note that the nonlinear Van der Pol oscillator looks like a damped harmonic oscillator with a damping "constant" that is really a function d = x^2 - r , that depends on the state x(t) itself. For small values of x(t), the damping constant is negative, the system gets excited rather than damped, and x(t) and x'(t) tend to increase in value. For large values of x(t), the damping constant is positive and solutions tend to decay toward the origin. It turns out that all initial conditions (except the fixed point x1=x2=0 at the origin) end up on a UNIQUE periodic orbit. The Van der Pol equation was originally invented to model oscillations in a vacuum tube (we are talking old technology here), and later applied by Van der Pol to beating of the heart. COMMENT: The van der Pol oscillator is just our two-dimensional cubic nullcline problem in disguise. We can write the first two pieces of the van der Pol equation as follows: x'' + r (x^2 - 1) x' = d/dt ( x' + r (x^3/3 - x) ) so if we define: v = x' + r (x^3/3 - x) , and relabel x to be the variable u, we have: u' = v - r (u^3/3 - u) , v' = x'' + r (x^2 -1) x' = - u . We then immediately expect to have oscillatory behavior because of the orientation of the nullclines. TANGENCY TO ORBITS AND NONEXISTENCE OF CYCLES, TORI: - Still another use of tangency is NONEXISTENCE, to demonstrate that periodic orbits or quasiperiodic orbits can't exist in certain parts of phase space. We can do this using the divergence theorem of multivariate calculus: Integral[ G(X) . (dl nhat) , curve ] = Integral[ Div[ G ] , area ] , = Integral[ D(G1)/d(x1) + D(G2)/d(x2), area ] , where dl is an infinitesimal length and nhat is the normal unit vector to the curve. In English, this equation says that the integral along some closed curve C with integrand given by the dot product of a vector field G(X) with the local infinitesimal vector dN that is normal to the curve at point X, is the same as the two-dimensional integral over the area enclosed by the curve C with integrand given by the divergence of the vector field G. - You should review your understanding about the famous divergence theorem (also known as Gauss's theorem) of multivariate calculus if you don't feel comfortable with this equation. A clear and physically motivated derivation of the divergence theorem (and the related Stoke's theorem) is given in Edward Purcell's classic book "Electricity and Magnetism", Berkeley Physics Course, Volume 2, Section 2.10. The physical interpretation of Gauss's theorem is perhaps clearer than its mathematical formulation: the integral over the surface of some region is the total flux of some quantity (e.g., mass or charge or points in phase space) leaving the surface covering some region. The divergence of a vector field div(F) = Trace[Jacobian[F]] is the total flux leaving some infinitesimal volume of space. So the flux through the surface must exactly balance the total flux through all infinitesimal volumes inside the surface. - You are probably most familiar with applications of the divergence theorem to 2d surfaces enclosing a 3d volume, rather than to a closed curve containing a planar area. It is rather easy to deduce the planar version of the theorem from the usual 3d version by taking the 2d vector field F(x,y) and just extending it everywhere above and below the plane, so that the vectors F(x,y) are still parallel to the x-y plane. Then invent a distorted cylinder of height h whose intersection with the plane is the given curve C, and whose top and bottom planes are parallel to the x-y plane. Applying the 3d divergence theorem to the volume contained by this cylinder, we get: Integral[ Surface, F(X) . dA ] = h Integral[ C , F(X) . dN ] since F(X) is perpendicular to the normal infinitesimal vector dA everywhere except on the sides of the cylinder that are perpendicular to the x-y plane. Similarly, for the 3d divergence part of the formula, we find: Integral[ Volume, d(Fx)/dx + d(Fy)/dy + d(Fz)/dz ] = h Integral[ Area, d(Fx)/dx + d(Fy)/dy ] , since Fz is zero by construction and the integral is the same on every plane parallel to the x-y plane. The height factor h now cancels from both sides, giving us the 2d version of the theorem. - Now assume that a periodic closed orbit C appears in some part of 2d phase space and integrate F(X).dN along the boundary of C, giving contour integral. Here dN is the vector normal to the curve with infinitesimal length dl, i.e., dN = N dl, where N is the unit vector normal to C at position X. The normal is, by definition, orthogonal to the tangent vector, and points in the direction exterior to the closed curve. (I'll let you worry about the subtlety of how you know whether some direction is exterior or interior to some long intricate curve.) The 2d version of the divergence theorem then gives: 0 = Integral[ F(X) . dN , Curve ] = Integral[ Div[ F(X) ] dx dy , Area ] = Integral[ ( d(Fx)/dx + d(Fy)/dy ) dx dy , area ] The first line follows from our assumption that the closed curve C is an orbit of the dynamical system dX/dt = F(X) and so F(X).dN = 0 at every point along the curve. This algebra leads to a contradiction if Div[F(x)] is nonzero and has the same sign everywhere inside the periodic orbit C and we have a cute little theorem: If Div[F] has constant sign in a simply connected region R then there can not be any periodic orbits in R. - As a simple example, consider the damped simple pendulum: x'' + b(x) x' + Sin[x] = 0 , with arbitrary positive damping b(x) > 0 that is an ARBITRARY positive function of the position x. This can be written in first-order autonomous form as x1' = x2 x2' = - b(x1) x2 - Sin[x1] , and has no known analytical solution. The flow has a constant and negative divergence everywhere: Div[ F[x] ] = D[ x2, x1] + D[ - b(x1) x2 - Sin[x1], x2 ] = - b(x1) < 0 Bendixson's theorem then proves that this equation has no periodic solutions anywhere in phase space. In fact, the same proof holds if Sin[x] is replaced by an arbitrary nonlinear function, f(x), and so is quite powerful. - This same conclusion does NOT hold for the DRIVEN damped simple pendulum x'' + b x' + Sin[x] = B Cos[w t] . In autonomous first-order form, x1' = x2 , x2' = - b x2 - Sin[x1] + B Cos[w x3] , x3' = 1, the divergence is still negative everywhere, but Bendixson's criterion doesn't work since this flow lives in a THREE-dimensional phase space (the t variable appears explicitly on the right side) and a limit cycle does not need to lie in a plane. And indeed, there are periodic orbits and they are twisted loops that do not lie in a plane. GENERALIZATION OF BENDIXSON'S CRITERION: NO QP2 FOR LORENZ EQS - Can we generalize Bendixson's criterion to higher dimensional phase spaces? We need to find orbits that trace out a CLOSED bounded surface in the limit of infinitely long time, and then apply the divergence theorem to F(X).dA over that closed surface. The surface integral vanishes since F(X) is tangent to orbits and hence orthogonal to the normal infinitesimal-area vector dA, and this gives a contradiction of existence if Div[F] is nonzero with constant sign inside the closed surface. - Since the orbit traces out the closed surface, the surface must have the property of having a FINITE TANGENT VECTOR everywhere. (If the tangent vector were not finite everywhere, it must be zero, in which case the orbit dies upon approaching some fixed point.) This is easily accomplished in 2d, e.g., a circle, but is harder in higher dimensions, e.g., there is a theorem that spheres in odd dimensions (3,5,...) can not have a finite tangent vector field everywhere. This is the so-called "cowlick theorem" since it explains why people must have a singularity (the cowlick) in their hair when trying to comb their hair flat on their essentially spherical head. However, people with four-dimensional heads can comb their hair perfectly flat without any topological singularities. (In two-dimensions, it is also easy to comb hair that is everywhere tangent to a circle.) - The simplest closed surfaces with finite tangent vector fields everywhere are the torus (donut), the 3-sphere in 4-space, and the Klein bottle, also in 4-space. I don't know the answer to the general question of finding all such closed surfaces in higher dimensions, don't even know if the answer is known. INTERESTING TIDBIT: The fact that the torus is the only bounded closed surface in three dimensions with a finite tangent vector field everywhere is the precise reason why most magnetic fusion confinement devices (e.g., tokamaks and stellarators) are donut shaped. The plasma particles flow at near relativistic speeds tangentially along magnetic field lines and so can be confined if these field lines trace out a torus over infinite time. The hard part of designing a magnetic fusion reactor is how to create a magnetic field out of currents so that the magnetic field lines trace out concentric tori. NONEXISTENCE OF QUASIPERIODIC SOLUTIONS IN THE LORENZ EQUATIONS - Let's apply these ideas and show that the Lorenz equations do not have any quasiperiodic solutions for the classical parameters sigma=10 and b=8/3. We consider again the Lorenz equations, x' = - s x + s y , y' = - y + r x - x z , z' = - b z + x y , The divergence of the Lorenz equations is seen to be: Div[ F ] = - (1 + b + s) , which is constant and negative for parameter values b > 0 and sigma > 0 for all values of r. The parameter b is a wavenumber and sigma is the Prandtl number so these physically should be positive. Thus the Lorenz equations can NOT have QP2 solutions for physical values of the parameters, a beautiful result considering how hard these equations are to work with analytically. This also shows that the Lorenz equations are seriously flawed as a model of convecting fluids, since Gollub and Benson's experiments and other experiments show that QP2 behavior is common. - CHALLENGE: I don't know whether there are quasiperiodic solutions to the Lorenz equations for unphysical parameter values, say for b and sigma such that 1 + b + s = 0, in which case the impossibility theorem no longer holds. Would someone like to explore this possibility? DEFINITION OF DISSIPATION: PHASE SPACE CONTRACTION - Now that you have had some experience in thinking about phase space and with orbits in phase space, we are finally in a position to define the important concepts of CONSERVATIVE and DISSIPATIVE nonlinear dynamical systems. The Bendixson criterion for the absence of periodic orbits in a 2d phase space and its generalization via the divergence theorem to higher dimensions suggests the importance of looking at the divergence of the vector field which can be expressed in several ways: Div[F] = Trace[ Jacobian[F] ] = sum of eigenvalues of Jacobian[F] at point X. This quantity turns out to have the meaning of the relative rate of contraction or expansion of infinitesimal volumes in phase space (areas in 2d, volumes in 3d, hypervolumes in N > 3 dimensions) along an orbit X'=F(X). More precisely, we will prove that: Div[ F ] = Tr[ Jacobian[F] ] = (1 / dV) d(dV)/dt , where dV(t) is an infinitesimal volume in phase space centered on the orbit at X(t). I postpone a proof for now so that we can focus on the implications of this fact. More specifically, you should think of every point in an infinitesimal volume dV as being an initial condition for the flow X'=F(X). Then as time evolves, each point moves to a nearby point in phase space along an orbit and the volume dV changes as it is stretched and squished. Div[F] gives the relative rate of change at a particular point in phase space. - We will make the following definitions: An orbit of a flow is DISSIPATIVE if phase space volumes decrease on average along the orbit. < Div[ F[X(t)] ] > < 0 , dissipative orbit An orbit of a dynamical system is CONSERVATIVE if phase space volumes are conserved on average along an orbit. < Div[ F[X(t)] ] > = 0 . conservative orbit An orbit of a dynamical system is EXPANSIVE if phase space volumes increase on average along an orbit. < Div[ F[X(t)] ] > > 0 . expansive orbit Here the notation <...> denotes a time average along the orbit: if A(X) is any function of points X in phase space, then < A(X(t)) > = (1 / T) Integral[ A(X(t)) dt, {t,t0,t1} ] , for a flow X(t) and < A(X[i]) > = (1 / N) Sum[ A(X[i]) , {i, 1, N} ] , for a map X[i]. In both cases, we ideally would try to take the limit of time going to infinity. Note that the definition involves orbits X(t), not the dynamical system dX/dt=F(X) itself. A given dynamical system can simultaneously have dissipative, conservative, and expanding (phase-volume increasing) orbits in different parts of its phase space. As an example, the time average of the divergence of the vector field at a fixed point: < Div[F](X*) > = Div[F](X*) , is the divergence itself since the fixed point doesn't change in time. So a fixed point X* is dissiptive, conservative, or expansive if Div[F](X*) is negative, zero, or positive respectively. It is easy to find systems with fixed points of different divergence simultaneously in phase space. Often "most" orbits are dissipative in which case people refer to the dynamical system as dissipative but this is a bit sloppy. It is also usually the case that a conservative system has a CONSTANT divergence, rather than a fluctuating divergence of positive and negative values of mean zero. We will generally not worry about expansive orbits such that > 0 is positive on the average. Such orbits correspond to unbounded dynamics and hence are rejected in most models of nature. CHALLENGE: Find a flow such that it is conservative (time-averaged divergence is zero) yet the divergence is itself not equal to zero. CHALLENGE: Consider a N-dimensional flow X'=F(X) with fixed point X* that is expansive, i.e., Div[J[F](X*)] > 0. Show that X* must be an unstable fixed point. So expansiveness implies instability for a fixed point. Does the converse hold: if a fixed point X* of a flow is unstable, is the fixed point also expansive? CHALLENGE: Does this hold for maps also? If X*=G(X*) is a fixed point of the N-dimensional map X[n+1]=G(X[n]), and if X* is expansive, |det[J[F](X*)]| > 1, is this fixed point unstable? CHALLENGE: Find an explicit dynamical system X'=F(X) with three distinct fixed points that are dissipative, conservative, and expansive respectively. CHALLENGE: If a periodic orbit X(t) of a flow is expansive, is it also unstable? EXAMPLES OF THE DISSIPATION DEFINITION - The damped harmonic oscillator x'' + d x' - g(x) = 0 for some arbitrary force function g(x) x1' = x2 , x2' = g(x1) - d x2 , has divergence Div[F] = - d , which is constant and negative if the damping constant d is positive. This oscillator is conservative if d=0 (phase space volumes conserved) and is expansive if d > 0. The function g(x) drops out completely in this definition. You can see that all orbits away from the origin converge towards the origin. So if you choose an area of initial conditions away from the origin, it gets squeezed to zero volume as all the initial conditions move in towards the origin x1=x2=0. From this physical example, we can show that dissipation is related to lack of energy conservation but this is misleading since for many dynamical systems (economics, biology, Lorenz equations), an energy can not be defined. The total energy E is given by the sum of the kinetic energy (1/2) (x')^2 and the sum of the potential energy: E = E[x,x'] = (1/2) (x')^2 + h(x) , where h(x) = - Integral[ g(x) dx ] , is the potential energy. But then by the chain rule, dE/dt = ( x'' + h'(x) ) x' = ( x'' - g(x) ) x' = - d (x')^2 < 0 steadily decreases towards E=0 until the fixed point x=x'=0 is reached. Note that if the divergence is negative (d < 0), energy steadily increases leading to diverging position and velocities. Energy arguments are not as useful as our divergence definition since the latter can always be calculated, and for many non-physics models, energy can not be defined (e.g., what is the kinetic energy of a fox population or of gross national product)? - THE VAN DER POL OSCILLATOR x'' + (x^2 - r) x' + x = 0 , or x1' = x2 , x2' = (r - x1^2) x2 - x1 . is a simple example for which Div[F] = r - x1^2 , takes on both positive and negative values with a mean that is negative and so is dissipative. In fact, one can show analytically for small positive r << 1 that x1(t) ~= 2 Sqrt[r] Sin(t) and so: < Div[F]> = < r - x1^2 > = < r - 4 r Sin(t)^2 > ~= - r < 0 , since < Sin(t)^2 > = 1/2. But you can see that the quantity "1 - 4 Sin(t)^2" can go positive when t is close to an integer multiple of 2Pi. Hint of proof that x1(t) = 2 Sqrt[r] Sin[t] to lowest order in the small quantity r: simply substitute this into the van der Pol equation and require that the Fourier coefficient of Sin[t] to lowest order in the small quantity r is zero. - Lotka-Volterra predator-prey equations for two species: x1' = a x1 - b x1 x2 prey equation , x2' = - c x2 + d x1 x2 predator equation , where x1(t) and x2(t) can be population or biomass. Historically, this was the first predator-prey equation to be written down although we know now, with our greater unerstanding of nonlinear dynamics, that this model is unphysical (structurally unstable). The divergence of the flow is given by: Div[F] = (a - b x2) + (-c + d x1) = a - c + d x1 - b x2 , and it is now not clear whether a given orbit is dissipative without further assumptions about the four parameters a, b, c, and d, or without knowing the time dependence of x1 and x2. The Lotka-Volterra equations do have two fixed point solutions: X = ( 0 , 0 ) extinction and X = ( c/d, a/b ) coexistence . Note that we must have c/d >= 0 and a/b >= 0 , in order for the nonzero fixed point to be biologically reasonable, with nonnegative populations. If we evaluate Div[F] at these fixed points, we get "a-c" and "0" respectively. So the extinction fixed point can be dissipative, conservative, or expansive depending on the relative sizes of the parameters a and c, while the coexistence fixed point is conservative. This latter is an unlikely coincidence and is unphysical and so this model is flawed; arbitrarily small perturbations of the equations (e.g., adding a small random function) will change the conservative nature of the coexistence fixed point. - The Lorenz equations have a strange attractor for the parameter values r=28, sigma=10, b=8/3, and for most initial conditions. For these values, the average dissipation is: < (1/dV) d(dV)/dt > = < Tr[Jacobian[F(X)]] > = - (1 + b + sigma ) approx - 13.7 , a constant value. All phase space volumes decay exponentially as V(0) Exp[-13.7 t]. In one unit of time, each volume is shrunk by a factor of about 1.2e-6, one millionth. This is a STRONGLY dissipative system, volumes collapse quite rapidly. DERIVATION OF PHASE-SPACE VOLUME CHANGES: FLOWS - Strogatz gives a global proof of the meaning of Div[F] by studying how a closed surface S in phase space is swept along by the flow F(X). Here I give a local proof by following how an infinitesimal volume changes volume under the flow. Consider a particular point in phase space, X(t), and a small hypercube of phase-space volume whose diagonal is X(t) + DX, where DX = ( dx[1], dx[2], ... dx[N] ) are infinitesimal increments added to each component of X(t). The small volume of this cube is simply the product of all sides, dV = Product[ dx[i] , {i,N} ] If each vertex of this cube is used as an initial condition, which evolves under the flow dX/dt=F(X), the increments and the volume evolve with time. One easily sees that: (1/dV) d(dV)/dt = Sum[ (1/dx[i]) d(dx[i])/dt , {i,N} ] , the relative change in the infinitesimal volume dV is the sum of the relative changes of the infinitesimals that make up the sides of the volume dV. - To find the relative change in volume, we simply need to find the relative change in the increment of each coordinate independently. Let us calculate the relative change in dx[1](t), the others follow by symmetry. Then we have the following approximate relations: at time t increment is dx[1] along x1 direction at time t+dt increment is difference between x1(t+dt) and (x1+dx[1])(t+dt). The flow equation dX/dt=F(X) implies dX = F(X) dt for small time intervals dt. We deduce: x1(t + dt) approx x1(t) + f1[X(t)] dt (x1+dx[1])(t+dt) approx x1(t) + dx[1](t) + f1[X(t) + dX1] dt We can use the Taylor series to lowest order to write: f1[X(t) + dX[1]] approx f1[X(t)] + (Df1/Dx1) dx[1] + h.o.t.s The increment at time t+dt is then: (x1+dx[1])(t+dt) - x1(t+dt) = dx[1](t) + (Df1/Dx1) dx[1] dt The CHANGE in the increment over the time dt is the difference between this increment at time t+dt and the original increment, dx[1](t): d(dx[1]) = (Df1/dx1) dx[1] dt , which proves that: (1 / dx[1]) d(dx[1])/dt = Df1/d1(X) . Adding the similar terms for the other N-1 components then proves the claimed result, (1 / dV) d(dV)/dt = Div[ F ] . This describes how an infinitesimal volume dV centered on a point X in phase space changes its size under the flow F[X].