#set document(title: "1.6 Autonomous equations", author: "Jiří Lebl")
#set page(width: 8.5in, height: auto, margin: 1in)
#import "@preview/cetz:0.5.2"
#set text(font: ("STIX Two Text", "Libertinus Serif", "New Computer Modern"), size: 10.5pt, lang: "en")
#show math.equation: set text(font: ("STIX Two Math", "New Computer Modern Math"))
#set par(justify: true, leading: 0.62em, spacing: 0.9em)
#set enum(spacing: 1.1em)   // room between list items so tall inline fractions don't collide
#set list(spacing: 1.1em)
#set table(stroke: 0.5pt + rgb("#c7ccd3"))

#let BLUE = rgb("#183B6F")     // brand navy — section bars + example/solution labels (white on navy 11.09:1)
#let ORANGE = rgb("#A94509")   // brand primary-700 — AA-safe deep orange for TEXT (5.93:1 on white; raw brand #F37021 is 2.94:1 and must never carry text)
#let RED = rgb("#DC2626")      // brand error-600
#let GREEN = rgb("#059669")    // brand success-600 (decoration only; small green text uses green-text #007942)

#show heading.where(level: 1): it => block(width: 100%, above: 0pt, below: 16pt,
  fill: gradient.linear(BLUE, rgb("#2C5AA0")), inset: (x: 14pt, y: 12pt), radius: 3pt,
  text(fill: white, weight: "bold", size: 19pt, it.body))
#show heading.where(level: 2): it => block(width: 100%, above: 18pt, below: 10pt,
  fill: BLUE, inset: (x: 10pt, y: 6pt), radius: 2pt,
  text(fill: white, weight: "bold", size: 12pt, it.body))
#show heading.where(level: 3): it => text(fill: ORANGE, weight: "bold", size: 12.5pt, it.body)
#show heading.where(level: 4): it => text(fill: BLUE, weight: "bold", size: 10.5pt, it.body)

#let examplebox(label, title, body) = block(width: 100%, breakable: true,
  fill: rgb("#EFF1F5"), stroke: 0.5pt + rgb("#CFDDF0"), radius: 4pt, inset: 10pt,
  above: 12pt, below: 12pt)[
  #block(below: 6pt)[#box(fill: BLUE, inset: (x: 6pt, y: 2pt), radius: 2pt,
    text(fill: white, weight: "bold", size: 8.5pt, label)) #h(0.4em) #strong[#title]]
  #body]

// rail = decorative left rule (raw brand token); labelcolor = AA-safe label text shade
#let notebox(label, rail, labelcolor, tint, body) = block(width: 100%, breakable: true, fill: tint,
    stroke: (left: 3pt + rail), inset: (left: 10pt, rest: 8pt), radius: (right: 4pt),
    above: 11pt, below: 11pt)[
    #text(fill: labelcolor, weight: "bold", size: 7.5pt, tracking: 0.5pt)[#upper(label)]
    #linebreak() #body]

#let solutionbox(body) = block(above: 4pt, below: 8pt)[
  #text(fill: BLUE, weight: "bold", size: 8.5pt)[Solution] #linebreak() #body]

#let figph(msg) = block(width: 100%, height: 60pt, fill: rgb("#f6f7f9"),
  stroke: (paint: rgb("#c7ccd3"), dash: "dashed"), radius: 4pt, inset: 10pt)[
  #align(center + horizon, text(fill: rgb("#889"), style: "italic", size: 9pt, msg))]

// Standardize inlined figure sizes: measure the natural CeTZ canvas, then scale to a
// consistent envelope (aspect-aware; see build_typst.py FIG_* constants). Unlike the
// print preamble, dimensions are FLOORED: in an editor a user can trim a figure to a
// degenerate 1-D shape (a bare line), and w/h or tw/w would then divide by zero.
#let _STD_W = 3.5
#let _WIDE_W = 5.6
#let _MAX_H = 3.4
#let _ASPECT_WIDE = 2.2
#let _UPSCALE_MAX = 1.15
#let stdfig(body) = context {
  let m = measure(body)
  let w = calc.max(m.width / 1in, 0.01)
  let h = calc.max(m.height / 1in, 0.01)
  let tw = if w / h > _ASPECT_WIDE { _WIDE_W } else { _STD_W }
  let s = calc.min(tw / w, _MAX_H / h, _UPSCALE_MAX)
  align(center, box(scale(x: s * 100%, y: s * 100%, reflow: true, body)))
}

#show figure: set block(breakable: false)
#set figure(gap: 8pt)
#show figure.caption: set text(size: 8.5pt, fill: rgb("#555"))


    
== 1.6#h(0.6em)Autonomous equations


    Let us consider general differential equation problems of the form


    #math.equation(block: true, alt: "the fraction d x over d t equals f open parenthesis x close parenthesis")[$frac(d x, d t) = f ( x )$]


    where the derivative of solutions depends only on #math.equation(block: false, alt: "x")[$x$] (the dependent variable). Such equations are called #strong[autonomous] equations. If we think of #math.equation(block: false, alt: "t")[$t$] as time, the naming comes from the fact that the equation is independent of time.


    Let us come back to the cooling coffee problem (see Example 1.3.3). Newton’s law of cooling says that


    #math.equation(block: true, alt: "the fraction d x over d t equals minus k open parenthesis x minus A close parenthesis")[$frac(d x, d t) = − k ( x − A )$]


    where #math.equation(block: false, alt: "x")[$x$] is the temperature, #math.equation(block: false, alt: "t")[$t$] is time, #math.equation(block: false, alt: "k")[$k$] is some constant, and #math.equation(block: false, alt: "A")[$A$] is the ambient temperature. See Figure #math.equation(block: false, alt: "1")[$1$] for an example with #math.equation(block: false, alt: "k equals 0.3")[$k = 0.3$] and #math.equation(block: false, alt: "A equals 5")[$A = 5$].


    Note the solution #math.equation(block: false, alt: "x equals A")[$x = A$] (in the figure #math.equation(block: false, alt: "x equals 5")[$x = 5$]). We call these constant solutions the equilibrium solutions. The points on the #math.equation(block: false, alt: "x")[$x$] axis where #math.equation(block: false, alt: "f open parenthesis X close parenthesis equals 0")[$f ( X ) = 0$] are called critical points. The point #math.equation(block: false, alt: "x equals A")[$x = A$] is a critical point. In fact, each critical point corresponds to an equilibrium solution. Note also, by looking at the graph, that the solution #math.equation(block: false, alt: "x equals A")[$x = A$] is “stable” in that small perturbations in #math.equation(block: false, alt: "x")[$x$] do not lead to substantially different solutions as #math.equation(block: false, alt: "t")[$t$] grows. If we change the initial condition a little bit, then as #math.equation(block: false, alt: "t → infinity")[$t → ∞$] we get #math.equation(block: false, alt: "x → A")[$x → A$]. We call such critical points stable. In this simple example it turns out that all solutions in fact go to #math.equation(block: false, alt: "A")[$A$] as #math.equation(block: false, alt: "t → infinity")[$t → ∞$]. If a critical point is not stable we would say it is unstable.


    #figure(figph[Slope field over time 0 to 20 and values -10 to 10, flat along the level 5. Five red solutions starting at 10, 5, 0, -5 and -10 all bend into the line at 5 by about time 13.], alt: "Slope field over time 0 to 20 and values -10 to 10, flat along the level 5. Five red solutions starting at 10, 5, 0, -5 and -10 all bend into the line at 5 by about time 13.", caption: [Figure #math.equation(block: false, alt: "1")[$1$]: The slope field and some solutions of #math.equation(block: false, alt: "x prime equals 0.3 open parenthesis 5 minus x close parenthesis")[$x^(′) = 0.3 ( 5 − x )$].])


    Let us consider the logistic equation


    #math.equation(block: true, alt: "the fraction d x over d t equals k x open parenthesis M minus x close parenthesis")[$frac(d x, d t) = k x ( M − x )$]


    for some positive #math.equation(block: false, alt: "k")[$k$] and #math.equation(block: false, alt: "M")[$M$]. This equation is commonly used to model population if we know the limiting population #math.equation(block: false, alt: "M")[$M$], that is the maximum sustainable population. The logistic equation leads to less catastrophic predictions on world population than #math.equation(block: false, alt: "x prime equals k x")[$x^(′) = k x$]. In the real world there is no such thing as negative population, but we will still consider negative #math.equation(block: false, alt: "x")[$x$] for the purposes of the math (see Figure #math.equation(block: false, alt: "2")[$2$] for an example).


    #figure(figph[Slope field over time 0 to 20 and values -5 to 10 with red equilibria at 5 and 0. Solutions above 0 head toward 5, those just below 0 plunge off the bottom, so 5 attracts and 0 repels.], alt: "Slope field over time 0 to 20 and values -5 to 10 with red equilibria at 5 and 0. Solutions above 0 head toward 5, those just below 0 plunge off the bottom, so 5 attracts and 0 repels.", caption: [Figure #math.equation(block: false, alt: "2")[$2$]: The slope field and some solutions of #math.equation(block: false, alt: "x prime equals 0.1 x open parenthesis 5 minus x close parenthesis")[$x^(′) = 0.1 x ( 5 − x )$].])


    Note two critical points, #math.equation(block: false, alt: "x equals 0")[$x = 0$] and #math.equation(block: false, alt: "x equals 5")[$x = 5$]. The critical point at #math.equation(block: false, alt: "x equals 5")[$x = 5$]. is stable. On the other hand the critical point at #math.equation(block: false, alt: "x equals 0")[$x = 0$]. is unstable.


    It is not really necessary to find the exact solutions to talk about the long term behavior of the solutions. For example, from the above slope field plot, we can easily see that


    #math.equation(block: true, alt: "5 if x open parenthesis 0 close parenthesis greater than 0 ,; 0 if x open parenthesis 0 close parenthesis equals 0 ,; DNE or minus infinity if x open parenthesis 0 close parenthesis less than 0 .")[$5 "if" x ( 0 ) > 0 , \ 0 "if" x ( 0 ) = 0 , \ "DNE or" − ∞ "if" x ( 0 ) < 0 .$]


    Where DNE means “does not exist.” From just looking at the slope field we cannot quite decide what happens if #math.equation(block: false, alt: "x open parenthesis 0 close parenthesis less than 0")[$x ( 0 ) < 0$]. It could be that the solution does not exist for #math.equation(block: false, alt: "t")[$t$] all the way to #math.equation(block: false, alt: "infinity")[$∞$]. Think of the equation #math.equation(block: false, alt: "x prime equals x squared")[$x^(′) = x^(2)$], we have seen that it only exists for some finite period of time. Same can happen here. In our example equation above it will actually turn out that the solution does not exist for all time, but to see that we would have to solve the equation. In any case, the solution does go to #math.equation(block: false, alt: "minus infinity")[$− ∞$], but it may get there rather quickly.


    If we are interested only in the long term behavior of the solution, we would be doing unnecessary work if we solved the equation exactly. We could draw the slope field, but it is easier to just look at the phase diagram or phase portrait, which is a simple way to visualize the behavior of autonomous equations. In this case there is one dependent variable #math.equation(block: false, alt: "x")[$x$]. We draw the #math.equation(block: false, alt: "x")[$x$]-axis, we mark all the critical points, and then we draw arrows in between. Since #math.equation(block: false, alt: "x")[$x$] is the dependent variable we draw the axis vertically, as it appears in the slope field diagrams above. If #math.equation(block: false, alt: "f open parenthesis x close parenthesis greater than 0")[$f ( x ) > 0$], we draw an up arrow. If #math.equation(block: false, alt: "f open parenthesis x close parenthesis less than 0")[$f ( x ) < 0$], we draw a down arrow. To figure this out, we could just plug in some #math.equation(block: false, alt: "x")[$x$] between the critical points, #math.equation(block: false, alt: "f open parenthesis x close parenthesis")[$f ( x )$] will have the same sign at all #math.equation(block: false, alt: "x")[$x$] between two critical points as long #math.equation(block: false, alt: "f open parenthesis x close parenthesis")[$f ( x )$] is continuous. For example, #math.equation(block: false, alt: "f open parenthesis 6 close parenthesis equals minus 0.6 less than 0")[$f ( 6 ) = − 0.6 < 0$], so #math.equation(block: false, alt: "f open parenthesis x close parenthesis less than 0")[$f ( x ) < 0$] for #math.equation(block: false, alt: "x greater than 5")[$x > 5$], and the arrow above #math.equation(block: false, alt: "x equals 5")[$x = 5$] is a down arrow. Next, #math.equation(block: false, alt: "f open parenthesis 1 close parenthesis equals 0.4 greater than 0")[$f ( 1 ) = 0.4 > 0$], so #math.equation(block: false, alt: "f open parenthesis x close parenthesis greater than 0")[$f ( x ) > 0$] whenever #math.equation(block: false, alt: "0 less than x less than 5")[$0 < x < 5$], and the arrow points up. Finally, #math.equation(block: false, alt: "f open parenthesis minus 1 close parenthesis equals minus 0.6 less than 0")[$f ( − 1 ) = − 0.6 < 0$] so #math.equation(block: false, alt: "f open parenthesis x close parenthesis less than 0")[$f ( x ) < 0$] when #math.equation(block: false, alt: "x less than 0")[$x < 0$], and the arrow points down.


    #figure(figph[Diagram of a vertical phase line with critical points at x = 5 and x = 0. Arrows on both sides of 5 point toward it and the arrow below 0 points away, so 5 is stable and 0 unstable.], alt: "Diagram of a vertical phase line with critical points at x = 5 and x = 0. Arrows on both sides of 5 point toward it and the arrow below 0 points away, so 5 is stable and 0 unstable.", caption: [Figure #math.equation(block: false, alt: "3")[$3$]])


    Armed with the phase diagram, it is easy to sketch the solutions approximately: As time #math.equation(block: false, alt: "t")[$t$] moves from left to right, the graph of a solution goes up if the arrow is up, and it goes down if the arrow is down.


    #notebox("Your Turn", rgb("#8a94a6"), rgb("#556666"), rgb("#f7f8fa"))[
Try sketching a few solutions simply from looking at the phase diagram. Check with the preceding graphs if you are getting the type of curves.


      Once we draw the phase diagram, we can easily classify critical points as stable or unstable.#math.equation(block: false, alt: "to the power 1")[$1$]


      #figure(figph[Diagram of four vertical phase lines with a critical point on each. Three, labelled unstable, have arrows pointing away, both up, or both down; the fourth, labelled stable, has both arrows aimed inward.], alt: "Diagram of four vertical phase lines with a critical point on each. Three, labelled unstable, have arrows pointing away, both up, or both down; the fourth, labelled stable, has both arrows aimed inward.", caption: [Figure #math.equation(block: false, alt: "4")[$4$]])
]


    Since any mathematical model we cook up will only be an approximation to the real world, unstable points are generally bad news.


    Let us think about the logistic equation with harvesting. Suppose an alien race really likes to eat humans. They keep a planet with humans on it and harvest the humans at a rate of #math.equation(block: false, alt: "h")[$h$] million humans per year. Suppose#math.equation(block: false, alt: "x")[$x$] is the number of humans in millions on the planet and #math.equation(block: false, alt: "t")[$t$] is time in years. Let #math.equation(block: false, alt: "M")[$M$] be the limiting population when no harvesting is done and #math.equation(block: false, alt: "k greater than 0")[$k > 0$] is some constant depending on how fast humans multiply. Our equation becomes


    #math.equation(block: true, alt: "the fraction d x over d t equals k x open parenthesis M minus x close parenthesis minus h")[$frac(d x, d t) = k x ( M − x ) − h$]


    We expand the right hand side and solve for critical points


    #math.equation(block: true, alt: "the fraction d x over d t equals minus k x squared plus k M x minus h")[$frac(d x, d t) = − k x^(2) + k M x − h$]


    Solving for the critical points #math.equation(block: false, alt: "A")[$A$] and #math.equation(block: false, alt: "B")[$B$] from the quadratic equations:


    #math.equation(block: true, alt: "A equals the fraction k M plus the square root of open parenthesis k M close parenthesis squared minus 4 h k over 2 k , B equals the fraction k M minus the square root of open parenthesis k M close parenthesis squared minus 4 h k over 2 k")[$A = frac(k M + sqrt(( k M )^(2) − 4 h k), 2 k) , #h(1em) B = frac(k M − sqrt(( k M )^(2) − 4 h k), 2 k)$]


    #notebox("Your Turn", rgb("#8a94a6"), rgb("#556666"), rgb("#f7f8fa"))[
Sketch a phase diagram for different possibilities. Note that these possibilities are #math.equation(block: false, alt: "A greater than B")[$A > B$], or #math.equation(block: false, alt: "A equals B")[$A = B$], or #math.equation(block: false, alt: "A")[$A$] and #math.equation(block: false, alt: "B")[$B$] both complex (i.e. no real solutions). Hint: Fix some simple #math.equation(block: false, alt: "k")[$k$] and #math.equation(block: false, alt: "M")[$M$] and then vary #math.equation(block: false, alt: "h")[$h$].
]


    For example, let #math.equation(block: false, alt: "M equals 8")[$M = 8$] and #math.equation(block: false, alt: "k equals 0.1")[$k = 0.1$]. When #math.equation(block: false, alt: "h equals 1")[$h = 1$], then #math.equation(block: false, alt: "A")[$A$] and #math.equation(block: false, alt: "B")[$B$] are distinct and positive. The slope field we get is in Figure #math.equation(block: false, alt: "5")[$5$]. As long as the population starts above #math.equation(block: false, alt: "B")[$B$], which is approximately 1.55 million, then the population will not die out. It will in fact tend towards #math.equation(block: false, alt: "A approximately equals 6.45")[$A ≈ 6.45$] million. If ever some catastrophe happens and the population drops below #math.equation(block: false, alt: "B")[$B$], humans will die out, and the fast food restaurant serving them will go out of business.


    
      Figures 5, 6 and 7 are three screenshots of one equation, x′ = 0.1x(8 − x) − h, at three harvest rates; this is all of them at once, with time along the horizontal axis and population up the vertical. Drag h up from 1 and watch the critical points A = 4 + √(16 − 10h) and B = 4 − √(16 − 10h) slide toward each other, collide at h = 1.6 where A = B = 4, then vanish — past which every curve in the frame plummets to zero no matter how well stocked the planet starts.
      {"functions":\[{"color":"\#1f4e79","domain":null,"expression":"0.1\*y\*(8-y) - h","expression2":"","id":"harvest-field","inequality":"lt","lineDash":\[\],"lineWidth":2,"mode":"slope-field","seed":\[0,7\],"seeds":\[{"x":0,"y":1},{"x":0,"y":2},{"x":0,"y":5},{"x":0,"y":8.5}\],"tMax":10,"tMin":-10,"variable":"x","visible":true}\],"grid":{"axisColor":"\#333333","gridColor":"\#cccccc","majorSpacing":1,"minorGridColor":"\#eeeeee","showAxes":true,"showLabels":true},"parameters":\[{"id":"p-h","label":"Harvest rate h","max":2.5,"min":0,"name":"h","step":0.05,"unit":"million/yr","value":1}\],"title":"Logistic harvesting: x' = 0.1x(8-x) - h","version":1,"viewport":{"centerX":7,"centerY":4.2,"scale":50}}
    
    #figure(figph[Slope field over time 0 to 20 and population 0 to 10 with red equilibria at about 6.4 and 1.6. Solutions above the lower line climb or fall toward 6.4; one just below it drops to zero, so 6.4 attracts and 1.6 repels.], alt: "Slope field over time 0 to 20 and population 0 to 10 with red equilibria at about 6.4 and 1.6. Solutions above the lower line climb or fall toward 6.4; one just below it drops to zero, so 6.4 attracts and 1.6 repels.", caption: [Figure #math.equation(block: false, alt: "5")[$5$]: Slope field and some solutions of #math.equation(block: false, alt: "x prime equals 0.1 x open parenthesis 8 minus x close parenthesis minus 1")[$x^(′) = 0.1 x ( 8 − x ) − 1$].])


    When #math.equation(block: false, alt: "h equals 1.6")[$h = 1.6$], then #math.equation(block: false, alt: "A equals B equals 4")[$A = B = 4$] and there is only #strong[one]critical point and it is #strong[unstable]. When the population starts above 4 million it will tend towards 4 million. If it ever drops below 4 million, humans will die out on the planet. This scenario is not one that we (as the human fast food proprietor) want to be in. A small perturbation of the equilibrium state and we are out of business; there is no room for error (see Figure #math.equation(block: false, alt: "6")[$6$]).


    #figure(figph[Slope field over time 0 to 20 and population 0 to 10, flat only along the red line at 4. A solution from above eases toward 4, while one starting near 3.5 plummets to zero at about time 15.5.], alt: "Slope field over time 0 to 20 and population 0 to 10, flat only along the red line at 4. A solution from above eases toward 4, while one starting near 3.5 plummets to zero at about time 15.5.", caption: [Figure #math.equation(block: false, alt: "6")[$6$]: The slope field and some solutions of #math.equation(block: false, alt: "x prime equals 0.1 x open parenthesis 8 minus x close parenthesis minus 1.6")[$x^(′) = 0.1 x ( 8 − x ) − 1.6$].])


    Finally if we are harvesting at 2 million humans per year, there are no critical points. The population will always plummet towards zero, no matter how well stocked the planet starts (see Figure #math.equation(block: false, alt: "7")[$7$]).


    #figure(figph[Slope field over time 0 to 20 and population 0 to 10 with every segment sloping downward, so no equilibrium. Three red solutions fall monotonically to zero, horizontal shifts of one another.], alt: "Slope field over time 0 to 20 and population 0 to 10 with every segment sloping downward, so no equilibrium. Three red solutions fall monotonically to zero, horizontal shifts of one another.", caption: [Figure #math.equation(block: false, alt: "7")[$7$]: Slope field and some solutions of #math.equation(block: false, alt: "x prime equals 0.1 x open parenthesis 8 minus x close parenthesis minus 2")[$x^(′) = 0.1 x ( 8 − x ) − 2$].])


    
=== Footnotes


    \[1\] Unstable points with one of the arrows pointing towards the critical point are sometimes called semistable.


    
=== References


    + Paul W. Berg and James L. McGregor, #link("http://books.google.com/books?id=EfJQAAAAMAAJ")[Elementary Partial Differential Equations], Holden-Day, San Francisco, CA, 1966.
+ William E. Boyce, Richard C. DiPrima, #link("http://books.google.com/books?id=nYWcQgAACAAJ")[Elementary Differential Equations and Boundary Value Problems], 9th edition, John Wiley & Sons Inc., New York, NY, 2008.
+ C.H. Edwards and D.E. Penney, #link("http://books.google.com/books?id=qi6ePwAACAAJ")[Differential Equations and Boundary Value Problems: Computing and Modeling], 4th edition, Prentice Hall, 2008.
+ Stanley J. Farlow, #link("http://books.google.com/books?id=_ozWAAAAMAAJ")[An Introduction to Differential Equations and Their Applications], McGraw-Hill, Inc., Princeton, NJ, 1994. (Published also by Dover Publications, 2006.)
+ E.L. Ince, #link("http://books.google.com/books?id=uYz-pqUD75gC")[Ordinary Differential Equations], Dover Publications, Inc., New York, NY, 1956.