# MATH 3600

## Introduction to Computational Mathematics

### Daniel D. Warner, Modified by Neil J. Calkin

Watch the video on Spirals, Fibonacci, and Being a Plant.

Is Ms. Vihart doing Mathematics?

What is the 100th Fibonacci number . . . the 1,000th . . . the 1,000,000 . . . the $10^{100}$?

Can we use the computer to calculate these numbers?

Can a computer count to $10^{100}$?

Are Ms. Vihart's rectangles related to the Golden Rectangle that was so prized by the ancient Greek architects?

If so, how?

What about drawing those spirals, can we get the computer to draw them?

Understanding how to effectively use a computer as a tool in doing mathematics is the subject of Computational Mathematics.

The first step is to learn how to program, and for this class we will use today's most powerful "platform" for programming in mathematics.

# Introduction to Sage

## WWW.SAGEMATH.ORG

Sage is a platform that brings together many powerful open source software packages for doing computational mathematics, and we talk to these packages using the Python programming language.

Python is a high-level scripting language that is particularly well suited for supporting the exploratory nature of Computational Mathematics.

## We start by making sure that everyone can access Sage.

1. Access the Sage Server at  https://sage.math.clemson.edu:34567
2. Set up an account and click on the "Published Button."  All worksheets for this course will start with MATH 3600.
3. You can easily download a copy of Sage from WWW.SAGEMATH.ORG, and install it on your personal computer, but be sure that you have high speed connection and at least 3 GB of space available on your hard drive.
4. You can also create an account at sagemath.org and run Sage on their cloud servers.

## And now a short demo...

# This is a simple introduction to Sage. # # The Sage notebook uses blocks of text, which are essentially Python statements that get executed. # Lines that start with the pound sign are comments and do not get echoed when the block is evaluated. # # A block is evaluated by clicking on the evaluate link or by typing shift-enter. # # Pressing tab will give auto-complete options, or function descriptions after a parens ( # A '?' after a command gives detailed information, e.g. 'plot?' tells you about plot # Sage can do simple calculations 3+2*6^2
 75 75
# Sage can use variables googol = 10^100

Note the slightly different spelling here: explore why it is spelled that way.

googol
 100000000000000000000000000000000000000000000000000000000000000000000000\ 00000000000000000000000000000 10000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
# and more elaborate calculations factorial(100)
 933262154439441526816992388562667004907159682643816214685929638952175999\ 932299156089414639761565182862536979208272237582511852109168640000000000\ 00000000000000 93326215443944152681699238856266700490715968264381621468592963895217599993229915608941463976156518286253697920827223758251185210916864000000000000000000000000
# Sage can factor integers factor(798)
 2 * 3 * 7 * 19 2 * 3 * 7 * 19
# The same operations can be applied to algebraic expressions factor(x^5-1)
 (x - 1)*(x^4 + x^3 + x^2 + x + 1) (x - 1)*(x^4 + x^3 + x^2 + x + 1)
p = x^5-1 p
 x^5 - 1 x^5 - 1
p^3
 (x^5 - 1)^3 (x^5 - 1)^3
# Sage can plot these algebraic expressions P1 = plot( p , -0.5, 0.5, rgbcolor=(1,0,0) ) P2 = plot( p^3, -0.5, 0.5, rgbcolor=(0,1,0) ) show(P1+P2)
factor(p^3)
 (x - 1)^3*(x^4 + x^3 + x^2 + x + 1)^3 (x - 1)^3*(x^4 + x^3 + x^2 + x + 1)^3
# Sage handles rational functions with ease rf = x / (x^2 - 1) rf
 x/(x^2 - 1) x/(x^2 - 1)
P3 = plot( rf , -0.5, 0.5, rgbcolor=(1,0,0) ) show(P3)
# What are x, p, and rf? print type(rf) # and what can we do with them?
  
# Partial fraction expansions rf.partial_fraction()
 1/2/(x - 1) + 1/2/(x + 1) 1/2/(x - 1) + 1/2/(x + 1)
# Sage can display results in the more familiar mathematical format show(rf.partial_fraction())

## Solving Algebraic Equations

# A quadratic example eq1 = x^2 - 2 == x show(eq1) solve([eq1],x)
# Cubic polynomials are just as easy eq2 = x^3 + 2*x^2 - 5*x - 6 == 0 show(eq2) solve([eq2],x)
# Sage can even find the ugly answers eq3 = x^3 - 2 == x show(eq3) sol = solve([eq3],x) show(sol)
for s in sol: show(s)
# A linear systems of equations var('x,y') eq4 = x-1 == y+1 eq5 = x+1 == 2*(y-1) eq4;eq5
sol = solve([eq4,eq5],x,y) sol
# A nonlinear system of equations var('x y') eq6 = (x-1)^2 + (y-1)^2 == 4; show(eq6) eq7 = 4*y - 3*x == -1; show(eq7) sol3 = solve([eq6,eq7],x,y,solution_dict=True) show(sol3)

## Now back to Fibonacci...

1 + 1 = 2

1 + 2 = 3

2 + 3 = 5

3 + 5 = 8

5 + 8 = 13

and so on . . .

# One approach # Use 3 variables and shift just like the example a = 1 b = 1 c = b + a print a,'+',b,'=',c
# Do it again and again a = b b = c c = b + a print a,'+',b,'=',c a = b b = c c = b + a print a,'+',b,'=',c
# Use a loop -- but make it stop a = 1 b = 1 c = b + a while c < 200: print a,'+',b,'=',c a = b b = c c = b + a

# Suppose we want to generate the first 25 Fibonacci numbers # We can put in a counter a = 1 b = 1 c = b + a count = 3 while count <= 25: print a,'+',b,'=',c a = b b = c c = b + a count = count + 1
# Same code but with some shorthand notation a = 1 b = 1 count = 3 while count <= 25: print a,'+',b,'=',(a+b) a, b = b, (a+b) count += 1
range?
range(2,25)
a = 1 b = 1 for count in range(2,25): c = b + a print a,'+',b,'=',c a = b b = c

Print the 1,000th Fibonacci number, but not all the others.

a = 1 b = 1 for count in range(2,1000): c = b + a a,b = b,c print a,'+',b,'=',c

How long do you think it would take to compute the 1,000,000th Fibonacci number?

Give it a shot!

How long do you think it would take to compute the $10^{100}$th Fibonacci number?

## Arrays and Lists

An array is an ordered homogeneous collection.  A list is an ordered inhomogeneous collection.  In low level languages like C when an array is set up it has a fixed size.  In higher level languages like Matlab, arrays are dynamic.  That is, they can grow as needed.  In Python lists are both inhomogeneous and dynamic.  They generalize dynamic arrays.

We will start by considering some array operations that are common in most languages, and then explore the additional features provided in Python.

# The range command gives us an array of integers a = range(15)
len(a)
a[8] = 2*a[7] a[9] = 3*a[7] a
# Generate an array of the first 25 Fibonacci numbers a = [0 for k in range(25)] print a a[0] = 1 a[1] = 1 for k in range(2,25): a[k] = a[k-1] + a[k-2] print a
# What can we learn about a? type(a)
a = [1,1] for k in range(2,25): a.append(a[k-1]+a[k-2]) a
# Entries in a list can be any object that Sage can work with. a = [0, -15, sqrt(2), 3/5, 'John Doe', [1,2,3]] a
for k in range(len(a)): print a[k], type(a[k])