import theme
theme.load_style()
This lecture by Tim Fuller is licensed under the Creative Commons Attribution 4.0 International License. All code examples are also licensed under the MIT license.
This book is intended for students in science/engineering/math and presumes an introductory knowledge of computer programming (be it whatever language). This book is not an introduction to computer programming.
For a complete introduction to IPython Notebooks, see the [IPython Notebook Tutorial](IPython Notebook Tutorial.ipynb).
As your read through each section, you will see cells labeled with
which are places meant for you to try the new concept introduced. Keep in mind, that these Notebooks are interactive, fee free to modify values, move things, test things out.
Python is an open source, general-purpose language, with hundreds of modules geared toward scientific computing. Some main features of Python are:
Some relative advantages of Python
Some relative disadvantages
There is no way a short tutorial can describe everything in Python (I've been using/developing Python for 10+ years and learn new stuff all the time). But, it will introduce enough concepts to complete the first homework assignment and provide a foundation for future learning.
Python is synonomous with the programming language, the interactive Python shell, and the interpreter under the hood that actually reads and processes the code. There are several ways to write/interact with Python: the interactive python shell, the IPython notebook, and by interpreting program files containing python statemtns. This is one of Python's relative strengths: it allows users to select the environment that best suits their needs.
The Python shell is a program that reads and executes Python statements as you enter them. It is opened by invoking the following at a command prompt
$ python
Interactive Python sessions are fantastic for quick calculations, or testing snippets of code.
Ipython Notebooks provide an interactive environment for Python similar to Maple's worksheets or Mathematica's notebooks. It is based on the IPython shell (a modified interactive Python shell). IPython notebook sessions are stored locally as json
files with a .ipynb
file extension and are interpreted graphically as a web session using IPython's notebook server. This file is an example of an IPython Notebook file.
The IPython Notebook server is launched locally by executing the following at the command prompt:
$ ipython notebook
This will open a new web browser (or tab in an existing window) in which notebooks can be created/edited/deleted. The fact that you are reading this means you have already launched the notebook server.
IPython notebooks consist of a collection of text cells and code cells. code
cells execute Python code in the notebook's kernel. Text cells are labeled markdown
, Raw NBConvert
, and Heading1-6
. Different text cell types display their contents differently. Looking at the toolbar above, can you tell what type of cell this is? (Markdown
).
See the full IPython Tutorial for more details.
Often, it is not convenient to work in an interactive or notebook environment. In those situations, Python statements can be saved in Python program files that are later interpreted and executed by the Python interpreter (this is not strictly true, but is true enough for now). Python program files (usually) have the file extension ".py
", e.g. baz.py
.
With the exception of comment lines, every line in a Python program file is considered to be a Python statement and statements are executed in order from the top of the file to the bottom. Any statement proceding a "#
" (until the end of the line) is considered a comment. Python program files can be executed individually by executing the following at a command prompt:
python filename.py
The strength of Python for scientific computing lies in the large number of add-on modules on which users can build applications. Many of the most important modules, namely
are not part of the standard Python distribution. Historically, users had to install each module (and their dependencies) individually. This painful process can, in large part, be bypassed by using one of the fantastic commercially available Python distributions. I recommend annoconda or enthought. Both offer free and paid versions. The free versions offer all the utility that we will need in this book.
The numpy
package is the foundation of most scientific computations performed in Python. It is implemented in C/Fortran, so performance is greatly improved over native Python data types. numpy
provides:
ndarray
(n dimensinal array) is numpy
's primary objectmatplotlib
is a 2D and 3D graphics library for generating scientific figures. matplotlib
:
sympy
provides a Computer Algebra System (CAS) for Python. sympy
is a regular Python module and integrates very well with IPython notebooks.
Since python is a dynamically interpreted language, a variable's type need not be declared - its type is inferred at the time of assignment. Variables are assigned a value by the assignment operateor "=
"
spam = 4
spam
An assignment statement assigns to the variable name on the left of the assignment operator the expression on the right.
Reminder: To execute the cell above move focus to the cell and press
Python supports many type
s of variables. The most commonly used in this course are str
, float
, int
, bool
, list
, tuple
, and dict
.
Strings are groups of one or more characters of type str
enclosed by single or double quotes
string = "foo"
type(string)
String concatenation is by the "+" operator
post = " bar"
string + post
Python's support for string parsing, manipulation, etc. is, in my opinion, one of its greatest strengths.
Python supports many types of numbers, the most common being reals and integers
Real numbers are of type float
a = 4.
a = float(4)
a = float("4.")
type(a)
Integers have type int
b = 4
b = int(4.)
type(b)
Be careful when converting strings to integers
b = int("4.")
The above code "raised" a ValueError
. Instead, first convert to float, then integer
b = int(float("4."))
b
Arithmetic is through the standard arithmetic operators
a = 3.
b = 2.
a + b # addition
a - b # subtraction
a * b # multiplication
a / b # division
1 / 2
a ** 2 # raising to power (note, do not use ^)
Python supports a boolean type bool
with two members True
and False
. Python also supports the related None
.
Lists are mutable container objects of type list
. Lists are composed of comma separated members enclosed by paired brackets
a = [1, 2, 3]
type(a)
Lists can be appended to and extended by other lists
a.append(4)
a.extend([5, 6, 7])
a
Important: In the cell above a is an object of type list. Put differently, a is an "instance" of the list class. The dot . following a above, is a special operator that accesses "methods" associated with a. Simplistically, methods are functions that are owned by a class and are accessible only by instances of the class through the dot operator. Classes will be covered more in depth later, it is sufficient for now to recognize the action of the dot operator.
List members need not be of the same type
b = [1, 2, "spam", True, [3, 4]]
b
first = b[0]
first
last = b[-1]
last
Slices of lists can be accessed
c = a[0:3]
c
Tuples are container objects like lists, but are immutable (not modifiable after creation) and use parenthesis in place of brackets
a = (1, 2, 3, 4)
Once created, immutable objects cannot be modified. For instance, attempting a number to the end of a
will result in an error
# a.append(4) # tuples are immutable, so this will result in an error
Like lists and tuples, dictionaries are yet another container object (type: dict
). Dictionaries consist of key:value
pairs enclosed by braces. Valid keys are any immutable object (other items can be used, but for now immutable objects will suffice). For example, in the cell below a dictionary d
is created with key a
and value 5
. A new key:value
pair is then added to d
.
d = {"a": 5}
d["b"] = "spam"
d
d["b"]
A KeyError
is raised when attempting to access a non-existent key
# d[12]
The dict
constructor can be used to create a dictionary from a list of tuples
d = dict([('a', 5), ('b', 'spam'), (0, (1,2,3))])
d
Try it!
Add key:value pairs "c": [1,2,3] and "d": (3,4,5) to d.
Unlike many other languages that explicitly define the scope of code blocks (opening/closing braces, begin/end keywords, etc.), Python implicitly defines code block scope through indentation
compound_statement ":"
code block
The Python language does not define a standard indentation depth, but the following advice should be followed:
Branching within a code is handled by if
, elif
, and else
clauses
if True:
print "I'm True!"
a = 4
if a < 3:
print "Less than 3"
elif a < 6:
print "Between 3 and 6"
else:
print "Greater than 6"
Notice that the indentation of the above code defines the scope of each block.
Python provides two forms of looping: for
and while
.
for
Loops¶for
loops in Python is one area that took me some time to become accustomed to, in comparison to C or Fortran. In C, Fortran, and many other languages, loops are constructed by incrementing an integer value from a starting value to an ending value. For example, in C
for(int i=0; i<10; i++){
std::cout << i << "\n";
}
In Python, for
loops are performed by iterating through an "iterable". Iterables are objects, such as lists
, tuples
, and dict
s, that support iteration. In python, the previous loop could be written:
for i in [0, 1, 2, 3, 4, 5, 6, 7, 8, 9]:
print i
or,
for i in range(10):
print i
Strings also support iteration. To loop through the values of a string,
for letter in "a very long string":
print letter
while
Loops¶Similar to for
loops, while loops iterate until a condition is met.
a = 1
while a < 5:
a += a # lhs += rhs is a shortcut for lhs = lhs + rhs
print a
A list comprehension is a concise way of creating a list
. Consider the following
a = []
for i in range(4):
a.append(i ** 2)
a
range(9)
The above code can be simplified to
a = [i ** 2 for i in range(4)]
a
List comprehensions consist of an opening bracket, expression, a for loop, followed by 1 or more for loops and/or if statements and are ended with a closing bracket.
Functions contain collections of related Python code. Functions have the following syntax
def function_name "(" [parameter_list] ")" ":"
suite
For example
def print_hello():
print "Hello, World!"
print_hello()
def print_hello_with_args(a, b):
print "{0} says hello to {1}!".format(a, b)
print_hello_with_args("Fred", "Sue")
"Mary {0} {3} {2} {1}".format("had", "lambs", "little", 5)
Functions can return 1 or more objects
def add(a, b):
c = a + b
return c
add(4, 5)
Another type of function is the lambda
function, or "anonymous" function. lambda
s allow the creation of functions on the fly
f = lambda x: x ** 2
f(4.)
lambda
functions are useful in many situations and we'll have a chance to use them often throughtout this book.
Python modules are importable files containing Python code. Importing a module in to the current file exposes that modules contents to the current file. The syntax for importing a module takes one of several forms, demonstrated below.
import math
math.pi
math.cos(0.)
Modules names can be aliased for convenience with the as
designator
import math as m
m.pi
Alternatively, all contents of a module can be imported in to the current namespace
from math import *
pi
However, this can lead to pollution of the current namespace and is not recommended - unless the module was designed to be imported in this way. For example, suppose you do
from spam import *
from eggs import *
baz = foo
from which module did foo
come from?
The Python standard library comes with many useful modules. In this course we will also use heavily the nonstandard numpy
, sympy
, and matplotlib
modules.
numpy
¶numpy
is a powerful Python module that is the defacto standard for scientific computing in Python. It is customary to import numpy
as
import numpy as np
The power of numpy
lies in the array object that it provides. numpy
arrays are similar to Python lists, but are statically typed and contain only one object type. Arrays can be instantiated in a number of ways,
a = np.array([1., 2., 3.])
a
numpy.linspace
functiona = np.linspace(0, 10, 5) # np.linspace takes (start, stop, number of elements)
a
numpy.arange
functiona = np.arange(0, 5, 1.5) # np.arange takes (start, stop, step size)
a
And still others.
numpy
¶numpy
provides many common operations from linear algebra to operate on n dimensional arrays. Consider the following arrays
# Define 1 and 2D arrays. 2D arrays are similar to matrices in matlab
a = np.array([3, 2, 6])
b = np.array([1, 1, 7])
M = np.array([[1, 5, 4], [5, 5, 3], [8, 2, 9]])
L = np.array([[3, 1, 5], [5, 7, 2], [1, 0, 3]])
np.transpose(M)
The *
operator operates element wise on arrays, so a * b
will not yield a scalar
a * b
dot
Product¶Scalar product of two vectors
np.dot(a, b)
Matrix/Vector multiplication is also through the dot
method
np.dot(M, a)
Think of the dot
method in terms of the indicial summation representation of matrix operations, e.g. matrix/vector multiplication
$\{c\} = [A]\{b\} \Rightarrow c_i = \sum_jA_{ij}b_{j}$
would be written
c = np.dot(A, b)
Or matrix multiplication,
$[C] = [A][B] \Rightarrow C_{ij} = \sum_m A_{im}B_{mj}$
is
C = np.dot(A, B)
Or the vector/vector (scalar) product
$c = \{v\}\cdot\{w\} = \sum_i v_i w_i$
is
c = np.dot(v, w)
np.linalg.inv(M)
# Determinant of a 2D array
np.linalg.det(M)
solve for $\{x\}$: $[M]\{x\} = \{b\}$
# Solution to system of equations (using inv is *very* inefficient)
x = np.dot(np.linalg.inv(M), b)
x
Better to use the solve
method
x = np.linalg.solve(M, b)
x
sympy
¶sympy
is a Computer Algebra System implemented in python. sympy
allows the creation of symbols and functions (among many other objects) and allows algebraic manipulation of those objects.
import sympy as sp
sp.init_printing() # allow pretty math output
Variables can be designated as Symbol
s and used in symbolic computation
x = sp.Symbol("x")
y, z = sp.symbols("y z")
f = x + y * sp.pi
f
Values can be substituted in place of symbols
f.subs({x: 2, y: 3})
Evaluate the value numerically
f.subs({x: 2, y: 3}).evalf()
sympy
¶Evaluate the following integral: $\int_{0}^{2} x \,dx$
expr = sp.integrate("x", ("x", 0, 2))
expr
Lets compare with numpy
's built in trapz
method
func = lambda x: x
xvals = np.linspace(0, 2, 5)
func_vals = np.array([func(x) for x in xvals])
a = np.trapz(func_vals, x=xvals)
a
expr.evalf() == a
Find $u(x)$ for $u'(x) + u(x) = x$, with $u(0) = 3$. Plot the solution on $x\in[0,10]$
x = sp.Symbol("x")
u = sp.Function("u")
# recast equation so that all terms appear on left hand side
ics = {u(0): 3}
de = sp.diff(u(x), x) + u(x) - x
gen_sol = sp.dsolve(de, u(x)).rhs
spec_sol = sp.dsolve(de, u(x), ics=ics, simplify=False)
gen_sol
spec_sol
Find the coefficients
coeffs = sp.solve([gen_sol.subs(x, 0) - 3], "C1")
sol = gen_sol.subs(coeffs)
sol
Now we plot with matplotlib
. matplotlib
offers plotting functionality very similar to Matlab
# "magic" command below allows inline printing of plots
import matplotlib.pyplot as plt
%matplotlib inline
dx = .25
xvals = np.arange(0, 10, dx)
sol_vals = [float(sol.subs({x: xval}).evalf()) for xval in xvals]
plt.plot(xvals, sol_vals)
Let's now find $u(x)$ numerically using finite differencing.
The forward difference operator is
$\frac{df}{dx} \approx \frac{f(x + \Delta{x}) - f(x)}{\Delta{x}}$
Substituting the finite difference relation for $u'$ gives
$u(x + \Delta{x}) = u(x)(1 - \Delta{x}) + x \Delta{x}$
Since $u(x=0)$ is known ($u(0) = 3$), we are now ready to solve for $u(x + \Delta{x})$ iteratively
# dx and xvals defined above
u = [3.]
for x in xvals[1:]:
u.append(u[-1] * (1 - dx) + x * dx)
u = np.array(u)
# Evaluate solution over interval [0, 10] and plot
plt.plot(xvals, sol_vals, label="Analytic")
plt.plot(xvals, u, "r+", label="Approximate")
plt.xlabel("x")
plt.ylabel("u(x)")
plt.legend(loc="best");