Sage Beginner's Guide: Unlock the full potential of Sage for simplifying and automating mathematical computing with this book and eBook

Much of the difficulty of higher mathematics actually lies in the extensive algebraic manipulations that are required to obtain a result. Sage can save you many hours, and many sheets of paper, by automating some tedious tasks in mathematics. We'll start with basic calculus. For example, let's compute the derivative of the following equation:

The following code defines the equation and computes the derivative:

var('x')
f(x) = (x^2 - 1) / (x^4 + 1)
show(f)
show(derivative(f, x))

The results will look like this:

The first line defines a symbolic variable x (Sage automatically assumes that x is always a symbolic variable, but we will define it in each example for clarity). We then defined a function as a quotient of polynomials. Taking the derivative of f(x) would normally require the use of the quotient rule, which can be very tedious to calculate. Sage computes the derivative effortlessly.

Now, we'll move on to integration, which can be one of the most daunting tasks in calculus. Let's compute the following indefinite integral symbolically:

The code to compute the integral is very simple:

f(x) = e^x * cos(x)
f_int(x) = integrate(f, x)
show(f_int)

The result is as follows:

To perform this integration by hand, integration by parts would have to be done twice, which could be quite time consuming. If we want to better understand the function we just defined, we can graph it with the following code:

f(x) = e^x * cos(x)
plot(f, (x, -2, 8))

Sage will produce the following plot:

Sage can also compute definite integrals symbolically:

To compute a definite integral, we simply have to tell Sage the limits of integration:

f(x) = sqrt(1 - x^2)
f_integral = integrate(f, (x, 0, 1))
show(f_integral)

The result is:

This would have required the use of a substitution if computed by hand.

There is actually a clever way to evaluate the integral from the previous problem without doing any calculus. If it isn't immediately apparent, plot the function f(x) from 0 to 1 and see if you recognize it. Note that the aspect ratio of the plot may not be square.

The partial fraction decomposition is another technique that Sage can do a lot faster than you. The solution to the following example covers two full pages in a calculus textbook —assuming that you don't make any mistakes in the algebra!

f(x) = (3 * x^4 + 4 * x^3 + 16 * x^2 + 20 * x + 9) / ((x + 2) * (x^2 + 3)^2)
g(x) = f.partial_fraction(x)
show(g)

The result is as follows:

We'll use partial fractions again when we talk about solving ordinary differential equations symbolically.

Linear algebra is one of the most fundamental tasks in numerical computing. Sage has many facilities for performing linear algebra, both numerical and symbolic. One fundamental operation is solving a system of linear equations:

Although this is a tedious problem to solve by hand, it only requires a few lines of code in Sage:

A = Matrix(QQ, [[0, -1, -1, 1], [1, 1, 1, 1], [2, 4, 1, -2],
    [3, 1, -2, 2]])
B = vector([0, 6, -1, 3])
A.solve_right(B)

The answer is as follows:

Notice that Sage provided an exact answer with integer values. When we created matrix A, the argument QQ specified that the matrix was to contain rational values. Therefore, the result contains only rational values (which all happen to be integers for this problem). Chapter 5 describes in detail how to do linear algebra with Sage.

Solving ordinary differential equations by hand can be time consuming. Although many differential equations can be handled with standard techniques such as the Laplace transform, other equations require special methods of solution. For example, let's try to solve the following equation:

The following code will solve the equation:

var('x, y, v')
y=function('y', x)
assume(v, 'integer')
f = desolve(x^2 * diff(y,x,2) + x*diff(y,x) + (x^2 - v^2) * y == 0,
    y, ivar=x)
show(f)

The answer is defined in terms of Bessel functions:

It turns out that the equation we solved is known as Bessel's equation. This example illustrates that Sage knows about special functions, such as Bessel and Legendre functions. It also shows that you can use the assume function to tell Sage to make specific assumptions when solving problems. In Chapter 7, we will explore Sage's powerful symbolic capabilities.

Sage does not restrict you to making plots in two dimensions. To demonstrate the 3D capabilities of Sage, we will create a parametric plot of a mathematical surface known as the "figure 8" immersion of the Klein bottle. You will need to have Java enabled in your web browser to see the 3D plot.

var('u,v')
r = 2.0
f_x = (r + cos(u / 2) * sin(v) - sin(u / 2) 
    * sin(2 * v)) * cos(u)
f_y = (r + cos(u / 2) * sin(v) - sin(u / 2)
    * sin(2 * v)) * sin(u)
f_z = sin(u / 2) * sin(v) + cos(u / 2) * sin(2 * v)
parametric_plot3d([f_x, f_y, f_z], (u, 0, 2 * pi), 
    (v, 0, 2 * pi), color="red")

In the Sage notebook interface, the 3D plot is fully interactive. Clicking and dragging with the mouse over the image changes the viewpoint. The scroll wheel zooms in and out, and right-clicking on the image brings up a menu with further options.

Sage can be used in conjunction with the LaTeX typesetting system to create publication-quality typeset mathematical expressions. In fact, all of the mathematical expressions in this chapter were typeset using Sage and exported as graphics. Chapter 10 explains how to use LaTeX and Sage together.

We introduced some new concepts in this example. On the first line, the statement import numpy allows us to access functions and classes from a module called NumPy. NumPy is based upon a fast, efficient array class, which we will use to store our data. We created a NumPy array and hard-coded the data values for OD, and created another array to store values of concentration (in practice, we would read these values from a file) We then defined a Python function called standard_curve, which we will use to convert optical density values to concentrations. We used the find_fit function to fit the slope and intercept parameters to the experimental data points. Finally, we plotted the data points with the scatter_plot function and the plotted the fitted line with the plot function. Note that we had to use a function called zip to combine the two NumPy arrays into a single list of points before we could plot them with scatter_plot. We'll learn all about Python functions in Chapter 4, and Chapter 8 will explain more about fitting routines and other numerical methods in Sage.

We defined one NumPy array of sample times, and another NumPy array of optical density values. As in the previous example, these values could easily be read from a file. We used the standard_curve function and the fitted parameter values from the previous example to convert the optical density to concentration. We then plotted the data points using the scatter_plot function.

We defined another Python function called gompertz to model the growth of bacteria in the presence of limited resources. Based on the data plot from the previous example, we estimated values for the parameters of the model to use an initial guess for the fitting routine. We used the find_fit function again to fit the model to the experimental data, and displayed the fitted values. Finally, we plotted the fitted model and the experimental data on the same axes.