Open-Source Math Software
Open-Source Math Software
Open-Source Math Software
Release 5.6
CONTENTS
Introduction 1.1 Installation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.2 Ways to Use Sage . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.3 Longterm Goals for Sage . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . A Guided Tour 2.1 Assignment, Equality, and Arithmetic 2.2 Getting Help . . . . . . . . . . . . . 2.3 Functions, Indentation, and Counting 2.4 Basic Algebra and Calculus . . . . . 2.5 Plotting . . . . . . . . . . . . . . . . 2.6 Some Common Issues with Functions 2.7 Basic Rings . . . . . . . . . . . . . . 2.8 Linear Algebra . . . . . . . . . . . . 2.9 Polynomials . . . . . . . . . . . . . 2.10 Parents, Conversion and Coercion . . 2.11 Finite Groups, Abelian Groups . . . . 2.12 Number Theory . . . . . . . . . . . 2.13 Some More Advanced Mathematics .
3 4 4 4 7 7 9 10 13 18 21 24 26 29 33 38 39 41 51 51 53 54 54 55 56 57 58 60 61 62 65 65 66 67 68
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
The Interactive Shell 3.1 Your Sage Session . . . . . . . . . . . 3.2 Logging Input and Output . . . . . . . 3.3 Paste Ignores Prompts . . . . . . . . . 3.4 Timing Commands . . . . . . . . . . . 3.5 Other IPython tricks . . . . . . . . . . 3.6 Errors and Exceptions . . . . . . . . . 3.7 Reverse Search and Tab Completion . . 3.8 Integrated Help System . . . . . . . . 3.9 Saving and Loading Individual Objects 3.10 Saving and Loading Complete Sessions 3.11 The Notebook Interface . . . . . . . . Interfaces 4.1 GP/PARI 4.2 GAP . . 4.3 Singular . 4.4 Maxima .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
. . . .
Sage, LaTeX and Friends 5.1 Overview . . . . . . . . . . . . . . . . . . . . . . . 5.2 Basic Use . . . . . . . . . . . . . . . . . . . . . . . 5.3 Customizing LaTeX Generation . . . . . . . . . . . 5.4 Customizing LaTeX Processing . . . . . . . . . . . 5.5 An Example: Combinatorial Graphs with tkz-graph . 5.6 A Fully Capable TeX Installation . . . . . . . . . . 5.7 External Programs . . . . . . . . . . . . . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
Programming 6.1 Loading and Attaching Sage les . . . . . . . . . . . . 6.2 Creating Compiled Code . . . . . . . . . . . . . . . . . 6.3 Standalone Python/Sage Scripts . . . . . . . . . . . . . 6.4 Data Types . . . . . . . . . . . . . . . . . . . . . . . . 6.5 Lists, Tuples, and Sequences . . . . . . . . . . . . . . . 6.6 Dictionaries . . . . . . . . . . . . . . . . . . . . . . . . 6.7 Sets . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6.8 Iterators . . . . . . . . . . . . . . . . . . . . . . . . . . 6.9 Loops, Functions, Control Statements, and Comparisons 6.10 Proling . . . . . . . . . . . . . . . . . . . . . . . . . Using SageTeX
. . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . .
7 8
Afterword 8.1 Why Python? . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8.2 I would like to contribute somehow. How can I? . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8.3 How do I reference Sage? . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Appendix 9.1 Arithmetical binary operator precedence . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
ii
Sage is free, open-source math software that supports research and teaching in algebra, geometry, number theory, cryptography, numerical computation, and related areas. Both the Sage development model and the technology in Sage itself are distinguished by an extremely strong emphasis on openness, community, cooperation, and collaboration: we are building the car, not reinventing the wheel. The overall goal of Sage is to create a viable, free, open-source alternative to Maple, Mathematica, Magma, and MATLAB. This tutorial is the best way to become familiar with Sage in only a few hours. You can read it in HTML or PDF versions, or from the Sage notebook (click Help, then click Tutorial to interactively work through the tutorial from within Sage). This work is licensed under a Creative Commons Attribution-Share Alike 3.0 License.
CONTENTS
CONTENTS
CHAPTER
ONE
INTRODUCTION
This tutorial should take at most 3-4 hours to fully work through. You can read it in HTML or PDF versions, or from the Sage notebook click Help, then click Tutorial to interactively work through the tutorial from within Sage. Though much of Sage is implemented using Python, no Python background is needed to read this tutorial. You will want to learn Python (a very fun language!) at some point, and there are many excellent free resources for doing so including [PyT] and [Dive]. If you just want to quickly try out Sage, this tutorial is the place to start. For example:
sage: 2 + 2 4 sage: factor(-2007) -1 * 3^2 * 223 sage: A = matrix(4,4, range(16)); A [ 0 1 2 3] [ 4 5 6 7] [ 8 9 10 11] [12 13 14 15] sage: factor(A.charpoly()) x^2 * (x^2 - 30*x - 80) sage: m = matrix(ZZ,2, range(4)) sage: m[0,0] = m[0,0] - 3 sage: m [-3 1] [ 2 3] sage: E = EllipticCurve([1,2,3,4,5]); sage: E Elliptic Curve defined by y^2 + x*y + 3*y = x^3 + 2*x^2 + 4*x + 5 over Rational Field sage: E.anlist(10) [0, 1, 1, 0, -1, -3, 0, -1, -3, -3, -3] sage: E.rank() 1 sage: k = 1/(sqrt(3)*I + 3/4 + sqrt(73)*5/9); k 1/(I*sqrt(3) + 5/9*sqrt(73) + 3/4) sage: N(k) 0.165495678130644 - 0.0521492082074256*I sage: N(k,30) # 30 "bits" 0.16549568 - 0.052149208*I sage: latex(k) \frac{1}{i \, \sqrt{3} + \frac{5}{9} \, \sqrt{73} + \frac{3}{4}}
1.1 Installation
If you do not have Sage installed on a computer and just want to try some commands, use online at http://www.sagenb.org. See the Sage Installation Guide in the documentation section of the main Sage webpage [SA] for instructions on installing Sage on your computer. Here we merely make a few comments. 1. The Sage download le comes with batteries included. In other words, although Sage uses Python, IPython, PARI, GAP, Singular, Maxima, NTL, GMP, and so on, you do not need to install them separately as they are included with the Sage distribution. However, to use certain Sage features, e.g., Macaulay or KASH, you must install the relevant optional package or at least have the relevant programs installed on your computer already. Macaulay and KASH are Sage packages (for a list of available optional packages, type sage -optional, or browse the Download page on the Sage website). 2. The pre-compiled binary version of Sage (found on the Sage web site) may be easier and quicker to install than the source code version. Just unpack the le and run sage. 3. If youd like to use the SageTeX package (which allows you to embed the results of Sage computations into a LaTeX le), you will need to make SageTeX known to your TeX distribution. To do this, see the section Make SageTeX known to TeX in the Sage installation guide (this link should take you to a local copy of the installation guide). Its quite easy; you just need to set an environment variable or copy a single le to a directory that TeX will search. The documentation for using SageTeX is located in $SAGE_ROOT/local/share/texmf/tex/generic/sagetex/, where $SAGE_ROOT refers to the directory where you installed Sage for example, /opt/sage-4.2.1.
Chapter 1. Introduction
in a paper you publish, you can rest assured that your readers will always have free access to Sage and all its source code, and you are even allowed to archive and re-distribute the version of Sage you used. Easy to compile: Sage should be easy to compile from source for Linux, OS X and Windows users. This provides more exibility for users to modify the system. Cooperation: Provide robust interfaces to most other computer algebra systems, including PARI, GAP, Singular, Maxima, KASH, Magma, Maple, and Mathematica. Sage is meant to unify and extend existing math software. Well documented: Tutorial, programming guide, reference manual, and how-to, with numerous examples and discussion of background mathematics. Extensible: Be able to dene new data types or derive from built-in types, and use code written in a range of languages. User friendly: It should be easy to understand what functionality is provided for a given object and to view documentation and source code. Also attain a high level of user support.
Chapter 1. Introduction
CHAPTER
TWO
A GUIDED TOUR
This section is a guided tour of some of what is available in Sage. For many more examples, see Sage Constructions, which is intended to answer the general question How do I construct ...?. See also the Sage Reference Manual, which has thousands more examples. Also note that you can interactively work through this tour in the Sage notebook by clicking the Help link. (If you are viewing the tutorial in the Sage notebook, press shift-enter to evaluate any input cell. You can even edit the input before pressing shift-enter. On some Macs you might have to press shift-return rather than shift-enter.)
4 * (10 // 4) + 10 % 4 == 10
The computation of an expression like 3^2*4 + 2%5 depends on the order in which the operations are applied; this is specied in the operator precedence table in Arithmetical binary operator precedence. Sage also provides many familiar mathematical functions; here are just a few examples:
sage: sqrt(3.4) 1.84390889145858 sage: sin(5.135) -0.912021158525540 sage: sin(pi/3) 1/2*sqrt(3)
As the last example shows, some mathematical expressions return exact values, rather than numerical approximations. To get a numerical approximation, use either the function n or the method n (and both of these have a longer name, numerical_approx, and the function N is the same as n)). These take optional arguments prec, which is the requested number of bits of precision, and digits, which is the requested number of decimal digits of precision; the default is 53 bits of precision.
sage: exp(2) e^2 sage: n(exp(2)) 7.38905609893065 sage: sqrt(pi).numerical_approx() 1.77245385090552 sage: sin(10).n(digits=5) -0.54402 sage: N(sin(10),digits=10) -0.5440211109 sage: numerical_approx(pi, prec=200) 3.1415926535897932384626433832795028841971693993751058209749
Python is dynamically typed, so the value referred to by each variable has a type associated with it, but a given variable may hold values of any Python type within a given scope:
sage: sage: <type sage: sage: <type sage: sage: <type a = 5 # a is an integer type(a) sage.rings.integer.Integer> a = 5/3 # now a is a rational number type(a) sage.rings.rational.Rational> a = hello # now a is a string type(a) str>
The C programming language, which is statically typed, is much different; a variable declared to hold an int can only hold an int in its scope. A potential source of confusion in Python is that an integer literal that begins with a zero is treated as an octal number, i.e., a number in base 8.
sage: 9 sage: 9 sage: sage: 11 011 8 + 1 n = 011 n.str(8)
The tangent function EXAMPLES: sage: tan(pi) 0 sage: tan(3.1415) -0.0000926535900581913 sage: tan(3.1415/4) 0.999953674278156 sage: tan(pi/4) 1 sage: tan(1/2) tan(1/2) sage: RR(tan(1/2)) 0.546302489843790 sage: log2? Type: <class sage.functions.constants.Log2> Definition: log2( [noargspec] ) Docstring: The natural logarithm of the real number 2. EXAMPLES: sage: log2 log2 sage: float(log2) 0.69314718055994529 sage: RR(log2) 0.693147180559945 sage: R = RealField(200); R Real Field with 200 bits of precision sage: R(log2) 0.69314718055994530941723212145817656807550013436025525412068 sage: l = (1-log2)/(1+log2); l (1 - log(2))/(log(2) + 1) sage: R(l) 0.18123221829928249948761381864650311423330609774776013488056 sage: maxima(log2) log(2) sage: maxima(log2).float() .6931471805599453 sage: gp(log2) 0.6931471805599453094172321215 # 32-bit 0.69314718055994530941723212145817656807 # 64-bit sage: sudoku?
Solve the 9x9 Sudoku puzzle defined by the matrix A. EXAMPLE: sage: A = matrix(ZZ,9,[5,0,0, 0,8,0, 0,4,9, 0,0,0, 5,0,0, 0,3,0, 0,6,7, 3,0,0, 0,0,1, 1,5,0, 0,0,0, 0,0,0, 0,0,0, 2,0,8, 0,0,0, 0,0,0, 0,0,0, 0,1,8, 7,0,0, 0,0,4, 1,5,0, 0,3,0, 0,0,2, 0,0,0, 4,9,0, 0,5,0, 0,0,3]) sage: A [5 0 0 0 8 0 0 4 9] [0 0 0 5 0 0 0 3 0] [0 6 7 3 0 0 0 0 1] [1 5 0 0 0 0 0 0 0] [0 0 0 2 0 8 0 0 0] [0 0 0 0 0 0 0 1 8] [7 0 0 0 0 4 1 5 0] [0 3 0 0 0 2 0 0 0] [4 9 0 0 5 0 0 0 3] sage: sudoku(A) [5 1 3 6 8 7 2 4 9] [8 4 9 5 2 1 6 3 7] [2 6 7 3 4 9 5 8 1] [1 5 8 4 6 3 9 7 2] [9 7 4 2 1 8 3 6 5] [3 2 6 7 9 5 4 1 8] [7 8 2 9 3 4 1 5 6] [6 3 5 1 7 2 8 9 4] [4 9 1 8 5 6 7 2 3]
Sage also provides Tab completion: type the rst few letters of a function and then hit the tab key. For example, if you type ta followed by TAB, Sage will print tachyon, tan, tanh, taylor. This provides a good way to nd the names of functions and other structures in Sage.
Note: Depending on which version of the tutorial you are viewing, you may see three dots ... on the second line of this example. Do not type them; they are just to emphasize that the code is indented. Whenever this is the case, press [Return/Enter] once at the end of the block to insert a blank line and conclude the function denition. You do not specify the types of any of the input arguments. You can specify multiple inputs, each of which may have an optional default value. For example, the function below defaults to divisor=2 if divisor is not specied.
10
You can also explicitly specify one or either of the inputs when calling the function; if you specify the inputs explicitly, you can give them in any order:
sage: is_divisible_by(6, divisor=5) False sage: is_divisible_by(divisor=2, number=6) True
In Python, blocks of code are not indicated by curly braces or begin and end blocks as in many other languages. Instead, blocks of code are indicated by indentation, which must match up exactly. For example, the following is a syntax error because the return statement is not indented the same amount as the other lines above it.
sage: def even(n): ... v = [] ... for i in range(3,n): ... if i % 2 == 0: ... v.append(i) ... return v Syntax Error: return v
Semicolons are not needed at the ends of lines; a line is in most cases ended by a newline. However, you can put multiple statements on one line, separated by semicolons:
sage: a = 5; b = a + 3; c = b^2; c 64
If you would like a single line of code to span multiple lines, use a terminating backslash:
sage: 2 + \ ... 3 5
In Sage, you count by iterating over a range of integers. For example, the rst line below is exactly like for(i=0; i<3; i++) in C++ or Java:
sage: for i in range(3): ... print i 0
11
1 2
The third argument controls the step, so the following is like for(i=1;i<6;i+=2).
sage: for i in range(1,6,2): ... print i 1 3 5
Often you will want to create a nice table to display numbers you have computed using Sage. One easy way to do this is to use string formatting. Below, we create three columns each of width exactly 6 and make a table of squares and cubes.
sage: for i in range(5): ... print %6s %6s %6s%(i, i^2, i^3) 0 0 0 1 1 1 2 4 8 3 9 27 4 16 64
The most basic data structure in Sage is the list, which is as the name suggests just a list of arbitrary objects. For example, the range command that we used creates a list:
sage: range(2,10) [2, 3, 4, 5, 6, 7, 8, 9]
Use len(v) to get the length of v, use v.append(obj) to append a new object to the end of v, and use del v[i] to delete the ith entry of v:
sage: len(v) 4 sage: v.append(1.5) sage: v [1, hello, 2/3, sin(x^3), 1.50000000000000] sage: del v[1] sage: v [1, 2/3, sin(x^3), 1.50000000000000]
12
Another important data structure is the dictionary (or associative array). This works like a list, except that it can be indexed with almost any object (the indices must be immutable):
sage: d = {hi:-2, sage: d[hi] -2 sage: d[e] pi 3/8:pi, e:pi}
You can also dene new data types using classes. Encapsulating mathematical objects with classes is a powerful technique that can help to simplify and organize your Sage programs. Below, we dene a class that represents the list of even positive integers up to n; it derives from the builtin type list.
sage: class Evens(list): ... def __init__(self, n): ... self.n = n ... list.__init__(self, range(2, n+1, 2)) ... def __repr__(self): ... return "Even positive numbers up to n."
The __init__ method is called to initialize the object when it is created; the __repr__ method prints the object out. We call the list constructor method in the second line of the __init__ method. We create an object of class Evens as follows:
sage: e = Evens(10) sage: e Even positive numbers up to n.
Note that e prints using the __repr__ method that we dened. To see the underlying list of numbers, use the list function:
sage: list(e) [2, 4, 6, 8, 10]
13
The following example of using Sage to solve a system of non-linear equations was provided by Jason Grout: rst, we solve the system symbolically:
sage: var(x y p q) (x, y, p, q) sage: eq1 = p+q==9 sage: eq2 = q*y+p*x==-6 sage: eq3 = q*y^2+p*x^2==24 sage: solve([eq1,eq2,eq3,p==1],p,q,x,y) [[p == 1, q == 8, x == -4/3*sqrt(10) - 2/3, y == 1/6*sqrt(2)*sqrt(5) - 2/3], [p == 1, q == 8, x == 4/3*sqrt(10) - 2/3, y == -1/6*sqrt(2)*sqrt(5) - 2/3]]
(The function n prints a numerical approximation, and the argument is the number of bits of precision.) Solving Equations Numerically Often times, solve will not be able to nd an exact solution to the equation or equations specied. When it fails, you can use find_root to nd a numerical solution. For example, solve does not return anything interesting for the following equation:
sage: theta = var(theta) sage: solve(cos(theta)==sin(theta), theta) [sin(theta) == cos(theta)]
On the other hand, we can use find_root to nd a solution to the above equation in the range 0 < < /2:
sage: phi = var(phi) sage: find_root(cos(phi)==sin(phi),0,pi/2) 0.785398163397448...
14
x sin(x2 ) dx and
dx
1 x2 1 :
This uses Sages interface to Maxima [Max], and so its output may be a bit different from other Sage output. In this case, this says that the general solution to the differential equation is x(t) = et (et + c). You can compute Laplace transforms also; the Laplace transform of t2 et sin(t) is computed as follows:
sage: s = var("s") sage: t = var("t") sage: f = t^2*exp(t) - sin(t) sage: f.laplace(t,s) 2/(s - 1)^3 - 1/(s^2 + 1)
Here is a more involved example. The displacement from equilibrium (respectively) for a coupled spring attached to a wall on the left
|------\/\/\/\/\---|mass1|----\/\/\/\/\/----|mass2| spring1 spring2
15
is modeled by the system of 2nd order differential equations m1 x1 + (k1 + k2 )x1 k2 x2 = 0 m2 x2 + k2 (x2 x1 ) = 0, where mi is the mass of object i, xi is the displacement from equilibrium of mass i, and ki is the spring constant for spring i. Example: Use Sage to solve the above problem with m1 = 2, m2 = 1, k1 = 4, k2 = 2, x1 (0) = 3, x1 (0) = 0, x2 (0) = 3, x2 (0) = 0. Solution: Take the Laplace transform of the rst equation (with the notation x = x1 , y = x2 ):
sage: de1 = maxima("2*diff(x(t),t, 2) + 6*x(t) - 2*y(t)") sage: lde1 = de1.laplace("t","s"); lde1 2*(-?%at(diff(x(t),t,1),t=0)+s^2*laplace(x(t),t,s)-x(0)*s)-2*laplace(y(t),t,s)+6*laplace(x(t),t,s
This is hard to read, but it says that 2x (0) + 2s2 X(s) 2sx(0) 2Y (s) + 6X(s) = 0
(where the Laplace transform of a lower case function like x(t) is the upper case function X(s)). Take the Laplace transform of the second equation:
sage: de2 = maxima("diff(y(t),t, 2) + 2*y(t) - 2*x(t)") sage: lde2 = de2.laplace("t","s"); lde2 -?%at(diff(y(t),t,1),t=0)+s^2*laplace(y(t),t,s)+2*laplace(y(t),t,s)-2*laplace(x(t),t,s)-y(0)*s
Plug in the initial conditions for x(0), x (0), y(0), and y (0), and solve the resulting two equations:
sage: var(s X Y) (s, X, Y) sage: eqns = [(2*s^2+6)*X-2*Y == 6*s, -2*X +(s^2+2)*Y == 3*s] sage: solve(eqns, X,Y) [[X == 3*(s^3 + 3*s)/(s^4 + 5*s^2 + 4), Y == 3*(s^3 + 5*s)/(s^4 + 5*s^2 + 4)]]
16
t = var(t) P = parametric_plot((cos(2*t) + 2*cos(t), 4*cos(t) - cos(2*t) ),\ (t, 0, 2*pi), rgbcolor=hue(0.9)) show(P)
For more on plotting, see Plotting. See section 5.5 of [NagleEtAl2004] for further information on differential equations.
we want to nd the approximate value of the solution at x = b with b > a. Recall from the denition of the derivative that y (x) y(x + h) y(x) , h
y(x+h)y(x) . h
where h > 0 is given and small. This and the DE together give f (x, y(x)) y(x + h) y(x) + h f (x, y(x)).
If we call hf (x, y(x)) the correction term (for lack of anything better), call y(x) the old value of y, and call y(x + h) the new value of y, then this approximation can be re-expressed as ynew yold + h f (x, yold ).
If we break the interval from a to b into n steps, so that h = in a table. x a a+h a + 2h ... b = a + nh y c c + hf (a, c) ... ??? hf (x, y) hf (a, c) ...
ba n ,
...
The goal is to ll out all the blanks of the table, one row at a time, until we reach the ??? entry, which is the Eulers method approximation for y(b). The idea for systems of ODEs is similar. Example: Numerically approximate z(t) at t = 1 using 4 steps of Eulers method, where z + tz + z = 0, z(0) = 1, z (0) = 0. 2.4. Basic Algebra and Calculus 17
We must reduce the 2nd order ODE down to a system of two rst order DEs (using x = z, y = z ) and apply Eulers method:
sage: t,x,y = PolynomialRing(RealField(10),3,"txy").gens() sage: f = y; g = -x - y * t sage: eulers_method_2x2(f,g, 0, 1, 0, 1/4, 1) t x h*f(t,x,y) y 0 1 0.00 0 1/4 1.0 -0.062 -0.25 1/2 0.94 -0.12 -0.48 3/4 0.82 -0.16 -0.66 1 0.65 -0.18 -0.74
Therefore, z(1) 0.65. We can also plot the points (x, y) to get an approximate picture of the curve. The function eulers_method_2x2_plot will do this; in order to use it, we need to dene functions f and g which takes one argument with three coordinates: (t, x, y).
sage: f = lambda z: z[2] # f(t,x,y) = y sage: g = lambda z: -sin(z[1]) # g(t,x,y) = -sin(x) sage: P = eulers_method_2x2_plot(f,g, 0.0, 0.75, 0.0, 0.1, 1.0)
At this point, P is storing two plots: P[0], the plot of x vs. t, and P[1], the plot of y vs. t. We can plot both of these as follows:
sage: show(P[0] + P[1])
At this point, Sage has only wrapped these functions for numerical use. For symbolic use, please use the Maxima interface directly, as in the following example:
sage: maxima.eval("f:bessel_y(v, w)") bessel_y(v,w) sage: maxima.eval("diff(f,w)") (bessel_y(v-1,w)-bessel_y(v+1,w))/2
2.5 Plotting
Sage can produce two-dimensional and three-dimensional plots. 18 Chapter 2. A Guided Tour
You can also create a circle by assigning it to a variable; this does not plot it:
sage: c = circle((0,0), 1, rgbcolor=(1,1,0))
Alternatively, evaluating c.save(filename.png) will save the plot to the given le. Now, these circles look more like ellipses because the axes are scaled differently. You can x this:
sage: c.show(aspect_ratio=1)
The command show(c, aspect_ratio=1) accomplishes the same thing, or you can save the picture using c.save(filename.png, aspect_ratio=1). Its easy to plot basic functions:
sage: plot(cos, (-5,5))
Once you specify a variable name, you can create parametric plots also:
sage: x = var(x) sage: parametric_plot((cos(x),sin(x)^3),(x,0,2*pi),rgbcolor=hue(0.6))
Its important to notice that the axes of the plots will only intersect if the origin is in the viewing range of the graph, and that with sufciently large values scientic notation may be used:
sage: plot(x^2,(x,300,500))
A good way to produce lled-in shapes is to produce a list of points (L in the example below) and then use the polygon command to plot the shape with boundary formed by those points. For example, here is a green deltoid:
sage: ... sage: sage: L = [[-1+cos(pi*i/100)*(1+cos(pi*i/100)),\ 2*sin(pi*i/100)*(1-cos(pi*i/100))] for i in range(200)] p = polygon(L, rgbcolor=(1/8,3/4,1/2)) p
2.5. Plotting
19
Type show(p, axes=false) to see this without any axes. You can add text to a plot:
sage: ... sage: sage: sage: L = [[6*cos(pi*i/100)+5*cos((6/2)*pi*i/100),\ 6*sin(pi*i/100)-5*sin((6/2)*pi*i/100)] for i in range(200)] p = polygon(L, rgbcolor=(1/8,1/4,1/2)) t = text("hypotrochoid", (5,4), rgbcolor=(1,0,0)) show(p+t)
Calculus teachers draw the following plot frequently on the board: not just one branch of arcsin but rather several of them: i.e., the plot of y = sin(x) for x between 2 and 2, ipped about the 45 degree line. The following Sage commands construct this:
sage: v = [(sin(x),x) for x in srange(-2*float(pi),2*float(pi),0.1)] sage: line(v)
Since the tangent function has a larger range than sine, if you use the same trick to plot the inverse tangent, you should change the minimum and maximum coordinates for the x-axis:
sage: v = [(tan(x),x) for x in srange(-2*float(pi),2*float(pi),0.01)] sage: show(line(v), xmin=-20, xmax=20)
Sage also computes polar plots, contour plots and vector eld plots (for special types of functions). Here is an example of a contour plot:
sage: f = lambda x,y: cos(x*y) sage: contour_plot(f, (-4, 4), (-4, 4))
Alternatively, you can use parametric_plot3d to graph a parametric surface where each of x, y, z is determined by a function of one or two variables (the parameters, typically u and v). The previous plot can be expressed parametrically as follows:
sage: sage: sage: sage: sage: u, v = var(u, v) f_x(u, v) = u f_y(u, v) = v f_z(u, v) = u^2 + v^2 parametric_plot3d([f_x, f_y, f_z], (u, -2, 2), (v, -2, 2))
The third way to plot a 3D surface in Sage is implicit_plot3d, which graphs a contour of a function like f (x, y, z) = 0 (this denes a set of points). We graph a sphere using the classical formula:
sage: x, y, z = var(x, y, z) sage: implicit_plot3d(x^2 + y^2 + z^2 - 4, (x,-2, 2), (y,-2, 2), (z,-2, 2))
20
u, v = var(u,v) fx = u*v fy = u fz = v^2 parametric_plot3d([fx, fy, fz], (u, -1, 1), (v, -1, 1), frame=False, color="yellow")
Cross cap:
sage: sage: sage: sage: sage: ... u, v = var(u,v) fx = (1+cos(v))*cos(u) fy = (1+cos(v))*sin(u) fz = -tanh((2/3)*(u-pi))*sin(v) parametric_plot3d([fx, fy, fz], (u, 0, 2*pi), (v, 0, 2*pi), frame=False, color="red")
Twisted torus:
sage: sage: sage: sage: sage: ... u, v = var(u,v) fx = (3+sin(v)+cos(u))*cos(2*v) fy = (3+sin(v)+cos(u))*sin(2*v) fz = sin(u)+2*cos(v) parametric_plot3d([fx, fy, fz], (u, 0, 2*pi), (v, 0, 2*pi), frame=False, color="red")
Lemniscate:
sage: x, y, z = var(x,y,z) sage: f(x, y, z) = 4*x^2 * (x^2 + y^2 + z^2 + z) + y^2 * (y^2 + z^2 - 1) sage: implicit_plot3d(f, (x, -0.5, 0.5), (y, -1, 1), (z, -1, 1))
In the last line, note the syntax. Using plot(f(z), 0, 2) instead will give an error, because z is a dummy variable in the denition of f and is not dened outside of that denition. Indeed, just f(z) returns an error. The following will work in this case, although in general there are issues and so it should probably be avoided (see item 4 below).
sage: var(z) # define z to be a variable z sage: f(z) z^2 sage: plot(f(z), 0, 2)
21
At this point, f(z) is a symbolic expression, the next item in our list. 2. Dene a callable symbolic expression. These can be plotted, differentiated, and integrated.
sage: g(x) = x^2 sage: g # g sends x to x^2 x |--> x^2 sage: g(3) 9 sage: Dg = g.derivative(); Dg x |--> 2*x sage: Dg(3) 6 sage: type(g) <type sage.symbolic.expression.Expression> sage: plot(g, 0, 2)
Note that while g is a callable symbolic expression, g(x) is a related, but different sort of object, which can also be plotted, differentated, etc., albeit with some issues: see item 5 below for an illustration.
sage: x^2 sage: <type sage: 2*x sage: g(x) type(g(x)) sage.symbolic.expression.Expression> g(x).derivative() plot(g(x), 0, 2)
3. Use a pre-dened Sage calculus function. These can be plotted, and with a little help, differentiated, and integrated.
sage: type(sin) <class sage.functions.trig.Function_sin> sage: plot(sin, 0, 2) sage: type(sin(x)) <type sage.symbolic.expression.Expression> sage: plot(sin(x), 0, 2)
Using f = sin(x) instead of sin works, but it is probably even better to use f(x) = sin(x) to dene a callable symbolic expression.
sage: S(x) = sin(x) sage: S.derivative() x |--> cos(x)
22
The issue: plot(h(x), 0, 4) plots the line y = x 2, not the multi-line function dened by h. The reason? In the command plot(h(x), 0, 4), rst h(x) is evaluated: this means plugging x into the function h, which means that x<2 is evaluated.
sage: type(x<2) <type sage.symbolic.expression.Expression>
When a symbolic equation is evaluated, as in the denition of h, if it is not obviously true, then it returns False. Thus h(x) evaluates to x-2, and this is the function that gets plotted. The solution: dont use plot(h(x), 0, 4); instead, use
sage: plot(h, 0, 4)
The problem: g(3), for example, returns an error, saying ValueError: the number of arguments must be less than or equal to 0.
sage: <type sage: <type type(f) sage.symbolic.expression.Expression> type(g) sage.symbolic.expression.Expression>
g is not a function, its a constant, so it has no variables associated to it, and you cant plug anything into it. The solution: there are several options. Dene f initially to be a symbolic expression.
sage: f(x) = x # instead of f = x sage: g = f.derivative() sage: g x |--> 1 sage: g(3) 1 sage: type(g) <type sage.symbolic.expression.Expression>
Or with f and g as dened originally, specify the variable for which you are substituting.
sage: f = x sage: g = f.derivative() sage: g 1
23
sage: g(x=3) 1
# instead of g(3)
Finally, heres one more way to tell the difference between the derivatives of f = x and f(x) = x
sage: sage: sage: () sage: (x,) sage: sage: sage: () sage: () f(x) = x g = f.derivative() g.variables() # the variables present in g g.arguments() # the arguments which can be plugged into g
As this example has been trying to illustrate, h accepts no arguments, and this is why h(3) returns an error.
Similar comments apply to matrices: the row-reduced form of a matrix can depend on the ring over which it is dened, as can its eigenvalues and eigenvectors. For more about constructing polynomials, see Polynomials, and for more about matrices, see Linear Algebra. The symbol I represents the square root of 1; i is a synonym for I. Of course, this is not a rational number:
24
Note: The above code may not work as expected if the variable i has been assigned a different value, for example, if it was used as a loop variable. If this is the case, type
sage: reset(i)
to get the original complex value of i. There is one subtlety in dening complex numbers: as mentioned above, the symbol i represents a square root of 1, but it is a formal or symbolic square root of 1. Calling CC(i) or CC.0 returns the complex square root of 1. Arithmetic involving different kinds of numbers is possible by so-called coercion, see Parents, Conversion and Coercion.
sage: i = CC(i) # floating point complex number sage: i == CC.0 True sage: a, b = 4/3, 2/3 sage: z = a + b*i sage: z 1.33333333333333 + 0.666666666666667*I sage: z.imag() # imaginary part 0.666666666666667 sage: z.real() == a # automatic coercion before comparison True sage: a + b 2 sage: 2*b == a True sage: parent(2/3) Rational Field sage: parent(4/2) Rational Field sage: 2/3 + 0.1 # automatic coercion before addition 0.766666666666667 sage: 0.1 + 2/3 # coercion rules are symmetric in SAGE 0.766666666666667
Here are more examples of basic rings in Sage. As noted above, the ring of rational numbers may be referred to using QQ, or also RationalField() (a eld is a ring in which the multiplication is commutative and in which every nonzero element has a reciprocal in that ring, so the rationals form a eld, but the integers dont):
sage: RationalField() Rational Field sage: QQ Rational Field sage: 1/2 in QQ True
The decimal number 1.2 is considered to be in QQ: decimal numbers which happen to also be rational can be coerced into the rational numbers (see Parents, Conversion and Coercion). The numbers and 2 are not rational, though:
sage: 1.2 in QQ True sage: pi in QQ
25
For use in higher mathematics, Sage also knows about other rings, such as nite elds, p-adic integers, the ring of algebraic numbers, polynomial rings, and matrix rings. Here are constructions of some of these:
sage: GF(3) Finite Field of size 3 sage: GF(27, a) # need to name the generator if not a prime field Finite Field in a of size 3^3 sage: Zp(5) 5-adic Ring with capped relative precision 20 sage: sqrt(3) in QQbar # algebraic closure of QQ True
Note that in Sage, the kernel of a matrix A is the left kernel, i.e. the space of vectors w such that wA = 0. Solving matrix equations is easy, using the method solve_right. Evaluating A.solve_right(Y) returns a matrix (or vector) X so that AX = Y :
sage: Y sage: X sage: X (-2, 1, sage: A (0, -4, = vector([0, -4, -1]) = A.solve_right(Y) 0) * X -1)
26
sage: A.solve_right(w) Traceback (most recent call last): ... ValueError: matrix equation has no solutions
Similarly, use A.solve_left(Y) to solve for X in XA = Y . Sage can also compute eigenvalues and eigenvectors:
sage: A = matrix([[0, 4], [-1, 0]]) sage: A.eigenvalues () [-2*I, 2*I] sage: B = matrix([[1, 3], [3, 1]]) sage: B.eigenvectors_left() [(4, [ (1, 1) ], 1), (-2, [ (1, -1) ], 1)]
(The syntax for the output of eigenvectors_left is a list of triples: (eigenvalue, eigenvector, multiplicity).) Eigenvalues and eigenvectors over QQ or RR can also be computed using Maxima (see Maxima below). As noted in Basic Rings, the ring over which a matrix is dened affects some of its properties. In the following, the rst argument to the matrix command tells Sage to view the matrix as a matrix of integers (the ZZ case), a matrix of rational numbers (QQ), or a matrix of reals (RR):
sage: AZ = matrix(ZZ, [[2,0], [0,1]]) sage: AQ = matrix(QQ, [[2,0], [0,1]]) sage: AR = matrix(RR, [[2,0], [0,1]]) sage: AZ.echelon_form() [2 0] [0 1] sage: AQ.echelon_form() [1 0] [0 1] sage: AR.echelon_form() [ 1.00000000000000 0.000000000000000] [0.000000000000000 1.00000000000000]
(To specify the space of 3 by 4 matrices, you would use MatrixSpace(QQ,3,4). If the number of columns is omitted, it defaults to the number of rows, so MatrixSpace(QQ,3) is a synonym for MatrixSpace(QQ,3,3).) The space of matrices has a basis which Sage stores as a list:
sage: B = M.basis() sage: len(B) 9 sage: B[1] [0 1 0]
27
[0 0 0] [0 0 0]
The basis of S used by Sage is obtained from the non-zero rows of the reduced row echelon form of the matrix of generators of S.
28
The multi-modular algorithm in Sage is good for square matrices (but not so good for non-square matrices):
sage: sage: sage: sage: sage: sage: M A E M A E = = = = = = MatrixSpace(QQ, 50, 100, sparse=True) M.random_element(density = 0.05) A.echelon_form() MatrixSpace(GF(2), 20, 40, sparse=True) M.random_element() A.echelon_form()
2.9 Polynomials
In this section we illustrate how to create and use polynomials in Sage.
This creates a polynomial ring and tells Sage to use (the string) t as the indeterminate when printing to the screen. However, this does not dene the symbol t for use in Sage, so you cannot use it to enter a polynomial (such as t2 + 1) belonging to R. An alternate way is
sage: S = QQ[t] sage: S == R True
This has the same issue regarding t. A third very convenient way is
sage: R.<t> = PolynomialRing(QQ)
or
sage: R.<t> = QQ[t]
or even
2.9. Polynomials
29
This has the additional side effect that it denes the variable t to be the indeterminate of the polynomial ring, so you can easily construct elements of R, as follows. (Note that the third way is very similar to the constructor notation in Magma, and just as in Magma it can be used for a wide range of objects.)
sage: poly = (t+1) * (t+2); poly t^2 + 3*t + 2 sage: poly in R True
Whatever method you use to dene a polynomial ring, you can recover the indeterminate as the 0th generator:
sage: R = PolynomialRing(QQ, t) sage: t = R.0 sage: t in R True
Note that a similar construction works with the complex numbers: the complex numbers can be viewed as being generated over the real numbers by the symbol i; thus we have the following:
sage: CC Complex Field with 53 bits of precision sage: CC.0 # 0th generator of CC 1.00000000000000*I
For polynomial rings, you can obtain both the ring and its generator, or just the generator, during ring creation as follows:
sage: sage: sage: sage: R, t = QQ[t].objgen() t = QQ[t].gen() R, t = objgen(QQ[t]) t = gen(QQ[t])
Notice that the factorization correctly takes into account and records the unit part. If you were to use, e.g., the R.cyclotomic_polynomial function a lot for some research project, in addition to citing Sage you should make an attempt to nd out what component of Sage is being used to actually compute the cyclotomic polynomial and cite that as well. In this case, if you type R.cyclotomic_polynomial?? to see the source code, youll quickly see a line f = pari.polcyclo(n) which means that PARI is being used for computation of the cyclotomic polynomial. Cite PARI in your work as well. 30 Chapter 2. A Guided Tour
Dividing two polynomials constructs an element of the fraction eld (which Sage creates automatically).
sage: x = QQ[x].0 sage: f = x^3 + 1; g = x^2 - 17 sage: h = f/g; h (x^3 + 1)/(x^2 - 17) sage: h.parent() Fraction Field of Univariate Polynomial Ring in x over Rational Field
Using Laurent series, one can compute series expansions in the fraction eld of QQ[x]:
sage: R.<x> = LaurentSeriesRing(QQ); R Laurent Series Ring in x over Rational Field sage: 1/(1-x) + O(x^10) 1 + x + x^2 + x^3 + x^4 + x^5 + x^6 + x^7 + x^8 + x^9 + O(x^10)
The ring is determined by the variable. Note that making another ring with variable called x does not return a different ring.
sage: sage: sage: True sage: True sage: True R = PolynomialRing(QQ, "x") T = PolynomialRing(QQ, "x") R == T R is T R.0 == T.0
Sage also has support for power series and Laurent series rings over any base ring. In the following example, we create an element of F7 [[T ]] and divide to create an element of F7 ((T )).
sage: R.<T> = PowerSeriesRing(GF(7)); R Power Series Ring in T over Finite Field of size 7 sage: f = T + 3*T^2 + T^3 + O(T^4) sage: f^3 T^3 + 2*T^4 + 2*T^5 + O(T^6) sage: 1/f T^-1 + 4 + T + O(T^2) sage: parent(1/f) Laurent Series Ring in T over Finite Field of size 7
You can also create power series rings using a double-brackets shorthand:
sage: GF(7)[[T]] Power Series Ring in T over Finite Field of size 7
2.9. Polynomials
31
Just as for dening univariate polynomial rings, there are alternative ways:
sage: GF(5)[z0, z1, z2] Multivariate Polynomial Ring in z0, z1, z2 over Finite Field of size 5 sage: R.<z0,z1,z2> = GF(5)[]; R Multivariate Polynomial Ring in z0, z1, z2 over Finite Field of size 5
Also, if you want the variable names to be single letters then you can use the following shorthand:
sage: PolynomialRing(GF(5), 3, xyz) Multivariate Polynomial Ring in x, y, z over Finite Field of size 5
You can also use more mathematical notation to construct a polynomial ring.
sage: R = GF(5)[x,y,z] sage: x,y,z = R.gens() sage: QQ[x] Univariate Polynomial Ring in x over Rational Field sage: QQ[x,y].gens() (x, y) sage: QQ[x].objgens() (Univariate Polynomial Ring in x over Rational Field, (x,))
Multivariate polynomials are implemented in Sage using Python dictionaries and the distributive representation of a polynomial. Sage makes some use of Singular [Si], e.g., for computation of gcds and Grbner basis of ideals.
sage: sage: sage: sage: x^2 R, (x, y) = PolynomialRing(RationalField(), 2, xy).objgens() f = (x^3 + 2*y^2*x)^2 g = x^2*y^2 f.gcd(g)
Next we create the ideal (f, g) generated by f and g, by simply multiplying (f,g) by R (we could also write ideal([f,g]) or ideal(f,g)).
sage: I = (f, g)*R; I Ideal (x^6 + 4*x^4*y^2 + 4*x^2*y^4, x^2*y^2) of Multivariate Polynomial Ring in x, y over Rational Field sage: B = I.groebner_basis(); B [x^6, x^2*y^2] sage: x^2 in I False
Incidentally, the Grbner basis above is not a list but an immutable sequence. This means that it has a universe, parent, and cannot be changed (which is good because changing the basis would break other routines that use the Grbner 32 Chapter 2. A Guided Tour
basis).
sage: B.parent() Category of sequences in Multivariate Polynomial Ring in x, y over Rational Field sage: B.universe() Multivariate Polynomial Ring in x, y over Rational Field sage: B[1] = x Traceback (most recent call last): ... ValueError: object is immutable; please change a copy instead.
Some (read: not as much as we would like) commutative algebra is available in Sage, implemented via Singular. For example, we can compute the primary decomposition and associated primes of I:
sage: I.primary_decomposition() [Ideal (x^2) of Multivariate Polynomial Ring in x, y over Rational Field, Ideal (y^2, x^6) of Multivariate Polynomial Ring in x, y over Rational Field] sage: I.associated_primes() [Ideal (x) of Multivariate Polynomial Ring in x, y over Rational Field, Ideal (y, x) of Multivariate Polynomial Ring in x, y over Rational Field]
2.10.1 Elements
If one wants to implement a ring in Python, a rst approximation is to create a class for the elements X of that ring and provide it with the required double underscore methods such as __add__, __sub__, __mul__, of course making sure that the ring axioms hold. As Python is a strongly typed (yet dynamically typed) language, one might, at least at rst, expect that one implements one Python class for each ring. After all, Python contains one type <int> for the integers, one type <float> for the reals, and so on. But that approach must soon fail: There are innitely many rings, and one can not implement innitely many classes. Instead, one may create a hierarchy of classes designed to implement elements of ubiquitous algebraic structures, such as groups, rings, skew elds, commutative rings, elds, algebras, and so on. But that means that elements of fairly different rings can have the same type.
sage: P.<x,y> = GF(3)[] sage: Q.<a,b> = GF(4,z)[] sage: type(x)==type(a) True
On the other hand, one could also have different Python classes providing different implementations of the same mathematical structure (e.g., dense matrices versus sparse matrices)
sage: P.<a> = PolynomialRing(ZZ) sage: Q.<b> = PolynomialRing(ZZ, sparse=True) sage: R.<c> = PolynomialRing(ZZ, implementation=NTL)
33
sage: type(a); type(b); type(c) <type sage.rings.polynomial.polynomial_integer_dense_flint.Polynomial_integer_dense_flint> <class sage.rings.polynomial.polynomial_element_generic.Polynomial_generic_sparse> <type sage.rings.polynomial.polynomial_integer_dense_ntl.Polynomial_integer_dense_ntl>
That poses two problems: On the one hand, if one has elements that are two instances of the same class, then one may expect that their __add__ method will allow to add them; but one does not want that, if the elements belong to very different rings. On the other hand, if one has elements belonging to different implementations of the same ring, then one wants to add them, but that is not straight forward if they belong to different Python classes. The solution to these problems is called coercion and will be explained below. However, it is essential that each element knows what it is element of. That is available by the method parent():
sage: a.parent(); b.parent(); c.parent() Univariate Polynomial Ring in a over Integer Ring Sparse Univariate Polynomial Ring in b over Integer Ring Univariate Polynomial Ring in c over Integer Ring (using NTL)
In algebra, objects sharing the same kind of algebraic structures are collected in so-called categories. So, there is a rough analogy between the class hierarchy in Sage and the hierarchy of categories. However, this analogy of Python classes and categories shouldnt be stressed too much. After all, mathematical categories are implemented in Sage as well:
sage: Rings() Category of rings sage: ZZ.category() Category of euclidean domains sage: ZZ.category().is_subcategory(Rings()) True sage: ZZ in Rings() True sage: ZZ in Fields() False sage: QQ in Fields() True
While Sages class hierarchy is centered at implementation details, Sages category framework is more centered on mathematical structure. It is possible to implement generic methods and tests independent of a specic implementation in the categories.
34
Parent structures in Sage are supposed to be unique Python objects. For example, once a polynomial ring over a certain base ring and with a certain list of generators is created, the result is cached:
sage: RR[x,y] is RR[x,y] True
While parents are unique, equal elements of a parent in Sage are not necessarily identical. This is in contrast to the behaviour of Python for some (albeit not all) integers:
sage: int(1) is int(1) # Python int True sage: int(-15) is int(-15) False sage: 1 is 1 # Sage Integer False
It is important to observe that elements of different rings are in general not distinguished by their type, but by their parent:
sage: a = GF(2)(1) sage: b = GF(5)(1) sage: type(a) is type(b) True sage: parent(a) Finite Field of size 2 sage: parent(b) Finite Field of size 5
Hence, from an algebraic point of view, the parent of an element is more important than its type.
35
Or If an element r1 of one ring R1 can somehow be interpreted in another ring R2, then all arithmetic operations involving r1 and any element of R2 are allowed. The multiplicative unit exists in all elds and many rings, and they should all be equal. Sage favours a compromise. If P1 and P2 are parent structures and p1 is an element of P1, then the user may explicitly ask for an interpretation of p1 in P2. This may not be meaningful in all cases or not be dened for all elements of P1, and it is up to the user to ensure that it makes sense. We refer to this as conversion:
sage: sage: sage: True sage: True a = GF(2)(1) b = GF(5)(1) GF(5)(a) == b GF(2)(b) == a
However, an implicit (or automatic) conversion will only happen if this can be done thoroughly and consistently. Mathematical rigour is essential at that point. Such an implicit conversion is called coercion. If coercion is dened, then it must coincide with conversion. Two conditions must be satised for a coercion to be dened: 1. A coercion from P1 to P2 must be given by a structure preserving map (e.g., a ring homomorphism). It does not sufce that some elements of P1 can be mapped to P2, and the map must respect the algebraic structure of P1. 2. The choice of these coercion maps must be consistent: If P3 is a third parent structure, then the composition of the chosen coercion from P1 to P2 with the coercion from P2 to P3 must coincide with the chosen coercion from P1 to P3. In particular, if there is a coercion from P1 to P2 and P2 to P1, the composition must be the identity map of P1. So, although it is possible to convert each element of GF(2) into GF(5), there is no coercion, since there is no ring homomorphism between GF(2) and GF(5). The second aspect - consistency - is a bit more difcult to explain. We illustrate it with multivariate polynomial rings. In applications, it certainly makes most sense to have name preserving coercions. So, we have:
sage: sage: sage: True sage: x sage: y R1.<x,y> = ZZ[] R2 = ZZ[y,x] R2.has_coerce_map_from(R1) R2(x) R2(y)
If there is no name preserving ring homomorphism, coercion is not dened. However, conversion may still be possible, namely by mapping ring generators according to their position in the list of generators:
sage: sage: False sage: z sage: x R3 = ZZ[z,x] R3.has_coerce_map_from(R1) R3(x) R3(y)
But such position preserving conversions do not qualify as coercion: By composing a name preserving map from ZZ[x,y] to ZZ[y,x] with a position preserving map from ZZ[y,x] to ZZ[a,b], a map would result that is neither name preserving nor position preserving, in violation to consistency. 36 Chapter 2. A Guided Tour
If there is a coercion, it will be used to compare elements of different rings or to do arithmetic. This is often convenient, but the user should be aware that extending the ==-relation across the borders of different parents may easily result in overdoing it. For example, while == is supposed to be an equivalence relation on the elements of one ring, this is not necessarily the case if different rings are involved. For example, 1 in ZZ and in a nite eld are considered equal, since there is a canonical coercion from the integers to any nite eld. However, in general there is no coercion between two different nite elds. Therefore we have
sage: True sage: True sage: False sage: True GF(5)(1) == 1 1 == GF(2)(1) GF(5)(1) == GF(2)(1) GF(5)(1) != GF(2)(1)
Similarly, we have
sage: R3(R1.1) == R3.1 True sage: R1.1 == R3.1 False sage: R1.1 != R3.1 True
Another consequence of the consistency condition is that coercions can only go from exact rings (e.g., the rationals QQ) to inexact rings (e.g., real numbers with a xed precision RR), but not the other way around. The reason is that the composition of the coercion from QQ to RR with a conversion from RR to QQ is supposed to be the identity on QQ. But this is impossible, since some distinct rational numbers may very well be treated equal in RR, as in the following example:
sage: RR(1/10^200+1/10^100) == RR(1/10^100) True sage: 1/10^200+1/10^100 == 1/10^100 False
When comparing elements of two parents P1 and P2, it is possible that there is no coercion between the two rings, but there is a canonical choice of a parent P3 so that both P1 and P2 coerce into P3. In this case, coercion will take place as well. A typical use case is the sum of a rational number and a polynomial with integer coefcients, yielding a polynomial with rational coefcients:
sage: P1.<x> = ZZ[] sage: p = 2*x+3 sage: q = 1/2 sage: parent(p) Univariate Polynomial Ring in x over Integer Ring sage: parent(p+q) Univariate Polynomial Ring in x over Rational Field
Note that in principle the result would also make sense in the fraction eld of ZZ[x]. However, Sage tries to choose a canonical common parent that seems to be most natural (QQ[x] in our example). If several potential common parents seem equally natural, Sage will not pick one of them at random, in order to have a reliable result. The mechanisms which that choice is based upon is explained in the tutorial worksheet No coercion into a common parent will take place in the following example:
sage: R.<x> = QQ[] sage: S.<y> = QQ[] sage: x+y
37
Traceback (most recent call last): ... TypeError: unsupported operand parent(s) for +: Univariate Polynomial Ring in x over Rational Fiel
The reason is that Sage would not choose one of the potential candidates QQ[x][y], QQ[y][x], QQ[x,y] or QQ[y,x], because all of these four pairwise different structures seem natural common parents, and there is no apparent canonical choice.
You can also obtain the character table (in LaTeX format) in Sage:
sage: G = PermutationGroup([[(1,2),(3,4)], [(1,2,3)]]) sage: latex(G.character_table()) \left(\begin{array}{rrrr} 1 & 1 & 1 & 1 \\ 1 & 1 & -\zeta_{3} - 1 & \zeta_{3} \\ 1 & 1 & \zeta_{3} & -\zeta_{3} - 1 \\ 3 & -1 & 0 & 0 \end{array}\right)
Sage also includes classical and matrix groups over nite elds:
sage: MS = MatrixSpace(GF(7), 2) sage: gens = [MS([[1,0],[-1,1]]),MS([[1,1],[0,1]])] sage: G = MatrixGroup(gens) sage: G.conjugacy_class_representatives() [ [1 0] [0 1], [0 1] [6 1], ... [6 0]
38
[0 6] ] sage: G = Sp(4,GF(7)) sage: G._gap_init_() Sp(4, 7) sage: G Symplectic Group of rank 2 over Finite Field of size 7 sage: G.random_element() # random output [5 5 5 1] [0 2 6 3] [5 0 1 0] [4 6 3 4] sage: G.order() 276595200
You can also compute using abelian groups (innite and nite):
sage: F = AbelianGroup(5, [5,5,7,8,9], names=abcde) sage: (a, b, c, d, e) = F.gens() sage: d * b**2 * c**3 b^2*c^3*d sage: F = AbelianGroup(3,[2]*3); F Multiplicative Abelian Group isomorphic to C2 x C2 x C2 sage: H = AbelianGroup([2,3], names="xy"); H Multiplicative Abelian Group isomorphic to C2 x C3 sage: AbelianGroup(5) Multiplicative Abelian Group isomorphic to Z x Z x Z x Z x Z sage: AbelianGroup(5).order() +Infinity
39
[22, 10, 6, 3, 2, 1, 1, 1] sage: next_prime(2005) 2011 sage: previous_prime(2005) 2003 sage: divisors(28); sum(divisors(28)); 2*28 [1, 2, 4, 7, 14, 28] 56 56
We next illustrate the extended Euclidean algorithm, Eulers -function, and the Chinese remainder theorem:
sage: d,u,v = xgcd(12,15) sage: d == u*12 + v*15 True sage: n = 2005 sage: inverse_mod(3,n) 1337 sage: 3 * 1337 4011 sage: prime_divisors(n) [5, 401] sage: phi = n*prod([1 - 1/p for p in prime_divisors(n)]); phi 1600 sage: euler_phi(n) 1600 sage: prime_to_m_part(n, 5) 401
40
sage: [kronecker(m,13) for m in range(1,13)] [1, -1, 1, 1, -1, -1, -1, -1, 1, 1, -1, 1] sage: n = 10000; sum([moebius(m) for m in range(1,n)]) -23 sage: Partitions(4).list() [[4], [3, 1], [2, 2], [2, 1, 1], [1, 1, 1, 1]]
Much work has been done implementing rings of integers in p-adic elds or number elds other than . The interested reader is invited to ask the experts on the sage-support Google group for further details. A number of related methods are already implemented in the NumberField class.
sage: R.<x> = PolynomialRing(QQ) sage: K = NumberField(x^3 + x^2 - 2*x + 8, a) sage: K.integral_basis() [1, 1/2*a^2 + 1/2*a, a^2] sage: K.galois_group(type="pari") Galois group PARI group [6, -1, 2, "S3"] of degree 3 of the Number Field in a with defining polynomial x^3 + x^2 - 2*x + 8 sage: K.polynomial_quotient_ring() Univariate Quotient Polynomial Ring in a over Rational Field with modulus x^3 + x^2 - 2*x + 8 sage: K.units() [3*a^2 + 13*a + 13] sage: K.discriminant() -503 sage: K.class_group() Class group of order 1 of Number Field in a with defining polynomial x^3 + x^2 - 2*x + 8 sage: K.class_number() 1
sage: x, y = AffineSpace(2, QQ, xy).gens() sage: C2 = Curve(x^2 + y^2 - 1) sage: C3 = Curve(x^3 + y^3 - 1) sage: D = C2 + C3 sage: D Affine Curve over Rational Field defined by x^5 + x^3*y^2 + x^2*y^3 + y^5 - x^3 - y^3 - x^2 - y^2 + 1 sage: D.irreducible_components() [ Closed subscheme of Affine Space of dimension 2 over Rational Field defined by: x^2 + y^2 - 1, Closed subscheme of Affine Space of dimension 2 over Rational Field defined by: x^3 + y^3 - 1 ]
We can also nd all points of intersection of the two curves by intersecting them and computing the irreducible components.
sage: V = C2.intersection(C3) sage: V.irreducible_components() [ Closed subscheme of Affine Space of dimension 2 over Rational Field defined by: y - 1, x, Closed subscheme of Affine Space of dimension 2 over Rational Field defined by: y, x - 1, Closed subscheme of Affine Space of dimension 2 over Rational Field defined by: x + y + 2, 2*y^2 + 4*y + 3 ]
Thus, e.g., (1, 0) and (0, 1) are on both curves (visibly clear), as are certain (quadratic) points whose y coordinates satisfy 2y 2 + 4y + 3 = 0. Sage can compute the toric ideal of the twisted cubic in projective 3 space:
sage: R.<a,b,c,d> = PolynomialRing(QQ, 4) sage: I = ideal(b^2-a*c, c^2-b*d, a*d-b*c) sage: F = I.groebner_fan(); F Groebner fan of the ideal: Ideal (b^2 - a*c, c^2 - b*d, -b*c + a*d) of Multivariate Polynomial Ring in a, b, c, d over Rational Field sage: F.reduced_groebner_bases () [[-c^2 + b*d, -b*c + a*d, -b^2 + a*c], [-c^2 + b*d, b^2 - a*c, -b*c + a*d], [-c^2 + b*d, b*c - a*d, b^2 - a*c, -c^3 + a*d^2], [c^3 - a*d^2, -c^2 + b*d, b*c - a*d, b^2 - a*c], [c^2 - b*d, -b*c + a*d, -b^2 + a*c], [c^2 - b*d, b*c - a*d, -b^2 + a*c, -b^3 + a^2*d], [c^2 - b*d, b*c - a*d, b^3 - a^2*d, -b^2 + a*c], [c^2 - b*d, b*c - a*d, b^2 - a*c]] sage: F.polyhedralfan() Polyhedral fan in 4 dimensions of dimension 4
42
where the ai s are coerced into the parent of a1 . If all the ai have parent Z, they are coerced into Q. EllipticCurve([a4 , a6 ]): Same as above, but a1 = a2 = a3 = 0. EllipticCurve(label): Returns the elliptic curve over from the Cremona database with the given (new!) Cremona label. The label is a string, such as "11a" or "37b2". The letter must be lower case (to distinguish it from the old labeling). EllipticCurve(j): Returns an elliptic curve with j-invariant j. EllipticCurve(R, [a1 , a2 , a3 , a4 , a6 ]): Create the elliptic curve over a ring R with given ai s as above. We illustrate each of these constructors:
sage: EllipticCurve([0,0,1,-1,0]) Elliptic Curve defined by y^2 + y = x^3 - x over Rational Field sage: EllipticCurve([GF(5)(0),0,1,-1,0]) Elliptic Curve defined by y^2 + y = x^3 + 4*x over Finite Field of size 5 sage: EllipticCurve([1,2]) Elliptic Curve defined by y^2
sage: EllipticCurve(37a) Elliptic Curve defined by y^2 + y = x^3 - x over Rational Field sage: EllipticCurve_from_j(1) Elliptic Curve defined by y^2 + x*y = x^3 + 36*x + 3455 over Rational Field sage: EllipticCurve(GF(5), [0,0,1,-1,0]) Elliptic Curve defined by y^2 + y = x^3 + 4*x over Finite Field of size 5
The pair (0, 0) is a point on the elliptic curve E dened by y 2 + y = x3 x. To create this point in Sage type E([0,0]). Sage can add points on such an elliptic curve (recall elliptic curves support an additive group structure where the point at innity is the zero element and three co-linear points on the curve add to zero):
sage: E = EllipticCurve([0,0,1,-1,0]) sage: E Elliptic Curve defined by y^2 + y = x^3 - x over Rational Field sage: P = E([0,0]) sage: P + P (1 : 0 : 1) sage: 10*P (161/16 : -2065/64 : 1) sage: 20*P (683916417/264517696 : -18784454671297/4302115807744 : 1)
43
sage: E.conductor() 37
The elliptic curves over the complex numbers are parameterized by the j-invariant. Sage computes j-invariant as follows:
sage: E = EllipticCurve([0,0,0,-4,2]); E Elliptic Curve defined by y^2 = x^3 - 4*x + 2 over Rational Field sage: E.conductor() 2368 sage: E.j_invariant() 110592/37
If we make a curve with the same j-invariant as that of E, it need not be isomorphic to E. In the following example, the curves are not isomorphic because their conductors are different.
sage: F = EllipticCurve_from_j(110592/37) sage: F.conductor() 37
We can compute the coefcients an of the L-series or modular form computation uses the PARI C-library:
n=0
sage: E = EllipticCurve([0,0,1,-1,0]) sage: print E.anlist(30) [0, 1, -2, -3, 2, -2, 6, -1, 0, 6, 4, -5, -6, -2, 2, 6, -4, 0, -12, 0, -4, 3, 10, 2, 0, -1, 4, -9, -2, 6, -12] sage: v = E.anlist(10000)
Elliptic curves can be constructed using their Cremona labels. This pre-loads the elliptic curve with information about its rank, Tamagawa numbers, regulator, etc.
sage: E = EllipticCurve("37b2") sage: E Elliptic Curve defined by y^2 + y = x^3 + x^2 - 1873*x - 31833 over Rational Field sage: E = EllipticCurve("389a") sage: E Elliptic Curve defined by y^2 + y = x^3 + x^2 - 2*x over Rational Field sage: E.rank() 2 sage: E = EllipticCurve("5077a") sage: E.rank() 3
We can also access the Cremona database directly. 44 Chapter 2. A Guided Tour
sage: db = sage.databases.cremona.CremonaDatabase() sage: db.curves(37) {a1: [[0, 0, 1, -1, 0], 1, 1], b1: [[0, 1, 1, -23, -50], 0, 3]} sage: db.allcurves(37) {a1: [[0, 0, 1, -1, 0], 1, 1], b1: [[0, 1, 1, -23, -50], 0, 3], b2: [[0, 1, 1, -1873, -31833], 0, 1], b3: [[0, 1, 1, -3, 1], 0, 3]}
The objects returned from the database are not of type EllipticCurve. They are elements of a database and have a couple of elds, and thats it. There is a small version of Cremonas database, which is distributed by default with Sage, and contains limited information about elliptic curves of conductor 10000. There is also a large optional version, which contains extensive data about all curves of conductor up to 120000 (as of October 2005). There is also a huge (2GB) optional database package for Sage that contains the hundreds of millions of elliptic curves in the Stein-Watkins database.
Having created the group, we next create an element and compute with it.
sage: G = DirichletGroup(21) sage: chi = G.1; chi Dirichlet character modulo 21 of conductor 7 mapping 8 |--> 1, 10 |--> zeta6 sage: chi.values() [0, 1, zeta6 - 1, 0, -zeta6, -zeta6 + 1, 0, 0, 1, 0, zeta6, -zeta6, 0, -1, 0, 0, zeta6 - 1, zeta6, 0, -zeta6 + 1, -1] sage: chi.conductor() 7 sage: chi.modulus() 21 sage: chi.order() 6 sage: chi(19) -zeta6 + 1 sage: chi(40) -zeta6 + 1
It is also possible to compute the action of the Galois group Gal(Q(N )/Q) on these characters, as well as the direct product decomposition corresponding to the factorization of the modulus.
sage: chi.galois_orbit() [Dirichlet character modulo 21 of conductor 7 mapping 8 |--> 1, 10 |--> zeta6,
45
Dirichlet character modulo 21 of conductor 7 mapping 8 |--> 1, 10 |--> -zeta6 + 1] sage: go = G.galois_orbits() sage: [len(orbit) for orbit in go] [1, 2, 2, 1, 1, 2, 2, 1] sage: [ Group 6 and Group 6 and ] G.decomposition() of Dirichlet characters of modulus 3 over Cyclotomic Field of order degree 2, of Dirichlet characters of modulus 7 over Cyclotomic Field of order degree 2
Next, we construct the group of Dirichlet characters mod 20, but with values in Q(i):
sage: sage: sage: Group K.<i> = NumberField(x^2+1) G = DirichletGroup(20,K) G of Dirichlet characters of modulus 20 over Number Field in i with defining polynomial x^2 + 1
In this example we create a Dirichlet character with values in a number eld. We explicitly specify the choice of root of unity by the third argument to DirichletGroup below.
sage: x = polygen(QQ, x) sage: K = NumberField(x^4 + 1, a); a = K.0 sage: b = K.gen(); a == b True sage: K Number Field in a with defining polynomial x^4 + 1 sage: G = DirichletGroup(5, K, a); G Group of Dirichlet characters of modulus 5 over Number Field in a with defining polynomial x^4 + 1 sage: chi = G.0; chi Dirichlet character modulo 5 of conductor 5 mapping 2 |--> a^2 sage: [(chi^i)(2) for i in range(4)] [1, a^2, -1, -a^2]
Here NumberField(x^4 + 1, a) tells Sage to use the symbol a in printing what K is (a Number Field in a with dening polynomial x4 + 1). The name a is undeclared at this point. Once a = K.0 (or equivalently a = K.gen()) is evaluated, the symbol a represents a root of the generating polynomial x4 + 1.
46
Next we illustrate computation of Hecke operators on a space of modular symbols of level 1 and weight 12.
sage: M = ModularSymbols(1,12) sage: M.basis() ([X^8*Y^2,(0,0)], [X^9*Y,(0,0)], [X^10,(0,0)]) sage: t2 = M.T(2) sage: t2 Hecke operator T_2 on Modular Symbols space of dimension 3 for Gamma_0(1) of weight 12 with sign 0 over Rational Field sage: t2.matrix() [ -24 0 0] [ 0 -24 0] [4860 0 2049] sage: f = t2.charpoly(x); f x^3 - 2001*x^2 - 97776*x - 1180224 sage: factor(f) (x - 2049) * (x + 24)^2 sage: M.T(11).charpoly(x).factor() (x - 285311670612) * (x - 534612)^2
47
Here is another example of how Sage can compute the action of Hecke operators on a space of modular forms.
sage: T = ModularForms(Gamma0(11),2) sage: T Modular Forms space of dimension 2 for Congruence Subgroup Gamma0(11) of weight 2 over Rational Field sage: T.degree() 2 sage: T.level() 11 sage: T.group() Congruence Subgroup Gamma0(11) sage: T.dimension() 2 sage: T.cuspidal_subspace() Cuspidal subspace of dimension 1 of Modular Forms space of dimension 2 for Congruence Subgroup Gamma0(11) of weight 2 over Rational Field sage: T.eisenstein_subspace() Eisenstein subspace of dimension 1 of Modular Forms space of dimension 2 for Congruence Subgroup Gamma0(11) of weight 2 over Rational Field sage: M = ModularSymbols(11); M Modular Symbols space of dimension 3 for Gamma_0(11) of weight 2 with sign 0 over Rational Field sage: M.weight() 2 sage: M.basis() ((1,0), (1,8), (1,9)) sage: M.sign() 0
Let Tp denote the usual Hecke operators (p prime). How do the Hecke operators T2 , T3 , T5 act on the space of modular symbols?
sage: M.T(2).matrix() [ 3 0 -1] [ 0 -2 0] [ 0 0 -2] sage: M.T(3).matrix() [ 4 0 -1]
48
49
50
CHAPTER
THREE
sage:
The wall time is the time that elapsed on the clock hanging from your wall. This is relevant, since CPU time does not track time used by subprocesses like GAP or Singular. (Avoid killing a Sage process with kill -9 from a terminal, since Sage might not kill child processes, e.g., Maple processes, or cleanup temporary les from $HOME/.sage/tmp.)
Here is an example:
sage: factor(100) _1 = 2^2 * 5^2 sage: kronecker_symbol(3,5)
51
_2 = -1 sage: %hist #This only works from the interactive shell, not the notebook. 1: factor(100) 2: kronecker_symbol(3,5) 3: %hist sage: _oh _4 = {1: 2^2 * 5^2, 2: -1} sage: _i1 _5 = factor(ZZ(100))\n sage: eval(_i1) _6 = 2^2 * 5^2 sage: %hist 1: factor(100) 2: kronecker_symbol(3,5) 3: %hist 4: _oh 5: _i1 6: eval(_i1) 7: %hist
We omit the output numbering in the rest of this tutorial and the other Sage documentation. You can also store a list of input from session in a macro for that session.
sage: E = EllipticCurve([1,2,3,4,5]) sage: M = ModularSymbols(37) sage: %hist 1: E = EllipticCurve([1,2,3,4,5]) 2: M = ModularSymbols(37) 3: %hist sage: %macro em 1-2 Macro em created. To execute, type its name (without quotes). sage: E Elliptic Curve defined by y^2 + x*y + 3*y = x^3 + 2*x^2 + 4*x + 5 over Rational Field sage: E = 5 sage: M = None sage: em Executing Macro... sage: E Elliptic Curve defined by y^2 + x*y + 3*y = x^3 + 2*x^2 + 4*x + 5 over Rational Field
When using the interactive shell, any UNIX shell command can be executed from Sage by prefacing it by an exclamation point !. For example,
sage: !ls auto example.sage glossary.tex t tmp tut.log tut.tex
returns the listing of the current directory. The PATH has the Sage bin directory at the front, so if you run gp, gap, singular, maxima, etc., you get the versions included with Sage.
sage: !gp Reading GPRC: /etc/gprc ...Done. GP/PARI CALCULATOR Version 2.2.11 (alpha) i686 running linux (ix86/GMP-4.1.4 kernel) 32-bit version
52
... sage: !singular SINGULAR A Computer Algebra System for Polynomial Computations 0< by: G.-M. Greuel, G. Pfister, H. Schoenemann FB Mathematik der Universitaet, D-67653 Kaiserslautern \ \ October 2005 / / Development version 3-0-1
If you use Sage in the Linux KDE terminal konsole then you can save your session as follows: after starting Sage in konsole, select settings, then history..., then set unlimited. When you are ready to save your session, select edit then save history as... and type in a name to save the text of your session to your computer. After saving this le, you could then load it into an editor, such as xemacs, and print it.
53
This means that 0.66 seconds total were taken, and the Wall time, i.e., the amount of time that elapsed on your wall clock, is also 0.66 seconds. If your computer is heavily loaded with other programs, the wall time may be much larger than the CPU time. Next we time exponentiation using the native Sage Integer type, which is implemented (in Cython) using the GMP library:
sage: %time a = 1938^99484 CPU times: user 0.04 s, sys: 0.00 s, total: 0.04 s Wall time: 0.04
GMP is better, but only slightly (as expected, since the version of PARI built for Sage uses GMP for integer arithmetic). You can also time a block of commands using the cputime command, as illustrated below:
sage: sage: sage: sage: sage: 0.64 t = cputime() a = int(1938)^int(99484) b = 1938^99484 c = pari(1938)^pari(99484) cputime(t)
sage: cputime? ... Return the time in CPU second since SAGE started, or with optional argument t, return the time since time t.
54
INPUT: t -- (optional) float, time in CPU seconds OUTPUT: float -- time in CPU seconds
The walltime command behaves just like the cputime command, except that it measures wall time. We can also compute the above power in some of the computer algebra systems that Sage includes. In each case we execute a trivial command in the system, in order to start up the server for that program. The most relevant time is the wall time. However, if there is a signicant difference between the wall time and the CPU time then this may indicate a performance issue worth looking into.
sage: time 1938^99484; CPU times: user 0.01 s, sys: 0.00 s, total: Wall time: 0.01 sage: gp(0) 0 sage: time g = gp(1938^99484) CPU times: user 0.00 s, sys: 0.00 s, total: Wall time: 0.04 sage: maxima(0) 0 sage: time g = maxima(1938^99484) CPU times: user 0.00 s, sys: 0.00 s, total: Wall time: 0.30 sage: kash(0) 0 sage: time g = kash(1938^99484) CPU times: user 0.00 s, sys: 0.00 s, total: Wall time: 0.04 sage: mathematica(0) 0 sage: time g = mathematica(1938^99484) CPU times: user 0.00 s, sys: 0.00 s, total: Wall time: 0.03 sage: maple(0) 0 sage: time g = maple(1938^99484) CPU times: user 0.00 s, sys: 0.00 s, total: Wall time: 0.11 sage: gap(0) 0 sage: time g = gap.eval(1938^99484;;) CPU times: user 0.00 s, sys: 0.00 s, total: Wall time: 1.02 0.01 s
0.00 s
0.00 s
0.00 s
0.00 s
0.00 s
0.00 s
Note that GAP and Maxima are the slowest in this test (this was run on the machine sage.math.washington.edu). Because of the pexpect interface overhead, it is perhaps unfair to compare these to Sage, which is the fastest.
(The comments not tested are here because the %bg syntax doesnt work well with Sages automatic testing facility. If you type this in yourself, it should work as written. This is of course most useful with commands which take a while to complete.)
sage: def quick(m): return 2*m sage: %bg quick(20) # not tested Starting job # 0 in a separate thread. sage: jobs.status() # not tested Completed jobs: 0 : quick(20) sage: jobs[0].result # the actual answer, not tested 40
Note that jobs run in the background dont use the Sage preparser see The Pre-Parser: Differences between Sage and Python for more information. One (perhaps awkward) way to get around this would be to run
sage: %bg eval(preparse(quick(20))) # not tested
It is safer and easier, though, to just use %bg on commands which dont require the preparser. You can use %edit (or %ed or ed) to open an editor, if you want to type in some complex code. Before you start Sage, make sure that the EDITOR environment variable is set to your favorite editor (by putting export EDITOR=/usr/bin/emacs or export EDITOR=/usr/bin/vim or something similar in the appropriate place, like a .profile le). From the Sage prompt, executing %edit will open up the named editor. Then within the editor you can dene a function:
def some_function(n): return n**2 + 3*n + 2
Save and quit from the editor. For the rest of your Sage session, you can then use some_function. If you want to modify it, type %edit some_function from the Sage prompt. If you have a computation and you want to modify its output for another use, perform the computation and type %rep: this will place the output from the previous command at the Sage prompt, ready for you to edit it.
sage: f(x) = cos(x) sage: f(x).derivative(x) -sin(x)
At this point, if you type %rep at the Sage prompt, you will get a new Sage prompt, followed by -sin(x), with the cursor at the end of the line. For more, type %quickref to get a quick reference guide to IPython. As of this writing (April 2011), Sage uses version 0.9.1 of IPython, and the documentation for its magic commands is available online.
56
sage: EllipticCurve([0,infinity]) -----------------------------------------------------------Traceback (most recent call last): ... TypeError: Unable to coerce Infinity (<class sage...Infinity>) to Rational
The interactive debugger is sometimes useful for understanding what went wrong. You can toggle it on or off using %pdb (the default is off). The prompt ipdb> appears if an exception is raised and the debugger is on. From within the debugger, you can print the state of any local variable, and move up and down the execution stack. For example,
sage: %pdb Automatic pdb calling has been turned ON sage: EllipticCurve([1,infinity]) --------------------------------------------------------------------------<type exceptions.TypeError> Traceback (most recent call last) ... ipdb>
l list n next p
tbreak u unalias up w
57
sage: V = QQ^3
Then it is easy to list all member functions for V using tab completion. Just type V., then type the [tab key] key on your keyboard:
sage: V.[tab key] V._VectorSpace_generic__base_field ... V.ambient_space V.base_field V.base_ring V.basis V.coordinates ... V.zero_vector
If you type the rst few letters of a function, then [tab key], you get only functions that begin as indicated.
sage: V.i[tab key] V.is_ambient V.is_dense V.is_full V.is_sparse
If you wonder what a particular function does, e.g., the coordinates function, type V.coordinates? for help or V.coordinates?? for the source code, as explained in the next section.
As shown above, the output tells you the type of the object, the le in which it is dened, and a useful description of the function with examples that you can paste into your current session. Almost all of these examples are regularly automatically tested to make sure they work and behave exactly as claimed.
58
Another feature that is very much in the spirit of the open source nature of Sage is that if f is a Python function, then typing f?? displays the source code that denes f. For example,
sage: V = QQ^3 sage: V.coordinates?? Type: instancemethod ... Source: def coordinates(self, v): """ Write $v$ in terms of the basis for self. ... """ return self.coordinate_vector(v).list()
This tells us that all the coordinates function does is call the coordinate_vector function and change the result into a list. What does the coordinate_vector function do?
sage: V = QQ^3 sage: V.coordinate_vector?? ... def coordinate_vector(self, v): ... return self.ambient_vector_space()(v)
The coordinate_vector function coerces its input into the ambient space, which has the effect of computing the vector of coefcients of v in terms of V . The space V is already ambient since its just Q3 . There is also a coordinate_vector function for subspaces, and its different. We create a subspace and see:
sage: V = QQ^3; W = V.span_of_basis([V.0, V.1]) sage: W.coordinate_vector?? ... def coordinate_vector(self, v): """ ... """ # First find the coordinates of v wrt echelon basis. w = self.echelon_coordinate_vector(v) # Next use transformation matrix from echelon basis to # user basis. T = self.echelon_to_user_matrix() return T.linear_combination_of_rows(w)
(If you think the implementation is inefcient, please sign up to help optimize linear algebra.) You may also type help(command_name) or help(class) for a manpage-like help le about a given class.
sage: help(VectorSpace) Help on class VectorSpace ... class VectorSpace(__builtin__.object) | Create a Vector Space. | | To create an ambient space over a field with given dimension | using the calling syntax ... : :
When you type q to exit the help system, your session appears just as it was. The help listing does not clutter up your session, unlike the output of function_name? sometimes does. Its particularly helpful to type 3.8. Integrated Help System 59
help(module_name). For example, vector spaces are dened in sage.modules.free_module, so type help(sage.modules.free_module) for documentation about that whole module. When viewing documentation using help, you can search by typing / and in reverse by typing ?.
You should now quit Sage and restart. Then you can get A back:
sage: sage: [ 15 [ 42 [ 69 A = load(A) A 18 21] 54 66] 90 111]
You can do the same with more complicated objects, e.g., elliptic curves. All data about the object that is cached is stored with the object. For example,
sage: sage: sage: sage: E = EllipticCurve(11a) v = E.anlist(100000) save(E, E) quit # takes a while
The saved version of E takes 153 kilobytes, since it stores the rst 100000 an with it.
~/tmp$ ls -l E.sobj -rw-r--r-- 1 was was 153500 2006-01-28 19:23 E.sobj ~/tmp$ sage [...] sage: E = load(E) sage: v = E.anlist(100000) # instant!
(In Python, saving and loading is accomplished using the cPickle module. In particular, a Sage object x can be saved via cPickle.dumps(x, 2). Note the 2!)
60
Sage cannot save and load individual objects created in some other computer algebra systems, e.g., GAP, Singular, Maxima, etc. They reload in a state marked invalid. In GAP, though many objects print in a form from which they can be reconstructed, many dont, so reconstructing from their print representation is purposely not allowed.
sage: a = gap(2) sage: a.save(a) sage: load(a) Traceback (most recent call last): ... ValueError: The session in which this object was defined is no longer running.
GP/PARI objects can be saved and loaded since their print representation is enough to reconstruct them.
sage: a = gp(2) sage: a.save(a) sage: load(a) 2
Saved objects can be re-loaded later on computers with different architectures or operating systems, e.g., you could save a huge matrix on 32-bit OS X and reload it on 64-bit Linux, nd the echelon form, then move it back. Also, in many cases you can even load objects into versions of Sage that are different than the versions they were saved in, as long as the code for that object isnt too different. All the attributes of the objects are saved, along with the class (but not source code) that denes the object. If that class no longer exists in a new version of Sage, then the object cant be reloaded in that newer version. But you could load it in an old version, get the objects dictionary (with x.__dict__), and save the dictionary, and load that into the newer version.
61
Next we save our session, which saves each of the above variables into a le. Then we view the le, which is about 3K in size.
sage: save_session(misc) Saving a Saving M Saving t Saving E sage: quit was@form:~/tmp$ ls -l misc.sobj -rw-r--r-- 1 was was 2979 2006-01-28 19:47 misc.sobj
Finally we restart Sage, dene an extra variable, and load our saved session.
sage: b = 19 sage: load_session(misc) Loading a Loading M Loading E Loading t
Each saved variable is again available. Moreover, the variable b was not overwritten.
sage: M Full Modular Symbols space for Gamma_0(37) of weight 2 with sign 0 and dimension 5 over Rational Field sage: E Elliptic Curve defined by y^2 + y = x^3 - x^2 - 10*x - 20 over Rational Field sage: b 19 sage: a 389
on the command line of Sage. This starts the Sage notebook and opens your default web browser to view it. The servers state les are stored in $HOME/.sage/sage\_notebook. Other options include:
sage: notebook("directory")
which starts a new notebook server using les in the given directory, instead of the default directory $HOME/.sage/sage_notebook. This can be useful if you want to have a collection of worksheets associated with a specic project, or run several separate notebook servers at the same time. When you start the notebook, it rst creates the following les in $HOME/.sage/sage_notebook:
62
(the notebook SAGE object file) (a directory containing SAGE objects) (a directory containing SAGE worksheets).
After creating the above les, the notebook starts a web server. A notebook is a collection of user accounts, each of which can have any number of worksheets. When you create a new worksheet, the data that denes it is stored in the worksheets/username/number directories. In each such directory there is a plain text le worksheet.txt - if anything ever happens to your worksheets, or Sage, or whatever, that human-readable le contains everything needed to reconstruct your worksheet. From within Sage, type notebook? for much more about how to start a notebook server. The following diagram illustrates the architecture of the Sage Notebook:
---------------------| | | | | firefox/safari | | | | javascript | | program | | | | | ---------------------| ^ | AJAX | V | ---------------------| | | sage | | web | ------------> | server | pexpect | | | | ----------------------
(Python processes)
For help on a Sage command, cmd, in the notebook browser box, type cmd? <shift-enter>).
For help on the keyboard shortcuts available in the notebook interface, click on the Help link.
63
64
CHAPTER
FOUR
INTERFACES
A central facet of Sage is that it supports computation with objects in many different computer algebra systems under one roof using a common interface and clean programming language. The console and interact methods of an interface do very different things. For example, using GAP as an example: 1. gap.console(): This opens the GAP console - it transfers control to GAP. Here Sage is serving as nothing more than a convenient program launcher, similar to the Linux bash shell. 2. gap.interact(): This is a convenient way to interact with a running GAP instance that may be full of Sage objects. You can import Sage objects into this GAP session (even from the interactive interface), etc.
4.1 GP/PARI
PARI is a compact, very mature, highly optimized C program whose primary focus is number theory. There are two very distinct interfaces that you can use in Sage: gp - the G o P ARI interpreter, and pari - the PARI C libraxry. For example, the following are two ways of doing the same thing. They look identical, but the output is actually different, and what happens behind the scenes is drastically different.
sage: gp(znprimroot(10007)) Mod(5, 10007) sage: pari(znprimroot(10007)) Mod(5, 10007)
In the rst case, a separate copy of the GP interpreter is started as a server, and the string znprimroot(10007) is sent to it, evaluated by GP, and the result is assigned to a variable in GP (which takes up space in the child GP processes memory that wont be freed). Then the value of that variable is displayed. In the second case, no separate program is started, and the string znprimroot(10007) is evaluated by a certain PARI C library function. The result is stored in a piece of memory on the Python heap, which is freed when the variable is no longer referenced. The objects have different types:
sage: type(gp(znprimroot(10007))) <class sage.interfaces.gp.GpElement> sage: type(pari(znprimroot(10007))) <type sage.libs.pari.gen.gen>
So which should you use? It depends on what youre doing. The GP interface can do absolutely anything you could do in the usual GP/PARI command line program, since it is running that program. In particular, you can load complicated PARI programs and run them. In contrast, the PARI interface (via the C library) is much more restrictive. First, not all
65
member functions have been implemented. Second, a lot of code, e.g., involving numerical integration, wont work via the PARI interface. That said, the PARI interface can be signicantly faster and more robust than the GP one. (If the GP interface runs out of memory evaluating a given input line, it will silently and automatically double the stack size and retry that input line. Thus your computation wont crash if you didnt correctly anticipate the amount of memory that would be needed. This is a nice trick the usual GP interpreter doesnt seem to provide. Regarding the PARI C library interface, it immediately copies each created object off of the PARI stack, hence the stack never grows. However, each object must not exceed 100MB in size, or the stack will overow when the object is being created. This extra copying does impose a slight performance penalty.) In summary, Sage uses the PARI C library to provide functionality similar to that provided by the GP/PARI interpreter, except with different sophisticated memory management and the Python programming language. First we create a PARI list from a Python list.
sage: v = pari([1,2,3,4,5]) sage: v [1, 2, 3, 4, 5] sage: type(v) <type sage.libs.pari.gen.gen>
Every PARI object is of type py_pari.gen. The PARI type of the underlying object can be obtained using the type member function.
sage: v.type() t_VEC
In PARI, to create an elliptic curve we enter ellinit([1,2,3,4,5]). Sage is similar, except that ellinit is a method that can be called on any PARI object, e.g., our t\_VEC v.
sage: e = v.ellinit() sage: e.type() t_VEC sage: pari(e)[:13] [1, 2, 3, 4, 5, 9, 11, 29, 35, -183, -3429, -10351, 6128487/10351]
Now that we have an elliptic curve object, we can compute some things about it.
sage: e.elltors() [1, [], []] sage: e.ellglobalred() [10351, [1, -1, 0, -1], 1] sage: f = e.ellchangecurve([1,-1,0,-1]) sage: f[:5] [1, -1, 0, 4, 3]
4.2 GAP
Sage comes with GAP 4.4.10 for computational discrete mathematics, especially group theory. Heres an example of GAPs IdGroup function, which uses the optional small groups database that has to be installed separately, as explained below.
sage: G = gap(Group((1,2,3)(4,5), (3,4))) sage: G Group( [ (1,2,3)(4,5), (3,4) ] ) sage: G.Center() Group( () )
66
Chapter 4. Interfaces
We can do the same computation in Sage without explicitly invoking the GAP interface as follows:
sage: G = PermutationGroup([[(1,2,3),(4,5)],[(3,4)]]) sage: G.center() Subgroup of (Permutation Group with generators [(3,4), (1,2,3)(4,5)]) generated by [()] sage: G.group_id() # requires optional database_gap package [120, 34] sage: n = G.order(); n 120
(For some GAP functionality, you should install two optional Sage packages. Type sage -optional for a list and choose the one that looks like gap\_packages-x.y.z, then type sage -i gap\_packages-x.y.z. Do the same for database\_gap-x.y.z. Some non-GPLd GAP packages may be installed by downloading them from the GAP web site [GAPkg], and unpacking them in $SAGE_ROOT/local/lib/gap-4.4.10/pkg. )
4.3 Singular
Singular provides a massive and mature library for Grbner bases, multivariate polynomial gcds, bases of RiemannRoch spaces of a plane curve, and factorizations, among other things. We illustrate multivariate polynomial factorization using the Sage interface to Singular (do not type the ...):
sage: R1 = singular.ring(0, (x,y), dp) sage: R1 // characteristic : 0 // number of vars : 2 // block 1 : ordering dp // : names x y // block 2 : ordering C sage: f = singular(9*y^8 - 9*x^2*y^7 - 18*x^3*y^6 - 18*x^5*y^6 + \ ... 9*x^6*y^4 + 18*x^7*y^5 + 36*x^8*y^4 + 9*x^10*y^4 - 18*x^11*y^2 - \ ... 9*x^12*y^3 - 18*x^13*y^2 + 9*x^16)
sage: f 9*x^16-18*x^13*y^2-9*x^12*y^3+9*x^10*y^4-18*x^11*y^2+36*x^8*y^4+18*x^7*y^5-18*x^5*y^6+9*x^6*y^4-18*x^ sage: f.parent() Singular sage: F = f.factorize(); F [1]: _[1]=9 _[2]=x^6-2*x^3*y^2-x^2*y^3+y^4 _[3]=-x^5+y^2 [2]: 1,1,2 sage: F[1][2] x^6-2*x^3*y^2-x^2*y^3+y^4
As with the GAP example in GAP, we can compute the above factorization without explicitly using the Singular interface (however, behind the scenes Sage uses the Singular interface for the actual computation). Do not type the ...:
4.3. Singular
67
x, y = QQ[x, y].gens() f = 9*y^8 - 9*x^2*y^7 - 18*x^3*y^6 - 18*x^5*y^6 + 9*x^6*y^4\ + 18*x^7*y^5 + 36*x^8*y^4 + 9*x^10*y^4 - 18*x^11*y^2 - 9*x^12*y^3\ - 18*x^13*y^2 + 9*x^16 factor(f) (-x^5 + y^2)^2 * (x^6 - 2*x^3*y^2 - x^2*y^3 + y^4)
4.4 Maxima
Maxima is included with Sage, as well as a Lisp implementation. The gnuplot package (which Maxima uses by default for plotting) is distributed as a Sage optional package. Among other things, Maxima does symbolic manipulation. Maxima can integrate and differentiate functions symbolically, solve 1st order ODEs, most linear 2nd order ODEs, and has implemented the Laplace transform method for linear ODEs of any degree. Maxima also knows about a wide range of special functions, has plotting capabilities via gnuplot, and has methods to solve and manipulate matrices (such as row reduction, eigenvalues and eigenvectors), and polynomial equations. We illustrate the Sage/Maxima interface by constructing the matrix whose i, j entry is i/j, for i, j = 1, . . . , 4.
sage: f = maxima.eval(ij_entry[i,j] := i/j) sage: A = maxima(genmatrix(ij_entry,4,4)); A matrix([1,1/2,1/3,1/4],[2,1,2/3,1/2],[3,3/2,1,3/4],[4,2,4/3,1]) sage: A.determinant() 0 sage: A.echelon() matrix([1,1/2,1/3,1/4],[0,0,0,0],[0,0,0,0],[0,0,0,0]) sage: A.eigenvalues() [[0,4],[3,1]] sage: A.eigenvectors() [[[0,4],[3,1]],[[[1,0,0,-4],[0,1,0,-2],[0,0,1,-4/3]],[[1,2,3,4]]]]
Finally, we give an example of using Sage to plot using openmath. Many of these were modied from the Maxima reference manual. A 2D plot of several functions (do not type the ...): 68 Chapter 4. Interfaces
A live 3D plot which you can move with your mouse (do not type the ...):
sage: ... sage: ... maxima.plot3d ("2^(-u^2 + v^2)", "[u, -3, 3]", "[v, -2, 2]",\ [plot_format, openmath]) # not tested maxima.plot3d("atan(-x^2 + y^3/4)", "[x, -4, 4]", "[y, -4, 4]",\ "[grid, 50, 50]",[plot_format, openmath]) # not tested
The next plot is the famous Mbius strip (do not type the ...):
sage: maxima.plot3d("[cos(x)*(3 + y*cos(x/2)), sin(x)*(3 + y*cos(x/2)),\ ... y*sin(x/2)]", "[x, -4, 4]", "[y, -4, 4]",\ ... [plot_format, openmath]) # not tested
The next plot is the famous Klein bottle (do not type the ...):
sage: maxima("expr_1: 5*cos(x)*(cos(x/2)*cos(y) + sin(x/2)*sin(2*y)+ 3.0)\ ... - 10.0") 5*cos(x)*(sin(x/2)*sin(2*y)+cos(x/2)*cos(y)+3.0)-10.0 sage: maxima("expr_2: -5*sin(x)*(cos(x/2)*cos(y) + sin(x/2)*sin(2*y)+ 3.0)") -5*sin(x)*(sin(x/2)*sin(2*y)+cos(x/2)*cos(y)+3.0) sage: maxima("expr_3: 5*(-sin(x/2)*cos(y) + cos(x/2)*sin(2*y))") 5*(cos(x/2)*sin(2*y)-sin(x/2)*cos(y)) sage: maxima.plot3d ("[expr_1, expr_2, expr_3]", "[x, -%pi, %pi]",\ ... "[y, -%pi, %pi]", "[grid, 40, 40]",\ ... [plot_format, openmath]) # not tested
4.4. Maxima
69
70
Chapter 4. Interfaces
CHAPTER
FIVE
5.1 Overview
It may be easiest to understand the various uses of LaTeX with a brief overview of the mechanics of the three principal methods employed by Sage. 1. Every object in Sage is required to have a LaTeX representation. You can access this representation by executing, in the notebook or at the sage command line, latex(foo) where foo is some object in Sage. The output is a string that should render a reasonably accurate representation of foo when used in TeXs math-mode (for example, when enclosed between a pair of single dollar signs). Some examples of this follow below. In this way, Sage can be used effectively for constructing portions of a LaTeX document: create or compute an object in Sage, print latex() of the object and cut/paste it into your document. 2. The notebook interface is congured to use MathJax to render mathematics cleanly in a web browser. MathJax is an open source JavaScript display engine for mathematics that works in all modern browsers. It is able to render a large, but not totally complete, subset of TeX. It has no support for things like complicated tables, sectioning or document management, as it is oriented towards accurately rendering snippets of TeX. Seemingly automatic rendering of math in the notebook is provided by converting the latex() representation of an object (as described above) into a form of HTML palatable to MathJax. Since MathJax uses its own scalable fonts, it is superior to other methods that rely on converting equations, or other snippets of TeX, into static inline images. 3. At the Sage command-line, or in the notebook when LaTeX code is more involved than MathJax can handle, a system-wide installation of LaTeX can be employed. Sage includes almost everything you need to build and use Sage, but a signicant exception is TeX itself. So in these situations you need to have TeX installed, along with some associated conversion utilities, to utilize the full power. Here we demonstrate some basic uses of the latex() function.
sage: var(z) z sage: latex(z^12) z^{12} sage: latex(integrate(z^4, z)) \frac{1}{5} \, z^{5} sage: latex(a string)
71
\verb|a|\phantom{\verb!x!}\verb|string| sage: latex(QQ) \Bold{Q} sage: latex(matrix(QQ, 2, 3, [[2,4,6],[-1,-1,-1]])) \left(\begin{array}{rrr} 2 & 4 & 6 \\ -1 & -1 & -1 \end{array}\right)
Basic MathJax functionality is largely automatic in the notebook, but we can partially demonstrate this support with the MathJax class, The eval function of this class converts a Sage object to its LaTeX representation and then wraps it in HTML that invokes the CSS math class, which then employs MathJax.
sage: from sage.misc.latex import MathJax sage: mj = MathJax() sage: var(z) z sage: mj(z^12) <html><script type="math/tex; mode=display">\newcommand{\Bold}[1]{\mathbf{#1}}z^{12}</script></html> sage: mj(QQ) <html><script type="math/tex; mode=display">\newcommand{\Bold}[1]{\mathbf{#1}}\Bold{Q}</script></html sage: mj(ZZ[x]) <html><script type="math/tex; mode=display">\newcommand{\Bold}[1]{\mathbf{#1}}\Bold{Z}[x]</script></h sage: mj(integrate(z^4, z)) <html><script type="math/tex; mode=display">\newcommand{\Bold}[1]{\mathbf{#1}}\frac{1}{5} \, z^{5}</s
72
The notebook has two other features for employing TeX. The rst is the Typeset button just above the rst cell of a worksheet, to the right of the four drop-down boxes. When checked, any subsequent evaluations of cells will result in output interpreted by MathJax, hence of a typeset quality. Note that this effect is not retroactive previously evaluated cells need to be re-evaluated. Essentially, checking the Typeset button is identical to wrapping the output of each cell in the view() command. A second feature of the notebook is entering TeX as part of annotating a worksheet. When the cursor is placed between cells of a worksheet so that a blue bar appears, then a shift-click will open a mini-word-processor, TinyMCE. This allows for the entry of text, using a WSIWYG editor to create HTML and CSS command for styled text. So it is possible to add formatted text as commentary within a worksheet. However, text between pairs of dollar signs, or pairs of double dollar signs is interpreted by MathJax as inline or display math (respectively).
The latex.vector_delimiters method works similarly. The way common rings and elds (integers, rational, reals, etc.) are typeset can be controlled by the latex.blackboard_bold method. These sets are by default typeset in bold, but may optionally be written in a double-struck fashion as sometimes done in written work. This is accomplished by redening the \Bold{} macro which is built-in to Sage.
sage: latex(QQ) \Bold{Q} sage: from sage.misc.latex import MathJax sage: mj=MathJax() sage: mj(QQ) <html><script type="math/tex; mode=display">\newcommand{\Bold}[1]{\mathbf{#1}}\Bold{Q}</script></html sage: latex.blackboard_bold(True) sage: mj(QQ)
73
It is possible to take advantage of the extensible nature of tex by adding in new macros and new packages. First, individual macros can be added so that they are used when MathJax interprets a snippet of TeX in the notebook.
sage: latex.extra_macros() sage: latex.add_macro("\\newcommand{\\foo}{bar}") sage: latex.extra_macros() \\newcommand{\\foo}{bar} sage: var(x y) (x, y) sage: latex(x+y) x + y sage: from sage.misc.latex import MathJax sage: mj=MathJax() sage: mj(x+y) <html><script type="math/tex; mode=display">\newcommand{\Bold}[1]{\mathbf{#1}}\newcommand{\foo}{bar}x
Additional macros added this way will also be used in the event that the system-wide version of TeX is called on something larger than MathJax can handle. The command latex_extra_preamble is used to build the preamble of a complete LaTeX document, so the following illustrates how this is accomplished. As usual note the need for the double-backslashes in the Python strings.
sage: latex.extra_macros() sage: latex.extra_preamble() sage: from sage.misc.latex import latex_extra_preamble sage: print latex_extra_preamble() \newcommand{\ZZ}{\Bold{Z}} ... \newcommand{\Bold}[1]{\mathbf{#1}} sage: latex.add_macro("\\newcommand{\\foo}{bar}") sage: print latex_extra_preamble() \newcommand{\ZZ}{\Bold{Z}} ... \newcommand{\Bold}[1]{\mathbf{#1}} \newcommand{\foo}{bar}
Again, for larger or more complicated LaTeX expressions, it is possible to add packages (or anything else) to the preamble of the LaTeX le. Anything may be incorporated into the preamble with the latex.add_to_preamble command, and the specialized command latex.add_package_to_preamble_if_available will rst check if a certain package is actually available before trying to add it to the preamble. Here we add the geometry package to the preamble and use it to set the size of the region on the page that TeX will use (effectively setting the margins). As usual, note the need for the double-backslashes in the Python strings.
sage: from sage.misc.latex import latex_extra_preamble sage: latex.extra_macros() sage: latex.extra_preamble() sage: latex.add_to_preamble(\\usepackage{geometry}) sage: latex.add_to_preamble(\\geometry{letterpaper,total={8in,10in}}) sage: latex.extra_preamble() \\usepackage{geometry}\\geometry{letterpaper,total={8in,10in}} sage: print latex_extra_preamble() \usepackage{geometry}\geometry{letterpaper,total={8in,10in}} \newcommand{\ZZ}{\Bold{Z}} ... \newcommand{\Bold}[1]{\mathbf{#1}}
74
A particular package may be added along with a check on its existence, as follows. As an example, we just illustrate an attempt to add to the preamble a package that presumably does not exist.
sage: latex.extra_preamble() sage: latex.extra_preamble() sage: latex.add_to_preamble(\\usepackage{foo-bar-unchecked}) sage: latex.extra_preamble() \\usepackage{foo-bar-unchecked} sage: latex.add_package_to_preamble_if_available(foo-bar-checked) sage: latex.extra_preamble() \\usepackage{foo-bar-unchecked}
Suppose a LaTeX expression is produced in the notebook with view() or while the Typeset button is checked, and then recognized as requiring the external LaTeX installation through the mathjax avoid list. Then the selected executable (as specied by latex.engine()) will process the LaTeX. However, instead of then spawning an external viewer (which is the command-line behavior), Sage will attempt to convert the result into a single, tightlycropped image, which is then inserted into the worksheet as the output of the cell. Just how this conversion proceeds depends on several factors mostly which executable you have specied as the engine and which conversion utilities are available on your system. Four useful converters that will cover all eventualities are dvips, ps2pdf, dvipng and from the ImageMagick suite, convert. The goal is to produce a PNG le as the output for inclusion back into the worksheet. When a LaTeX expression can be converted successfully to a dvi by the latex engine, then dvipng should accomplish the conversion. If the LaTeX expression and chosen engine creates a dvi with specials that dvipng cannot handle, then dvips will create a PostScript le. Such a PostScript le,
75
or a PDF le created by an engine such as pdflatex, is then processed into a PNG with the convert utility. The presence of two of these converters can be tested with the have_dvipng() and have_convert() routines. These conversions are done automatically if you have the necessary converters installed; if not, then an error message is printed telling you whats missing and where to download it. For a concrete example of how complicated LaTeX expressions can be processed, see the example in the next section (An Example: Combinatorial Graphs with tkz-graph) for using the LaTeX tkz-graph package to produce high-quality renderings of combinatorial graphs. For other examples, there are some pre-packaged test cases. To use these, it is necessary to import the sage.misc.latex.latex_examples object, which is an instance of the sage.misc.latex.LatexExamples class, as illustrated below. This class currently has examples of commutative diagrams, combinatorial graphs, knot theory and pstricks, which respectively exercise the following packages: xy, tkz-graph, xypic, pstricks. After the import, use tab-completion on latex_examples to see the pre-packaged examples. Calling each example will give you back some explanation about what is required to make the example render properly. To actually see the examples, it is necessary to use view() (once the preamble, engine, etc are all set properly).
sage: from sage.misc.latex import latex_examples sage: latex_examples.diagram() LaTeX example for testing display of a commutative diagram produced by xypic. To use, try to view this object -- it wont work. Now try latex.add_to_preamble("\\usepackage[matrix,arrow,curve,cmtip]{xy}"), and try viewing again -- it should work in the command line but not from the notebook. In the notebook, run latex.add_to_mathjax_avoid_list("xymatrix") and try again -- you should get a picture (a part of the diagram arising from a filtered chain complex).
76
At this point, a command like view(graphs.CompleteGraph(4)) should produce a graphic version of the graph pasted into the notebook, having used pdflatex to process tkz-graph commands to realize the graph. Note that there is a variety of options to affect how a graph is rendered in LaTeX via tkz-graph, which is again outside the scope of this section, see the section of the Reference manual titled LaTeX Options for Graphs for instructions and details.
77
78
CHAPTER
SIX
PROGRAMMING
6.1 Loading and Attaching Sage les
Next we illustrate how to load programs written in a separate le into Sage. Create a le called example.sage with the following content:
print "Hello World" print 2^3
You can read in and execute example.sage le using the load command.
sage: load "example.sage" Hello World 8
You can also attach a Sage le to a running session using the attach command:
sage: attach "example.sage" Hello World 8
Now if you change example.sage and enter one blank line into Sage (i.e., hit return), then the contents of example.sage will be automatically reloaded into Sage. In particular, attach automatically reloads a le whenever it changes, which is handy when debugging code, whereas load only loads a le once. When Sage loads example.sage it converts it to Python, which is then executed by the Python interpreter. This conversion is minimal; it mainly involves wrapping integer literals in Integer() oating point literals in RealNumber(), replacing ^s by **s, and replacing e.g., R.2 by R.gen(2)}. The converted version of example.sage is contained in the same directory as example.sage and is called example.sage.py. This le contains the following code:
print "Hello World" print Integer(2)**Integer(3)
Integer literals are wrapped and the ^ is replaced by a **. (In Python ^ means exclusive or and ** means exponentiation.) This preparsing is implemented in sage/misc/interpreter.py.) You can paste multi-line indented code into Sage as long as there are newlines to make new blocks (this is not necessary in les). However, the best way to enter such code into Sage is to save it to a le and use attach, as described above.
79
Here the trailing L indicates a Python long integer (see The Pre-Parser: Differences between Sage and Python). Note that Sage will recompile factorial.spyx if you quit and restart Sage. The compiled shared object library is stored under $HOME/.sage/temp/hostname/pid/spyx. These les are deleted when you exit Sage. NO Sage preparsing is applied to spyx les, e.g., 1/3 will result in 0 in a spyx le instead of the rational number 1/3. If foo is a function in the Sage library, to use it from a spyx le import sage.all and use sage.all.foo.
import sage.all def foo(n): return sage.all.factorial(n)
80
Chapter 6. Programming
If an additional library foo is needed to compile the C code generated from a Cython le, add the line clib foo to the Cython source. Similarly, an additional C le bar can be included in the compilation with the declaration cfile bar.
In order to use this script, your SAGE_ROOT must be in your PATH. If the above script is called factor, here is an example usage:
bash $ 2 * 17 bash $ (2*x ./factor 2006 * 59 ./factor "32*x^5-1" 1) * (16*x^4 + 8*x^3 + 4*x^2 + 2*x + 1)
81
Only certain functions can be called on V. In other math software systems, these would be called using the functional notation foo(V,...). In Sage, certain functions are attached to the type (or class) of V, and are called using an object-oriented syntax like in Java or C++, e.g., V.foo(...). This helps keep the global namespace from being polluted with tens of thousands of functions, and means that many different functions with different behavior can be named foo, without having to use type-checking of arguments (or case statements) to decide which to call. Also, if you reuse the name of a function, that function is still available (e.g., if you call something zeta, then want to compute the value of the Riemann-Zeta function at 0.5, you can still type s=.5; s.zeta()).
sage: zeta = -1 sage: s=.5; s.zeta() -1.46035450880959
In some very common cases, the usual functional notation is also supported for convenience and because mathematical expressions might look confusing using object-oriented notation. Here are some examples.
sage: n = 2; n.sqrt() sqrt(2) sage: sqrt(2) sqrt(2) sage: V = VectorSpace(QQ,2) sage: V.basis() [ (1, 0), (0, 1) ] sage: basis(V) [ (1, 0), (0, 1) ] sage: M = MatrixSpace(GF(7), 2); M Full MatrixSpace of 2 by 2 dense matrices over Finite Field of size 7 sage: A = M([1,2,3,4]); A [1 2] [3 4] sage: A.charpoly(x) x^2 + 2*x + 5 sage: charpoly(A, x) x^2 + 2*x + 5
To list all member functions for A, use tab completion. Just type A., then type the [tab] key on your keyboard, as explained in Reverse Search and Tab Completion.
82
Chapter 6. Programming
(When indexing into a list, it is OK if the index is not a Python int!) A Sage Integer (or Rational, or anything with an __index__ method) will work just ne.
sage: sage: 3 sage: sage: 3 sage: 3 v = [1,2,3] v[2] n = 2 v[n] v[int(n)] # SAGE Integer # Perfectly OK! # Also OK.
The range function creates a list of Python ints (not Sage Integers):
sage: range(1, 15) [1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14]
For more about how to create lists using list comprehensions, see [PyT]. List slicing is a wonderful feature. If L is a list, then L[m:n] returns the sublist of L obtained by starting at the mth element and stopping at the (n 1)st element, as illustrated below.
sage: L = [factor(n) for n in range(1, 20)] sage: L[4:9] [5, 2 * 3, 7, 2^3, 3^2] sage: print L[:4] [1, 2, 3, 2^2] sage: L[14:4] [] sage: L[14:] [3 * 5, 2^4, 17, 2 * 3^2, 19]
Tuples are similar to lists, except they are immutable, meaning once they are created they cant be changed.
sage: v = (1,2,3,4); v (1, 2, 3, 4) sage: type(v) <type tuple> sage: v[1] = 5 Traceback (most recent call last):
83
Sequences are a third list-oriented Sage type. Unlike lists and tuples, Sequence is not a built-in Python type. By default, a sequence is mutable, but using the Sequence class method set_immutable, it can be set to be immutable, as the following example illustrates. All elements of a sequence have a common parent, called the sequences universe.
sage: v = Sequence([1,2,3,4/5]) sage: v [1, 2, 3, 4/5] sage: type(v) <class sage.structure.sequence.Sequence_generic> sage: type(v[1]) <type sage.rings.rational.Rational> sage: v.universe() Rational Field sage: v.is_immutable() False sage: v.set_immutable() sage: v[0] = 3 Traceback (most recent call last): ... ValueError: object is immutable; please change a copy instead.
Sequences derive from lists and can be used anywhere a list can be used:
sage: v = Sequence([1,2,3,4/5]) sage: isinstance(v, list) True sage: list(v) [1, 2, 3, 4/5] sage: type(list(v)) <type list>
As another example, basis for vector spaces are immutable sequences, since its important that you dont change them.
sage: V = QQ^3; B = V.basis(); B [ (1, 0, 0), (0, 1, 0), (0, 0, 1) ] sage: type(B) <class sage.structure.sequence.Sequence_generic> sage: B[0] = B[1] Traceback (most recent call last): ... ValueError: object is immutable; please change a copy instead. sage: B.universe() Vector space of dimension 3 over Rational Field
6.6 Dictionaries
A dictionary (also sometimes called an associative array) is a mapping from hashable objects (e.g., strings, numbers, and tuples of such; see the Python documentation http://docs.python.org/tut/node7.html and http://docs.python.org/lib/typesmapping.html for details) to arbitrary objects.
84
Chapter 6. Programming
sage: d = {1:5, sage:17, ZZ:GF(7)} sage: type(d) <type dict> sage: d.keys() [1, sage, Integer Ring] sage: d[sage] 17 sage: d[ZZ] Finite Field of size 7 sage: d[1] 5
The third key illustrates that the indexes of a dictionary can be complicated, e.g., the ring of integers. You can turn the above dictionary into a list with the same data:
sage: d.items() [(1, 5), (sage, 17), (Integer Ring, Finite Field of size 7)]
6.7 Sets
Python has a built-in set type. The main feature it offers is very fast lookup of whether an element is in the set or not, along with standard set-theoretic operations.
sage: X = set([1,19,a]); sage: X set([a, 1, 19]) sage: Y set([1, 2/3]) sage: a in X True sage: a in Y False sage: X.intersection(Y) set([1]) Y = set([1,1,1, 2/3])
Sage also has its own set type that is (in some cases) implemented using the built-in Python set type, but has a little bit of extra Sage-related functionality. Create a Sage set using Set(...). For example,
sage: X = Set([1,19,a]); Y = Set([1,1,1, 2/3]) sage: X {a, 1, 19} sage: Y {1, 2/3} sage: X.intersection(Y) {1} sage: print latex(Y) \left\{1, \frac{2}{3}\right\} sage: Set(ZZ) Set of elements of Integer Ring
6.7. Sets
85
6.8 Iterators
Iterators are a recent addition to Python that are particularly useful in mathematics applications. Here are several examples; see [PyT] for more details. We make an iterator over the squares of the nonnegative integers up to 10000000.
sage: sage: 0 sage: 1 sage: 4 v = (n^2 for n in xrange(10000000)) v.next() v.next() v.next()
We create an iterate over the primes of the form 4p + 1 with p also prime, and look at the rst few values.
sage: w = (4*p + 1 for p in Primes() if is_prime(4*p+1)) sage: w # in the next line, 0xb0853d6c is a random 0x number <generator object at 0xb0853d6c> sage: w.next() 13 sage: w.next() 29 sage: w.next() 53
Certain rings, e.g., nite elds and the integers have iterators associated to them:
sage: [x for x in GF(7)] [0, 1, 2, 3, 4, 5, 6] sage: W = ((x,y) for x in ZZ for y in ZZ) sage: W.next() (0, 0) sage: W.next() (0, 1) sage: W.next() (0, -1)
Note the colon at the end of the for statement (there is no do or od as in GAP or Maple), and the indentation before the body of the loop, namely print(i). This indentation is important. In Sage, the indentation is automatically put in for you when you hit enter after a :, as illustrated below.
86
Chapter 6. Programming
The symbol = is used for assignment. The symbol == is used to check for equality:
sage: for i in range(15): ... if gcd(i,15) == 1: ... print(i) 1 2 4 7 8 11 13 14
Keep in mind how indentation determines the block structure for if, for, and while statements:
sage: def legendre(a,p): ... is_sqr_modp=-1 ... for i in range(p): ... if a % p == i^2 % p: ... is_sqr_modp=1 ... return is_sqr_modp sage: legendre(2,7) 1 sage: legendre(3,7) -1
Of course this is not an efcient implementation of the Legendre symbol! It is meant to illustrate various aspects of Python/Sage programming. The function {kronecker}, which comes with Sage, computes the Legendre symbol efciently via a C-library call to PARI. Finally, we note that comparisons, such as ==, !=, <=, >=, >, <, between numbers will automatically convert both numbers into the same type if possible:
sage: 2 < 3.1; 3.1 <= 1 True False sage: 2/3 < 3/2; 3/2 < 3/1 True True
Almost any two objects may be compared; there is no assumption that the objects are equipped with a total ordering.
sage: 2 < CC(3.1,1) True sage: 5 < VectorSpace(QQ,3) True
87
When comparing objects of different types in Sage, in most cases Sage tries to nd a canonical coercion of both objects to a common parent (see Parents, Conversion and Coercion for more details). If successful, the comparison is performed between the coerced objects; if not successful, the objects are considered not equal. For testing whether two variables reference the same object use is. As we see in this example, the Python int 1 is unique, but the Sage Integer 1 is not:
sage: False sage: True sage: False sage: True 1 is 2/2 int(1) is int(2)/int(2) 1 is 1 1 == 2/2
In the following two lines, the rst equality is False because there is no canonical morphism Q F5 , hence no canonical way to compare the 1 in F5 to the 1 Q. In contrast, there is a canonical map Z F5 , hence the second comparison is True. Note also that the order doesnt matter.
sage: GF(5)(1) == QQ(1); QQ(1) == GF(5)(1) False False sage: GF(5)(1) == ZZ(1); ZZ(1) == GF(5)(1) True True sage: ZZ(1) == QQ(1) True
WARNING: Comparison in Sage is more restrictive than in Magma, which declares the 1 F5 equal to 1 Q.
sage: magma(GF(5)!1 eq Rationals()!1) true # optional magma required
6.10 Proling
Section Author: Martin Albrecht (malb@informatik.uni-bremen.de) Premature optimization is the root of all evil. - Donald Knuth Sometimes it is useful to check for bottlenecks in code to understand which parts take the most computational time; this can give a good idea of which parts to optimize. Python and therefore Sage offers several prolingas this process is calledoptions. The simplest to use is the prun command in the interactive shell. It returns a summary describing which functions took how much computational time. To prole (the currently slow! - as of version 1.0) matrix multiplication over nite elds, for example, do:
sage: k,a = GF(2**8, a).objgen() sage: A = Matrix(k,10,10,[k.random_element() for _ in range(10*10)]) sage: %prun B = A*A 32893 function calls in 1.100 CPU seconds
88
Chapter 6. Programming
Ordered by: internal time ncalls tottime percall cumtime percall filename:lineno(function) 12127 0.160 0.000 0.160 0.000 :0(isinstance) 2000 0.150 0.000 0.280 0.000 matrix.py:2235(__getitem__) 1000 0.120 0.000 0.370 0.000 finite_field_element.py:392(__mul__) 1903 0.120 0.000 0.200 0.000 finite_field_element.py:47(__init__) 1900 0.090 0.000 0.220 0.000 finite_field_element.py:376(__compat) 900 0.080 0.000 0.260 0.000 finite_field_element.py:380(__add__) 1 0.070 0.070 1.100 1.100 matrix.py:864(__mul__) 2105 0.070 0.000 0.070 0.000 matrix.py:282(ncols) ...
Here ncalls is the number of calls, tottime is the total time spent in the given function (and excluding time made in calls to sub-functions), percall is the quotient of tottime divided by ncalls. cumtime is the total time spent in this and all sub-functions (i.e., from invocation until exit), percall is the quotient of cumtime divided by primitive calls, and filename:lineno(function) provides the respective data of each function. The rule of thumb here is: The higher the function in that listing, the more expensive it is. Thus it is more interesting for optimization. As usual, prun? provides details on how to use the proler and understand the output. The proling data may be written to an object as well to allow closer examination:
sage: %prun -r A*A sage: stats = _ sage: stats?
Note: entering stats = prun -r A\*A displays a syntax error message because prun is an IPython shell command, not a regular function. For a nice graphical representation of proling data, you can use the hotshot proler, a small script called hotshot2cachetree and the program kcachegrind (Unix only). The same example with the hotshot proler:
sage: sage: sage: sage: sage: k,a = GF(2**8, a).objgen() A = Matrix(k,10,10,[k.random_element() for _ in range(10*10)]) import hotshot filename = "pythongrind.prof" prof = hotshot.Profile(filename, lineevents=1)
This results in a le pythongrind.prof in the current working directory. It can now be converted to the cachegrind format for visualization. On a system shell, type
hotshot2calltree -o cachegrind.out.42 pythongrind.prof
The output le cachegrind.out.42 can now be examined with kcachegrind. Please note that the naming convention cachegrind.out.XX needs to be obeyed.
6.10. Proling
89
90
Chapter 6. Programming
CHAPTER
SEVEN
USING SAGETEX
The SageTeX package allows you to embed the results of Sage computations into a LaTeX document. It comes standard with Sage. To use it, you will need to install it into your local TeX system; here install means copying a single le. See Installation in this tutorial and the Make SageTeX known to TeX section of the Sage installation guide (this link should take you to a local copy of the installation guide) for more information on doing that. Here is a very brief example of using SageTeX. The full documentation can be found in SAGE_ROOT/local/share/texmf/tex/generic/sagetex, where SAGE_ROOT is the directory where your Sage installation is located. That directory contains the documentation, an example le, and some possibly useful Python scripts. To see how SageTeX works, follow the directions for installing SageTeX (in Installation) and copy the following text into a le named, say, st_example.tex: Warning: The text below will have several errors about unknown control sequences if you are viewing this in the live help. Use the static version to see the correct text.
\documentclass{article} \usepackage{sagetex} \begin{document} Using Sage\TeX, one can use Sage to compute things and put them into your \LaTeX{} document. For example, there are $\sage{number_of_partitions(1269)}$ integer partitions of $1269$. You dont need to compute the number yourself, or even cut and paste it from somewhere. Heres some Sage code: \begin{sageblock} f(x) = exp(x) * sin(2*x) \end{sageblock} The second derivative of $f$ is \[ \frac{\mathrm{d}^{2}}{\mathrm{d}x^{2}} \sage{f(x)} = \sage{diff(f, x, 2)(x)}. \] Heres a plot of $f$ from $-1$ to $1$: \sageplot{plot(f, -1, 1)}
91
\end{document}
Run LaTeX on st_example.tex as usual. Note that LaTeX will have some complaints, which will include:
Package sagetex Warning: Graphics file sage-plots-for-st_example.tex/plot-0.eps on page 1 does not exist. Plot command is on input line 25. Package sagetex Warning: There were undefined Sage formulas and/or plots. Run Sage on st_example.sage, and then run LaTeX on st_example.tex again.
Notice that, in addition to the usual collection of les produced by LaTeX, there is a le called st_example.sage. That is a Sage script produced when you run LaTeX on st_example.tex. The warning message told you to run Sage on st_example.sage, so take its advice and do that. It will tell you to run LaTeX on st_example.tex again, but before you do that, notice that a new le has been created: st_example.sout. That le contains the results of Sages computations, in a format that LaTeX can use to insert into your text. A new directory containing an EPS le of your plot has also been created. Run LaTeX again and youll see that everything that Sage computed and plotted is now included in your document. The different macros used above should be pretty easy to understand. A sageblock environment typesets your code verbatim and also executes the code when you run Sage. When you do \sage{foo}, the result put into your document is whatever you get from running latex(foo) inside Sage. Plot commands are a bit more complicated, but in their simplest form, \sageplot{foo} inserts the image you get from doing foo.save(filename.eps). In general, the mantra is: run LaTeX on your .tex le; run Sage on the generated .sage le; run LaTeX again. You can omit running Sage if you havent changed around any Sage commands in your document. Theres a lot more to SageTeX, and since both Sage and LaTeX are complex, powerful tools, its a good idea to read the documentation for SageTeX, which is in SAGE_ROOT/local/share/texmf/tex/generic/sagetex.
92
CHAPTER
EIGHT
AFTERWORD
8.1 Why Python?
8.1.1 Advantages of Python
The primary implementation language of Sage is Python (see [Py]), though code that must be fast is implemented in a compiled language. Python has several advantages: Object saving is well-supported in Python. There is extensive support in Python for saving (nearly) arbitrary objects to disk les or a database. Excellent support for documentation of functions and packages in the source code, including automatic extraction of documentation and automatic testing of all examples. The examples are automatically tested regularly and guaranteed to work as indicated. Memory management: Python now has a well thought out and robust memory manager and garbage collector that correctly deals with circular references, and allows for local variables in les. Python has many packages available now that might be of great interest to users of Sage: numerical analysis and linear algebra, 2D and 3D visualization, networking (for distributed computations and servers, e.g., via twisted), database support, etc. Portability: Python is easy to compile from source on most platforms in minutes. Exception handling: Python has a sophisticated and well thought out system of exception handling, whereby programs gracefully recover even if errors occur in code they call. Debugger: Python includes a debugger, so when code fails for some reason, the user can access an extensive stack trace, inspect the state of all relevant variables, and move up and down the stack. Proler: There is a Python proler, which runs code and creates a report detailing how many times and for how long each function was called. A Language: Instead of writing a new language for mathematics as was done for Magma, Maple, Mathematica, Matlab, GP/PARI, GAP, Macaulay 2, Simath, etc., we use the Python language, which is a popular computer language that is being actively developed and optimized by hundreds of skilled software engineers. Python is a major open-source success story with a mature development process (see [PyDev]).
93
This use of ^ may appear odd, and it is inefcient for pure math research, since the exclusive or function is rarely used. For convenience, Sage pre-parses all command lines before passing them to Python, replacing instances of ^ that are not in strings with **:
sage: 2^8 256 sage: 3^2 9 sage: "3^2" 3^2
Integer division: The Python expression 2/3 does not behave the way mathematicians might expect. In Python, if m and n are ints, then m/n is also an int, namely the quotient of m divided by n. Therefore 2/3=0. There has been talk in the Python community about changing Python so 2/3 returns the oating point number 0.6666..., and making 2//3 return 0. We deal with this in the Sage interpreter, by wrapping integer literals in Integer( ) and making division a constructor for rational numbers. For example:
sage: 2/3 2/3 sage: (2/3).parent() Rational Field sage: 2//3 0 sage: int(2)/int(3) 0
Long integers: Python has native support for arbitrary precision integers, in addition to C-ints. These are signicantly slower than what GMP provides, and have the property that they print with an L at the end to distinguish them from ints (and this wont change any time soon). Sage implements arbitrary precision integers using the GMP C-library, and these print without an L. Rather than modifying the Python interpreter (as some people have done for internal projects), we use the Python language exactly as is, and write a pre-parser for IPython so that the command line behavior of IPython is what a mathematician expects. This means any existing Python code can be used in Sage. However, one must still obey the standard Python rules when writing packages that will be imported into Sage. (To install a Python library, for example that you have found on the Internet, follow the directions, but run sage -python instead of python. Very often this means typing sage -python setup.py install.)
94
Chapter 8. Afterword
in your bibliography (replacing 4.3 with the version of Sage you used). Moreover, please attempt to track down what components of Sage are used for your computation, e.g., PARI?, GAP?, Singular? Maxima? and also cite those systems. If you are in doubt about what software your computation uses, feel free to ask on the sage-devel Google group. See Univariate Polynomials for further discussion of this point.
If you happen to have just read straight through this tutorial, and have some sense of how long it took you, please let us know on the sage-devel Google group. Have fun with Sage!
95
96
Chapter 8. Afterword
CHAPTER
NINE
APPENDIX
9.1 Arithmetical binary operator precedence
What is 3^2*4 + 2%5? The value (38) is determined by this operator precedence table. The table below is based on the table in 5.14 of the Python Language Reference Manual by G. Rossum and F. Drake. the operations are listed here in increasing order of precedence. Operators or and not in, not in is, is not >, <=, >, >=, ==, != +, *, /, % **, ^ Description boolean or boolean and boolean not membership identity test comparison addition, subtraction multiplication, division, remainder exponentiation
Therefore, to compute 3^2*4 + 2%5, Sage brackets the computation this way: ((3^2)*4) + (2%5). Thus, rst compute 3^2, which is 9, then compute both (3^2)*4 and 2%5, and nally add these.
97
98
Chapter 9. Appendix
CHAPTER
TEN
BIBLIOGRAPHY
99
100
CHAPTER
ELEVEN
101
102
BIBLIOGRAPHY
[Cyt] Cython, http://www.cython.org. [Dive] Dive into Python, Freely available online at http://diveintopython.org. [GAP] The GAP Group, GAP - Groups, Algorithms, and Programming, Version 4.4; 2005, http://www.gapsystem.org [GAPkg] GAP Packages, http://www.gap-system.org/Packages/packages.html [GP] PARI/GP http://pari.math.u-bordeaux.fr/. [Ip] The IPython shell http://ipython.scipy.org. [Jmol] Jmol: an open-source Java viewer for chemical structures in 3D http://www.jmol.org/. [Mag] Magma http://magma.maths.usyd.edu.au/magma/. [Max] Maxima http://maxima.sf.net/ [NagleEtAl2004] Nagle, Saff, and Snider. Fundamentals of Differential Equations. 6th edition, Addison-Wesley, 2004. [Py] The Python language http://www.python.org/ Reference Manual http://docs.python.org/ref/ref.html. [PyDev] Guido, Some Guys, and http://www.python.org/dev/dev_intro.html. a Mailing List: How Python is Developed,
[Pyr] Pyrex, http://www.cosc.canterbury.ac.nz/~greg/python/Pyrex/. [PyT] The Python Tutorial http://www.python.org/. [SA] Sage web site http://www.sagemath.org/. [Si] G.-M. Greuel, G. Pster, and H. Schnemann. Singular 3.0. A Computer Algebra System for Polynomial Computations. Center for Computer Algebra, University of Kaiserslautern (2005). http://www.singular.uni-kl.de. [SJ] William Stein, David Joyner, Sage: System for Algebra and Geometry Experimentation, Comm. Computer Algebra {39}(2005)61-64.
103
104
Bibliography
INDEX
E
EDITOR, 56 environment variable EDITOR, 56
105