Wednesday, September 30th, 2026

Last class, we continued to look at plotting using matplotlib.pyplot. In particular, we had looked at how to generate a list of equally-spaced points spanning a finite interval, and how to plug those points into some function so that it could be plotted.

NumPy

The NumPy (Numerical Python) module contains many useful tools for numerical calculations in Python. We typically import the module and assign the name np.

import numpy as np

#help(np)

The basic building blocks in NumPy are arrays, which in many ways behave like lists. We can use the np.array function to convert a list to an array.

my_array = np.array([0,1,2,3,4,5,6])
my_array

array([0, 1, 2, 3, 4, 5, 6])

Just like with lists, we can acccess elements of an array by index using square brackets.

print(my_array[0])

0

print(my_array[3])

3

print(my_array[-1])

6

We can also access portions of an array with slicing, similar to lists.

print(my_array[:4])

[0 1 2 3]

print(my_array[::2])

[0 2 4 6]

print(my_array[1:-1:3])

[1 4]

Unlike lists, NumPy arrays are built to perform arithmetic operations on an element-by-element basis. For example, compare what happens when we multiply a list by an integer and what happens when we multiply an array by an integer.

my_list = [1,2,3,4]
my_array = np.array([1,2,3,4])

2*my_list

[1, 2, 3, 4, 1, 2, 3, 4]

2*my_array

array([2, 4, 6, 8])

In a similar way:

• Adding/subtracting/multiplying/dividing an array by an integer/float adds/subtracts/multiplies/divides each element by the integer/float.
• Exponentiating/modular dividing an array by an integer/float is performed on each element.
• If we have two arrays of the same shape, we can add/subtract/multiply/divide/exponentiate/modular divide one by the other. The operation will be performed element-by-element (that is, the first elements from each array will be added/subtracted/multiplied/etc., the second elements from each array will be...).

my_array + 10

array([11, 12, 13, 14])

my_array - 5

array([-4, -3, -2, -1])

my_array / 2

array([0.5, 1. , 1.5, 2. ])

2**my_array

array([ 2,  4,  8, 16])

my_array**2

array([ 1,  4,  9, 16])

my_array % 2

array([1, 0, 1, 0])

Let's return to the earlier exercise of plotting $y(t)$ using a densely-packed list (or array) of . The np.linspace function can be used to easliy generate an array of evenly spaced points over some interval. The basic syntax

np.linspace(a,b,N)

generates an array of $N$ evenly spaced points over the interval $[a, b]$.

np.linspace(0,2,5)

array([0. , 0.5, 1. , 1.5, 2. ])

np.linspace(0,10,11)

array([ 0.,  1.,  2.,  3.,  4.,  5.,  6.,  7.,  8.,  9., 10.])

points = np.linspace(-5,5,1001)
print(points[:10])

[-5.   -4.99 -4.98 -4.97 -4.96 -4.95 -4.94 -4.93 -4.92 -4.91]

Exercise: Use NumPy and plt.plot to plot the previously described trajectories of objects in freefall on Earth and Mars.

def y(t,y0,v0=0,g=-9.8):
    return (g/2)*t**2 + v0*t + y0

N = 1000
a = 0
b = 16
y0 = 1000

t_array = np.linspace(a,b,N+1)
y_Earth = y(t_array, y0)

import matplotlib.pyplot as plt

plt.plot(t_array, y_Earth, 'b--', label='Earth')
plt.xlabel('Time $t$ (in seconds)')
plt.ylabel('Height $y(t)$ (in meters)')
plt.title('Objects in free-fall')
plt.legend()

<matplotlib.legend.Legend at 0x1cbbc5b5d30>

Python loops versus NumPy

Python is a relatively easy language to code with (in comparison to C/C++, for example). On the other hand, it is a somewhat slow language. For example, Python is much slower to complete a for loop iteration than C/C++. The backbone of NumpPy is written in C/C++, which comes with speed advantages. The module essentially gives us the best of both worlds: the ease of use of Python with the speed of C/C++.

Let's compare how long it takes to generate $1,000,000$ equally spaced $t$-values over the interval $[0,2\pi]$ and compute both $x=\cos(t)$ and $y=\sin(t)$ using:

• for loops to construct the lists of $t$-, $x$-, and $y$-values.
• list comprehensions to construct the lists of $t$-, $x$-, and $y$-values.
• np.linspace, np.sin, and np.cos to generate the arrays of $t$-, $x$-, and $y$-values.

Reminder: We can use the time function from the time module to time our code.

import time

Recall, earlier we discussed how to divide an interval $[a,b]$ into $N$ equally spaced sub-intervals with endpoints $t_0, t_1, t_2, \dots, t_N$, where $$t_i = a + i\Delta t \qquad \text{for} \qquad i = 0, 1, ..., N, \qquad \text{with} \qquad \Delta t = \frac{b-a}{N}.$$

Let's time how long it takes to build our list of $t$-values, $x$-values, and $y$-values using for loops only (without list comprehension).

from math import sin, cos, pi

t0 = time.time()

a = 0
b = 2*pi
N = 10**6
dt = (b - a)/N

t_list = []
for i in range(N+1):
    t_list.append(a + i*dt)

x_list = []
for t in t_list:
    x_list.append(cos(t))

y_list = []
for t in t_list:
    y_list.append(sin(t))

t1 = time.time()

print('Elapsed time:', t1-t0)

Elapsed time: 0.6989145278930664

Let's see if using list comprehensions brings any speed improvements.

t0 = time.time()

a = 0
b = 2*pi
N = 10**6
dt = (b - a)/N

t_list = [(a + i*dt) for i in range(N+1)]
x_list = [cos(t) for t in t_list]
y_list = [sin(t) for t in t_list]

t1 = time.time()

print('Elapsed time:', t1-t0)

Elapsed time: 0.5027945041656494

Finally, let's see if using NumPy arrays brings any speed improvements.

t0 = time.time()

a = 0
b = 2*np.pi
N = 10**6

t_array = np.linspace(a,b,N+1)
x_array = np.cos(t_array)
y_array = np.sin(t_array)

t1 = time.time()

print('Elapsed time:', t1-t0)

Elapsed time: 0.03254246711730957

plt.plot(t_array, x_array)
plt.plot(t_array, y_array)

[<matplotlib.lines.Line2D at 0x1cbc1461bd0>]

When possible, it is usually good practice to try to take advantage of NumPy for computations.

Generating Numpy arrays

We've already seen how we can convert lists to arrays and how we can use the np.linspace function to create arrays. There are many other ways to build arrays. For example:

• np.zeros can generate arrays full of 0s.
• np.ones can generate arrays full of 1s.
• np.arange works just like the normal range function, except that it returns an array instead.

np.zeros(10)

array([0., 0., 0., 0., 0., 0., 0., 0., 0., 0.])

np.ones(5)

array([1., 1., 1., 1., 1.])

np.arange(10)

array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9])

np.arange(4,21)

array([ 4,  5,  6,  7,  8,  9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20])

np.arange(4,21,6)

array([ 4, 10, 16])

Note: the np.zeros function and np.ones function produce arrays filled with floats, while the np.arange function returns an array full of integers.

In general, NumPy arrays can only be filled with one datatype (unlike lists which can mix-and-match). We can check what datatype an array holds using the .dtype attribute:

my_array = np.linspace(0,1,100)
print(my_array.dtype)

float64

my_array = np.arange(11)
print(my_array.dtype)

int64

my_array = np.zeros(100)
print(my_array.dtype)

float64

We can change the datatype when defining arrays using np.zeros or np.ones by using the dtype keyword argument:

np.zeros(10, dtype=int)

array([0, 0, 0, 0, 0, 0, 0, 0, 0, 0])

np.zeros(10, dtype=bool)

array([False, False, False, False, False, False, False, False, False,
       False])

np.ones(5, dtype=bool)

array([ True,  True,  True,  True,  True])

We can also convert an existing array to a new datatype using the .astype method:

my_zeros = np.zeros(10)
print(my_zeros.dtype)

float64

my_zeros.astype(int)

array([0, 0, 0, 0, 0, 0, 0, 0, 0, 0])

my_zeros.astype(bool)

array([False, False, False, False, False, False, False, False, False,
       False])

Note: when converting to a Boolean data type, any non-zero integer/float is considered True. Only 0 is considered False.

Slicing with NumPy

Just like with lists, we use slicing to access portions of a NumPy array. Consider the following experiment.

First, we define a list and an array, each containing the integers 0, ..., 19.

N = 20

my_list = [i for i in range(N)]
my_array = np.arange(N)

print(my_list)
print(my_array)

[0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19]
[ 0  1  2  3  4  5  6  7  8  9 10 11 12 13 14 15 16 17 18 19]

Now, let's take slices of the list and array that starts at index 1 and takes every other element.

list_slice = my_list[1::2]
array_slice = my_array[1::2]

Before continuing, let's see what the slices look like.

print(list_slice)
print(array_slice)

[1, 3, 5, 7, 9, 11, 13, 15, 17, 19]
[ 1  3  5  7  9 11 13 15 17 19]

Small note: lists and NumPy arrays are displayed slighty differently when printed. Namely, when printing NumPy arrays, the commas separating the elements are omitted.

What happens if we modify these slices? For example, let's change the first entry of the list slice to be 99.

list_slice[0] = 99
print(list_slice)

[99, 3, 5, 7, 9, 11, 13, 15, 17, 19]

print(my_list)

[0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19]

We can see that, with the list, the changes to the slice do not propogate back to the original list. What happens if we perform the same experiment with the NumPy array?

array_slice[0] = 99
print(array_slice)

[99  3  5  7  9 11 13 15 17 19]

print(my_array)

[ 0 99  2  3  4  5  6  7  8  9 10 11 12 13 14 15 16 17 18 19]

With the NumPy array, the changes to the slice do propogate back to the original list.

What's happening here?

• When slicing a list, we obtain a new list object that is unattached to the original list. Changes to one do not affect the other.
• When slicing an array, we obtain a "view" of the original array. Changes to the slice affect the original array, and vice-versa.

If we want to obtain a slice of an array that is unattached to the original array, we can use the .copy method to sever the connection.

Back to Project 2: Pythagorean triples

In-class exercises start

Last class, we started to look at putting code together to find and visualize Pythagorean triples.

Exercise: Write a function is_square that takes in an integer n and returns True if n is a perfect square and False otherwise.

Exercise: Write a function is_pdouble that takes in a double of integers $(a,b)$ and returns True if the double is a Pythagorean double, and False otherwise.

Exercise: Write a function get_ptriples that takes in an integer $n$ and returns a list of all Pythagorean triples $(a,b,c)$ such that $1 \leq a,b \leq n$.

In-class exercises stop

Exercise: Use your get_ptriples function to generate all Pythagorean triples with $1 \leq a,b \leq 50$. Then, create a figure that plots $a$ versus $b$.

You will need to generate separate lists that contain just the $a$ parts and just the $b$ parts of each Pythagorean triple.