Up One Level SAGE Tutorial Contents Index
Previous: Up: SAGE Tutorial Next:
Subsections

2. A Guided Tour

This section is a guided tour of some of what is available in SAGE 1.0. For more examples, see the SAGE documentation ``SAGE constructions'', which is intended to answer the general question ``How do I construct ...?''.


2.1 Basic, and not-so-basic, Rings

We illustrate some basic rings in SAGE. For example, the field $ \mathbf{Q}$ of rational numbers may be referred to using either RationalField() or QQ: 

sage: RationalField()
Rational Field
sage: QQ
Rational Field
sage: 1/2 in QQ
True
The decimal number $ 1.2$ is considered in $ \mathbf{Q}$ , since there is a coercion map from the reals to the rationals: 
sage: 1.2 in QQ
True
However, the following doesn't, since there is no coercion: 
sage: I = ComplexField().0
sage: I in QQ
False
Also, of course, the symbolic constant $ \pi$ is not in $ \mathbf{Q}$ : 
sage: pi in QQ
False

If you use QQ as a variable, you can still fetch the rational numbers using the command RationalField(). By the way, some other pre-defined SAGE rings include the integers ZZ, the real numbers RR, the complex numbers CC (which uses I (or i), as usual, for the square root of $ -1$ ). We discuss polynomial rings in Section 2.2.

Do not redefine Integer or RealNumber unless you really know what you are doing. They are used by the SAGE interpreter to wrap integer and real literals. For example, if you type Integer = int, then integer literals will behave as they usually do in Python, so e.g., 4/3 evaluates to the Python int 1. For example 

sage: 4/3
4/3
sage: parent(_)
Rational Field
sage: prev = Integer
sage: Integer = int
sage: 4/3
1
sage: parent(_)
<type 'int'>
sage: Integer = prev
sage: 4/3
4/3

Now we illustrate some arithmetic involving various numbers. 

sage: a, b = 4/3, 2/3
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.76666666666666661
sage: 0.1 + 2/3       # coercion rules are symmetric in SAGE
0.76666666666666661
sage: z = a + b*I
sage: z
1.3333333333333333 + 0.66666666666666663*I
sage: z.real() == a        # automatic coercion before comparision
True
sage: QQ(11.1)
111/10

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: a = 5
sage: type(a)
<type 'integer.Integer'>
sage: a = 5/3
sage: type(a)
<type 'rational.Rational'>
sage: a = 'hello'
sage: type(a)
<type '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.

The field of p-adic numbers is implements as well: 

sage: K = Qp(11); K.prec(10)
sage: a = K(211/17); a
4 + 4*11 + 11^2 + 7*11^3 + 9*11^5 + 5*11^6 + 4*11^7 + 8*11^8 + 7*11^9 +
O(11^10)
sage: a.denominator()
1
sage: b = K(3211/11^2); b
10*11^-2 + 5*11^-1 + 4 + 2*11 + O(11^Infinity)
sage: b.denominator()
121

Rings of integers in p-adic fields or number fields other than $ \mathbf{Q}$ have not yet been implemented. However,a number of related methods are already implemented in the NumberField class. 

sage: x = PolynomialRing(QQ).gen()
sage: K = NumberField(x^3 + x^2 - 2*x + 8, 'a')
sage: K.integral_basis()
[1, a, 1/2*a^2 + 1/2*a]
sage: K.galois_group()                # requires optional GAP database package
Transitive group number 2 of degree 3
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()
Abelian group on 0 generators () with invariants []
sage: K.class_number()
1


2.2 Polynomials

In this section we illustrate how to create and use polynomials in SAGE.


2.2.1 Univariate Polynomials

There are three ways to create polynomial rings. 

sage: R = PolynomialRing(QQ, 'x')
sage: R
Univariate Polynomial Ring in x over Rational Field

An alternate way is 

sage: S = QQ['x']
sage: S == R
True

A third very convenient way is 

sage: R.<x> = PolynomialRing(QQ)
or 
sage: R.<x> = QQ['x']
This has the additional side effect that it defines the variable $ x$ to be the indeterminate of the polynomial ring. (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.)

The indeterminate of the polynomial ring is the 0 th generator: 

sage: R = PolynomialRing(QQ, 'x')
sage: x = R.0
sage: x in R
True

Alternatively, you can obtain both the ring and its generator, or just the generator during ring creation as follows: 

sage: R, x = QQ['x'].objgen()
sage: x    = QQ['x'].gen()
sage: R, x = objgen(QQ['x'])
sage: x    = gen(QQ['x'])

Finally we do some arithmetic in $ \mathbf{Q}[x]$ . 

sage: R, x = QQ['x'].objgen()
sage: f = 2*x^7 + 3*x^2 - 15/19
sage: f^2
4*x^14 + 12*x^9 - 60/19*x^7 + 9*x^4 - 90/19*x^2 + 225/361
sage: cyclo = R.cyclotomic_polynomial(7); cyclo
x^6 + x^5 + x^4 + x^3 + x^2 + x + 1
sage: g = 7 * cyclo * x^5 * (x^5 + 10*x + 2)
sage: g
7*x^16 + 7*x^15 + 7*x^14 + 7*x^13 + 77*x^12 + 91*x^11 + 91*x^10 + 84*x^9 
       + 84*x^8 + 84*x^7 + 84*x^6 + 14*x^5
sage: F = factor(g); F
(7) * x^5 * (x^5 + 10*x + 2) * (x^6 + x^5 + x^4 + x^3 + x^2 + x + 1)
sage: F.unit()
7
sage: list(F)
[(x, 5), (x^5 + 10*x + 2, 1), (x^6 + x^5 + x^4 + x^3 + x^2 + x + 1, 1)]
Notice that that the factorization correctly takes into account and records the unit part, unlike many other programs (e.g., PARI, Magma).

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 find 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, you'll 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.

Dividing two polynomials constructs an element of the fraction field. 

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 field of QQ[x]: 

sage: R = LaurentSeriesRing(QQ, 'x'); R
      Laurent Series Ring in x over Rational Field
sage: x = R.gen()
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)

If we name the variable differently, we obtain a different univariate polynomial ring. 

sage: R.<x> = PolynomialRing(QQ)
sage: S.<y> = PolynomialRing(QQ)
sage: x == y
False
sage: R == S
False
sage: R(y)
x
sage: R(y^2 - 17)
x^2 - 17

The ring is determined by the variable. Note that making another ring with variable called x does return a different ring. 

sage: R = PolynomialRing(QQ, "x")
sage: T = PolynomialRing(QQ, "x")
sage: R == T
True      
sage: R is T
False
sage: R.0 == T.0
True

SAGE also has support for power series and Laurent series rings over any base ring. In the following example we create an element of $ \mathbf{F}_7[[T]]$ and divide to create an element of $ \mathbf{F}_7((T))$ . 

sage: R = PowerSeriesRing(GF(7), 'T'); R
Power Series Ring in T over Finite Field of size 7
sage: T = R.0
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.2.2 Multivariate Polynomials

To work with polynomials of several variables, we declare the polynomial ring and variables first, in one of two ways. 

sage: R = MPolynomialRing(GF(5),3,"z")
sage: R
Polynomial Ring in z0, z1, z2 over Finite Field of size 5

Just as for univariate polynomials, there is an alternative more compact notation: 

sage: GF(5)['z0, z1, z2']
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: MPolynomialRing(GF(5), 3, 'xyz')
Polynomial Ring in x, y, z over Finite Field of size 5

Next lets do some arithmetic. 

sage: z = GF(5)['z0, z1, z2'].gens()
sage: z
(z0, z1, z2)
sage: (z[0]+z[1]+z[2])^2
z2^2 + 2*z1*z2 + z1^2 + 2*z0*z2 + 2*z0*z1 + z0^2
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 the Python dictionaries and the ``distributive representation'' of a polynomial. SAGE makes some use of Singular [Si], e.g., for computation of gcd's and Gröbner basis of ideals. 

sage: R, (x, y) = PolynomialRing(RationalField(), 2, 'xy').objgens()
sage: f = (x^3 + 2*y^2*x)^2
sage: g = x^2*y^2
sage: f.gcd(g)
x^2

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^2*y^2, 4*x^2*y^4 + 4*x^4*y^2 + x^6) of Polynomial Ring in x, y over Rational Field
sage: B = I.groebner_basis(); B
[x^2*y^2, x^6]
sage: x^2 in I
False

Incidentally, the Groebner basis above is not just 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 Groebner basis). 

sage: B.parent()
Category of sequences in Polynomial Ring in x, y over Rational Field
sage: B.universe()
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$ : 

ssage: I.primary_decomposition()
[Ideal (x^2) of Polynomial Ring in x, y over Rational Field,
 Ideal (y^2, 4*x^2*y^4 + 4*x^4*y^2 + x^6) of Polynomial Ring in x, y over Rational Field]
sage: I.associated_primes()
[Ideal (x) of Polynomial Ring in x, y over Rational Field,
 Ideal (y, x) of Polynomial Ring in x, y over Rational Field]

2.3 Number Theory

SAGE has extensive functionality for number theory. For example, we can do arithmetic in $ \mathbf{Z}/N\mathbf{Z}$ as follows: 

sage: R = IntegerModRing(97)
sage: a = R(2) / R(3)
sage: a
33
sage: a.rational_reconstruction()
2/3
sage: b = R(47)
sage: b^20052005
50
sage: b.modulus()
97
sage: b.is_square()
True

SAGE contains standard number theoretic functions. For example, 

sage: gcd(515,2005)
5
sage: factor(2005)
5 * 401
sage: c = factorial(25); c
15511210043330985984000000
sage: [valuation(c,p) for p in prime_range(2,23)]
[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
Perfect!

SAGE's sigma(n,k) function adds up the $ k$ th powers of the divisors of $ n$ (note the order of $ n$ and $ k$ !): 

sage: sigma(28,0); sigma(28,1); sigma(28,2)
6
56
1050

We next illustrate the extended Euclidean algorithm, the Euler's $ \phi$ -function, and the Chinese remainder theorem: 

sage: d,u,v = xgcd(12,15)
sage: d == u*12 + v*15
True
sage: inverse_mod(3,2005)
1337
sage: 3 * 1337
4011
sage: n = 2005
sage: prime_divisors(n)
[5, 401]
sage: phi = n*prod([1 - 1/p for p in prime_divisors(n)]); phi
1600
sage: euler_phi(2005)
1600
sage: prime_to_m_part(n, 5)
401

We next verify something about the $ 3n+1$ problem. 

sage: n = 2005
sage: for i in range(1000):
          n = 3*odd_part(n) + 1
          if odd_part(n)==1:
              print i
              break
38

Finally we illustrate the Chinese remainder theorem. 

sage: x = crt(2, 1, 3, 5); x   
-4
sage: x % 3  # x mod 3 = 2
2
sage: x % 5  # x mod 5 = 1
1
sage: [binomial(13,m) for m in range(14)]
[1, 13, 78, 286, 715, 1287, 1716, 1716, 1287, 715, 286, 78, 13, 1]
sage: [binomial(13,m)%2 for m in range(14)]
[1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1]
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: list(partitions(4))
[(1, 1, 1, 1), (1, 1, 2), (2, 2), (1, 3), (4,)]


2.3.1 Dirichlet Characters

A Dirichlet character is the extension of a homomorphism

$\displaystyle (\mathbf{Z}/N\mathbf{Z})^* \to R^*, $

for some ring $ R$ , to the map $ \mathbf{Z}/N\mathbf{Z}\to R$ obtained by sending those $ x\in\mathbf{Z}/N\mathbf{Z}$ with $ \gcd(N,x)>1$ to 0 . 

sage: G = DirichletGroup(21)
sage: list(G)
[[1, 1], [-1, 1], [1, zeta6], [-1, zeta6], [1, zeta6 - 1], 
 [-1, zeta6 - 1], [1, -1], [-1, -1], [1, -zeta6], [-1, -zeta6], 
 [1, -zeta6 + 1], [-1, -zeta6 + 1]]
sage: G.gens()
([-1, 1], [1, zeta6])
sage: len(G)
12

Having created the group, we next create an element and compute with it. 

sage: chi = G.1; chi
[1, 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

It is also possible to compute the action of the Galois group $ {\rm Gal}(\mathbf{Q}(\zeta_n)/\mathbf{Q})$ on these characters, as well as the direct product decomposition corresponding to the factorization of the modulus. 

sage: G.galois_orbits()
    [
    [[1, 1]],
    [[-1, 1]],
    [[1, zeta6], [1, -zeta6 + 1]],
    [[-1, zeta6], [-1, -zeta6 + 1]],
    [[1, zeta6 - 1], [1, -zeta6]],
    [[-1, zeta6 - 1], [-1, -zeta6]],
    [[1, -1]],
    [[-1, -1]]
    ]
sage: G.decomposition()
    [
    Group of Dirichlet characters of modulus 3 over Cyclotomic Field of order 6 and degree 2,
    Group of Dirichlet characters of modulus 7 over Cyclotomic Field of order 6 and degree 2
    ]

Next, we construct the group of Dirichlet character mod 20, but with values in $ \mathbf{Q}(i)$ : 

sage: G = DirichletGroup(20)
sage: G.list()
[[1, 1], [-1, 1], [1, zeta4], [-1, zeta4], [1, -1], 
    [-1, -1], [1, -zeta4], [-1, -zeta4]]

We next compute several invariants of G: 

sage: G.gens()
([-1, 1], [1, zeta4])
sage: G.unit_gens()
[11, 17]
sage: G.zeta()
zeta4
sage: G.zeta_order()
4

In this example we create a Dirichlet character with values in a number field. We explicitly specify the choice of root of unity by the third argument to DirichletGroup below. 

sage: x = PolynomialRing(QQ).gen()
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: G.list()
[[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 defining polynomial $ x^4 + 1$ ''). The name ``a'' is undeclared at this point. Once a = K.0 (or equivalently a = K.gen()) is typed, the symbol ``a'' represents a root of the generating polynomial, $ x^4 + 1$ .

2.4 Linear Algebra

SAGE provides standard linear algebra commands, e.g., characteristic polynomial, echelon form, trace, decomposition, etc., of a matrix.

We create the space Mat$ _{3\times 3}(\mathbf{Q})$ : 

sage: M = MatrixSpace(QQ,3)
sage: M
Full MatrixSpace of 3 by 3 dense matrices over Rational Field

The space of matrices has a basis: 

sage: B = M.basis()
sage: len(B)
9
sage: B[1]
[0 1 0]
[0 0 0]
[0 0 0]

We create a matrix as an element of M. 

sage: A = M(range(9)); A
[0 1 2]
[3 4 5]
[6 7 8]

Next we compute its reduced row echelon form and kernel. 

sage: A.echelon_form()
[ 1  0 -1]
[ 0  1  2]
[ 0  0  0]
sage: A^20
[ 2466392619654627540480  3181394780427730516992 3896396941200833493504] 
[ 7571070245559489518592  9765907978125369019392 11960745710691248520192]
[12675747871464351496704 16350421175823007521792 20025094480181663546880]
sage: A.kernel()
Vector space of degree 3 and dimension 1 over Rational Field
Basis matrix:
[ 1 -2  1]

Eigenvalues and eigenvectors over $ \mathbf{Q}$ or $ \mathbf{R}$ can be computed using maxima (see section 4.4 below).

Next we illustrate computation of matrices defined over finite fields: 

sage: M = MatrixSpace(GF(2),4,8)
sage: A = M([1,1,0,0,1,1,1,1,0,1,0,0,1,0,1,1,0,0,1,0,1,1,0,1,0,0,1,1,1,1,1,0])
sage: A
[1 1 0 0 1 1 1 1]
[0 1 0 0 1 0 1 1]
[0 0 1 0 1 1 0 1]
[0 0 1 1 1 1 1 0]
sage: rows = A.rows()
sage: A.columns()
[(1, 0, 0, 0), (1, 1, 0, 0), (0, 0, 1, 1), (0, 0, 0, 1), 
 (1, 1, 1, 1), (1, 0, 1, 1), (1, 1, 0, 1), (1, 1, 1, 0)]
sage: rows
[(1, 1, 0, 0, 1, 1, 1, 1), (0, 1, 0, 0, 1, 0, 1, 1), 
 (0, 0, 1, 0, 1, 1, 0, 1), (0, 0, 1, 1, 1, 1, 1, 0)]

We make the subspace over $ \mathbf{F}_2$ spanned by the above rows. 

sage: V = VectorSpace(GF(2),8)
sage: S = V.subspace(rows)
sage: S
Vector space of degree 8 and dimension 4 over Finite Field of size 2
Basis matrix:
[1 0 0 0 0 1 0 0]
[0 1 0 0 1 0 1 1]
[0 0 1 0 1 1 0 1]
[0 0 0 1 0 0 1 1]
sage: A.echelon_form()
[1 0 0 0 0 1 0 0]
[0 1 0 0 1 0 1 1]
[0 0 1 0 1 1 0 1]
[0 0 0 1 0 0 1 1]

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$ .

2.4.1 Sparse Linear Algebra

SAGE has support for sparse linear algebra over PID's. 

sage: M = MatrixSpace(QQ, 100, sparse=True)
sage: A = M.random_element(prob = 0.05)
sage: E = A.echelon_form()

The multi-modular algorithm in SAGE is good for square matrices (but not so good for non-square matrices): 

sage: M = MatrixSpace(QQ, 50, 100, sparse=True)
sage: A = M.random_element(prob = 0.05)
sage: E = A.echelon_form()                  
sage: M = MatrixSpace(GF(2), 20, 40, sparse=True)
sage: A = M.random_element()
sage: E = A.echelon_form()

Note that Python is case sensitive: 

sage: M = MatrixSpace(QQ, 10,10, Sparse=True)
Traceback (most recent call last):
...
TypeError: MatrixSpace() got an unexpected keyword argument 'Sparse'

SAGE can compute eigenvalues and eigenvectors: 

sage: MS = MatrixSpace(GF(7),2,2)
sage: g = MS([[5, 1], [4, 1]])
sage: g.eigenvalues()
[4, 2]
sage: g.eigenvectors()
[(1, 5), (1, 1)]


2.4.2 Numerical Linear Algebra

SAGE includes Numeric, which is a standard Python package for numerical linear algebra. If you have the appropriate numerical libraries installed on your computer when you built SAGE, then Numeric will use them for highly optimized matrix algorithms. To use Numeric, type import Numeric and proceed as described in the Numeric documentation (type help(Numeric)). Also, if $ A$ is a SAGE matrix, you can obtain the corresponding Numeric array as follows. 

sage: import Numeric
sage: A = Matrix(QQ,3,3,range(9))
sage: N = A.numeric_array(); N
    [[ 0., 1., 2.,]
     [ 3., 4., 5.,]
     [ 6., 7., 8.,]]
sage: Numeric.matrixmultiply(N,N)
    [[  15.,  18.,  21.,]
     [  42.,  54.,  66.,]
     [  69.,  90., 111.,]]
sage: A*A
[ 15  18  21]
[ 42  54  66]
[ 69  90 111]


2.5 Finite Groups

SAGE has some support for computing with permutation groups, finte classical grous (such as $ SU(n,q)$ ), finite matrix groups (with your own generators), and abelian groups (even infinite ones). Much of this is implemented using the interface to GAP.

For example, to create a permutation group, give a list of generators, as in the following example. 

sage: G = PermutationGroup(['(1,2,3)(4,5)', '(3,4)'])
sage: G
Permutation Group with generators [(1,2,3)(4,5), (3,4)]
sage: G.order()
120
sage: G.is_abelian()
False
sage: G.derived_series()           # random-ish output
[Permutation Group with generators [(1,2,3)(4,5), (3,4)],
 Permutation Group with generators [(1,5)(3,4), (1,5)(2,4), (1,3,5)]]
sage: G.center()
Permutation Group with generators [()]
sage: G.random()
(1,5,3)(2,4)
sage: print latex(G)
\langle (1,2,3)(4,5), (3,4) \rangle

You can also obtain the character table (in a semi-legal latex format even) 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&-zeta3 - 1&zeta3\\
1&1&zeta3&-zeta3 - 1\\
3&-1&0&0
\end{array}\right)

Though this latex array needs some editing, it is useful when trying to create the table for a technical paper.

Also implemented are classical and matrix groups over finite fields: 

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], [0 1]
    [6 3], [0 1]
    [6 2], [0 1]
    [6 6], [0 1]
    [6 4], [0 1]
    [6 5], [0 3]
    [2 2], [0 3]
    [2 5], [0 1]
    [6 0], [6 0]
    [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()
[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 (infinite and finite): 

sage: A=AbelianGroup(5,[3, 5, 5, 7, 8], names="abcde")
sage: a,b,c,d,e=A.gens()
sage: b1 = a^3*b*c*d^2*e^5
sage: b2 = a^2*b*c^2*d^3*e^3
sage: b3 = a^7*b^3*c^5*d^4*e^4
sage: b4 = a^3*b^2*c^2*d^3*e^5
sage: b5 = a^2*b^4*c^2*d^4*e^5
sage: e.word_problem([b1,b2,b3,b4,b5],display=False)
[[b^2*c^2*d^3*e^5, 245]]
sage: (b^2*c^2*d^3*e^5)^245
e
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


2.6 Elliptic Curves

Elliptic curves functionality includes most of the elliptic curve functionality of PARI, access to the data in Cremona's online tables (requires optional database package), the functionality of mwrank, i.e., $ 2$ -descents with computation of the full Mordell-Weil group, the SEA algorithm, computation of all isogenies, much new code for curves over $ \mathbf{Q}$ , and some of Denis Simon's algebraic descent software.

The command EllipticCurve for creating an elliptic curve has many forms:

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  = x^3 + x + 2 over Rational Field

sage: EllipticCurve('37a')
Elliptic Curve defined by y^2 + y = x^3 - x over Rational Field

sage: EllipticCurve(5)
Elliptic Curve defined by y^2 + x*y  = x^3 + 36/1723*x + 1/1723 
            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$ defined by $ y^2 + y = x^3 - 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 infinity 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)
sage: E.conductor()
37

The elliptic curves over the complex numbers are parameterized by the $ j$ -invariant. SAGE computes $ j$ -invariants as follows: 

sage: E = EllipticCurve([0,0,1,-1,0]); E
Elliptic Curve defined by y^2 + y = x^3 - x over Rational Field
sage: E.j_invariant()
110592/37
If we make a curve with $ j$ -invariant the same 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(110592/37)
sage: factor(F.conductor())
2^6 * 37

However, the twist of $ F$ by $ 2$ gives an isomorphic curve. 

sage: G = F.quadratic_twist(2); G
Elliptic Curve defined by y^2 + y = x^3 - x over Rational Field
sage: G.conductor()
37
sage: G.j_invariant()
110592/37

We can compute the coefficients $ a_n$ of the $ L$ -series or modular form $ \sum_{n=0}^\infty a_nq^n$ attached to the elliptic curve. This computation uses the PARI C-library: 

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)

It only takes a second to compute all $ a_n$ for $ n\leq 10^5$ : 

sage: time v = E.anlist(100000)
CPU times: user 0.98 s, sys: 0.06 s, total: 1.04 s
Wall time: 1.06

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. 

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 fields, and that's it. There is a small version of Cremona's database, which is distributed by default with SAGE, and contains limited information about elliptic curves of conductor $ \leq 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.

2.7 Plotting

The ``Constructions'' SAGE documentation has some examples of using SAGE for plotting, as do sections 2.8.4 and 4.4 below. We shall give some other examples here of using matplotlib. To view any one of these, after entering the commands below for the picture you want, type p.save("<path>/my_plot.png") and view the plot in a graphics viewer such as GIMP.

Here's a yellow circle: 

sage: L = [[cos(pi*i/100),sin(pi*i/100)] for i in range(200)]
sage: p = polygon(L, rgbcolor=(1,1,0))

A green deltoid: 

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)]
sage: p = polygon(L, rgbcolor=(1/8,3/4,1/2))

A blue figure 8: 

sage: L = [[2*cos(pi*i/100)*sqrt(1-sin(pi*i/100)^2),2*sin(pi*i/100)*sqrt(1-sin(pi*i/100)^2)] for i in range(200)]
sage: p = polygon(L, rgbcolor=(1/8,1/4,1/2))

A blue hypotrochoid: 

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)]
sage: p = polygon(L, rgbcolor=(1/8,1/4,1/2))

A purple epicycloid: 

sage: m = 9; b = 1
sage: L = [[m*cos(pi*i/100)+b*cos((m/b)*pi*i/100),m*sin(pi*i/100)-b*sin((m/b)*pi*i/100)] for i in range(200)]
sage: p = polygon(L, rgbcolor=(7/8,1/4,3/4))

A blue 8-leaved petal: 

sage: L = [[sin(5*pi*i/100)^2*cos(pi*i/100)^3,sin(5*pi*i/100)^2*sin(pi*i/100)] for i in range(200)]
sage: p = polygon(L, rgbcolor=(1/3,1/2,3/5))

You can also add text to a plot: 

L = [[cos(pi*i/100)^3,sin(pi*i/100)] for i in range(200)]
p = line(L, rgbcolor=(1/4,1/8,3/4))
t = text("a bulb", (-1.7, 0.5))
x = text("x axis", (2,-0.2))
y = text("y axis", (0.6,1.3))
g = p+t+x+y
view(g, xmin=-1.5, xmax=2, ymin=-1, ymax=1.3)

2.8 Calculus

The ``Constructions'' SAGE documentation has some examples of using SAGE for calculus computations, such as integration, differentiation, and Laplace transforms. In this chapter, we present a few other examples.

2.8.1 Functions

SAGE allows one to construct piecewise-defined functions. To define

\begin{displaymath} f(x) = \left\{ \begin{array}{ll} 1,& 0<x<1,\\ 1-x, & 1<x<2,\\ 2x, & 2<x<3,\\ 10x-x^2, 3<x<10, \end{array}\right. \end{displaymath}

type 

sage: x = PolynomialRing(RationalField()).gen()
sage: f1 = x^0
sage: f2 = 1-x
sage: f3 = 2*x
sage: f4 = 10*x-x^2
sage: f = Piecewise([[(0,1),f1],[(1,2),f2],[(2,3),f3],[(3,10),f4]])
sage: f
Piecewise defined function with 4 parts, [[(0, 1), 1], [(1, 2), -x + 1], [(2, 3), 2*x], [(3, 10), -x^2 + 10*x]]

By convention, we assume this takes the average value of the jumps at each of the inner midpoints.

To compute critical points and function values, 

sage: f.critical_points()
[5.0]
sage: f(5)
      25
sage: f(1/2)
      1
sage: f(1)
      1/2
sage: f(0)
      1
sage: f(10)
      0
Several other methods are available for these functions, such as laplace transforms, Fourier series, and plotting (not only of the function itself but even the plotting of the partial sums of Fourier series of the function). See the section on ``piecewise functions'' in the SAGE reference manual for examples.

2.8.2 Elementary functions

Call a univariate function an "elementary function" if it can be written as a sum of functions of the form "polynomial times and exponential times a sine or a cosine".

The set E of elementary functions is an algebra over RR. If D is differentiation and A = RR[D] is the polynomial ring in D over RR, let us define a smooth function f to be finite if the vector space A(f) is finite dimensional.

Theorem: E is the algebra of all finite functions. 

sage: R = ElementaryFunctionRing(QQ,"t"); R
ElementaryFunctionRing over Rational Field in t
sage: t = R.polygen(); t
t
sage: f = exponential(2,t); f
Elementary function (1)exp(2*t)
sage: f.diff()
Elementary function (2)exp(2*t)
sage: f.int([])
Elementary function (1/2)exp(2*t)
sage: f.latex()
'(1)e^{2t}\\cos(0t)'
sage: f(1)
7.3890560989306504
sage: f.laplace_transform("s")
'1/(s - 2)'
sage: f^2
Elementary function (1)exp(4*t)

## The example below shows how to solve x'' - x = sin(2t):

sage: DR = PolynomialRing(QQ,"D")
sage: D = DR.gen()
sage: Phi = D^2 - 1
sage: R = ElementaryFunctionRing(QQ,"t")
sage: t = R.polygen()
sage: g = ElementaryFunction([(1*t^0,0,0,2)])
sage: g.desolve(Phi,"x")
"x(t) = e^t*(5*x'(0)) + 5*x(0) + 2)/10 - e^-t*(5*x'(0)) - 5*x(0) + 2)/10 - sin(2*t)/5"

2.8.3 Differentiation, integration, etc

To compute $ \frac{d^4\sin(x^2)}{dx^4}$ : 

sage: maxima('sin(x^2)').diff('x',4)             
16*x^4*sin(x^2) - 12*sin(x^2) - 48*x^2*cos(x^2)

To compute, $ \frac{\partial (x^2+17y^2)}{\partial x}$ , $ \frac{\partial (x^2+17y^2)}{\partial y}$ : 

sage: f = maxima('x^2 + 17*y^2')                 
sage: f.diff('x')
2*x
sage: f.diff('y')                                
34*y

To compute $ \int x\sin(x^2) dx$ , $ \int_0^1 \frac{x}{x^2+1} dx$ : 

sage: maxima('x*sin(x^2)').integrate('x')
      -cos(x^2)/2
sage: maxima('x/(x^2+1)').integral('x', 0, 1)
      log(2)/2

To compute the partial fraction decomposition of $ \frac{1}{x^2-1}$ : 

sage: f = maxima('1/((1+x)*(x-1))')            
sage: f.partial_fraction_decomposition('x')    
1/(2*(x - 1)) - 1/(2*(x + 1))
sage: f.partial_fraction_decomposition('x').display2d() 
                     1           1
                 --------- - ---------
                 2 (x - 1)   2 (x + 1)


2.8.4 Systems of DEs using Laplace transforms

In this section, we provide a few details which are useful to teaching a lower level ordinary differential equatoins course using SAGE.

The displacement from equilibrium (respectively) for a coupled spring attached to a wall on the left 

|------\/\/\/\/\---|mass1|----\/\/\/\/\/----|mass2|
         spring1               spring2
is modeled by the system of 2nd order ODEs

$\displaystyle m_1x_1'' + (k_1+k_2)x_1 - k_2x_2 = 0, \ m_2x_2''+ k_2(x_2-x_1) = 0, $

where $ x_1$ denotes the displacement from equilibrium of mass 1, denoted $ m_1$ , $ x_2$ denotes the displacement from equilibrium of mass 2, denoted $ m_2$ , and $ k_1$ , $ k_2$ are the respective spring constants.

Example: Use SAGE to solve the above problem with $ m_1 = 2$ , $ m_2 = 1$ , $ k_1 = 4$ , $ k_2 = 2$ , $ x_1(0) = 3$ , $ x_1'(0) = 0$ , $ x_2(0) = 3$ , $ x_2'(0) = 0$ .

Soln: Take Laplace transforms of the first DE (for simplicity of notation, let $ x = x_1$ , $ y = x_2$ ): 

sage: _ = maxima.eval("x2(t) := diff(x(t),t, 2)")
sage: maxima("laplace(2*x2(t)+6*x(t)-2*y(t),t,s)")
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 says $ -2x_1'(0) + 2s^2*X_1(s) - 2sx_1(0) - 2X_2(s) + 2X_1(s) = 0$ (where the Laplace transform of a lower case function is the upper case function). Take Laplace transforms of the second DE: 

sage: _ = maxima.eval("y2(t) := diff(y(t),t, 2)")
sage: maxima("laplace(y2(t)+2*y(t)-2*x(t),t,s)")
-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
This says $ s^2X_2(s) + 2X_2(s) - 2X_1(s) - 3s = 0$ . Solve these two equations: 

sage: eqns = ["(2*s^2+6)*X-2*Y=6*s","-2*X +(s^2+2)*Y = 3*s"] 
sage: vars = ["X","Y"]
sage: maxima.solve_linear(eqns, vars)
[X = (3*s^3 + 9*s)/(s^4 + 5*s^2 + 4),Y = (3*s^3 + 15*s)/(s^4 + 5*s^2 + 4)]
This says $ X_1(s) = (3s^3 + 9s)/(s^4 + 5s^2 + 4)$ , $ X_2(s) = (3s^3 + 15s)/(s^4 + 5s^2 + 4)$ . Take inverse Laplace transforms to get the answer: 

sage: maxima("ilt((3*s^3 + 9*s)/(s^4 + 5*s^2 + 4),s,t)")
cos(2*t) + 2*cos(t)
sage: maxima("ilt((3*s^3 + 15*s)/(s^4 + 5*s^2 + 4),s,t)")
4*cos(t) - cos(2*t)
Therefore, $ x_1(t) = \cos(2t) + 2\cos(t)$ , $ x_2(t) = 4\cos(t) - \cos(2t)$ . This can be plotted parametrically using 

maxima.plot2d_parametric(["cos(2*t) + 2*cos(t)","4*cos(t) - cos(2*t)"], "t",[0,1])
and individually using 

maxima.plot2d('cos(2*x) + 2*cos(x)','[x,0,1]')
maxima.plot2d('4*cos(x) - cos(2*x)','[x,0,1]')

REFERENCES: Nagle, Saff, Snider, Fundamentals of DEs, 6th ed, Addison-Wesley, 2004. (see §5.5).

2.8.5 Euler's method for systems of DEs

Finally, we show how the files in the examples/calculus subdirectory of the main SAGE_HOME directory can be used. (And please feel free to contribute your own by emailing them to the SAGE Forum or to William Stein.)

In the next example, we will illustrate Euler's method for $ 2nd$ order ODE's. We first recall the basic idea. The goal is to find an approximate solution to the problem

$\displaystyle y'=f(x,y), y(a)=c,$ (2.1)

where $ f(x,y)$ is some given function. We shall try to approximate the value of the solution at $ x=b$ , where $ b>a$ is given.

Recall from the definition of the derivative that

$\displaystyle y'(x)\cong \frac{y(x+h)-y(x)}{h}, $

$ h>0$ is a given and small. This and the DE together give $ f(x,y(x))\cong \frac{y(x+h)-y(x)}{h}$ . Now solve for $ y(x+h)$ :

$\displaystyle y(x+h)\cong y(x)+h\cdot f(x,y(x)). $

If we call $ h\cdot f(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

$\displaystyle y_{new}=y_{old}+h\cdot f(x,y_{old}). $

Tabular idea: Let $ n>0$ be an integer, which we call the step size. This is related to the increment by

$\displaystyle h=\frac{b-a}{n}. $

This can be expressed simplest using a table.

$ x$ $ y$ $ hf(x,y)$
$ a$ $ c$ $ hf(a,c)$
$ a+h$ $ c+hf(a,c)$ &vellip#vdots;
$ a+2h$ &vellip#vdots;  
&vellip#vdots;    
$ b$ ??? xxx

The goal is to fill out all the blanks of the table but the xxx entry and find the ??? entry, which is the Euler's 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 Euler's method, where $ z''+tz'+z=0$ , $ z(0)=1$ , $ z'(0)=0$ .

First, you must attach the appropriate file, so we type 

sage: attach os.environ['SAGE_ROOT'] + '/examples/calculus/eulers_method.sage'

Now one must reduce the 2nd order ODE doewn to a system of two first order DEs (using $ x=z$ , $ y=z'$ ) and apply Euler's method: 

sage: t,y1,y2 = PolynomialRing(RealField(10),3).gens()
sage: f = y2; g = -y1 - y2 * t
sage: eulers_method_2x2(f,g, 0, 1, 0, 1/4, 1)
         t                    x                h*f(t,x,y)                    y           h*g(t,x,y)
         0                    1                   0.00000                    0             -0.25000
       1/4               1.0000                 -0.062500             -0.25000             -0.23438
       1/2              0.93750                  -0.11719             -0.46875             -0.17578
       3/4              0.82031                  -0.15381             -0.61523            -0.089722
         1              0.66602                  -0.16650             -0.66602              0.00000
Therefore, $ z(1)\approx 0.75$ .

Other files useful for illustrating basic aspects of differential euqations os.environ['SAGE_ROOT'] '/examples/calculus/desolvers.sage'+ and os.environ['SAGE_ROOT'] '/examples/calculus/field_plot2d.sage'+. They have docstrings with many examples of the syntax.

2.8.6 Special functions

Several orthogonal polynomials and special functions are implemented, using both pari and maxima. These are documented in the appropriate section (``Orthogonal polynomials'') of the reference manual. 

sage: x = PolynomialRing(QQ, 'x').gen()
sage: chebyshev_U(2,x)
4*x^2 - 1

The special functions in SAGE are also documented in the appropriate section (``Special functions'') of the reference manual. 

sage: bessel_I(1,1,"pari",500)
0.56515910399248502720769602760986330732889962162109200948029448947925564096437113409266499776681441006467788605552630267685763768491717981204113120812118
sage: bessel_I(1,1)
0.56515910399248503
sage: bessel_I(2,1.1,"maxima")  # last few digits are random
0.16708949925104899

Maxima has fixed accuracy, whereas those functions implemented using pari have higher accuracy.

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) - v*bessel_y(v,w)/w'

2.9 Algebraic Geometry

You can define arbitrary algebraic varieties in SAGE, but sometimes nontrivial functionality is limited to rings over $ \mathbf{Q}$ or finite fields. For example, we compute the union of two affine plane curves, then recover the curves as the irreducible components of the union. 

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 1 - y^2 - y^3 + y^5 - x^2 
                 + x^2*y^3 - x^3 + x^3*y^2 + x^5
sage: D.irreducible_components()
[Closed subscheme of Affine Space of dimension 2 over Rational Field defined by:
  -1 + y^2 + x^2, 
Closed subscheme of Affine Space of dimension 2 over Rational Field defined by:
  -1 + y^3 + x^3]

We can also find 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:
      -1 + y
      x, Closed subscheme of Affine Space of dimension 2 over Rational Field defined by:
      2 + y + x
      3 + 4*y + 2*y^2]
Thus, e.g., $ (1,0)$ and $ (0,1)$ are on both curves (visibly clear), as are certain (quadratic) points whose $ y$ coordinate 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 (c^2 - b*d, -1*b*c + a*d, b^2 - a*c) of Polynomial Ring in a, b, c, d over Rational Field
sage: F.reduced_groebner_bases ()
 [[-1*c^2 + b*d, -1*b*c + a*d, -1*b^2 + a*c],
  [c^2 - b*d, -1*b*c + a*d, -1*b^2 + a*c],
  [c^2 - b*d, b*c - a*d, -1*b^2 + a*c, -1*b^3 + a^2*d],
  [c^2 - b*d, b*c - a*d, b^3 - a^2*d, -1*b^2 + a*c],
  [c^2 - b*d, b*c - a*d, b^2 - a*c],
  [-1*c^2 + b*d, b^2 - a*c, -1*b*c + a*d],
  [-1*c^2 + b*d, b*c - a*d, b^2 - a*c, -1*c^3 + a*d^2],
  [c^3 - a*d^2, -1*c^2 + b*d, b*c - a*d, b^2 - a*c]]
sage: F.fvector()
 (1, 8, 8)


2.10 Modular Forms

SAGE can do some computations related to modular forms, including dimensions, computing spaces of modular symbols, Hecke operators, and decompositions.

There are several functions available for computing dimensions of spaces of modular forms. For example, 

sage: dimension_cusp_forms(Gamma0(11),2)
1
sage: dimension_cusp_forms(Gamma0(1),12)
1
sage: dimension_cusp_forms(Gamma1(389),2)
6112

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(); f
x^3 - 2001*x^2 - 97776*x - 1180224
sage: factor(f)
(x - 2049) * (x + 24)^2
sage: M.T(11).charpoly().factor()
(x - 285311670612) * (x - 534612)^2

We can also create spaces for $ \Gamma_0(N)$ and $ \Gamma_1(N)$ . 

sage: ModularSymbols(11,2)
Modular Symbols space of dimension 3 for Gamma_0(11) of weight 2 with sign 0 over Rational Field
sage: ModularSymbols(Gamma1(11),2)
Modular Symbols space of dimension 11 for Gamma_1(11) of weight 2 with sign 0 and over Rational Field
Let's compute some characteristic polynomials and $ q$ -expansions. 
sage: M = ModularSymbols(Gamma1(11),2)
sage: M.T(2).charpoly()
x^11 - 8*x^10 + 20*x^9 + 10*x^8 - 145*x^7 + 229*x^6 + 58*x^5 
           - 360*x^4 + 70*x^3 - 515*x^2 + 1804*x - 1452
sage: M.T(2).charpoly().factor()
(x - 3) * (x + 2)^2 * (x^4 - 7*x^3 + 19*x^2 - 23*x + 11) 
        * (x^4 - 2*x^3 + 4*x^2 + 2*x + 11)
sage: S = M.cuspidal_submodule()
sage: S.T(2).matrix()
[-2  0]
[ 0 -2]
sage: S.q_expansion_basis(10)
[
    q - 2*q^2 - q^3 + 2*q^4 + q^5 + 2*q^6 - 2*q^7 - 2*q^9 + O(q^10)
]

We can even compute spaces of modular symbols with character. 

sage: G = DirichletGroup(13)
sage: e = G.0^2
sage: M = ModularSymbols(e,2); M
Modular Symbols space of dimension 4 and level 13, weight 2, character [zeta6], sign 0, over Cyclotomic Field of order 6 and degree 2
sage: M.T(2).charpoly().factor()
(x + -2*zeta6 - 1) * (x + -zeta6 - 2) * (x + zeta6 + 1)^2
sage: S = M.cuspidal_submodule(); S
Modular Symbols subspace of dimension 2 of Modular Symbols space of dimension 4 and level 13, weight 2, character [zeta6], sign 0, over Cyclotomic Field of order 6 and degree 2
sage: S.T(2).charpoly().factor()
(x + zeta6 + 1)^2
sage: S.q_expansion_basis(10)
[
q + (-zeta6 - 1)*q^2 + (2*zeta6 - 2)*q^3 + zeta6*q^4 + (-2*zeta6 + 1)*q^5 + (-2*zeta6 + 4)*q^6 + (2*zeta6 - 1)*q^8 + -zeta6*q^9 + O(q^10)
]

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: T.new_subspace()
 Modular Forms subspace of dimension 2 of Modular Forms space of dimension 2 for Congruence Subgroup Gamma0(11) of weight 2 over Rational Field
sage: M = ModularSymbols(11)
sage: 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 $ T_p$ denote the usual Hecke operators ($ p$ prime). How do the Hecke operators $ T_2$ , $ T_3$ , $ T_5$ 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]
[ 0 -1  0]
[ 0  0 -1]
sage: M.T(5).matrix()
[ 6  0 -1]
[ 0  1  0]
[ 0  0  1]

Up One Level SAGE Tutorial Contents Index
Previous: Up: SAGE Tutorial Next:
Release 1.4.1, documentation updated on October 18, 2006.
See About this document... for information on suggesting changes.