math_salad.sagews

A sampler of Sage mathematics

Sage offers tools for learning and for doing research in many different areas of mathematics. The purpose of this module is to provide a very brief glimpse into some these areas, including:

1. Multivariate calculus
2. Linear algebra
3. Ordinary differential equations
4. Abstract algebra/ Groups

Functions of 2 or more variables

There are several ways to define them. The easiest and most intuitive way is shown in the example below. it defines a function, evaluates it at a point, and takes its partial derivatives.

Example 1: Let $f(x,y,z)=(x^2+y^2)z$. Evaluate the function at some point in its domain, and find its 1st partial derivatives.

# A function of 3 variables: 
f(x,y,z) = (x^2 + y^2) * z

# Can evaluate the function at any point:
f(1,0,2)

# There are several ways to get partial derivatives:
f.diff(x)
f.derivative(x)
derivative(f,x)
derivative(f,z)
# Note that the answer is returned in function notation, 
# with inputs to the left and outputs to the right of the arrow.

# It is easy to compute the gradient of f:
f.gradient()

2 (x, y, z) |--> 2*x*z (x, y, z) |--> 2*x*z (x, y, z) |--> 2*x*z (x, y, z) |--> x^2 + y^2 (x, y, z) |--> (2*x*z, 2*y*z, x^2 + y^2)

Vector functions of 2 or more variables

When the dependent variable of a function is also a vector, we have the more general case. Here are some examples

${\bf f}(x,y) = \left[\begin{array}{c} x^2 + y^2 \\ 3 x y \end{array}\right ]$,       ${\bf g}(x,y,z) = \left[\begin{array}{c} x^2 y \\ x^2 z \\ y z \end{array}\right ]$

Such functions offer more interesting and fun variations for applying familiar calculus concepts! For example, taking their derivatives can be approached in the form of a "divergence" or a "curl," with very different conceptual meaning and results.

The Sage segment below defines each of these functions, evaluates them at a point, and computes their divergence and curl.

# Here's how to define vector functions of multiple variables. 
# Notice, it is an intuitive generalization of the previous syntax.
#
x, y, z = var('x, y, z')
f = vector( (x^2 + y^2, 3*x*y))
g = vector( (x^2 * y, x^2 * z, y*z) )

# You can evaluate the functions at any specific point by 
# similarly generalizing the single variable syntax.
# (I've used the print command just to show the output clearly):
print "f(0,2)=",  f(x=0, y=2)
print "g(1,2,3)=", g(x=1, y=2, z=3)

# Find divergence and curl. Note, the syntax requires 
# specifying the names of the variables in order.  (Again, 
# the print command is just to show the output clearly):
print "div(f)=", f.div([x,y])
print "curl(f)=", f.curl([x,y])
print "div(g)=", g.div([x,y,z])
print "curl(g)=", g.curl([x,y,z])

f(0,2)= (4, 0) g(1,2,3)= (2, 3, 6) div(f)= 5*x curl(f)= y div(g)= 2*x*y + y curl(g)= (-x^2 + z, 0, -x^2 + 2*x*z)

Some linear algebra

Sage offers powerful tools for working with vectors and matrices. Here is an example showing how to define them, and how to do basic arithmetic operations:

# Example showing how to work with matrices & vectors
# Define matrix by rows in brackets [] separated by commas
# Define vector by giving its column entries in [] 
A = matrix([ [1, 4, 1], [0, 1, 3], [0, -2, -5] ])
x = vector([0, 1, 2])
y = vector([6, 7, -12])

# Add the two vectors:
print "x+y =", x+y

# Multiply the matrix with vector "x":
print "A*x =", A*x

x+y = (6, 8, -10) A*x = (6, 7, -12)

Here is how to solve a matrix-vector system:

# To solve the matrix-vector system: A*x = y
answer = A\y
print "x =", answer

x = (0, 1, 2)

To find the eigenvalues and eigenvectors of a matrix:

# To find eigenvalues of matrix
A.eigenvalues()
A.eigenvectors_right()

[1, -3.732050807568878?, -0.2679491924311227?] [(1, [ (1, 0, 0) ], 1), (-3.732050807568878?, [(1, -1.953254218877944?, 3.080966067942896?)], 1), (-0.2679491924311227?, [(1, -0.3544380888143643?, 0.1498031628263346?)], 1)]

It is equally easy to transpose, invert, or raise a matrix to any power:

# To transpose a matrix:
transpose(A)
# To raise matrix to n_th power (e.g., n=6):
A^6
# To invert a matrix:
B = A^(-1)
print "A^(-1) =", B

[ 1 0 0] [ 4 1 -2] [ 1 3 -5] [ 1 510 1200] [ 0 -989 -2340] [ 0 1560 3691] A^(-1) = [ 1 18 11] [ 0 -5 -3] [ 0 2 1]

Ordinary differential equations

ODEs come in many flavors. Some have closed-form analytical solutions, while most can only be solved numerically. In addition, sometimes we have only single equations, and other times we have systems of ODEs. Sage provides tools for working with all these situations.

The following example shows how to find the general solution of the 1st order differential equation $x^\prime = 3/x$. Note that the key step is to rewrite the equation with $0$ on one side: $x^\prime - 3/x = 0$.

# Example showing how to solve the 1st ord. ODE: x' = 3/x
t = var('t')
x = function('x')(t)
DE = diff(x,t) - 3/x
desolve(DE, [x,t])

1/6*x(t)^2 == _C + t

It is also possible to solve ODEs that contain arbitrary parameters, as shown next:

# Example showing how to solve a 1st ord. ODE with parameters: y' = -ay(1-by)
t = var('t')
y = function('y')(t)
a, b = var('a  b')
DE = diff(y,t) + a*y*(1 - b*y)
desolve(DE, [y,t])

(log(b*y(t) - 1) - log(y(t)))/a == _C + t

Often, we want to plot direction fields for an ODE. Here is how to do it:

# Example showing how to plot direction fields for 1st ord. ODE: y' = 2y+4-t
t,y = var('t,y')
plot_vector_field((1,2*y+4-t), (t,0,3), (y,-4,0))

Sage can solve systems of linear 1st order ODEs in analytical form. Here are two examples -- one shows how to find the general solution, the other shows how to find a particular solution:

# Example showing how to find general solution of an ODE system.
# System is 1st order, linear: x'=2x+2y; y'=-x+4y

t = var('t')
x = function('x')(t)
y = function('y')(t)
de1 = diff(x,t) - 2*x - 2*y == 0
de2 = diff(y,t) + x - 4*y == 0
desolve_system ([de1, de2], [x, y])

[x(t) == -((x(0) - 2*y(0))*sin(t) - cos(t)*x(0))*e^(3*t), y(t) == -((x(0) - y(0))*sin(t) - cos(t)*y(0))*e^(3*t)]

# Example showing how to find particular solution.
# initial value problem is: x'=2x+2y; y'=-x+4y; x(0)=-2, y(0)=1

t = var('t')
x = function('x')(t)
y = function('y')(t)
de1 = diff(x,t) - 2*x - 2*y == 0
de2 = diff(y,t) + x - 4*y == 0
solution = desolve_system ([de1, de2], [x, y], ics=[0, -2, 1])
solution

[x(t) == -2*(cos(t) - 2*sin(t))*e^(3*t), y(t) == (cos(t) + 3*sin(t))*e^(3*t)]

It is also possible to solve ODE systems numerically, as shown below:

# Example showing how to solve the same system numerically
# Systems is: x'=2x+2y; y'=-x+4y; x(0)=-2, y(0)=1

#from sage.calculus.desolvers import desolve_system_rk4
x, y, t = var('x y t')
de1 = 2*x + 2*y
de2 = -x + 4*y
P = desolve_system_rk4 ([de1, de2], [x, y], ics=[0, -2, 1], ivar=t, end_points=[0,2] )
Q = [ [i,j] for i, j, k in P]  # Data for plot of x vs t
xp = line(Q, color='green', legend_label='$x$', legend_color='green')
Q = [ [i,k] for i, j, k in P]  # Data for plot of y vs t
yp = line(Q, axes_labels=['$t$', '$x,y$'], legend_label='$y$', legend_color='blue')
Q = [ [j,k] for i, j, k in P]  # Data for plot of y vs x
show(xp+yp)
line(Q, axes_labels=['$x$', '$y$'])

A taste of abstract algebra/ Groups

The group of all permutations of $n$ elements is a symmetric group. For example, the symmetric group of the set $\{1,2,3,4\}$ can be constructed using the command:

G = SymmetricGroup([1,2,3,4])
print G

Symmetric group of order 4! as a permutation group

It is now possible to find the number of elements in the group, and to list all of them:

G.cardinality()  # Number of elements in G.
G.list()         # List of elements in G.

24 [(), (1,2), (1,2,3,4), (1,3)(2,4), (1,3,4), (2,3,4), (1,4,3,2), (1,3,4,2), (1,3,2,4), (1,4,2,3), (1,2,4,3), (2,4,3), (1,4,3), (1,4)(2,3), (1,4,2), (1,3,2), (1,3), (3,4), (2,4), (1,4), (2,3), (1,2)(3,4), (1,2,3), (1,2,4)]

s = G.random_element(); s  # Pick a random element and print it.
t = G.random_element(); t  # Another random element.
s*t
t*s

(1,4,2,3) (1,4) (2,3,4) (1,2,3)