3848 days ago by SC2011



Introduction to Python programming within Sage

Neil Calkin, Holly Hirst, and Dan Warner



Mission Statement

Create a viable free open source alternative to Magma, Maple, Mathematica, and Matlab.

What is Sage?

  • Sage = Python + Math
  • A unified self-contained distribution of open source free mathematical software.
  • Nearly a half million lines of new Python (and Cython) code that implements new capabilities and algorithms.
  • A "cloud" application like GMail or Google Docs:  (currently about 30,000 users);  Of course Sage also runs on your desktop and supercomputer.



Looking a little more carefully at how Sage works with Python.

First, Sage makes extensive use of classes

# What is pi? type(pi) 

# How do we tell Sage to compute pi to 4493 digits if we don't remember the command? # The tab key after "." of "(" generates a help display which can be used to # autocomplete the statement. pi. 
# Sage can solve equations for one variable in terms of others. # var('x b c') show(solve([x^2 + b*x + c == 0],x)) 

This example illustrated the fact that Python variables by default are "programming variables" and must be initialized. 

To create a symbolic variable we must declare the variable with the var command. 

The variable x is a special case for convenience.

# Sage can even find the ugly answers eq3 = x^3 - 2 == x show(eq3) sol = solve([eq3],x) for z in sol: show(z) 

This example illustrated the for loop in Python. 

First, the for statement which walks through a list of objects. It ends with the colon. 

Second, the block inside the loop is indented.  There are no begin ... end constructs in Python. 

Good programming style strongly suggests that you indent subordinate blocks of code. 

Python insists that you indent and rewards you with shorter programs.

# Sage can solve linear systems of algebraic equations # Note the introduction of the symbolic variable y var('y') eq4 = x-1 == y+1 eq5 = x+1 == 2*(y-1) eq4;eq5 

This example underscored the fact that the single "=" is assignment and the double "==" is equality.

It also showed that multiple statements can be on the same line when separated by ";"

sol2 = solve([eq4,eq5],x,y) sol2 

In this example, the solution is returned in a list.  Lists are the principal data structure workhorse.

# Sage can also return the solutions in a dictionary # This provides access by name sol2 = solve([eq4,eq5],x,y,solution_dict=True) sd = sol2[0] print sd (sd[x].n(digits=2),sd[y].n(digits=3)) 

In this example we used the dictionary, which is the most sophisticated built-in data structure in Python.


Python Data Structures

The key data structures are: Lists, Tuples, and Dictionaries

Creating Lists I: [Square brackets]

L = [3, Permutation([5,1,4,2,3]), 17, 17, 3, 51] L 

Exercise: Create the list [63, 12, -10, 'a', 12], assign it to the variable L, and print the list.


Exercise: Create the empty list. (You will often need to do this.)


Creating Lists II: range

The range function provides an easy way to construct a list of integers. Here is the documentation of the range function.


range([start,] stop[, step]) -> list of integers

Return a list containing an arithmetic progression of integers.
range(i, j) returns [i, i+1, i+2, ..., j-1]; start (!) defaults to 0.
When step is given, it specifies the increment (or decrement).
For example, range(4) returns [0, 1, 2, 3].  The end point is omitted!
These are exactly the valid indices for a list of 4 elements.

Exercise: Use range to construct the list $[1,2,\ldots,50]$.


Exercise: Use range to construct the list of even numbers between 1 and 100 (including 100).


Exercise: The step argument for the range command can be negative. Use range to construct the list $[10, 7, 4, 1, -2]$.


Creating Lists III: list comprehensions

List comprehensions provide a concise way to create lists from other lists (or other data types). This technique is very similar to the mathematical notation for describing a set.

Example. We already know how to create the list $[1, 2, \ldots, 16]$:


Using a list comprehension, we can create the list of even numbers as follows:

E = [x for x in range(1,101) if x%2==0] E 

Using a list comprehension, we can also create the list $[1^2, 2^2, 3^2, ..., 16^2]$ as follows:

[i^2 for i in range(1,17)]
[i^2 for i in range(1,17)] 
sum([i^2 for i in range(1,17)]) 


The sum of the squares of the first ten natural numbers is,

1^(2) + 2^(2) + ... + 10^(2) = 385

The square of the sum of the first ten natural numbers is,

(1 + 2 + ... + 10)^(2) = 55^(2) = 3025

Hence the difference between the sum of the squares of the first ten natural numbers and the square of the sum is 3025 − 385 = 2640.

Find the difference between the sum of the squares of the first one hundred natural numbers and the square of the sum.

s1 = sum( 
s2 = ( 
s2 - s1 

Filtering lists with a list comprehension

A list can be filtered using a list comprehension.

Example: To create a list of the squares of the prime numbers between 1 and 100, we use a list comprehension as follows.

[p^2 for p in [1,2,..,100] if is_prime(p)] 
[p^2 for p in range(1,101) if is_prime(p)] 

Exercise: Use a list comprehension to list all the natural numbers below 20 that are multiples of 3 or 5.


  1. To get the remainder of 7 divided by 3 use 7%3.
  2. To test for equality use two equal signs (==); for example, 3 == 7.

Nested list comprehensions

List comprehensions can be nested!


[(x,y) for x in range(5) for y in range(3)] 
[[i^j for j in range(1,4)] for i in range(6)] 
matrix([[i^j for j in range(1,4)] for i in range(6)]) 


Pythagorean triple is a triple $(x,y,z)$ of positive integers satisfying $x^2+y^2=z^2$. 

The Pythagorean triples whose components are at most $10$ are :  $[(3, 4, 5), (4, 3, 5), (6, 8, 10), (8, 6, 10)]$.

Using a filtered list comprehension, construct the list of Pythagorean triples whose components are at most $50$.


Accessing individual elements of lists

To access an element of the list, use the syntax L[i], where i is the index of the item.


  1. Construct the list L = [1,2,3,4,3,5,6]. What is L[3]?
  2. What is L[1]?
  3. What is the index of the first element of L?
  4. What is L[-1]? What is L[-2]?
  5. What is L.index(2)? What is L.index(3)?


L = [1,2,..,6];L 

Modifying lists: changing an element in a list

To change the item in position i of a list L:

L = ["a", 4, 1, 8] L 
L[2] = 0 L 
L[1:3]=[1,2,-1] L 

Modifying lists: append and extend

To append an object to a list:

L = ["a", 4, 1, 8] L 
L.append(17) L 

To extend a list by another list:

L1 = [1,2,3] L2 = [7,8,9,0] 
print L1 print L2 
L1.extend(L2) L1 

Modifying lists: reverse, sort, ...

L = [4,2,5,1,3] L 
L.reverse() L 
L.sort() L 
L = [3,1,6,4] 

Concatenating Lists

To concatenate two lists, add them (+). This is not a commutative operation....

L1 = [1,2,3] L2 = [7,8,9,0] print L1 print L2 
L1 + L2 

Slicing Lists

You can slice a list using the syntax L[start : stop : step]. This will return a sublist of L.

Exercise: Below are some examples of slicing lists. Try to guess what the output will be before evaluating the cell.

L = range(20) L 
L[3:5] = [-2,'a'];L 

Exercise: The following function combines a loop with the some of the list operations above. What does the function do?

def f(number_of_iterations): L = [1] for n in range(2, number_of_iterations): L = [sum(L[:i]) for i in range(n-1, -1, -1)] return numerical_approx(2*L[0]*len(L)/sum(L), digits=50) 


A tuple is an immutable list. That is, it cannot be changed once it is created.

The syntax for creating a tuple is to the use parentheses. 

t = (3, 5, [3,1], (17,[2,3],17), 4) t 

We can create a tuple from a list, or vice-versa.


Tuples behave like lists in many respects:

Operation Syntax for lists Syntax for tuples
Accessing a letter
list[3] tuple[3]
Concatenation list1 + list2 tuple1 + tuple2
Slicing list[3:17:2] tuple[3:17:2]
A reversed copy
list[::-1] tuple[::-1]
Length len(list) len(tuple)

Trying to modify a tuple will fail.

t = (5, 'a', 6/5) t 
t[1] = 'b' 


"Tuple-comprehension" does not exist.  The syntax produces something called a generator.  A generator allows you to process a sequence of items one at a time.  Each item is created when it is needed, and then forgotten.  This can be very efficient if we only need to use each item once.

(i^2 for i in range(5)) 
g = (i^2 for i in range(5)) 
[x for x in g] 

g is now empty.

[x for x in g] 

A nice 'pythonic' trick is to use generators as the argument to functions.  We do not need double parentheses for this.

sum( i^2 for i in srange(100001) ) 


A dictionary is another built-in data type. Unlike lists, which are indexed by a range of numbers, dictionaries are "indexed" by keys, which can be any immutable object. Strings and numbers can always be keys (because they are immutable). Dictionaries are sometimes called "associative arrays" in other programming languages.

There are several ways to define dictionaries. One method is to use braces, {}, with comma-separated entries given in the form key:value.

d = {3:17, "key":[4,1,5,2,3], (3,1,2):"goo", 3/2 : 17} d 

Dictionaries behave as lists and tuples for several important operations.

Operation Syntax for lists Syntax for dictionaries
Accessing elements
L[3] D["key"]
Length len(L) len(D)
Modifying L[3] = 17 D["key"] = 17
Deleting items del L[3] del D["key"]

A dictionary can have the same value multiple times, but each key can only appear once and must be immutable.

d = {3 : 14, 4: 14} d 
d = {3: 13, 3:14} d 
d = {[1,2,3] : 12} d 
d = {(1,2,3) : 12} d 

Another way to add items to a dictionary is with the update() method:

d = {} d 
d.update( {10 : 'newvalue', 20: 'newervalue', 3: 14, 'a':[1,2,3], 15:[]} ) d 
d[15].append(7) d 

We can iterate through the keys, or values, or both, of a dictionary.

d = {10 : 'newvalue', 20: 'newervalue', 3: 14, 'a':[1,2,3]} 
[key for key in d] 
[key for key in d.iterkeys()] 
[value for value in d.itervalues()] 
[(key, value) for key, value in d.iteritems()] 

Exercise:  Consider a digraph whose vertices are {1, 2, 3} and whose directed edges are {(1,2), (1,3), (2,3)} 

Create a dictionary whose keys are the vertices of the above directed graph, and whose values are the lists of the vertices that it points to. For instance, the vertex 1 points to the vertices 2 and 3, so the dictionary will look like:

d = { ..., 1:[2,3], ... }

Then try

g = DiGraph(d)



An important class that is intimately related to dictionaries is the set.

greek = set(['alpha','beta','gamma','delta','epsilon','zeta','beta']) greek 
set(['epsilon', 'beta', 'delta', 'alpha', 'zeta', 'gamma'])
set(['epsilon', 'beta', 'delta', 'alpha', 'zeta', 'gamma'])
greek.add('theta') greek.remove('delta') greek 
set(['epsilon', 'beta', 'alpha', 'zeta', 'theta', 'gamma'])
set(['epsilon', 'beta', 'alpha', 'zeta', 'theta', 'gamma'])
greek.add('delta') greek.add('alpha') greek 
set(['epsilon', 'beta', 'delta', 'alpha', 'zeta', 'theta', 'gamma'])
set(['epsilon', 'beta', 'delta', 'alpha', 'zeta', 'theta', 'gamma'])

Using sage types: The srange command

Example:  Construct a $3 \times 3$ matrix whose $(i,j)$ entry is the rational number $\frac{i}{j}$.

matrix([[ i/j for j in range(1,4)] for i in range(1,4)]) 
[1 0 0]
[2 1 0]
[3 1 1]
[1 0 0]
[2 1 0]
[3 1 1]
matrix([[ i/j for j in srange(1,4)] for i in srange(1,4)]) 
[  1 1/2 1/3]
[  2   1 2/3]
[  3 3/2   1]
[  1 1/2 1/3]
[  2   1 2/3]
[  3 3/2   1]
[0.500000000000000, 0.700000000000000, 0.900000000000000,
1.10000000000000, 1.30000000000000, 1.50000000000000, 1.70000000000000,
1.90000000000000, 2.10000000000000, 2.30000000000000, 2.50000000000000,
2.70000000000000, 2.90000000000000, 3.10000000000000, 3.30000000000000]
[0.500000000000000, 0.700000000000000, 0.900000000000000, 1.10000000000000, 1.30000000000000, 1.50000000000000, 1.70000000000000, 1.90000000000000, 2.10000000000000, 2.30000000000000, 2.50000000000000, 2.70000000000000, 2.90000000000000, 3.10000000000000, 3.30000000000000]
<type 'sage.rings.integer.Integer'>
<type 'sage.rings.integer.Integer'>
<type 'int'>
<type 'sage.rings.integer.Integer'>
<type 'int'>
<type 'sage.rings.integer.Integer'>

Modifying lists has consequences!

Try to predict the results of the following commands.

a = [1, 2, 3] L = [a, a, a] L 
a.append(4) L 

Now try these:

a = [1, 2, 3] L = [a, a, a] L 
a = [1, 2, 3, 5] L 
L[0].append(4) L 

You can use the command deepcopy to avoid these issues.

a = [1,2,3] L = [deepcopy(a), deepcopy(a)] L 
a.append(4) L 

The same issues occur with dictionaries.

d = {1:'a', 2:'b', 3:'c'} dd = d d.update( { 4:'d' } ) dd