3

I am putting all of my code in one block for easy copy and paste.

sys = {u (1 - u + a v), r v (1 - v + b u)};
sys1 = sys /. {a -> 2, b -> 3, r -> 1};
equilibria = Solve[sys1 == {0, 0}];
equilibria = equilibria[[{2, 3, 4}]];
sol[{u0_, v0_}] := 
  NDSolveValue[{{u'[t], 
      v'[t]} == (sys1 /. {u -> u[t], v -> v[t]}), {u[0], v[0]} == {u0,
       v0},
    WhenEvent[Abs[u[t] - 1] > 1, "StopIntegration"],
    WhenEvent[Abs[v[t] - 1] > 1, "StopIntegration"]}, {u[t], 
    v[t]}, {t, -4, 4},
   "ExtrapolationHandler" -> {Indeterminate &, 
     "WarningMessage" -> False}];
initconds = {{0.1, 1.5}, {1.5, 0.1}, {0.1, 0.5}, {0.5, 0.1}};
toplot = Map[sol, initconds];
plt = ParametricPlot[toplot, {t, -4, 4}, PlotStyle -> Blue];
Show[vp, cp, plt, 
 Graphics[{Red, PointSize[Large], Point[{u, v}] /. equilibria,
   Text[Style["(0,0) nodal source", Black, Background -> White], {0.1,
      0.1}],
   Text[Style["(1,0) saddle", Black, Background -> White], {1, 0.1}],
   Text[Style["(0,1) saddle", Black, Background -> White], {0.1, 
     1.1}]}]]

The resulting figure follows.

enter image description here

I've only tested these concepts one time, on this one system. I am wondering what folks think. Will I get in trouble in general?

Grisha's Suggestion:

Clear["Global`*"];
sys = {u (1 - u + a v), r v (1 - v + b u)};
sys1 = sys /. {a -> 2, b -> 3, r -> 1};
vp = VectorPlot[sys1, {u, -0.1, 2}, {v, -0.1, 2},
   VectorScale -> {0.045, 0.9, None},
   VectorStyle -> {GrayLevel[0.8]},
   Axes -> True, AxesLabel -> {x, y}, AxesOrigin -> {0, 0}];
cp = ContourPlot[sys1, {u, -0.1, 2}, {v, -0.1, 2},
   ContourStyle -> {Orange, Green}];
equilibria = Solve[sys1 == {0, 0}];
equilibria = equilibria[[{2, 3, 4}]];
initconds = {{0.1, 1.5}, {1.5, 0.1}, {0.1, 0.5}, {0.5, 0.1}};
sp = StreamPlot[sys1, {u, -0.5, 2}, {v, -0.5, 2}, StreamScale -> None,
    StreamPoints -> initconds,
   StreamStyle -> {Blue, "Line"}];
Show[vp, cp, sp, 
 Graphics[{Red, PointSize[Large], Point[{u, v}] /. equilibria,
   Text[Style["(0,0) nodal source", Black, Background -> White], {0.1,
      0.1}],
   Text[Style["(1,0) saddle", Black, Background -> White], {1, 0.1}],
   Text[Style["(0,1) saddle", Black, Background -> White], {0.1, 
     1.1}]}]]

enter image description here

One can see that Grisha's suggestion worked nicely here.

Case where ab<1

Clear["Global`*"];
sys = {u (1 - u + a v), r v (1 - v + b u)};
sys2 = sys /. {a -> 1/3, b -> 1/4, r -> 1};
vp2 = VectorPlot[sys2, {u, -0.1, 2}, {v, -0.1, 2},
   VectorScale -> {0.045, 0.9, None},
   VectorStyle -> {GrayLevel[0.8]},
   Axes -> True, AxesLabel -> {x, y}, AxesOrigin -> {0, 0}];
cp2 = ContourPlot[sys2, {u, -0.1, 2}, {v, -0.1, 2},
   ContourStyle -> {Orange, Green}];
equilibria = Solve[sys2 == {0, 0}];
initconds2 = {{0.1, 1.5}, {1.5, 0.1}, {0.1, 0.5}, {0.5, 0.1}, {0.7, 
    2}, {2, 2}, {2, 0.7}};
sp2 = StreamPlot[sys2, {u, -0.1, 2}, {v, -0.1, 2}, 
   StreamScale -> None, StreamPoints -> initconds2,
   StreamStyle -> {Blue, "Line"}];
Show[vp2, cp2, sp2, 
 Graphics[{Red, PointSize[Large], Point[{u, v}] /. equilibria,
   Text[Style["(0,0) nodal source", Black, Background -> White], {0.1,
      0.1}],
   Text[Style["(1,0) saddle", Black, Background -> White], {1, 0.1}],
   Text[Style["(0,1) saddle", Black, Background -> White], {0.1, 1.1}],
   Text[Style["(16/11,15/11) nodal sink", Black, 
     Background -> White], {1.5, 1.5}]}]]

Produces the following image.

enter image description here

Notice that it failed to use three initial conditions. Changing the initial condition {0.7, 2} to {0.7,1.9} fixes the problem.

Clear["Global`*"];
sys = {u (1 - u + a v), r v (1 - v + b u)};
sys2 = sys /. {a -> 1/3, b -> 1/4, r -> 1};
vp2 = VectorPlot[sys2, {u, -0.1, 2}, {v, -0.1, 2},
   VectorScale -> {0.045, 0.9, None},
   VectorStyle -> {GrayLevel[0.8]},
   Axes -> True, AxesLabel -> {x, y}, AxesOrigin -> {0, 0}];
cp2 = ContourPlot[sys2, {u, -0.1, 2}, {v, -0.1, 2},
   ContourStyle -> {Orange, Green}];
equilibria = Solve[sys2 == {0, 0}];
initconds2 = {{0.1, 1.5}, {1.5, 0.1}, {0.1, 0.5}, {0.5, 0.1}, {0.7, 
    1.9}, {2, 2}, {2, 0.7}};
sp2 = StreamPlot[sys2, {u, -0.1, 2}, {v, -0.1, 2}, 
   StreamScale -> None, StreamPoints -> initconds2,
   StreamStyle -> {Blue, "Line"}];
Show[vp2, cp2, sp2, 
 Graphics[{Red, PointSize[Large], Point[{u, v}] /. equilibria,
   Text[Style["(0,0) nodal source", Black, Background -> White], {0.1,
      0.1}],
   Text[Style["(1,0) saddle", Black, Background -> White], {1, 0.1}],
   Text[Style["(0,1) saddle", Black, Background -> White], {0.1, 1.1}],
   Text[Style["(16/11,15/11) nodal sink", Black, 
     Background -> White], {1.5, 1.5}]}]]

Giving the following image.

enter image description here

Not bad, not bad at all. Curious why I had to make the initial condition point change I made.

However, note that now all the stream plots don't push very close to the nodal sink. Is there a way you can force them to move a bit closer?

My Technique Again My technique seemed to work again.

Clear["Global`*"];
sys = {u (1 - u + a v), r v (1 - v + b u)};
sys2 = sys /. {a -> 1/3, b -> 1/4, r -> 1};
vp2 = VectorPlot[sys2, {u, -0.1, 2}, {v, -0.1, 2},
   VectorScale -> {0.045, 0.9, None},
   VectorStyle -> {GrayLevel[0.8]},
   Axes -> True, AxesLabel -> {x, y}, AxesOrigin -> {0, 0}];
cp2 = ContourPlot[sys2, {u, -0.1, 2}, {v, -0.1, 2},
   ContourStyle -> {Orange, Green}];
equilibria = Solve[sys2 == {0, 0}];
sol[{u0_, v0_}] := 
  NDSolveValue[{{u'[t], 
      v'[t]} == (sys2 /. {u -> u[t], v -> v[t]}), {u[0], v[0]} == {u0,
       v0},
    WhenEvent[Abs[u[t] - 1] > 1, "StopIntegration"],
    WhenEvent[Abs[v[t] - 1] > 1, "StopIntegration"]}, {u[t], 
    v[t]}, {t, -4, 4},
   "ExtrapolationHandler" -> {Indeterminate &, 
     "WarningMessage" -> False}];
initconds2 = {{0.1, 1.5}, {1.5, 0.1}, {0.1, 0.5}, {0.5, 0.1}, {0.7, 
    1.9}, {2, 2}, {2, 0.7}};
toplot2 = Map[sol, initconds2];
plt2 = ParametricPlot[toplot2, {t, -4, 4}, PlotStyle -> Blue];
Show[vp2, cp2, plt2, 
 Graphics[{Red, PointSize[Large], Point[{u, v}] /. equilibria,
   Text[Style["(0,0) nodal source", Black, Background -> White], {0.1,
      0.1}],
   Text[Style["(1,0) saddle", Black, Background -> White], {1, 0.1}],
   Text[Style["(0,1) saddle", Black, Background -> White], {0.1, 1.1}],
   Text[Style["(16/11,15/11) nodal sink", Black, 
     Background -> White], {1.5, 1.5}]}]]

Producing this image:

enter image description here

However, this time I got a couple warning messages:

enter image description here

I'd love to hear suggestions from folks.

David
  • 14,883
  • 4
  • 44
  • 117
  • @MichaelE2 What do you think? Is this a good idea for my students? I am amazed that I found that ExtrapolationHandler thing. – David Apr 22 '15 at 06:43
  • 1
    Why do you not use StreamPlot[sys1, {u, -0.5, 2}, {v, -0.5, 2}, StreamScale -> Full, StreamPoints -> initconds, StreamStyle -> {Blue, "Line"}] for the selected initial conditions? – Grisha Kirilin Apr 22 '15 at 07:34
  • Why not plot in 4 quadrants to see arrows on all sides of saddle points? – Narasimham Apr 22 '15 at 10:29
  • @GrishaKirilin Used your suggestion and ran into a small problem. And I am looking for a few more suggestions. See the edits to my post. Thanks for your idea. It is really a good one that I might be able to easily use with my student. I will note however, that StreamStyle -> {Blue, "Line"} doesn't really work. Try StreamStyle -> {Blue, "ArrowArrow"}. Think we've seen these errors before. – David Apr 22 '15 at 20:47
  • (1) I often use Grisha's approach, but as you've discovered, it does many things automatically. I find it hard to control in some cases. It stops, deletes or moves streamlines that get too close to one another. (2) The messages occur because the initial points {2,2} and {2,0.7} start on the boundary of the WhenEvent conditions; in such a case, the WhenEvent is ignored and NDSolve tries to integrate to infinity. WhenEvent only detects crossing a boundary. If you change each 2 to 1.99, no messages are generated. (3) I don't quite understand the purpose of the ContourPlot. – Michael E2 Apr 22 '15 at 21:48
  • @MichaelE2 When you have have a system of two differential equations, x'=f(x,y) and y'=g(x,y), the values of x and y that make f(x,y)=0 is called the nullcline for the first equation, and the values of x and y that make g(x,y)=0 is called the nullcline of the second equations. I have my students draw them in two different colors. Where the orange and green nullclines intersect, you have an equilibrium point. Along the nullcline of the first equation, x'=0, so there is no right and left motion. For the y'=0 nullcline, there is no up and down motion. This helps create phase portraits by hand. – David Apr 22 '15 at 22:47
  • @MichaelE2 For an example, see here and here. – David Apr 22 '15 at 22:54
  • I see. For some reason my mind had gotten stuck on separatrices, both because of one of your other questions and because when I picture a saddle etc., I think of the flow in the regions partitioned by the separatrices. That's what I kept wanting to be drawn (even though I know they're not even straight lines in this case). – Michael E2 Apr 23 '15 at 02:29
  • I don't understand what the question is. –  Apr 23 '15 at 03:32

0 Answers0