Loops

for loops

A for loop allows us to execute a block of code multiple times with some parameters updated each time through the loop. A for loop begins with the for statement:

for item in iterable:
    # code block indented 4 spaces

The main points to observe are:

  • for and in keywords
  • iterable is a sequence object such as a list, tuple or range
  • item is a variable which takes each value in iterable
  • end for statement with a colon :
  • code block indented 4 spaces which executes once for each value in iterable

For example, let's print $n^2$ for $n$ from 0 to 5:

for n in [0,1,2,3,4,5]:
    square = n**2
    print(n,'squared is',square)
print('The for loop is complete!')
0 squared is 0
1 squared is 1
2 squared is 4
3 squared is 9
4 squared is 16
5 squared is 25
The for loop is complete!

Copy and paste this code and any of the examples below into the Python visualizer to see each step in a for loop!

while loops

What if we want to execute a block of code multiple times but we don't know exactly how many times? We can't write a for loop because this requires us to set the length of the loop in advance. This is a situation when a while loop is useful.

The following example illustrates a while loop:

n = 5
while n > 0:
    print(n)
    n = n - 1
5
4
3
2
1

The main points to observe are:

  • while keyword
  • a logical expression followed by a colon :
  • loop executes its code block if the logical expression evaluates to True
  • update the variable in the logical expression each time through the loop
  • BEWARE! If the logical expression always evaluates to True, then you get an infinite loop!

We prefer for loops over while loops because of the last point. A for loop will never result in an infinite loop. If a loop can be constructed with for or while, we'll always choose for.

Constructing Sequences

There are several ways to construct a sequence of values and to save them as a Python list. We have already seen Python's list comprehension syntax. There is also the append list method described below.

Sequences by a Formula

If a sequence is given by a formula then we can use a list comprehension to construct it. For example, the sequence of squares from 1 to 100 can be constructed using a list comprehension:

squares = [d**2 for d in range(1,11)]
print(squares)
[1, 4, 9, 16, 25, 36, 49, 64, 81, 100]

However, we can achieve the same result with a for loop and the append method for lists:

# Intialize an empty list
squares = []
for d in range(1,11):
    # Append the next square to the list
    squares.append(d**2)
print(squares)
[1, 4, 9, 16, 25, 36, 49, 64, 81, 100]

In fact, the two examples above are equivalent. The purpose of list comprehensions is to simplify and compress the syntax into a one-line construction.

Recursive Sequences

We can only use a list comprehension to construct a sequence when the sequence values are defined by a formula. But what if we want to construct a sequence where the next value depends on previous values? This is called a recursive sequence.

For example, consider the Fibonacci sequence:

$$ x_1 = 1, x_2 = 1, x_3 = 2, x_4 = 3, x_5 = 5, ... $$

where

$$ x_{n} = x_{n-1} + x_{n-2} $$

We can't use a list comprehension to build the list of Fibonacci numbers, and so we must use a for loop with the append method instead. For example, the first 15 Fibonacci numbers are:

fibonacci_numbers = [1,1]
for n in range(2,15):
    fibonacci_n = fibonacci_numbers[n-1] + fibonacci_numbers[n-2]
    fibonacci_numbers.append(fibonacci_n)
    print(fibonacci_numbers)
[1, 1, 2]
[1, 1, 2, 3]
[1, 1, 2, 3, 5]
[1, 1, 2, 3, 5, 8]
[1, 1, 2, 3, 5, 8, 13]
[1, 1, 2, 3, 5, 8, 13, 21]
[1, 1, 2, 3, 5, 8, 13, 21, 34]
[1, 1, 2, 3, 5, 8, 13, 21, 34, 55]
[1, 1, 2, 3, 5, 8, 13, 21, 34, 55, 89]
[1, 1, 2, 3, 5, 8, 13, 21, 34, 55, 89, 144]
[1, 1, 2, 3, 5, 8, 13, 21, 34, 55, 89, 144, 233]
[1, 1, 2, 3, 5, 8, 13, 21, 34, 55, 89, 144, 233, 377]
[1, 1, 2, 3, 5, 8, 13, 21, 34, 55, 89, 144, 233, 377, 610]

Computing Sums

Suppose we want to compute the sum of a sequence of numbers $x_0$, $x_1$, $x_2$, $x_3$, $\dots$, $x_n$. There are at least two approaches:

  1. Compute the entire sequence, store it as a list $[x_0,x_1,x_2,\dots,x_n]$ and then use the built-in function sum.
  2. Initialize a variable with value 0 (and name it result for example), create and add each element in the sequence to result one at a time.

The advantage of the second approach is that we don't need to store all the values at once. For example, here are two ways to write a function which computes the sum of squares.

For the first approach, use a list comprehension:

def sum_of_squares_1(N):
    "Compute the sum of squares 1**2 + 2**2 + ... + N**2."
    return sum([n**2 for n in range(1,N + 1)])
sum_of_squares_1(4)
30

For the second approach, use a for loop with the initialize-and-update construction:

def sum_of_squares_2(N):
    "Compute the sum of squares 1**2 + 2**2 + ... + N**2."
    # Initialize the output value to 0
    result = 0
    for n in range(1,N + 1):
        # Update the result by adding the next term
        result = result + n**2
    return result
sum_of_squares_2(4)
30

Again, both methods yield the same result however the second uses less memory!

Computing Products

There is no built-in function to compute products of sequences therefore we'll use an initialize-and-update construction similar to the example above for computing sums.

Write a function called factorial which takes a positive integer $N$ and return the factorial $N!$.

def factorial(N):
    "Compute N!=N(N-1) ... (2)(1) for N >= 1."
    # Initialize the output variable to 1
    product = 1
    for n in range(2,N + 1):
        # Update the output variable
        product = product * n
    return product

Let's test our function for input values for which we know the result:

factorial(2)
2
factorial(5)
120

We can use our function to approximate $e$ using the Taylor series for $e^x$:

$$ e^x = \sum_{k=0}^{\infty} \frac{x^k}{k!} $$

For example, let's compute the 100th partial sum of the series with $x=1$:

sum([1/factorial(k) for k in range(0,101)])
2.7182818284590455

Solutions of Polynomials

We can use nested for loops to search for integer solutions of equations. For example, suppose we would like to find all representations of a positive integer $N$ as a sum of two squares. In other words, we want to find all integer solutions $(x,y)$ of the equation:

$$ x^2 + y^2 = N $$

Write a function called reps_sum_squares which takes an integer $N$ and finds all representations of $N$ as a sum of squares $x^2 + y^2 = N$ for $0 \leq x \leq y$. The function returns the representations as a list of tuples. For example, if $N = 50$ then $1^2 + 7^2 = 50$ and $5^2 + 5^2 = 50$ and the function returns the list [(1, 7),(5, 5)].

Let's outline our approach before we write any code:

  1. Given $x \leq y$, the largest possible value for $x$ is $\sqrt{\frac{N}{2}}$
  2. For $x \leq \sqrt{\frac{N}{2}}$, the pair $(x,y)$ is a solution if $N - x^2$ is a square
  3. Define a helper function called is_square to test if an integer is square
def is_square(n):
    "Determine if the integer n is a square."
    if round(n**0.5)**2 == n:
        return True
    else:
        return False

def reps_sum_squares(N):
    '''Find all representations of N as a sum of squares x**2 + y**2 = N.

    Input:
        N : integer
    Output:
        List of tuples (x,y) of positive integers such that x**2 + y**2 = N
    Example:
        >>> reps_sum_squares(1105)
        [(4, 33), (9, 32), (12, 31), (23, 24)]
    '''
    reps = []
    if is_square(N/2):
        # If N/2 is a square, search up to x = (N/2)**0.5
        max_x = round((N/2)**0.5)
    else:
        # If N/2 is not a square, search up to x = floor((N/2)**0.5)
        max_x = int((N/2)**0.5)
    for x in range(0,max_x + 1):
        y_squared = N - x**2
        if is_square(y_squared):
            y = round(y_squared**0.5)
            reps.append((x,y))
    return reps
reps_sum_squares(1105)
[(4, 33), (9, 32), (12, 31), (23, 24)]

What is the smallest integer which can be expressed as the sum of squares in 5 different ways?

N = 1105
num_reps = 4
while num_reps != 5:
    N = N + 1
    reps = reps_sum_squares(N)
    num_reps = len(reps)
print(N,':',reps_sum_squares(N))
4225 : [(0, 65), (16, 63), (25, 60), (33, 56), (39, 52)]

Examples

Prime Numbers

A positive integer is prime if it is divisible only by 1 and itself. Write a function called is_prime which takes an input parameter n and returns True or False depending on whether n is prime or not.

Let's outline our approach before we write any code:

  1. $d$ divides $n$ if there is no remainder of $n$ divided by $d$
  2. Use the modulus operator % to compute the remainder
  3. If $d$ divides $n$ then $n = d q$ for some integer $q$ and either $d \leq \sqrt{n}$ or $q \leq \sqrt{n}$ (and not both), therefore we need only test if $d$ divides $n$ for integers $d \leq \sqrt{n}$
def is_prime(n):
    "Determine whether or not n is a prime number."
    if n <= 1:
        return False
    # Test if d divides n for d <= n**0.5
    for d in range(2,round(n**0.5) + 1):
        if n % d == 0:
            # n is divisible by d and so n is not prime
            return False
    # If we exit the for loop, then n is not divisible by any d
    # and therefore n is prime
    return True

Let's test our function on the first 30 numbers:

for n in range(0,31):
    if is_prime(n):
        print(n,'is prime!')
2 is prime!
3 is prime!
5 is prime!
7 is prime!
11 is prime!
13 is prime!
17 is prime!
19 is prime!
23 is prime!
29 is prime!

Our function works! Let's find all the primes between 20,000 and 20,100.

for n in range(20000,20100):
    if is_prime(n):
        print(n,'is prime!')
20011 is prime!
20021 is prime!
20023 is prime!
20029 is prime!
20047 is prime!
20051 is prime!
20063 is prime!
20071 is prime!
20089 is prime!

Divisors

Let's write a function called divisors which takes a positive integer $N$ and returns the list of positive integers which divide $N$.

def divisors(N):
    "Return the list of divisors of N."
    # Initialize the list of divisors (which always includes 1)
    divisor_list = [1]
    # Check division by d for d <= N/2
    for d in range(2,N // 2 + 1):
        if N % d == 0:
            divisor_list.append(d)
    divisor_list.append(N)
    return divisor_list

Let's test our function:

divisors(10)
[1, 2, 5, 10]
divisors(100)
[1, 2, 4, 5, 10, 20, 25, 50, 100]
divisors(59)
[1, 59]

Collatz Conjecture

Let $a$ be a positive integer and consider the recursive sequence where $x_0 = a$ and

$$ x_{n+1} = \left\{ \begin{array}{cl} x_n/2 & \text{if } x_n \text{ is even} \\ 3x_n+1 & \text{if } x_n \text{ is odd} \end{array} \right. $$

The Collatz conjecture states that this sequence will always reach 1. For example, if $a = 10$ then $x_0 = 10$, $x_1 = 5$, $x_2 = 16$, $x_3 = 8$, $x_4 = 4$, $x_5 = 2$ and $x_6 = 1$.

Write a function called collatz which takes one input parameter a and returns the sequence of integers defined above and ending with the first occurrence $x_n=1$.

def collatz(a):
    "Compute the Collatz sequence starting at a and ending at 1."
    # Initialize list with first value a
    sequence = [a]
    # Compute values until we reach 1
    while sequence[-1] > 1:
        # Check if the last element in the list is even
        if sequence[-1] % 2 == 0:
            # Compute and append the new value
            sequence.append(sequence[-1] // 2)
        else:
            # Compute and append the new value
            sequence.append(3*sequence[-1] + 1)
    return sequence

Let's test our function:

print(collatz(10))
[10, 5, 16, 8, 4, 2, 1]
collatz(22)
[22, 11, 34, 17, 52, 26, 13, 40, 20, 10, 5, 16, 8, 4, 2, 1]

The Collatz conjecture is quite amazing. No matter where we start, the sequence always terminates at zero!

collatz(393570)
[393570, 196785, 590356, 295178, 147589, 442768, 221384, 110692, 55346, 27673, 83020, 41510, 20755, 62266, 31133, 93400, 46700, 23350, 11675, 35026, 17513, 52540, 26270, 13135, 39406, 19703, 59110, 29555, 88666, 44333, 133000, 66500, 33250, 16625, 49876, 24938, 12469, 37408, 18704, 9352, 4676, 2338, 1169, 3508, 1754, 877, 2632, 1316, 658, 329, 988, 494, 247, 742, 371, 1114, 557, 1672, 836, 418, 209, 628, 314, 157, 472, 236, 118, 59, 178, 89, 268, 134, 67, 202, 101, 304, 152, 76, 38, 19, 58, 29, 88, 44, 22, 11, 34, 17, 52, 26, 13, 40, 20, 10, 5, 16, 8, 4, 2, 1]

Which $a < 1000$ produces the longest sequence?

max_length = 1
a_max = 1
for a in range(1,1001):
    seq_length = len(collatz(a))
    if seq_length > max_length:
        max_length = seq_length
        a_max = a
print('Longest sequence begins with a =',a_max,'and has length',max_length)
Longest sequence begins with a = 871 and has length 179

Exercises

  1. Fermat's theorem on the sum of two squares states that every prime number $p$ of the form $4k+1$ can be expressed as the sum of two squares. For example, $5 = 2^2 + 1^2$ and $13 = 3^2 + 2^2$. Find the smallest prime greater than $2018$ of the form $4k+1$ and write it as a sum of squares. (Use the functions is_prime and reps_sum_squares from this section.)

  2. What is the smallest prime number which can be represented as a sum of squares in 2 different ways?

  3. What is the smallest integer which can be represented as a sum of squares in 3 different ways?

  4. Write a function called primes_between which takes two integer inputs $a$ and $b$ and returns the list of primes in the closed interval $[a,b]$.

  5. Write a function called primes_d_mod_N which takes four integer inputs $a$, $b$, $d$ and $N$ and returns the list of primes in the closed interval $[a,b]$ which are congruent to $d$ mod $N$ (this means that the prime has remainder $d$ after division by $N$). This kind of list is called primes in an arithmetic progression.

  6. Write a function called reciprocal_recursion which takes three positive integers $x_0$, $x_1$ and $N$ and returns the sequence $[x_0,x_1,x_2,\dots,x_N]$ where

    $$ x_n = \frac{1}{x_{n-1}} + \frac{1}{x_{n-2}} $$

  7. Write a function called root_sequence which takes input parameters $a$ and $N$, both positive integers, and returns the $N$th term $x_N$ in the sequence:

    $$ \begin{align} x_0 &= a \\ x_n &= 1 + \sqrt{x_{n-1}} \end{align} $$

    Does the sequence converge to different values for different starting values $a$?

  8. Write a function called fib_less_than which takes one input $N$ and returns the list of Fibonacci numbers less than $N$.

  9. Write a function called fibonacci_primes which takes an input parameter $N$ and returns the list of Fibonacci numbers less than $N$ which are also prime numbers.

  10. Let $w(N)$ be the number of ways $N$ can be expressed as a sum of two squares $x^2 + y^2 = N$ with $1 \leq x \leq y$. Then

    $$ \lim_{N \to \infty} \frac{1}{N} \sum_{n=1}^{N} w(n) = \frac{\pi}{8} $$

    Compute the left side of the formula for $N=100$ and compare the result to $\pi / 8$.

  11. A list of positive integers $[a,b,c]$ (with $1 \leq a < b$) are a Pythagorean triple if $a^2 + b^2 = c^2$. Write a function called py_triples which takes an input parameter $N$ and returns the list of Pythagorean triples [a,b,c] with $c \leq N$.