4420-Using Sage for Calculus.sagews

[removed]

Using Sage for Calculus

April 24, 2011

Contents

1 Algebra
 1.1 Numbers
 1.2 Variables and Expressions
 1.3 Substitution
 1.4 Solving Equations
2 Plotting
 2.1 Functions of one variable
 2.2 Functions of two variables
 2.3 Polar coordinates; parametric functions
3 Calculus
 3.1 Limits
 3.2 Differentiation
 3.3 Implicit Differentiation
 3.4 Integration
4 Getting Help

Introduction

Sage is a powerful mathematical calculator, often called a computer algebra system (CAS) or symbolic mathematics program, because unlike most calculators it can do algebra and calculus symbolically. For example, it can expand the product ${(x + 1)}^{30}$ into the normal representation of this polynomial, and it can factor the result back into the original form. It can compute the derivative or antiderivative of $3x∕\sqrt{{x}^{2 } + 2x + 47}$. It can also serve as a numerical calculator, but again unlike most calculators it can do integer and rational number arithmetic exactly, no matter how many digits are involved. If it is undesirable or impossible to express a result exactly, Sage can do numerical approximation. Finally, Sage can plot curves and surfaces in two and three dimensions.

To learn how to use some of Sage’s features, you should work through this tutorial step by step. As you go, you may alter the Sage commands and type new ones to experiment; you will learn better and faster if you do this. However, some commands depend on the results of previous commands, so before you continue you may want to restore the original calculation. A good alternative is to experiment in a new command window, called a cell in Sage. If you hover your mouse just a bit above an existing cell you will see a blue bar appear; click on this bar and a new cell appears. The blocks of plain text between command cells are themselves cells, and you can open new command cells above them in the same way.

Sage remembers values from one cell to the next, but if you don’t execute any commands for a while these values are forgotten. You may have to return to an earlier cell and re-execute it to refresh a value.

1 Algebra

1.1 Numbers

To evaluate the contents of a command cell, click in the box; the outline should turn blue. Then hold down the shift key and press return. This takes some getting used to; if you press just return, Sage will simply start a new line. When this happens, just backspace and then press shift-return. If you prefer, you can click on the word Evaluate below the box.

Try evaluating the cell below. Then hover above the box and open a new cell. Then hover between the original cell and the following text (“You may type more than one command…”) and open a new cell there as well. Try evaluating some simple arithmetic expressions using “$+$”, “$−$”, “$∗$”, and “$∕$”.

123456789*987654321

You may type more than one command in a window by placing semicolons between them.

2+3; 34*98-45

The ordinary precedence rules apply to mathematical operations: multiplication and division are done first, then addition and subtraction. Multiplication and division have the same precedence, as do addition and subtraction; a sequence of mixed multiplications and divisions is evaluated left to right, and likewise for additions and subtractions. Make sure you understand these:

24/3/2; 24/3*2; 24/4+1

You may of course use parentheses to alter the order of operations. Sometimes using extra parentheses is a good idea because it makes it easier to read and understand an expression.

Sage will do arithmetic with rational numbers without converting to decimal form:

2/5 + 1/3; 23/45-167/345*2/3

Sage knows standard mathematical notation, including two symbols for exponentiation. The underscore character refers to the previous result.

2**5; 3^27; factor(_)

You can force Sage to approximate values using decimal notation.

12345/98765; n(_)

Sage understands $\pi$, $e$, and $i$ (the square root of $− 1$).

3*pi; n(_); n(e); i^2

Sage’s normal output can be a bit hard to read when it is complicated, though it becomes easier with practice, and it can be useful for copying results and pasting them in to new commands. At the top of the worksheet is a checkbox labeled Typeset. Check the box and then evaluate the previous cell; the output will be more familiar. Alternately, you may use the jsmath function to typeset a single result. Uncheck the Typeset box at the top of the page and then execute the next box.

sqrt(2); jsmath(sqrt(2))

The function $\mathop{ exp}\nolimits (x)$ is an alternate form of ${e}^{x}$.

e^2; exp(2)

1.2 Variables and Expressions

Sage can manipulate expressions involving variables as well as numbers. Variables can be assigned values, either actual numeric values or entire expressions.

f=(x+5)^2; f; expand(f); factor(_)

For readability you can also put different commands on different lines; just press return without the shift at the end of the line. Note that only the last result will be displayed in this case, so generally you use this when the other lines produce no output.

f=(x+5)^2
expand(f)
factor(_)

f=(x+5)^2
expand(f); factor(_)

We normally think of $f = {(x + 5)}^{2}$ as a function, but it is not officially a function.

f(5)

So that gives the result you expect, but only after a warning. Here’s how to make $f$ an official function:

f(x) = (x+5)^2; f(5)

You should get in the habit of defining functions this way. You can see the value of $f$ by evaluating it:

f

This notation indicates that Sage thinks of $f$ as a function, not merely an expression. Compare to this:

g=(x+5)^2; g

If you have a mere expression that you would like to be a function, you have a couple of options. You can redefine it, but the expression you have might be the result of many calculations, not a simple definition. To turn the previous $g$ into a function you can do this:

h(x) = g; h

There’s no need to use a different name if you don’t need the original one:

g(x) = g; g

Here’s something a little more realistic:

expand((x+4)^2 * (2*x+1)); g(x) = _; g; g(3)

You may need to apply a function to a list of input values; here is something a bit easier than doing each one separately:

L = [1,2,3,10]; map(g,L)

In ordinary mathematical notation you are used to leaving out some multiplication symbols; Sage does not allow this:

2x

Type an asterisk between the 2 and the $x$ and evaluate that cell again.

We saw above how to expand and factor polynomials. If a polynomial has integer or rational coefficients, the “factor” command tries to factor it into polynomials with integer coefficients, perhaps times a rational number:

expand((x+2)*(x-1/3)); factor(_)

If Sage can’t do this, it leaves the polynomial alone, even if it has a “nice” factorization:

factor(x^2-2); factor(x^2-x-1)

1.3 Substitution

We use “variables” in a variety of ways, to stand for “true variables”, for functions, or for unspecified constants. We might say, for example, $f(x) = A{(x + 5)}^{2}$, using all of these. If we try this, there is already a problem:

f(x)=A*(x+5)^2

Although Sage understands that $x$ is a variable, we must tell it explicitly about others:

A=var("A"); f(x)=A*(x+5)^2

We can manipulate this as you would expect:

expand(f); factor(_)

Notice that strangely the result of expand is a function but the result of factor is not!

We may then want to assign the constant $A$ a value. We can do this temporarily like this:

f(A=2); f(3,A=2); f; A

So we see that the definition of $f$ has not changed, it still contains the symbol $A$. This happens even if we give $A$ a permanent value:

A=2; f; A

This is probably a bad idea, since it is at least confusing. If you would like a new function with $A$ replaced by 2 you can do something like this:

g(x) = f(A=2); g; g(5)

Let’s restore $A$ to its status as a variable with no value:

A; A=var("A"); A

You can subsitute expressions for variables, not just values.

g(x) = f(A=x^3-1); g; h(x) = f(2*x+3,A=x^3-1); h

Make sure you understand how $h$ is defined. This substitution method works quite generally on expressions that may not have names:

(A^2+2)(A=3*x+1)

1.4 Solving Equations

Sage will solve a limited class of equations exactly and a broader class approximately. The solve command will give exact answers if possible, including complex answers.

solve(x^2-x-10==10,x);

If you specify just a function, Sage assumes you want to set it equal to zero and solve.

solutions=solve(g(x),x); solutions

Even when Sage can solve the equation exactly, you may want to see numerical approximations. You can use the n function on each solution, or with a bit of programming do them all at once. First we generate the solutions in a different form:

solutions=solve(g(x),x,solution_dict=True); solutions;

Then we can extract just the numerical values:

newsolutions = [ s[x] for s in solutions]; newsolutions

Then use the map function to apply n to each number on the list:

map(n,newsolutions)

Or we could combine the last two steps:

[ n(s[x]) for s in solutions]

If Sage can’t solve an equation it will simply reprint the equation. You can approximate a solution, but you have to limit the search to a range, which you might discover by graphing the function or by substituting some values into the function. For example, using $f(x) = {x}^{5} + 4{x}^{3} − 3{x}^{2} + 10$, you might see that $f(0) = 10$ and $f(−2) = −66$, so there must be a root between $0$ and $− 2$.

solve(x^5+4*x^3-3*x^2+10,x); find_root(x^5+4*x^3-3*x^2+10,-2,0)

find_root(tan(x)==x+1,0,pi/2);

Sage will also solve systems of equations; note the square brackets forming a list of equations.

y=var("y"); solve([3*x+5*y==6,4*x-6*y==-8],x,y)

A bit more involved:

h1=3*x-6*y-5; h2=-4*x*y+7; solve([h1,h2],x,y); jsmath(_)

2 Plotting

2.1 Functions of one variable

Sage can do two and three dimensional graphing. There are many ways to adjust the plots that Sage produces; we’ll cover just the basics here.

It’s easy to plot ordinary one-variable functions.

plot(g,-5,5);

The vertical scale makes it hard to see what’s going on in that graph; you can limit the vertical scale for a better look, using “ymin” and “ymax”.

plot(g,-5,5,ymax=10);

Some functions have vertical asymptotes, but Sage will just connect the dots:

plot(1/(x-1),-5,5,ymin=-10,ymax=10);

You can get rid of the spurious bit in the middle:

plot(1/(x-1),-5,5,ymin=-10,ymax=10,exclude=[1]);

You can also put a list of excluded values in the brackets; try plotting $f(x) = 1∕({x}^{2} − 1)$ between $− 5$ and $5$.

If you right-click on a plot you can download it as a png image.

You can put multiple plots on the same axes:

plot([g,1/(x-1)],-5,5,ymin=-110,ymax=20,exclude=[1]);

It’s often better to do the plots separately but display them together:

p1=plot(g,-5,5,ymin=-110,ymax=20,color='red')
p2=plot(1/(x-1),-5,5,ymin=-110,ymax=20,exclude=[1])
p1+p2

Sage normally uses different scales on the two axes; you can force it to use the same scales to see the “true shape”.

(p1+p2).show(aspect_ratio=1)

Here’s a better example: the line doesn’t look like it has slope 1 in the first graph.

p1=plot(x+1,-5,5); p1; p1.show(aspect_ratio=1)

2.2 Functions of two variables

Sage will plot functions of two variables in 3D. You can rotate the resulting graph with the mouse: click and hold, then drag. You can change the appearance from a menu that appears when you right-click on the graph. Note especially the stereographic option—try the cross-eyed version. Click on the “Get Image” link for a view that you can download.

y=var("y"); r=sqrt(x^2+y^2); f=sin(r)/r
plot3d(f,(x,-10,10),(y,-10,10))

2.3 Polar coordinates; parametric functions

You can do plots in polar coordinates.

polar_plot(1+2*sin(x),-pi,pi);

Sage will do parametric plots in two or three dimensions.

parametric_plot3d([cos(x),sin(x),cos(6*x)],(x,0,2*pi))

Plots in polar coordinates are really special examples of parametric plots.

polar_plot(1+cos(x),0,2*pi)

parametric_plot([cos(x)*(1+cos(x)),sin(x)*(1+cos(x))],(x,0,2*pi))

You can also do parametric surface plots:

(u,v)=var("u v")
parametric_plot3d([v*cos(2*u),v*sin(2*u),u],(u,0,2*pi),(v,0,1))

You may notice that the edges of this surface appear more polygonal than round; you can change the number of points in either direction on the plot. The $u$ direction is the one we want to increase; the $v$ direction is along radial lines, so a very small value is fine.

(u,v)=var("u v")
parametric_plot3d([v*cos(2*u),v*sin(2*u),u],(u,0,2*pi),(v,0,1),plot_points=[200,5])

Here’s a fancier plot: a sphere inside a cylinder, and the cylinder is partially transparent.

(u,v)=var("u v")
s=parametric_plot3d([sin(v)*cos(u),sin(v)*sin(u),cos(v)],(u,0,2*pi),(v,0,pi),color='red',aspect_ratio=1)
c=parametric_plot3d([v,cos(u),sin(u)],(u,0,2*pi),(v,-1,2),opacity=0.6)
s+c

Note that Sage has a bit of trouble where the cylinder and the sphere meet. You can fix this by cheating a little: make the cylinder just slightly bigger.

(u,v)=var("u v")
s=parametric_plot3d([sin(v)*cos(u),sin(v)*sin(u),cos(v)],(u,0,2*pi),(v,0,pi),color='red',aspect_ratio=1)
c=parametric_plot3d([v,1.02*cos(u),1.02*sin(u)],(u,0,2*pi),(v,-1,2),opacity=0.6)
s+c

You can plot quite complicated surfaces using multiple parametric plots together. The intersection of three cylinders of the same size at right angles to one another forms a solid. Note that we have plotted each color as four surfaces when in fact one would do; try plotting just the red surface with a single plot, using the range $(u, 0, 2π)$ and see what happens; can you explain what has gone wrong? You might try adding the original cylinders to this plot, colored appropriately, and partially transparent. First we set up some useful functions.

(x,y,u,t,a,b)=var("x y u t a b")
minimum(a,b)=((a+b)-abs(b-a))/2
                                                                                                
                                                                                                
x=cos(u); y=sin(u)
v=t*minimum(sqrt(1-x^2),sqrt(1-y^2))

Now notice that each of the colors is essentially the same, with the coordinates plotted against different axes for the different colors.

redsurf1=parametric_plot3d([x,v,y],(u,0,pi/2),(t,-1,1),color='red',aspect_ratio=(1,1,1))
redsurf2=parametric_plot3d([x,v,y],(u,pi/2,pi),(t,-1,1),color='red',aspect_ratio=(1,1,1))
redsurf3=parametric_plot3d([x,v,y],(u,pi,3*pi/2),(t,-1,1),color='red',aspect_ratio=(1,1,1))
redsurf4=parametric_plot3d([x,v,y],(u,3*pi/2,2*pi),(t,-1,1),color='red',aspect_ratio=(1,1,1))
bluesurf1=parametric_plot3d([x,y,v],(u,0,pi/2),(t,-1,1),color='blue',aspect_ratio=(1,1,1))
bluesurf2=parametric_plot3d([x,y,v],(u,pi/2,pi),(t,-1,1),color='blue',aspect_ratio=(1,1,1))
bluesurf3=parametric_plot3d([x,y,v],(u,pi,3*pi/2),(t,-1,1),color='blue',aspect_ratio=(1,1,1))
bluesurf4=parametric_plot3d([x,y,v],(u,3*pi/2,2*pi),(t,-1,1),color='blue',aspect_ratio=(1,1,1))
yellowsurf1=parametric_plot3d([v,x,y],(u,0,pi/2),(t,-1,1),color='yellow',aspect_ratio=(1,1,1))
yellowsurf2=parametric_plot3d([v,x,y],(u,pi/2,pi),(t,-1,1),color='yellow',aspect_ratio=(1,1,1))
yellowsurf3=parametric_plot3d([v,x,y],(u,pi,3*pi/2),(t,-1,1),color='yellow',aspect_ratio=(1,1,1))
yellowsurf4=parametric_plot3d([v,x,y],(u,3*pi/2,2*pi),(t,-1,1),color='yellow',aspect_ratio=(1,1,1))
redsurf1+redsurf2+redsurf3+redsurf4+bluesurf1+bluesurf2+bluesurf3+bluesurf4+yellowsurf1+yellowsurf2+yellowsurf3+yellowsurf4

3 Calculus

3.1 Limits

Sage can compute many typical limits as encountered in calculus. Here’s an easy one, for a continuous function.

f(x)=x^2+4*x+30
limit(f,x=3); f(3)

And everyone’s favorite non-trivial limit:

g(x)=sin(x)/x
limit(g,x=0); limit(g,x=infinity)

You can specify one-sided limits in two ways.

g(x)=x/abs(x)
limit(g,x=0); limit(g,x=0,dir='-'); limit(g,x=0,dir='right')

3.2 Differentiation

Here’s a typical easy calculus problem: find a derivative, find critical values, find the $y$-coordinates for each.

f(x)=x^3+4*x^2+3*x+4
fp=diff(f)
critvals=solve(fp,x,solution_dict=True)
f; fp; critvals; jsmath([expand(f(s[x])) for s in critvals])

3.3 Implicit Differentiation

Sage has no built-in implicit differentiation command, but we can do it using differentiation, in essentially the way we do it by hand. We tell Sage that $y$ is a function of $x$, then proceed essentially as usual. Here is the standard example computing $y' = −x∕y$ for $y$ defined by ${x}^{2} + {y}^{2} = 1$. Note that we convert this to ${x}^{2} + {y}^{2} − 1 = 0$, which guarantees that the derivative of the right side is zero. (In this case the derivative of the right side, 1, was already zero, but in general it might be more complicated.)

y=function('y',x)
temp=diff(x^2+y^2-1)
solve (temp,diff(y))

3.4 Integration

The “integral” command computes an antiderivative. Here, we take the derivative to check the antiderivative.

f(x)=(2*x-3)*sqrt(x-3)
g(x)=integral(f,x)
g; diff(g); factor(_)

If Sage can’t handle an integral, it does nothing.

integral(tan(x^2),x)

We could of course do definite integration by substituting and subtracting, but the integral command will also do it.

integral(f(x),x,5,10); g(10)-g(5)

Sage can do some improper integrals.

f(x)=1/sqrt(x); integral(f(x),x,0,1); integral(f(x),x,1,infinity)

f(x)=1/x^2; integral(f(x),x,1,infinity)

The function ${e}^{−{x}^{2}}$ doesn’t have a “nice” antiderivative, though it can be expressed in terms of the “erf” function. For most values of $x$, erf(x) can only be approximated, but an improper integral turns out to be quite nice. (The function ${e}^{−{x}^{2}}$ is a Gaussian function, a “bell curve”.)

f(x)=exp(-x^2);
integral(f,x); integral(f(x),x,0,1); integral(f(x),x,-infinity,infinity)

Sometimes the result is mysterious.

f(x)=sin(x); integral(0,infinity)

Sage can’t find antiderivatives of some functions. In such cases you can approximate the value of a definite integral using the $n$ function.

integral(tan(x^2),x,0,1); n(_)

4 Getting Help

If you click on “Help” at the top of the Sage window, you will find a number of options. There is a tutorial that is longer and more complete than this introduction, and there is a comprehensive reference manual. If you know the command you want help with you can enter “commandname?” in a command box and execute it. Google can also find helpful pages; generally you should search for “sagemath” plus some description of the topic you’re interested in.

plot3d?