Sage Education Day 1 Intro

95 days ago by kcrisman

Sage: Free and Open-Source Mathematics Software

December 5th, Clay Math Institute - Sage Education Day 1

Speaker: Karl-Dieter Crisman, Gordon College

Thanks to everyone for coming!  I hope this talk will have something for everyone, but we will start fairly basic, and try to focus on Sage in the classroom from a beginner's viewpoint.  Some of you may have seen some of this before, but hopefully that's okay; repetition makes the master!

2+2; 2^2; 2*2 
       

Let's start.  We usually expect to be able to do calculations like these.  But what about calculations like these?

integral(3*x^2/sqrt(25*x^2-3)); factorial(21); primitive_root(991); find_root(cos(x)==sin(x),0,pi/2) 
       

Unless we had a specific reason for them doing these by hand, most of us would be happy to be able to use these things in some other problem.  And doing computations like this (symbolic and numerical) is what comprehensive mathematics software does - often, very well indeed.

But there are two potential barriers to our using them.  Since you are here, presumably you don't need to be convinced that Sage can help overcome these barriers.

You can even try it out on the first one of these - if you wish, create an account, click on "Published", and this talk should be there along with a few other things.  You can even ignore me for the rest of the talk.  Warning: if you all do it at once, you might have to wait a bit for it to load up.

Okay, let's get right into the mechanics - always an issue for computer programs.  You have two options once you type in a command - either press Shift and Enter together (Shift-Enter), or click on the "evaluate" link you'll see.

2+2 
       

It's as easy as that; type in the command and Shift-Enter.  This is very similar to other systems you may be familiar with, and makes it easy to experiment.  There are tutorials and references, both at the Sage website and strewn elsewhere.  The "Help" button at the top has some links to these.

Getting back to the worksheet; if I did not get 4 as the answer above, but instead got the following for my answer, I probably did something wrong:

\zeta(s)=\sum_{n=1}^{\infty}\frac{1}{n^s}=\prod_p \left(\frac{1}{1-p^{-s}}\right)
 

Now, I did this partly in honor of being at the Clay Math Institute!  But how does one get that to show up?  Let's explore that.  Move your mouse between two math cells (that's what we call these) until you see the little blue line.  Press shift, and click on the line simultaneously; now you have a WYSIWYG text editor called TinyMCE available - and it understands TeX through the magic of jsmath!  

So you can type "x" or "$x$", and the result should be different.  Note that depending on the browser, you may have to install some extra fonts for full effect, but it should still look better than nothing.  The above formula is "$$\zeta(s)=\sum_{n=1}^{\infty}\frac{1}{n^s}=\prod_p \left(\frac{1}{1-p^{-s}}\right)$$", by the way.

Not everyone will want to use this, but it is great for putting comments in homework or sharing your attempts at a problem with other people in a study group.  Incidentally, the "Share" and "Publish" buttons at the top allow very easy collaboration; I haven't used them as much for my own needs, but others have done so very successfully.

A = Matrix([[1,2,3],[3,2,1],[1,1,1]]) A 
       

Back to mathematics.  Here, our example is from linear algebra.  Notice how making two lines is easy, so we can separate commands.  

Above, in the first line we have assigned a specific matrix to the variable A, and the second line tells Sage to represent A for us - in this case, by nicely typesetting it, if I have clicked the "Typeset" button at the very tippy-top of the page. 

Matrix([[1,2,3],[3,2,1],[1,1,1]]) 
       

Flexibility like this is very nice to have.  What else can I do?

w = vector([1,1,-4]); w; w*A; A*w 
       

Note that it's important to assign most things to a name, in order to obviate the tedium of retyping everything.  Just don't reassign things like \pi:

pi; pi=3; pi 
       

Luckily, I can fix such a mistake.

reset('pi'); pi 
       

Sage has support for a wide variety of things, and often has a couple different ways to get the same result.  Here, we are solving the linear system

Y = vector([0,-4,-1]) X = A.solve_right(Y) A \ Y; X; A*X 
       

It will also tell me if there isn't a solution, which of course can happen!

w; A.solve_right(w) 
       

Plenty of other things are available, too.

rank(A); det(A) 
       

Okay, this is nice, but what about finding out what else Sage can do?  Here, there is a very nice set of help features.  What can we do with a matrix, for instance?  Type the dot/period after your object's name, and then [tab] to see what options you have:

A. 
       
A.echelon_form();A.determinant();A.rank() 
       

Naturally, the only options that come up are ones you can actually do to a matrix.  

Also note the "()" syntax - it betrays Sage's foundation in the Python programming language.  You don't need to know Python to use Sage, though it doesn't hurt to learn just a tad of it to be able to do cooler things.  

The other nice feature in Sage for help is access to the documentation of each function or object.  This is done simply by appending "?" to your request.

A.rook_vector 
       

Sage can do the standard things you expect math software to do:

var('x,y') plot3d(sin(pi*(x^2+y^2))/2,(x,-2,2),(y,-2,2)) 
       

Note that one can make most Sage things interactive (this one already is, but only moving it around) in a very simple way.  Here, plot_points is the attribute you think it is:

var('x,y') @interact def _(p=[10,50,100,200,500]): show(plot3d(sin(pi*(x^2+y^2))/2,(x,-2,2),(y,-2,2),plot_points=p)) 
       
var('x,y') @interact def _(p=slider(50,500,50,50)): show(plot3d(sin(pi*(x^2+y^2))/2,(x,-2,2),(y,-2,2),plot_points=p)) 
       

Here is another typical thing you might want to do.

y = var('y') P=plot_slope_field(2-y,(x,0,3),(y,0,20)) y = function('y',x) # declare y to be a function of x f = desolve(diff(y,x) + y - 2, y, ics=[0,7]) # solve the DE, with Initial ConditionS of 2,2 Q=plot(f,0,3) html('$f=%s$'%(latex(f),)); show(P+Q) # I can use HTML in my outcomes 
       

One annoying, but crucial, point in the above is that you must use the syntax "var('y')" to define all variables other than "x".  In the long run, this is good, because it helps avoid doing silly things like the reassignment I did above, but in the short run it is annoying, I apologize.  

The only other such thing I can think of off the top of my head is that Sage doesn't like it if you pretend a symbolic expression is "callable" without explicitly saying so.  In one variable this seems pointless, though I suppose f(3) could be interpreted as f\cdot 3:

f=x^2; f(3) 
       
f(x)=x^2; f(3) 
       

But it makes lots of sense to require this in more than one variable! 

g=x^2*y h=y*x^2 g(2,1); h(2,1) 
       
g(x,y)=x^2*y h(y,x)=y*x^2 g(2,1);h(2,1) 
       

Going back to the example, I can make it interactive too, though in order to avoid having 'y' mean two different things, I have to do a little more work.

y = var('y') @interact def _(g=input_box(default = 1-y),a=2,b=2,auto_update=False): P=plot_slope_field(g,(x,0,3),(y,0,20)) yfun = function('yfun',x) f = desolve(diff(yfun,x) - g(y=yfun), yfun, ics=[a,b]) Q=plot(f,0,3) show(P+Q) 
       

There are several ways to access statistics.  The one with by far the most undergraduate resources is the R implementation of the S language; it is a standard.   One has to download a few extra things to do this particular example, though not just to use R; for numpy or Maxima stats, they are completely built in.  I have to turn the typesetting off to use this properly, since jsmath doesn't understand it.  Hopefully someone here today will fix both of these issues :-)

r.library("MASS") r.data('Cars93') r.mean('Cars93$MPG.highway') 
       

One can also actually evaluate R commands directly:

%r mean(Cars93$MPG.city) mean(Cars93$MPG.highway) xtabs( ~ Origin + Type, data=Cars93) 
       

There is eye candy too.  Some is not so serious:

plot([x^i for i in [1..10]],(x,0,1),fill = dict((i,[i+1]) for i in [0..9])) 
       

Other stuff can be.  I think William was the first person to do this example, and I know nothing about image processing, but I love it.  This picture is off the CMI website.

import pylab A_image = pylab.mean(pylab.imread(DATA + 'CMI.png'), 2) @interact def svd_image(i=(5,(1..30)),display_axes=True): u,s,v = pylab.linalg.svd(A_image) A = sum(s[j]*pylab.outer(u[0:,j],v[j,0:]) for j in range(i)) g = graphics_array([matrix_plot(A),matrix_plot(A_image)]) show(g,axes=display_axes, figsize=(8,3)) html('<h2>Compressed using %s eigenvalues</h2>'%i) 
       

Since I like teaching number theory, I'll include this example.  How many theorems can you see in this data?

@interact def power_table_plot(p=(7,prime_range(50))): P=matrix_plot(matrix(p-1,[mod(a,p)^b for a in range(1,p) for b in srange(p)]),cmap='jet') show(P) 
       

Dynamical systems/calculus labs are possible too.  Notice that I've hidden all the relevant programming I needed to do (and/or crib off of the Sage interact Wiki):

%auto %hide def my_eulers_method_plot(f,x0,y0,h,x1): n=int((1.0)*(x1-x0)/h) x00=x0; y00=y0 x01=x0; y01=y0 P=point((x00,y00),rgbcolor=hue(1)) # red Q=Graphics() # default is blue if f.number_of_arguments()==1: try: f = f._fast_float_(f.args()[0]) except AttributeError: pass for i in range(n+1): y01 = y00+h*f(x00) x01 = x00+h P=P+point((x01,y01),rgbcolor=hue(1)) Q=Q+line([(x00,y00),(x01,y01)]) x00=x01 y00=y01 if f.number_of_arguments()==2: try: f = f._fast_float_(f.args()[0],f.args()[1]) except AttributeError: pass for i in range(n+1): y01 = y00+h*f(x00,y00) x01 = x00+h P=P+point((x01,y01),rgbcolor=hue(1)) Q=Q+line([(x00,y00),(x01,y01)]) x00=x01 y00=y01 return(P+Q) var('t,x,y') def euler_logistic_plot(parameter,y_start,end=15,step=1): function(x,y)=parameter*y*(1-y) my_eulers_method_plot(function,0,y_start,step,end).show(ymin=0) 
       
@interact def _(start=slider(0,1,.05,.4),parameter=RR(2),end=slider(10,100,1,14),step=slider(.1,1,.1,1)): print 'starting population percent = %.1f%%' %float(100*start) print 'parameter =', parameter, '(Warning: do not let this get above 3)' print 'length of time =', end+1 print 'step size = %.1f' %float(step) euler_logistic_plot(parameter,start,end,step) 
       

Well, I could go on, but I really just wanted to give you a taste of how Sage is and isn't like similar systems; I assume you have your own experience and needs.  And naturally, there are some things that are easier to do in one system than another, and some that aren't possible in any of them (solving the Riemann Hypothesis, for instance).  

Finally, I do want to point out is the variety of pedagogical techniques one can use with this (or any other) software - but from the convenience of students' own rooms.  For instance:

Thank you!  Please feel free to ask questions, though we may delay some to talks or discussions later today.  And don't forget to check out the website:

http://www.sagemath.org