A distinguishing property of *programming* languages is that the programmer can create their own *functions*. Creating a *function* is like teaching the computer a new trick. Typically a function will receive some data as *input*, will perform an *algorithm* involving the input data, and will *output* data when the algorithm terminates.

In this part, we explore Python functions. We also explore control statements, which allow a program to behave in different ways for different inputs. We also introduce the *while loop*, a loop whose repetition can be more carefully controlled than a for loop. As an application of these techniques, we implement the Euclidean algorithm as a Python function in a few ways, to effectively find the GCD of integers and solve linear Diophantine equations. This complements Chapter 1 of An Illustrated Theory of Numbers.

A *function* in Python is a construction which takes input data, performs some actions, and outputs data. It is best to start with a few examples and break down the code. Here is a function `square`

. Run the code as usual by pressing *shift-Enter* when the code block is selected.

In [ ]:

```
def square(x):
answer = x * x
return answer
```

When you run the code block, you probably didn't see anything happen. But you have effectively taught your computer a new trick, increasing the vocabulary of commands it understands through the Python interpreter. Now you can use the `square`

command as you wish.

In [ ]:

```
square(12)
```

In [ ]:

```
square(1.5)
```

Let's break down the syntax of the *function declaration*, line by line.

```
def square(x):
answer = x * x
return answer
```

The first line begins with the Python reserved word `def`

. (So don't use `def`

as a variable name!). The word `def`

stands for "define" and it defines a function called `square`

. After the function name `square`

comes parentheses `(x)`

containing the **argument** `x`

. The *arguments* or *parameters* of a function refer to the input data. Even if your function has no arguments, you need parentheses. The argument `x`

is used to name whatever number is input into the `square`

function.

At the end of the function declaration line is a colon `:`

and the following two lines are indented. As in the case of for loops, the colon and indentation are signals of *scope*. Everything on the indented lines is considered the *scope of the function* and is carried out when the function is used later.

The second line `answer = x * x`

is the beginning of the scope of the function. It declares a variable `answer`

and sets the value to be `x * x`

. So if the argument `x`

is 12, then `answer`

will be set to 144. The variable `answer`

, being declared within the scope of the function, will not be accessible outside the scope of the function. It is called a **local variable**.

The last line `return answer`

contains the Python reserved word `return`

, which terminates the function and outputs the value of the variable `answer`

. So when you apply the function with the command `square(1.5)`

, the number `1.5`

is the argument `x`

, and `answer`

is `2.25`

, and that number `2.25`

becomes the output.

A function does not have to return a value. Some functions might just provide some information. Here is a function which displays the result of division with remainder as a sentence with addition and multiplication.

In [ ]:

```
def display_divmod(a,b):
quotient = a // b # Integer division
remainder = a % b #
print "%d = %d (%d) + %d"%(a,quotient,b,remainder)
```

In [ ]:

```
display_divmod(23,5)
```

Notice that this function has no `return`

line. The function terminates automatically at the end of its scope.

The function also uses a slick Python string manipulation called **string substitution**. This has changed between Python 2.x and 3.x, and as usual we are following Python 2.x syntax.

String substitution allows you to place variable numbers in a string by putting placeholders like `%d`

within the string, and then placing the list of variables after the string (as a tuple). Here are some examples to get used to the syntax, since we will use this often in this lesson.

In [ ]:

```
print "My favorite number is %d"%(17) # The % symbol substitutes 17 for %d.
```

In [ ]:

```
print "%d + %d = %d"%(13,12,13+12)
```

The placeholder `%d`

is meant to hold the place of an integer. For floating point numbers, use `%f`

instead. This will impose a bit of rounding before printing. One can modify the placeholder to change the amount of digits displayed. Here is the official reference for placeholder syntax and string substitutions. We will only use the most basic features.

In [ ]:

```
print "pi is approximately %f"%(3.14159265) # This uses a default precision for display.
```

In [ ]:

```
print "pi is approximately %.3f"%(3.14159265) # This sets 3 digits after the decimal point.
```

What are the signals of scope in Python?

How is the symbol

`%`

used in Python? Can you think of three ways?Write a function called area_circle, which takes one argument radius. The function should return the area of the circle, as a floating point number. Then add one line to the function, using string substitutions, so that it additionally prints a helpful sentence of the form "The area of a circle of radius 1.0 is 3.14159." (depending on the radius and the area it computes).

Can you think of a reason you might want to have a function with

*no*arguments?

In [ ]:

```
# Use this space to work on the Exercises.
# Remember that you can add a new cell above/below by clicking to the left of a cell,
# (the cell will have a blue bar at the left) and then pressing "a" or "b" on the keyboard.
```

It is important for a computer program to behave differently under different circumstances. The simplest control statements, `if`

and its relative `else`

, can be used to tell Python to carry out different actions depending on the value of a boolean variable. The following function exhibits the syntax.

In [ ]:

```
def is_even(n):
if n%2 == 0:
print "%d is even."%(n)
return True
else:
print "%d is odd."%(n)
return False
```

In [ ]:

```
is_even(17)
```

In [ ]:

```
is_even(1000)
```

The broad syntax of the function should be familiar. We have created a function called `is_even`

with one argument called `n`

. The body of the function uses the **control statement** `if n%2 == 0:`

. Recall that `n%2`

gives the remainder after dividing `n`

by `2`

. Thus `n%2`

is 0 or 1, depending on whether `n`

is even or odd. Therefore the **boolean** `n%2 == 0`

is `True`

if `n`

is even, and `False`

if `n`

is odd.

The next two lines (the first `print`

and `return`

statements) are within the **scope** of the `if <boolean>:`

statement, as indicated by the colon and the indentation. The `if <boolean>:`

statement tells the Python interpreter to perform the statements within the scope if the boolean is `True`

, and to ignore the statements within the scope if the boolean is `False`

.

Putting it together, we can analyze the code.

```
if n%2 == 0:
print "%d is even."%(n)
return True
```

If `n`

is even, then the Python interpreter will print a sentence of the form `n is even`

. Then the interpreter will return (output) the value `True`

and the function will terminate. If `n`

is odd, the Python interpreter will ignore the two lines of scope.

Often we don't just want Python to *do nothing* when a condition is not satisfied. In the case above, we would rather Python tell us that the number is odd. The `else:`

control statement tells Python what to do in case the `if <boolean>:`

control statement receives a `False`

boolean. We analyze the code

```
else:
print "%d is odd."%(n)
return False
```

The `print`

and `return`

commands are within the scope of the `else:`

control statement. So when the `if`

statement receives a false signal (the number `n`

is odd), the program prints a sentence of the form `n is odd.`

and then returns the value `False`

and terminates the function.

The function `is_even`

is a verbose, or "talkative" sort of function. Such a function is sometimes useful in an interactive setting, where the programmer wants to understand everything that's going on. But if the function had to be called a million times, the screen would fill with printed sentences! In practice, an efficient and silent function `is_even`

might look like the following.

In [ ]:

```
def is_even(n):
return (n%2 == 0)
```

In [ ]:

```
is_even(17)
```

A `for`

loop and an `if`

control statement, used together, allow us to carry out a **brute force** search. We can search for factors in order to check whether a number is prime. Or we can look for solutions to an equation until we find one.

One thing to note: the function below begins with a block of text between a triple-quote (three single-quotes when typing). That text is called a **docstring** and it is meant to document what the function does. Writing clear docstrings becomes more important as you write longer programs, collaborate with other programmers, and when you want to return months or years later to use a program again. There are different style conventions for docstrings; for example, here are Google's docstring conventions. We take a less formal approach.

In [ ]:

```
def is_prime(n):
'''
Checks whether the argument n is a prime number.
Uses a brute force search for factors between 1 and n.
'''
for j in range(2,n): # the list of numbers 2,3,...,n-1.
if n%j == 0: # is n divisible by j?
print "%d is a factor of %d."%(j,n)
return False
return True
```

An important note: the `return`

keyword **terminates** the function. So as soon as a factor is found, the function terminates and outputs `False`

. If no factor is found, then the function execution survives past the loop, and the line `return True`

is executed to terminate the function.

In [ ]:

```
is_prime(91)
```

In [ ]:

```
is_prime(101)
```

Try the `is_prime`

function on bigger numbers -- try numbers with 4 digits, 5 digits, 6 digits. Where does it start to slow down? Do you get any errors when the numbers are large? Make sure to save your work first, just in case this crashes your computer!

In [ ]:

```
# Experiment with is_prime here.
```

There are two limiting factors, which we study in more detail in the next lesson. These are **time** and **space** (your computer's memory space). As the loop of `is_prime`

goes on and on, it might take your computer a long time! If each step of the loop takes only a nanosecond (1 billionth of a second), the loop would take about a second when executing `is_prime(1000000001)`

. If you tried `is_prime`

on a much larger number, like `is_prime(2**101 - 1)`

, the loop would take longer than the lifetime of the Earth.

But if you try such a big number, you'll hit a more immediate problem with *space*. Namely, the `range(2,n)`

attempts to *store* the entire list of numbers `[2,3,4,...,n-1]`

in the memory of your computer. Your computer has some (4 or 8 or 16, perhaps) gigabytes of memory (RAM). A gigabyte is a billion bytes, and a byte is enough memory to store a number between 0 and 255. (More detail about this later!). So a gigabyte will not even hold a billion numbers, and your computer will probably run out of memory if you attempt to try `is_prime(n)`

when `n`

is close to a billion.

We will address both of these problems in the next lesson, to some extent. The memory (*space*) problem is easier, and we will be able to adapt the `is_prime(n)`

function to deal with `n`

in the trillions. To go beyond this, to work with `n`

of size $10^{100}$ or even $10^{1000}$, we will find completely different methods in a later lesson.

Create a function

`my_abs(x)`

which outputs the absolute value of the argument`x`

. (Note that Python already has a built-in`abs(x)`

function).Modify the

`is_prime`

function so that it prints a message`Number too big`

and returns`None`

if the input argument is bigger than one million. (Note that`None`

is a Python reserved word. You can use the one-line statement`return None`

.)Write a Python function

`thrarity`

which takes an argument`n`

, and outputs the string`threeven`

if`n`

is a multiple of three, or`throdd`

is`n`

is one more than a multiple of three, or`thrugly`

if`n`

is one less than a multiple of three. Example:`thrarity(31)`

should output`throdd`

and`thrarity(44)`

should output`thrugly`

. Hint: study the`if`

/`elif`

syntax at the official Python tutorialWrite a Python function

`sum_of_squares(n)`

which finds and prints a pair of natural numbers $x$, $y$, such that $x^2 + y^2 = n$. The function should use a brute force search.

In [ ]:

```
# Use this space for your solutions to the questions.
```

We *almost* have all the tools we need to implement the Euclidean algorithm. The last tool we will need is the **while loop**. We have seen the *for loop* already, for iterating over a range of numbers. The Euclidean algorithm involves repetition, but there is no way to know in advance how many steps it will take. The while loop allows us to repeat a process as long as a boolean value (sometimes called a **flag**) is True. The following countdown example illustrates the structure of a while loop.

In [ ]:

```
def countdown(n):
current_value = n
while current_value > 0: # The condition (current_value > 0) is checked before every instance of the scope!
print current_value
current_value = current_value - 1
```

In [ ]:

```
countdown(10)
```

The while loop syntax begins with `while <boolean>:`

and the following indented lines comprise the scope of the loop. If the boolean is `True`

, then the scope of the loop is executed. If the boolean is `True`

again afterwards, then the scope of the loop is executed again. And again and again and so on.

This can be a **dangerous process**! For example, what would happen if you made a little typo and the last line of the while loop read `current_value = current_value + 1`

? The numbers would increase and increase... and the boolean `current_value > 0`

would **always** be `True`

. Therefore the loop would never end. Bigger and bigger numbers would scroll down your computer screen.

You might panic under such a circumstance, and maybe turn your computer off to stop the loop. Here is some advice for when your computer gets stuck in such a neverending loop:

- Back up your work often. When you're programming, make sure everything else is saved just in case.
- Save your programming work (use "Save and checkpoint" under the "File" menu) often, especially before running a cell with a loop for the first time.
- If you
*do*get stuck in a neverending loop, click on "Kernel... Interrupt". This will often unstick the loop and allow you to pick up where you left off. - On a Mac, you might try a "Force Quit" of the Python process, using the Activity Manager.

Now, if you're feeling brave, save your work, change the while loop so that it never ends, and try to recover where you left off. But be aware that this could cause your computer to freeze or behave erratically, crashing your browser, etc. Don't panic... it won't break your computer permanently.

The neverending loop causes two problems here. One is with your computer processor, which will be essentially spinning its wheels. This is called busy waiting, and your computer will essentially be busy waiting forever. The other problem is that your loop is printing more and more lines of text into the notebook. This could easily crash your web browser, which is trying to store and display zillions of lines of numbers. So be ready for problems!

The **Euclidean Algorithm** is a process of repeated division with remainder. Beginning with two integers `a`

(dividend) and `b`

(divisor), one computes quotient `q`

and remainder `q`

to express `a = qb + r`

. Then `b`

becomes the dividend and `r`

becomes the divisor, and one repeats. The repetition continues, and the **last nonzero** remainder is the greatest common divisor of `a`

and `b`

. It is the subject of Chapter 1 of An Illustrated Theory of Numbers.

We implement the Euclidean algorithm in a few variations. The first will be a verbose version, to show the user what happens at every step. We use a while loop to take care of the repetition.

In [ ]:

```
def Euclidean_algorithm(a,b):
dividend = a
divisor = b
while divisor != 0: # Recall that != means "is not equal to".
quotient = dividend // divisor
remainder = dividend % divisor
print "%d = %d (%d) + %d"%(dividend, quotient, divisor, remainder)
dividend = divisor
divisor = remainder
```

In [ ]:

```
Euclidean_algorithm(133, 58)
```

In [ ]:

```
Euclidean_algorithm(1312331323, 58123123)
```

This is excellent if we want to know every step of the Euclidean algorithm. If we just want to know the GCD of two numbers, we can be less verbose. We carefully return the last nonzero remainder after the while loop is concluded. This last nonzero remainder becomes the divisor when the remainder becomes zero, and then it would become the dividend in the next (unprinted) line. That is why we return the (absolute value) of the dividend after the loop is concluded. You might insert a line at the end of the loop, like `print dividend, divisor, remainder`

to help you track the variables.

In [ ]:

```
def GCD(a,b):
dividend = a # The first dividend is a.
divisor = b # The first divisor is b.
while divisor != 0: # Recall that != means "not equal to".
quotient = dividend // divisor
remainder = dividend % divisor
dividend = divisor
divisor = remainder
return abs(dividend) # abs() is used, since we like our GCDs to be positive.
```

Note that the `return dividend`

statement occurs *after* the scope of the while loop. So as soon as the *divisor* variable equals zero, the funtion `GCD`

returns the *dividend* variable and terminates.

In [ ]:

```
GCD(111,27)
```

In [ ]:

```
GCD(111,-27)
```

We can refine our code in a few ways. First, note that the `quotient`

variable is never used! It was nice in the verbose version of the Euclidean algorithm, but plays no role in finding the GCD. Our refined code reads

```
def GCD(a,b):
dividend = a
divisor = b
while divisor != 0: # Recall that != means "not equal to".
remainder = dividend % divisor
dividend = divisor
divisor = remainder
return abs(dividend)
```

Now there are two slick Python tricks we can use to shorten the code. The first is called **multiple assignment**. It is possible to set the values of two variables in a single line of code, with a syntax like below.

In [ ]:

```
x,y = 2,3 # Sets x to 2 and y to 3.
```

This is particular useful for self-referential assignments, because as for ordinary assignment, the right side is evaluated first and then bound to the variables on the left side. For example, after the line above, try the line below. Use print statements to see what the values of the variables are afterwards!

In [ ]:

```
x,y = y,x # Guess what this does!
```

In [ ]:

```
print "x =",x
print "y =",y
```

Now we can use multiple assignment to turn three lines of code into one line of code. For the `remainder`

variable is only used temporarily before its value is given to the `divisor`

variable. Using multiple assignment, the three lines

```
remainder = dividend % divisor
dividend = divisor
divisor = remainder
```

can be written in one line,

```
dividend, divisor = divisor, dividend % divisor # Evaluations on the right occur before any assignments!
```

Our newly shortened GCD function looks like this.

```
def GCD(a,b):
dividend = a
divisor = b
while divisor != 0: # Recall that != means "not equal to".
dividend, divisor = divisor, dividend % divisor
return abs(dividend)
```

The next trick involves the while loop. The usual syntax has the form `while <boolean>:`

. But if `while`

is followed by a numerical type, e.g. `while <int>:`

, then the scope of the while loop will execute as long as the number is nonzero! Therefore, the line

```
while divisor != 0:
```

can be replaced by the shorter line

```
while divisor:
```

This is truly a trick. It probably won't speed anything up, and it does not make your program easier to read for beginners. So use it if you prefer communicating with experienced Python programmers! Here is the whole function again.

```
def GCD(a,b):
dividend = a
divisor = b
while divisor: # Executes the scope if divisor is nonzero.
dividend, divisor = divisor, dividend % divisor
return abs(dividend)
```

The next optimization is a bit more dangerous for beginners, but it works here. In general, it can be dangerous to operate directly on the arguments to a function. But in this setting, it is safe, and makes no real difference to the Python interpreter. Instead of creating new variables called `dividend`

and `divisor`

, one can manipulate `a`

and `b`

directly within the function. If you do this, the GCD function can be shortened to the following.

In [ ]:

```
def GCD(a,b):
while b: # Recall that != means "not equal to".
a, b = b, a % b
return abs(a)
```

In [ ]:

```
# Try it out. Try it on some big numbers and see how quick it runs!
```

This code is essentially optimal, if one wishes to execute the Euclidean algorithm to find the GCD of two integers. It almost matches the GCD code in a standard Python library. It might be slightly faster than our original code -- but there is a tradeoff here between execution speed and readability of code. In this and the following lessons, we often optimize enough for everyday purposes, but not so much that readability is lost.

Modify the

`is_prime`

function by using a while loop instead of`for j in range(2,n):`

. Hint: the function should contain the lines`j = 2`

and`while j < n:`

and`j = j + 1`

in various places. Why might this be an improvement from the for loop?Modify the

`Euclidean_algorithm`

function to create a function which returns the*number of steps*that the Euclidean algorithm requires, i.e., the number of divisions-with-remainder.Create a function which carries out division with minimal remainder. In other words, given integers a,b, the function expresses a = q(b) + r, where r is a

*positive or negative*integer of magnitude bounded by b/2. Use such a function to create a new Euclidean algorithm function which uses minimal remainder.What

`GCD(a,b)`

function do you think strikes the best balance between efficiency and readability?

In [ ]:

```
# Use this space to work on the exercises.
```

In Chapter 1 of An Illustrated Theory of Numbers, we not only used the Euclidean algorithm to find the GCD of two integers, but also to solve the linear Diophantine equation $ax + by = c$. On paper, this required us to perform the Euclidean algorithm, then "work backwards" to carefully solve a series of linear equations. This process is repetetive and error-prone... perfect for a computer.

So here we develop a function `solve_LDE(a,b,c)`

which will describe all integer solutions $x,y$ to the equation $ax + by = c$.

The idea of the algorithm is to keep track of "hops" and "skips" throughout the Euclidean algorithm. A general step in the Euclidean algorithm looks like `u = q(v) + r`

. The remainder can then be expressed by the formula `r = u - q(v)`

. If `u`

and `v`

can be built from hops and skips, then `r`

can be built from hops and skips. How many? Just tally the hops and skips to find:

`r_hops = u_hops - q (v_hops)` and `r_skips = u_skips - q (v_skips)`.

This sort of tallying is what makes the algorithm below work. The function below does not introduce any new programming concepts, but it assembles many ideas together.

In [ ]:

```
def hop_and_skip(a,b):
'''
Takes two integer arguments a,b, and prints a sentence of the form
GCD(a,b) = x(a) + y(b). The method is the Euclidean algorithm,
tallying hops (units of a) and skips (units of b) along the way.
'''
u = a # We use u instead of dividend.
v = b # We use v instead of divisor.
u_hops, u_skips = 1,0 # u is built from one hop (a) and no skips, for now.
v_hops, v_skips = 0,1 # v is built from no hops and one skip (b), for now.
while v != 0: # We could just write "while v:"
q = u // v # q stands for quotient.
r = u % v # r stands for remainder. So u = q(v) + r.
r_hops = u_hops - q * v_hops # Tally hops
r_skips = u_skips - q * v_skips # Tally skips
u,v = v,r # The new dividend,divisor is the old divisor,remainder.
u_hops, v_hops = v_hops, r_hops # The new u_hops, v_hops is the old v_hops, r_hops
u_skips, v_skips = v_skips, r_skips # The new u_skips, v_skips is the old v_skips, r_skips
print "%d = %d(%d) + %d(%d)"%(u,u_hops,a,u_skips,b)
```

In [ ]:

```
hop_and_skip(102,45)
```

Try out the hop_and_skip code on some integers of your choice. Does it behave correctly? Check the results, using Python as a calculator. Does it run quickly for large integers?

In [ ]:

```
# Experimentation space here.
```

To conclude this lesson, we put everything together to create a long(ish) function to solve linear Diophantine equations. We want this function to be smart enough to respond when an equation has no solutions, and to describe *all* solutions when they exist.

The first part of the function `solve_LDE`

is the same as the hop and skip function above. But rather than expressing $GCD(a,b)$ as $ax + by$, the function uses the GCD to determine the existence and the general form of a solution to $ax + by = c$. The formula for the general form comes from An Illustrated Theory of Numbers, Chapter 1, Corollary 1.25.

In [ ]:

```
def solve_LDE(a,b,c):
'''
Describes all of the solutions to the linear Diophantine equation
ax + by = c. There are either no solutions or infinitely many solutions.
Prints a description of the solution set, and returns None if there are no solutions
or returns a single solution if one exists.
'''
u = a # We use u instead of dividend.
v = b # We use v instead of divisor.
u_hops, u_skips = 1,0 # u is built from one hop (a) and no skips.
v_hops, v_skips = 0,1 # v is built from no hops and one skip (b).
while v != 0: # We could just write while v:
q = u // v # q stands for quotient.
r = u % v # r stands for remainder. So u = q(v) + r.
r_hops = u_hops - q * v_hops # Tally hops
r_skips = u_skips - q * v_skips # Tally skips
u,v = v,r # The new dividend,divisor is the old divisor,remainder.
u_hops, v_hops = v_hops, r_hops # The new u_hops, v_hops is the old v_hops, r_hops
u_skips, v_skips = v_skips, r_skips # The new u_skips, v_skips is the old v_skips, r_skips
g = u # The variable g now describes the GCD of a and b.
if c%g == 0: # When GCD(a,b) divides c...
d = c/g
x = d * u_hops
y = d * u_skips # Now ax + by = c is a specific solution!
print "%d x + %d y = %d if and only if "%(a, b, c)
print "x = %d + %d n and y = %d - %d n, for some integer n."%(x,b/g,y,-a/g)
return x,y
else: # When GCD(a,b) does not divide c...
print "There are no solutions to %d x + %d y = %d,"%(a,b,c)
print "because GCD(%d, %d) = %d, which does not divide %d."%(a,b,g,c)
return None
```

In [ ]:

```
solve_LDE(102,45,3)
```

In [ ]:

```
solve_LDE(72,100,17)
```

Solve problems 4-7 of Chapter 1 of An Illustrated Theory of Numbers using the

`solve_LDE`

function.Write an

`LCM`

function, using the previous`GCD`

function and the GCD-LCM product formula, Theorem 1.23 of An Illustrated Theory of Numbers.Sometimes it is important to find not the integer solutions, but the

*positive*integer solutions to a Diophantine equation. Modify the`solve_LDE`

function to create a`solve_LDE_positive(a,b,c)`

function. The output of the function should be all pairs of*positive*integers $x$, $y$, such that $ax + by = c$ (if any pairs exist), and a helpful message if no pairs exist (and a`return None`

should be used in this case).

In [ ]:

```
# Use this space to work on the exercises.
```