Plotting Direction Fields, Phase Portraits, and Contour Maps

[Alexander Jankowski]

PREFACE

In this tutorial, you will learn how to create a graph like the following one, which is the phase portrait for the system
\begin{equation*}
\left\{ \begin{array}{l} x'(t) = x\ln{(y^2)}, \\
y'(t) = y\ln{(x^2)}.
\end{array} \right.
\end{equation*}

Note that the right-hand expressions are not continuously differentiable at $(0,0)$ and there is no linearization at this point, which explains that while stable, the point $(0,0)$ does not look (after close examination) like a standard node.

In order to create such plots, you need computational software and the ability to write expressions that the software can interpret. There are three popular programs for this purpose: Waterloo Maple, MATLAB, and Wolfram Mathematica. However, I should caution you that these programs are expensive, even if you are buying a student license. There are a lot of open-source alternatives: Mathics, Sage, GNU Octave (similar to MATLAB), and Scilab, to name a few. Each software has its advantages and disadvantages. For example, it is commonly agreed that  Mathematica is good for symbolic manipulation of expressions (e.g. integration, solving differential equations, etc.) whereas MATLAB is more appropriate for working with large sets of data. This is a decision that I will leave to you. This guide is concerned with Mathematica.

1 DIRECTION FIELDS

A direction field for a first-order linear ordinary differential equation can be created using the VectorPlot function. For our purposes, the basic form is

Code: [Select]

VectorPlot[{1, <function>}, {x, <xmin>, <xmax>}, {y, <ymin>, <ymax>}]
Here, <function> is the differential equation in the form $y'(x) = f(x,y)$. By default, the direction field that you get will be a coloured vector field. If you want to get a classic direction field, you'll need to add the option VectorStyle -> Arrowheads[0]. Often, the colour unnecessarily increases image size, and if you intend the plot your solution on the direction field, then it is better to keep the direction field in black and the solution in colour. To do this, you can add the VectorStyle -> Black option. To see how these options are used, see the following example for the differential equation $y' = xy$.

Code: [Select]

VectorPlot[{1, x*y}, {x, -1, 1}, {y, -1, 1}, VectorStyle -> {Arrowheads[0], Black}]

2 PHASE PORTRAITS

A phase portrait for a two-dimensional autonomous system of ordinary differential equations can be created using the StreamPlot function. The generic form of what you will type will look like this:

Code: [Select]

StreamPlot[{<function1>, <function2>}, {x, <xmin>, <xmax>}, {y, <ymin>, <ymax>}]
Here, <function1> and <function2> correspond to $x'(t)$ and $y'(t)$, respectively. As with VectorPlot, you can add some options after these arguments. I often use only StreamStyle -> Black since the default colour is unnecessary. To exemplify this function, I will show you how to plot the phase portrait for case $(d)$ of the system in the third semester-end challenge.

Code: [Select]

StreamPlot[{-y + x*(-(Sin[Pi*Sqrt[x^2 + y^2]])^2), x + y*(-(Sin[Pi*Sqrt[x^2 + y^2]])^2)}, {x, -5, 5}, {y, -5, 5}, StreamStyle -> Black]

You can see that the integral curves on this phase portrait are not very long. To increase the length, you can make use of the StreamScale -> <size> option, where <size> is a real number in $[0,1]$. The length that you specify depends on the phase plane and on your intentions. Applying the StreamScale option not only scales the length of the curves but also the arrowheads. This is not what we want, and we can add the Arrowheads[<size>] option to avoid this. It is usually sufficient to use a size of Tiny or Small. In the following example, you can see a phase portrait of the very first system in this tutorial. The main difference is that the integral curves are much longer.

Code: [Select]

StreamPlot[{x*Log[y^2], y*Log[x^2]}, {x, -1, 1}, {y, -1, 1}, StreamScale -> 1, StreamStyle -> {Arrowheads[Small], Black}]

3 CONTOUR MAPS

A contour map for a system is made with the ContourPlot function. Its form is not very different from the above functions:

Code: [Select]

ContourPlot[{<implicitfunction>}, {x, <xmin>, <xmax>}, {y, <ymin>, <ymax>}]
Note that <implicitfunction> is simply a function in the form $H(x,y)=C$. Thus, ContourPlot is useful only if we can find $H$. For example, it was possible in the Easter challenge to determine $H$ for system $(a)$ to be $H(x,y) = \cos{x}+\cos{y}$. We can plot this with the expression

Code: [Select]

ContourPlot[{Cos[x] + Cos[y]}, {x, -2Pi, 2Pi}, {y, -2Pi, 2Pi}, ColorFunction -> "BlueGreenYellow"]

Note that, out of personal preference, I have added the option ColorFunction -> "BlueGreenYellow".

FURTHER READING

The most important advice that I can give anyone who is willing to pursue this kind of work is to read the program documentation. The documentation for Mathematica is located below, but you can also access it from within the program. You can learn a lot by simply looking at the examples that are given.

http://reference.wolfram.com/mathematica/guide/Mathematica.html

Sometimes, reading documentation is not good enough because you may not directly make a connection between a problem and a function. Then, you should search for tutorials. Many nice tutorials can be found on university websites. Here is an example.

http://www.math.mtu.edu/~msgocken/pdebook2/mathtut2.pdf

Last, but not least, is a technical comparison of Maple, MATHLAB, and Mathematica.

http://amath.colorado.edu/computing/mmm/

LINKS

Maplesoft's Maple
http://www.maplesoft.com/products/maple/

MathWorks' MATLAB
http://www.mathworks.com/products/matlab/

Wolfram Mathematica
http://www.wolfram.com/mathematica/

Mathics
http://www.mathics.org/

Sage
http://www.sagemath.org/

Scilab
http://www.scilab.org/

GNU Octave
http://www.gnu.org/software/octave/

Gnuplot
http://www.gnuplot.info

LOG

06.04.13: First version published.
09.04.13: Updated Section 2 with StreamScale option and new example.