Solutions#
Let’s see if your functions function.
1. There are several infinite series that converge to something proportional to \(\pi\). For example, the Basel series is defined as an infinite sum of reciprocal square numbers:
a) Write a function to approximate \(\pi\) by evaluating the Basel series up to a given term \(n\). Test your function for \(n = 1000\).
import math
def basel(n):
"""
Evaluate the Basel series (1/1^2 + 1/2^2 + 1/3^3 + ...) up to term n.
Args:
n (int): The term up to which the series should be evaluated.
Returns:
(float): The series evaluated up to term n.
"""
cumulative_sum = 0
for k in range(1, n + 1):
cumulative_sum += 1 / k ** 2
return math.sqrt(6 * cumulative_sum)
pi = basel(1000)
print(f'pi is approximately {pi}')
pi is approximately 3.1406380562059946
b) Another infinite series that converges to something proportional to \(\pi\) is the Leibniz formula:
Write a function to approximate \(\pi\) by evaluating the Leibniz formula up to a given term \(n\). Test your function for \(n = 1000\).
def leibniz(n):
"""
Evaluate the Leibniz formula (1 - 1/3 + 1/5 - 1/7 + 1/9 - ...) up to term n.
Args:
n (int): The term up to which the series should be evaluated.
Returns:
(float): The series evaluated up to term n.
"""
cumulative_sum = 0
for k in range(0, n + 1):
cumulative_sum += (-1) ** k / (2 * k + 1)
return 4 * cumulative_sum
pi = leibniz(1000)
print(f'pi is approximately {pi}')
pi is approximately 3.1425916543395442
c) The accuracy of an approximate value of \(\pi\) is simply measured by calculating the absolute difference from the exact value:
Write a function to evaluate the above expression and use this to determine which of your previous functions is more accurate when evaluated up to the \(50\)th term.
def error(pi_approx):
"""
Calculate the absolute difference between the exact value of pi and an approximation.
Args:
pi_approx (float): The approximation for pi.
Returns:
(float): The absolute difference between the exact and approximate values of pi.
"""
return abs(math.pi - pi_approx)
pi_basel = basel(50)
pi_leibniz = leibniz(50)
print(f'pi - pi_approx (basel) = {error(pi_basel)}')
print(f'pi - pi_approx (leibniz) = {error(pi_leibniz)}')
pi - pi_approx (basel) = 0.01896613065606667
pi - pi_approx (leibniz) = 0.01960595939725751
2. The distance between two points \(i\) and \(j\) in two-dimensional space can be calculated with the following equation:
where \(x_{i}\) is the \(x\)-coordinate of point \(i\), \(x_{j}\) is the \(x\)-coordinate of point \(j\) etc.
a) Write a function to calculate the distance between two points, where each point is a list of \(x\), \(y\) and \(z\) coordinates:
example_point_i = [0.1, 0.2]
example_point_j = [0.7, 0.5]
Test your function by calculating the distance between the example points above.
import math
def distance(point_i, point_j):
"""
Calculate the distance between two points in two-dimensional space.
Args:
point_i (list): The x and y coordinates of the first point.
point_j (list): The x and y coordinates of the second point.
Returns:
(float): The distance between points i and j.
"""
return math.sqrt((point_i[0] - point_j[0]) ** 2 + (point_i[1] - point_j[1]) ** 2)
example_point_i = [0.1, 0.2]
example_point_j = [0.7, 0.5]
distance(example_point_i, example_point_j)
0.6708203932499369
b) The atoms in a molecule of BF\(_{3}\) are located at the following coordinates:
bf3_coordinates = [[0.00, 0.00],
[0.00, 1.30],
[1.13, -0.65],
[-1.13, -0.65]]
This is a nested list (a list of lists), where each sublist contains the coordinates associated with one atom.
Using loops, calculate and print the distances between all pairs of atoms in BF\(_{3}\). Use the function you wrote in a) to calculate each distance.
bf3_coordinates = [[0.00, 0.00],
[0.00, 1.30],
[1.13, -0.65],
[-1.13, -0.65]]
for atom_i in bf3_coordinates:
for atom_j in bf3_coordinates:
print(distance(atom_i, atom_j))
0.0
1.3
1.3036103712382776
1.3036103712382776
1.3
0.0
2.2537524265100637
2.2537524265100637
1.3036103712382776
2.2537524265100637
0.0
2.26
1.3036103712382776
2.2537524265100637
2.26
0.0
c) Rather than using the print function to display all of your calculated distances, modify the code you wrote for b) so that the distances are instead used to build a distance matrix:
where \(d_{12}\) is the distance between atom \(1\) and atom \(2\), \(d_{14}\) is the distance between atoms \(1\) and \(4\) etc.
In Python, your final distance matrix should be a nested list:
[[d_11, d_12, d_13, d_14],
[d_21, d_22, d_23, d_24],
[d_31, d_32, d_33, d_34],
[d_41, d_42, d_43, d_44]]
You can build up the matrix one row at a time using the loops you wrote in b).
D = []
for atom_i in bf3_coordinates:
row = []
for atom_j in bf3_coordinates:
row.append(distance(atom_i, atom_j))
D.append(row)
# Loop over sublists so that it prints nicely.
for sublist in D:
print(sublist)
[0.0, 1.3, 1.3036103712382776, 1.3036103712382776]
[1.3, 0.0, 2.2537524265100637, 2.2537524265100637]
[1.3036103712382776, 2.2537524265100637, 0.0, 2.26]
[1.3036103712382776, 2.2537524265100637, 2.26, 0.0]
3. Consider the following blocks of code. Each one of them will cause an error when run. Taking each block one-by-one:
Predict which line of code will cause the error and why this line is problematic.
Run the block of code in a Jupyter notebook and see whether or not your were correct.
Modify the code to solve the problem and thus remove the error.
Block 1:
def calculate_product(sequence):
"""
Calculate the product of a sequence of numbers (n_1 * n_2 * n_3 * ...).
Args:
sequence (iterable): The sequence used to calculate the product.
Returns:
(float): The product.
"""
product = sequence[0]
for number in sequence[1:]:
product *= number
return product
product_of_sequence = calculate_product(5.12, 6.79, 2.22, 0.19, 13.6)
print(f'The cumulative product = {product_of_sequence}')
def calculate_product(sequence):
"""
Calculate the product of a sequence of numbers (n_1 * n_2 * n_3 * ...).
Args:
sequence (iterable): The sequence used to calculate the product.
Returns:
(float): The product.
"""
product = sequence[0]
for number in sequence[1:]:
product *= number
return product
product_of_sequence = calculate_product([5.12, 6.79, 2.22, 0.19, 13.6])
print(f'The cumulative product = {product_of_sequence}')
The cumulative product = 199.427579904
Block 2:
Note
This block also contains a problem which, whilst it will not cause an error, will lead to the wrong result. See if you can spot and fix this problem in addition to removing the error.
def check_palindrome(convert_to_smalls=False, word):
"""
Check if a word is a palindrome.
Args:
word (str): The word to be checked.
Returns:
(bool): Whether or not the word is a palindrome.
"""
if convert_to_smalls:
word = word.lower() # Convert to small-case.
if word == word[::-1]:
check = True
else:
check = False
example_word = 'Civic'
if check_palindrome(example_word, convert_to_smalls=True):
print(f'{example_word} is a palindrome.')
else:
print(f'{example_word} is not a palindrome.')
def check_palindrome(word, convert_to_smalls=False):
"""
Check if a word is a palindrome.
Args:
word (str): The word to be checked.
Returns:
(bool): Whether or not the word is a palindrome.
"""
if convert_to_smalls:
word = word.lower() # Convert to small-case.
if word == word[::-1]:
check = True
else:
check = False
return check
example_word = 'Civic'
if check_palindrome(example_word, convert_to_smalls=True):
print(f'{example_word} is a palindrome.')
else:
print(f'{example_word} is not a palindrome.')
Civic is a palindrome.
Block 3:
import scipy
def harmonic_energy(v, nu):
"""
Calculate the energy of a quantum harmonic oscillator.
Args:
v (int): Quantum number labelling each energy eigenstate.
nu (float): Vibrational frequency of the oscillator.
Returns:
(float): The energy of the oscillator.
"""
energy = (v + 0.5) * scipy.constants.hbar * nu
return energy
frequency_HCl = 5.40e14
harmonic_energy(0, frequency_HCl)
print(f'The zero-point energy of HCl = {energy} J')
import scipy
def harmonic_energy(v, nu):
"""
Calculate the energy of a quantum harmonic oscillator.
Args:
v (int): Quantum number labelling each energy eigenstate.
nu (float): Vibrational frequency of the oscillator.
Returns:
(float): The energy of the oscillator.
"""
energy = (v + 0.5) * scipy.constants.hbar * nu
return energy
frequency_HCl = 5.40e14
energy_HCl = harmonic_energy(0, frequency_HCl)
print(f'The zero-point energy of HCl = {energy_HCl} J')
The zero-point energy of HCl = 2.8473439076446224e-20 J