Homework 10

Student's Name: PLEASE INSERT YOUR NAME HERE

Directions: Add work to this notebook to solve the problems below. Some other cells have been included which are important for testing your work and should give you some feeback.

Check your work.

Your notebook should run without printing errors and without user input. If this is not the case, points may be deducted from your grade.

Problem Sources:

  • LL: Programming for Computations - Python by Svein Linge and Hans Petter Langtangen, 2nd edition.
  • L: A Primer on Scientific Programming with Python by Hans Petter Langtangen, 2nd edition.
  • TAK: Applied Scientific Computing With Python by Peter R. Turner, Thomas Arildsen, and Kathleen Kavanagh.
In [ ]:
# Standard imports
import numpy as np
import matplotlib.pyplot as plt
import math as m
from mpmath import mp, iv
from scipy import linalg

from numpy import random

import random as random_number

1. Coins that land on their edge

According to an authorative youtube video a nickel has a $1/6000$ chance of landing on its edge. Write a function nickel_flip() that takes no inputs and returns one of three strings, "Edge", "Heads" and "Tails". The function should return "Edge" with probability $1/6000$ and return "Heads" and Tails" with equal probability.

In [ ]:

2. Video game character's position

A video game character lives in a $2$-dimensional integer lattice. That is, his coordinates are always a pair of integers $(x, y)$. He can move in the four compass directions one unit at a time.

Write a class Position whose initialization method takes two parameters (in addition to self), the initial x-coordinate x0 and the initial y-coordinate y0. Both should be integers. If p is a Position object, it should have the following methods:

  • p.north() should increase the character's y-coordinate by one.
  • p.east() should increase the character's x-coordinate by one.
  • p.south() should decrease the character's y-coordinate by one.
  • p.west() should decrease the character's x-coordinate by one.
  • p.beam(xb, yb) should move the character directly to position (xb, yb)
  • p.x() should return the characters x-coordinate.
  • p.y() should return the characters y-coordinate.
  • p.coordinates() should return the character's coordinates as a pair.
In [ ]:
In [ ]:
# Test
p = Position(3, 5) # Start at position (3,5)
assert p.x() == 3
assert p.y() == 5

p.south()
assert p.x() == 3
assert p.y() == 4

p.west()
p.south()
assert p.coordinates() == (2, 3)

p.beam(-2, 3)
p.north()
p.east()
assert p.coordinates() == (-1, 4)

3. Keeping track of statistical data

Write a class named StatisticsTracker which keeps track of statistical quantities related to real number measurements. No data should be passed to the initializer of the class. A StatisticalTracker object st should have the following methods available:

  • st.number() should return the number of measurements received so far. (This should be a non-negative integer.)
  • st.mean() should return the mean of the measurements received so far.
  • st.variance() should return the variance of the measurements received so far.
  • st.new_measurement(x) should update the internal variables of the class so that the quantities returned by the above methods are correct when adding the float x to the prior received data set. This method should return nothing.

Constraint: You are not to store all the measurements received. Instead just store quantities necessary to keep track of the statistical quantities. Some hints are below.

To keep track of the number of measurements, your class should devote a variable to keeping count of the number of measurements received.

Recall that the mean of measurements $x_0, \ldots, x_{n-1}$ is $$\mu = \frac{1}{n} \sum_{i=0}^{n-1} x_i.$$ The number $n$ is the total number of measurements. So, it also makes sense to keep track of the sum of the measurements.

The variance should be computed by assuming all measurements made are equally likely. Thus the variance is given by the formula $$Var = \frac{1}{n} \sum_{i=0}^{n-1} (x_i - \mu)^2.$$ Since the mean $\mu$ will change as more data is received, it is best to treat $\mu$ as an unknown variable. Expanding the above formula we see: $$Var = \frac{1}{n} \left(\sum_{i=0}^{n-1} x_i^2\right) - \frac{2 \mu}{n}\left(\sum_{i=0}^{n-1} x_i\right) + \mu^2.$$ Thus it also makes sense to keep track of the sums of the squares of the measurements.

In [ ]:
In [ ]:
# Test
st = StatisticsTracker()
assert st.number() == 0, "The number of tests should be 0."

st.new_measurement(5.0) # Add one measurement of 5.0
assert st.number() == 1, "The number of tests should be 1."
assert st.mean() == 5.0, "The mean should be 5."
assert st.variance() == 0.0, "The variance should be 0."

st.new_measurement(1.0) # Add a second measurement of 1.0
assert st.number() == 2, "The number of tests should be 2."
assert st.mean() == 3.0, "The mean should be 3."
assert st.variance() == 4.0, "The variance should be 2."

st.new_measurement(7.0)  # Add a third measurement of 7.0
st.new_measurement(-4.0) # Add a fourth measurement of -4.0
assert st.number() == 4, "The number of tests should be 2."
assert st.mean() == 2.25, "The mean should be 2.25."
assert st.variance() == 17.6875, "The variance should be 17.6875."

4. Taylor Series

Write a class TaylorSeries so that the initializer of a TaylorSeries object takes two parameters, a coefficient function c mapping degrees (non-negative integers) to coefficients (floating point real numbers) and a center x0. The initialized object should represent the Taylor Series (expressed as a function of $x$): $$\sum_{d=0}^{+\infty} c_d \cdot (x - x_0)^d.$$ (Here we have used $c_d$ to represent c(d) and $x_0$ for x0.)

A TaylorSeries object ts should have the following methods:

  • ts.center() should return the center $x_0$.
  • ts.coefficient(d) should return the coefficient of the degree $d$ term. That is, this function should return $c_d$ in the expression above.
  • ts.partial_sum(x,N) should return the finite sum $$\sum_{d=0}^{N} c_d \cdot (x - x_0)^d.$$
  • ts.derivative() should return the derivative of this Taylor Series as a TaylorSeries. The derivative of the series is $$\sum_{d=0}^{+\infty} (d+1) c_{d+1} \cdot (x - x_0)^d.$$

Remark: To define a derivative you should first define a coefficient function for the derivative series using currying. Then you can pass that new function to the TaylorSeries initializer.

In [ ]:
In [ ]:
# Test

# Construct the exponential function
c = lambda d: 1.0/m.factorial(d)
ts = TaylorSeries(c, 0)

assert ts.partial_sum(0, 10) == 1.0
assert abs( ts.partial_sum(1.0, 15) - m.e ) < 10**-10, \
    "ts.partial_sum(1.0, 15) should be close to e."
assert abs( ts.partial_sum(-1.0, 15) - 1/m.e ) < 10**-10, \
    "ts.partial_sum(-1.0, 15) should be close to 1/e."

# The derivative should also be the exponential function
der = ts.derivative()

for d in range(10):
    assert abs( der.coefficient(d) * m.factorial(d) - 1) < 10**-8, \
        "The coefficient of der of degree {} is wrong.".format(d)

# Same checks as above.
assert der.partial_sum(0, 10) == 1.0
assert abs( der.partial_sum(1.0, 15) - m.e ) < 10**-10, \
    "der.partial_sum(1.0, 15) should be close to e."
assert abs( der.partial_sum(-1.0, 15) - 1/m.e ) < 10**-10, \
    "der.partial_sum(-1.0, 15) should be close to 1/e."

5. Random coefficients

We'd like to construct a random function using Taylor series. To define a function in terms of Taylor series you just need to determine the coefficients. We'll choose the coefficients $c_d$ at uniformly from the interval $[\frac{-1}{d!}, \frac{1}{d!}]$. (This seems like a weird choice but guarantees that every derivative $f^{(n)}(x_0)$ is taken at random from $[-1,1]$.

The following is a first attempt:

In [ ]:
def random_coef(d):
    bound = 1 / m.factorial(d)
    return 2*bound*np.random.random_sample() - bound

Unfortunately it doesn't work as the following plot demonstrates.

In [ ]:
ts = TaylorSeries(random_coef, 0)
x = np.linspace(-5,5,1000)
y = [ts.partial_sum(val, 20) for val in x]
plt.plot(x,y)

What is happening is that every time the TaylorSeries object needs a coefficient, it calls random_coef. Each time this occurs it generates a new random number. But the coefficicient of degree d should be a constant. Observe:

In [ ]:
c3 = random_coef(3)
print("The coeficient of degree 3 is {}.".format(c3))
c3 = random_coef(3)
print("The coeficient of degree 3 is {}.".format(c3))

To fix this, we will create a function class that caches (stores) coefficients that it has already selected.

Construct a class RandomCoefficients whose initializer takes no input. An object rc of type RandomCoefficients should have a __call__ method so that rc(d) returns a random number between $\frac{-1}{d!}$ and $\frac{1}{d!}$. The call method should return the same random number whenever called with the same value of d.

Suggestion: You might want to have a RandomCoefficients store a dictionary. When rc(d) is called for the first time with input d, the random number should be computed and stored as the image of d under the dictionary. Then when rc(d) is called a second time, the number should be restored from the dictionary and returned.

In [ ]:
In [ ]:
# Test

rc = RandomCoefficients()
val = rc(0)
assert -1 <= val <= 1
assert val == rc(0), "The return value of rc(0) has changed."

val = rc(2)
assert -1/2 <= val <= 1/2
assert val == rc(2), "The return value of rc(2) has changed."

assert val != rc(1), "The return value of rc(1) should differ from rc(2)."

6. Random Functions

Use objects created by the RandomCoefficients class to create some random functions. Experiment with plotting these functions.

In [ ]: