Download Array Computing and Curve Plotting: Visualizing Functions using Arrays and Vectorization and more Study Guides, Projects, Research Engineering in PDF only on Docsity!
Ch.5: Array computing and curve plotting
Hans Petter Langtangen^1 ,^2
Simula Research Laboratory^1
University of Oslo, Dept. of Informatics^2
Aug 21, 2016
Goal: learn to visualize functions
We need to learn about a new object: array
Curves y = f (x) are visualized by drawing straight lines
between points along the curve
Meed to store the coordinates of the points along the curve in
lists or arrays x and y
Arrays ≈ lists, but computationally much more e cient
To compute the y coordinates (in an array) we need to learn
about array computations or vectorization
Array computations are useful for much more than plotting
curves!
The minimal need-to-know about vectors
Vectors are known from high school mathematics, e.g.,
point (x, y ) in the plane, point (x, y , z) in space
In general, a vector v is an n-tuple of numbers:
v = (v 0 ,... , vn− 1 )
Vectors can be represented by lists: vi is stored as v[i],
but we shall use arrays instead
Vectors and arrays are key concepts in this chapter. It takes separate math courses to understand what vectors and arrays really are, but in this course we only need a small subset of the complete story. A learning strategy may be to just start using vectors/arrays in programs and later, if necessary, go back to the more mathematical details in the rst part of Ch. 5.
The minimal need-to-know about arrays
Arrays are a generalization of vectors where we can have multiple
indices: Ai,j , Ai,j,k
Example: table of numbers, one index for the row, one for the
column
A =
A 0 , 0 · · · A 0 ,n− 1
Am− 1 , 0 · · · Am− 1 ,n− 1
The no of indices in an array is the rank or number of
dimensions
Vector = one-dimensional array, or rank 1 array
In Python code, we use Numerical Python arrays instead of
nested lists to represent mathematical arrays (because this is
computationally more e cient)
Storing (x,y) points on a curve in lists
Collect points on a function curve y = f (x) in lists:
def f(x): ... return x**
...
n = 5 # no of points dx = 1.0/(n-1) # x spacing in [0,1] xlist = [i*dx for i in range(n)] ylist = [f(x) for x in xlist]
pairs = [[x, y] for x, y in zip(xlist, ylist)]
Turn lists into Numerical Python (NumPy) arrays:
import numpy as np # module for arrays x = np.array(xlist) # turn list xlist into array y = np.array(ylist)
Make arrays directly (instead of lists)
The pro drops lists and makes NumPy arrays directly:
n = 5 # number of points x = np.linspace(0, 1, n) # n points in [0, 1] y = np.zeros(n) # n zeros (float data type) for i in xrange(n): ... y[i] = f(x[i]) ...
Note:
xrange is like range but faster (esp. for large n - xrange
does not explicitly build a list of integers, xrange just lets you
loop over the values)
Entire arrays must be made by numpy (np) functions
Arrays are not as exible as list, but computational much
more e cient
List elements can be any Python objects
Array elements can only be of one object type
Arrays are very e cient to store in memory and compute with
if the element type is float, int, or complex
Rule: use arrays for sequences of numbers!
We can work with entire arrays at once - instead of one
element at a time
Compute the sine of an array:
from math import sin
for i in xrange(len(x)): y[i] = sin(x[i])
However, if x is array, y can be computed by
y = np.sin(x) # x: array, y: array
The loop is now inside np.sin and implemented in very e cient C
code.
Operating on entire arrays at once is called vectorization
Vectorization gives:
shorter, more readable code, closer to the mathematics
much faster code
Use %timeit in IPython to measure the speed-up for y = sin xe−x^ :
In [1]: n = 100000
In [2]: import numpy as np
In [3]: x = np.linspace(0, 2*np.pi, n+1)
In [4]: y = np.zeros(len(x))
In [5]: %timeit for i in xrange(len(x)):
y[i] = np.sin(x[i])*np.exp(-x[i]) 1 loops, best of 3: 247 ms per loop
In [6]: %timeit y = np.sin(x)*np.exp(-x) 100 loops, best of 3: 4.77 ms per loop
In [7]: 247/4. Out[7]: 51.781970649895186 # vectorization: 50x speed-up!
A function f(x) written for a number x usually works for
array x too
from numpy import sin, exp, linspace
def f(x): return x3 + sin(x)exp(-3x)
x = 1.2 # float object y = f(x) # y is float
x = linspace(0, 3, 10001) # 10000 intervals in [0,3] y = f(x) # y is array
Note: math is for numbers and numpy for arrays
import math, numpy x = numpy.linspace(0, 1, 11) math.sin(x[3])
math.sin(x) ... TypeError: only length-1 arrays can be converted to Python scalars numpy.sin(x) array([ 0. , 0.09983, 0.19866, 0.29552, 0.38941, 0.47942, 0.56464, 0.64421, 0.71735, 0.78332, 0.84147])
Array arithmetics is broken down to a series of unary/binary
array operations
Consider y = f(x), where f returns x**3 +
sin(x)exp(-3x)
f(x) leads to the following set of vectorized
sub-computations:
1 r1 = x**
for i in range(len(x)): r1[i] = x[i]**
(but with loop in C)
2 r2 = sin(x) (computed elementwise in C)
3 r3 = -3*x
4 r4 = exp(r3)
5 r5 = r3*r
6 r6 = r1 + r
7 y = r
Note: this is the same set of operations as you would do with
a calculator when x is a number
Easyviz (imported from scitools.std) allows a more
compact Pythonic syntax for plotting curves
Use keyword arguments instead of separate function calls:
plot(t, y, xlabel='t', ylabel='y', legend='t^2*exp(-t^2)', axis=[0, 3, -0.05, 0.6], title='My First Easyviz Demo', savefig='tmp1.png', show=True) # display on the screen (default)
(This cannot be done with Matplotlib.)
0.0 0.5 1.0 1.5t 2.0 2.5 3.
y
My First Matplotlib Demo t^2*exp(-t^2)
Plotting several curves in one plot
Plot t^2 e−t
2
and t^4 e−t
2
in the same plot:
from scitools.std import * # curve plotting + array computing
def f1(t): return t2exp(-t*2)
def f2(t): return t*2f1(t)
t = linspace(0, 3, 51) y1 = f1(t) y2 = f2(t)
plot(t, y1) hold('on') # continue plotting in the same plot plot(t, y2)
xlabel('t') ylabel('y') legend('t^2exp(-t^2)', 't^4exp(-t^2)') title('Plotting two curves in the same plot') savefig('tmp2.png')
Alternative, more compact plot command
plot(t, y1, t, y2, xlabel='t', ylabel='y', legend=('t^2exp(-t^2)', 't^4exp(-t^2)'), title='Plotting two curves in the same plot', savefig='tmp2.pdf')
equivalent to
plot(t, y1) hold('on') plot(t, y2)
xlabel('t') ylabel('y') legend('t^2exp(-t^2)', 't^4exp(-t^2)') title('Plotting two curves in the same plot') savefig('tmp2.pdf')
The resulting plot with two curves
0.0 0.5 1.0 1.5 2.0 2.5 3. t
y
Plotting two curves in the same plot
t^2exp(-t^2) t^4exp(-t^2)
Controlling line styles
When plotting multiple curves in the same plot, the dierent lines
(normally) look dierent. We can control the line type and color, if
desired:
plot(t, y1, 'r-') # red (r) line (-) hold('on') plot(t, y2, 'bo') # blue (b) circles (o)
or
plot(t, y1, 'r-', t, y2, 'bo')
Documentation of colors and line styles: see the book, Ch. 5, or
Unix> pydoc scitools.easyviz
Quick plotting with minimal typing
A lazy pro would do this:
t = linspace(0, 3, 51) plot(t, t2exp(-t2), t, t4exp(-t**2))
Plot function given on the command line
Task: plot function given on the command line
Terminal> python plotf.py expression xmin xmax Terminal> python plotf.py "exp(-0.2x)sin(2pix)" 0 4*pi
Should plot e−^0.^2 x^ sin( 2 πx), x ∈ [ 0 , 4 π]. plotf.py should work for
any mathematical expression.
Solution
Complete program:
from scitools.std import *
or alternatively
from numpy import * from matplotlib.pyplot import *
formula = sys.argv[1] xmin = eval(sys.argv[2]) xmax = eval(sys.argv[3])
x = linspace(xmin, xmax, 101) y = eval(formula) plot(x, y, title=formula)
Let's make a movie/animation
0
1
2
-6 -4 -2 0 2 4 6
s=0. s= s=
The Gaussian/bell function
f (x; m, s) =
s
exp
[
x − m
s
) 2 ]
m is the location of the peak
s is a measure of the width
of the function
Make a movie (animation)
of how f (x; m, s) changes
shape as s goes from 2 to
0
1
2
-6 -4 -2 0 2 4 6
s=0.2 s= s=
Movies are made from a (large) set of individual plots
Goal: make a movie showing how f (x) varies in shape as s
decreases
Idea: put many plots (for dierent s values) together
(exactly as a cartoon movie)
How to program: loop over s values, call plot for each s and
make hardcopy, combine all hardcopies to a movie
Very important: x the y axis! Otherwise, the y axis always
adapts to the peak of the function and the visual impression
gets completely wrong
The complete code for making the animation
from scitools.std import * import time
def f(x, m, s): return (1.0/(sqrt(2pi)s))exp(-0.5((x-m)/s)**2)
m = 0; s_start = 2; s_stop = 0. s_values = linspace(s_start, s_stop, 30)
x = linspace(m -3s_start, m + 3s_start, 1000)
f is max for x=m (smaller s gives larger max value)
max_f = f(m, m, s_stop)
Show the movie on the screen
and make hardcopies of frames simultaneously
import time frame_counter = 0
for s in s_values: y = f(x, m, s) plot(x, y, axis=[x[0], x[-1], -0.1, max_f], xlabel='x', ylabel='f', legend='s=%4.2f' % s, savefig='tmp_%04d.png' % frame_counter) frame_counter += 1 #time.sleep(0.2) # pause to control movie speed
Let's try to plot a discontinuous function
The Heaviside function is frequently used in science and
engineering:
H(x) =
0 , x < 0
1 , x ≥ 0
Python implementation:
def H(x): return (0 if x < 0 else 1)
4 3 2 1 0 1 2 3 4
Plotting the Heaviside function: rst try
Standard approach:
x = linspace(-10, 10, 5) # few points (simple curve) y = H(x) plot(x, y)
First problem: ValueError error in H(x) from if x < 0
Let us debug in an interactive shell:
x = linspace(-10,10,5) x array([-10., -5., 0., 5., 10.]) b = x < 0 b array([ True, True, False, False, False], dtype=bool) bool(b) # evaluate b in a boolean context ... ValueError: The truth value of an array with more than one element is ambiguous. Use a.any() or a.all()
if x < 0 does not work if x is array
Remedy 1: use a loop over x values
def H_loop(x): r = zeros(len(x)) # or r = x.copy() for i in xrange(len(x)): r[i] = H(x[i]) return r
n = 5 x = linspace(-5, 5, n+1) y = H_loop(x)
Downside: much to write, slow code if n is large
if x < 0 does not work if x is array
Remedy 2: use vectorize
from numpy import vectorize
Automatic vectorization of function H
Hv = vectorize(H)
Hv(x) works with array x
Downside: The resulting function is as slow as Remedy 1
if x < 0 does not work if x is array
Remedy 3: code the if test dierently
def Hv(x): return where(x < 0, 0.0, 1.0)
More generally:
def f(x): if condition: x = else: x = return x
def f_vectorized(x): def f_vectorized(x): x1 = x2 = r = np.where(condition, x1, x2) return r
Back to plotting the Heaviside function
With a vectorized Hv(x) function we can plot in the standard way
x = linspace(-10, 10, 5) # linspace(-10, 10, 50) y = Hv(x) plot(x, y, axis=[x[0], x[-1], -0.1, 1.1])
4 3 2 1 0 1 2 3 4
How to make the function look discontinuous in the plot?
Newbie: use a lot of x points; the curve gets steeper
Pro: plot just two horizontal line segments
one from x = −10 to x = 0, y = 0; and one from x = 0 to
x = 10, y = 1
plot([-10, 0, 0, 10], [0, 0, 1, 1], axis=[x[0], x[-1], -0.1, 1.1])
Draws straight lines between (− 10 , 0 ), ( 0 , 0 ), ( 0 , 1 ), ( 10 , 1 )
The nal plot of the discontinuous Heaviside function
4 3 2 1 0 1 2 3 4
4 3 2 1 0 1 2 3 4
Removing the vertical jump from the plot
Question
Some will argue and say that at high school they would draw H(x)
as two horizontal lines without the vertical line at x = 0, illustrating
the jump. How can we plot such a curve?
Some functions are challenging to visualize
Plot f (x) = sin( 1 /x)
def f(x): return sin(1.0/x)
x1 = linspace(-1, 1, 10) # use 10 points x2 = linspace(-1, 1, 1000) # use 1000 points plot(x1, f(x1), label='%d points' % len(x)) plot(x2, f(x2), label='%d points' % len(x))
Plot based on 10 points
-0.
-0.
-0.
-0.
0
1
-1 -0.5 0 0.5 1
Plot based on 1000 points
-0.
-0.
-0.
-0.
0
1
-1 -0.5 0 0.5 1
Generalized array indexing
Recall slicing: a[f:t:i], where the slice f:t:i implies a set of
indices (from, to, increment).
Any integer list or array can be used to indicate a set of indices:
a = linspace(1, 8, 8) a array([ 1., 2., 3., 4., 5., 6., 7., 8.]) a[[1,6,7]] = 10 a array([ 1., 10., 3., 4., 5., 6., 10., 10.]) a[range(2,8,3)] = -2 # same as a[2:8:3] = - a array([ 1., 10., -2., 4., 5., -2., 10., 10.])
Generalized array indexing with boolean expressions
a < 0 [False, False, True, False, False, True, False, False]
a[a < 0] # pick out all negative elements array([-2., -2.])
a[a < 0] = a.max() # if a[i]<10, set a[i]= a array([ 1., 10., 10., 4., 5., 10., 10., 10.])
Two-dimensional arrays; math intro
When we have a table of numbers,
(called matrix by mathematicians) it is natural to use a
two-dimensional array Ai,j with one index for the rows and one for
the columns:
A =
A 0 , 0 · · · A 0 ,n− 1
Am− 1 , 0 · · · Am− 1 ,n− 1
Two-dimensional arrays; Python code
Making and lling a two-dimensional NumPy array goes like this:
A = zeros((3,4)) # 3x4 table of numbers A[0,0] = - A[1,0] = 1 A[2,0] = 10 A[0,1] = - ... A[2,3] = -
can also write (as for nested lists)
A[2][3] = -
From nested list to two-dimensional array
Let us make a table of numbers in a nested list:
Cdegrees = [-30 + i10 for i in range(3)] Fdegrees = [9./5C + 32 for C in Cdegrees] table = [[C, F] for C, F in zip(Cdegrees, Fdegrees)] print table [[-30, -22.0], [-20, -4.0], [-10, 14.0]]
Turn into NumPy array:
table2 = array(table) print table
a[range(2,8,3)] = -2 # same as a[2:8:3] = - >>> a array([ 1., 10., -2., 4., 5., -2., 10., 10.]) ## Generalized array indexing with boolean expressions >>> a < 0 [False, False, True, False, False, True, False, False] >>> a[a < 0] # pick out all negative elements array([-2., -2.]) >>> a[a < 0] = a.max() # if a[i]<10, set a[i]= >>> a array([ 1., 10., 10., 4., 5., 10., 10., 10.]) ## Two-dimensional arrays; math intro ## When we have a table of numbers, ## (called matrix by mathematicians) it is natural to use a ## two-dimensional array Ai,j with one index for the rows and one for ## the columns: ## A = ## A 0 , 0 · · · A 0 ,n− 1 ## Am− 1 , 0 · · · Am− 1 ,n− 1 ## Two-dimensional arrays; Python code ## Making and lling a two-dimensional NumPy array goes like this: A = zeros((3,4)) # 3x4 table of numbers A[0,0] = - A[1,0] = 1 A[2,0] = 10 A[0,1] = - ... A[2,3] = - # can also write (as for nested lists) A[2][3] = - ## From nested list to two-dimensional array ## Let us make a table of numbers in a nested list: >>> Cdegrees = [-30 + i10 for i in range(3)] >>> Fdegrees = [9./5C + 32 for C in Cdegrees] >>> table = [[C, F] for C, F in zip(Cdegrees, Fdegrees)] >>> print table [[-30, -22.0], [-20, -4.0], [-10, 14.0]] ## Turn into NumPy array: >>> table2 = array(table) >>> print table [[-30. -22.] [-20. -4.] [-10. 14.]]
How to loop over two-dimensional arrays
table2.shape # see the number of elements in each dir. (3, 2) # 3 rows, 2 columns
A for loop over all array elements:
for i in range(table2.shape[0]): ... for j in range(table2.shape[1]): ... print 'table2[%d,%d] = %g' % (i, j, table2[i,j])
a array([ 1., 10., 10., 4., 5., 10., 10., 10.]) ## Two-dimensional arrays; math intro ## When we have a table of numbers, ## (called matrix by mathematicians) it is natural to use a ## two-dimensional array Ai,j with one index for the rows and one for ## the columns: ## A = ## A 0 , 0 · · · A 0 ,n− 1 ## Am− 1 , 0 · · · Am− 1 ,n− 1 ## Two-dimensional arrays; Python code ## Making and lling a two-dimensional NumPy array goes like this: A = zeros((3,4)) # 3x4 table of numbers A[0,0] = - A[1,0] = 1 A[2,0] = 10 A[0,1] = - ... A[2,3] = - # can also write (as for nested lists) A[2][3] = - ## From nested list to two-dimensional array ## Let us make a table of numbers in a nested list: >>> Cdegrees = [-30 + i10 for i in range(3)] >>> Fdegrees = [9./5C + 32 for C in Cdegrees] >>> table = [[C, F] for C, F in zip(Cdegrees, Fdegrees)] >>> print table [[-30, -22.0], [-20, -4.0], [-10, 14.0]] ## Turn into NumPy array: >>> table2 = array(table) >>> print table [[-30. -22.] [-20. -4.] [-10. 14.]] ## How to loop over two-dimensional arrays >>> table2.shape # see the number of elements in each dir. (3, 2) # 3 rows, 2 columns ## A for loop over all array elements: >>> for i in range(table2.shape[0]): ... for j in range(table2.shape[1]): ... print 'table2[%d,%d] = %g' % (i, j, table2[i,j]) ... table2[0,0] = - table2[0,1] = - ... table2[2,1] = 14
Alternative single loop over all elements:
for index_tuple, value in np.ndenumerate(table2): ... print 'index %s has value %g' % \
... (index_tuple, table2[index_tuple]) ... index (0,0) has value - index (0,1) has value - ... index (2,1) has value 14
type(index_tuple) <type 'tuple'>
How to take slices of two-dimensional arrays
Rule: can use slices start:stop:inc for each index
table2[0:table2.shape[0], 1] # 2nd column (index 1) array([-22., -4., 14.])
table2[0:, 1] # same array([-22., -4., 14.])
table2[:, 1] # same array([-22., -4., 14.])
t = linspace(1, 30, 30).reshape(5, 6) t[1:-1:2, 2:] array([[ 9., 10., 11., 12.], Cdegrees = [-30 + i10 for i in range(3)] >>> Fdegrees = [9./5C + 32 for C in Cdegrees] >>> table = [[C, F] for C, F in zip(Cdegrees, Fdegrees)] >>> print table [[-30, -22.0], [-20, -4.0], [-10, 14.0]] ## Turn into NumPy array: >>> table2 = array(table) >>> print table [[-30. -22.] [-20. -4.] [-10. 14.]] ## How to loop over two-dimensional arrays >>> table2.shape # see the number of elements in each dir. (3, 2) # 3 rows, 2 columns ## A for loop over all array elements: >>> for i in range(table2.shape[0]): ... for j in range(table2.shape[1]): ... print 'table2[%d,%d] = %g' % (i, j, table2[i,j]) ... table2[0,0] = - table2[0,1] = - ... table2[2,1] = 14 ## Alternative single loop over all elements: >>> for index_tuple, value in np.ndenumerate(table2): ... print 'index %s has value %g' %
... (index_tuple, table2[index_tuple]) ... index (0,0) has value - index (0,1) has value - ... index (2,1) has value 14 >>> type(index_tuple) <type 'tuple'> ## How to take slices of two-dimensional arrays ## Rule: can use slices start:stop:inc for each index table2[0:table2.shape[0], 1] # 2nd column (index 1) array([-22., -4., 14.]) >>> table2[0:, 1] # same array([-22., -4., 14.]) >>> table2[:, 1] # same array([-22., -4., 14.]) >>> t = linspace(1, 30, 30).reshape(5, 6) >>> t[1:-1:2, 2:] array([[ 9., 10., 11., 12.], [ 21., 22., 23., 24.]]) t array([[ 1., 2., 3., 4., 5., 6.], [ 7., 8., 9., 10., 11., 12.], table = [[C, F] for C, F in zip(Cdegrees, Fdegrees)] >>> print table [[-30, -22.0], [-20, -4.0], [-10, 14.0]] ## Turn into NumPy array: >>> table2 = array(table) >>> print table [[-30. -22.] [-20. -4.] [-10. 14.]] ## How to loop over two-dimensional arrays >>> table2.shape # see the number of elements in each dir. (3, 2) # 3 rows, 2 columns ## A for loop over all array elements: >>> for i in range(table2.shape[0]): ... for j in range(table2.shape[1]): ... print 'table2[%d,%d] = %g' % (i, j, table2[i,j]) ... table2[0,0] = - table2[0,1] = - ... table2[2,1] = 14 ## Alternative single loop over all elements: >>> for index_tuple, value in np.ndenumerate(table2): ... print 'index %s has value %g' %
... (index_tuple, table2[index_tuple]) ... index (0,0) has value - index (0,1) has value - ... index (2,1) has value 14 >>> type(index_tuple) <type 'tuple'> ## How to take slices of two-dimensional arrays ## Rule: can use slices start:stop:inc for each index table2[0:table2.shape[0], 1] # 2nd column (index 1) array([-22., -4., 14.]) >>> table2[0:, 1] # same array([-22., -4., 14.]) >>> table2[:, 1] # same array([-22., -4., 14.]) >>> t = linspace(1, 30, 30).reshape(5, 6) >>> t[1:-1:2, 2:] array([[ 9., 10., 11., 12.], [ 21., 22., 23., 24.]]) >>> t array([[ 1., 2., 3., 4., 5., 6.], [ 7., 8., 9., 10., 11., 12.], [ 13., 14., 15., 16., 17., 18.], [ 19., 20., 21., 22., 23., 24.], [ 25., 26., 27., 28., 29., 30.]])
Time for a question
Problem:
Given
t array([[ 1., 2., 3., 4., 5., 6.], [ 7., 8., 9., 10., 11., 12.], [ 13., 14., 15., 16., 17., 18.], [ 19., 20., 21., 22., 23., 24.], [ 25., 26., 27., 28., 29., 30.]])
What will t[1:-1:2, 2:] be?
Solution:
Slice 1:-1:2 for rst index results in
[ 7., 8., 9., 10., 11., 12.]
[ 19., 20., 21., 22., 23., 24.]
Slice 2: for the second index then gives
[ 9., 10., 11., 12.]
[ 21., 22., 23., 24.]
Summary of vectors and arrays
Vector/array computing: apply a mathematical expression to
every element in the vector/array (no loops in Python)
Ex: sin(x4)exp(-x*2), x can be array or scalar
for array the i'th element becomes
sin(x[i]4)exp(-x[i]*2)
Vectorization: make scalar mathematical computation valid for
vectors/arrays
Pure mathematical expressions require no extra vectorization
Mathematical formulas involving if tests require manual work
for vectorization:
scalar_result = expression1 if condition else expression vector_result = where(condition, expression1, expression2)
Summary of plotting y = f (x) curves
Curve plotting (uni ed syntax for Matplotlib and SciTools):
from matplotlib.pyplot import * #from scitools.std import *
plot(x, y) # simplest command
plot(t1, y1, 'r', # curve 1, red line t2, y2, 'b', # curve 2, blue line t3, y3, 'o') # curve 3, circles at data points axis([t1[0], t1[-1], -1.1, 1.1]) legend(['model 1', 'model 2', 'measurements']) xlabel('time'); ylabel('force') savefig('myframe_%04d.png' % plot_counter)
Note: straight lines are drawn between each data point
Alternativ plotting of y = f (x) curves
Single SciTools plot command with keyword arguments:
from scitools.std import *
plot(t1, y1, 'r', # curve 1, red line t2, y2, 'b', # curve 2, blue line t3, y3, 'o', # curve 3, circles at data points axis=[t1[0], t1[-1], -1.1, 1.1], legend=('model 1', 'model 2', 'measurements'), xlabel='time', ylabel='force', savefig='myframe_%04d.png' % plot_counter)
Summary of making animations
Make a hardcopy of each plot frame (PNG or PDF format)
Use avconv or ffmpeg to make movie
Terminal> avconv -r 5 -i tmp_%04d.png -vcodec flv movie.flv
