Homework 3

Student's Name: PLEASE INSERT YOUR NAME HERE (DOUBLE CLICK THIS BOX TO EDIT.)

Directions: Add work to this notebook to solve the problems below.

Check your work. I give some tests that should be passed if you get it correct.

Your notebook should run without printing errors and without user input. If this is not the case, points will 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 [ ]:
import numpy as np
import matplotlib.pyplot as plt
import math as m

1. Geometric Series

A geometric series is a series with a constant ratio $r$ between successive terms. Let $c$ be the first term in the series. Then the geometric series with ratio $r$ and first term $c$ is $$c + cr + cr^2 + cr^3+\ldots=\sum_{k=0}^\infty c r^k.$$

Write a function geometric_partial_sum(c, r, N) that returns the sum of the first $N$ terms of this series. That is, you should return $$\sum_{k=0}^{N-1} c r^k.$$ You should work with floats. The test included is insufficient. Test your code!

In [ ]:
# Tests
x = geometric_partial_sum(1, 0.5, 2)
assert abs(x - (1+ 1/2))<2**-50, "This series sum (1, 0.5, 2) is incorrect."
# Tests
x = geometric_partial_sum(5, 2, 3)
assert abs(x - (5+10+20))<2**-50, "This series sum (5, 2, 3) is incorrect."

2. Number of terms needed

It is well known that if $|r|<1$, then the geometric series converges and satisfies $$\frac{c}{1-r} = \sum_{k=0}^\infty c r^k.$$ You probably learned this in calculus.

Write a function geometric_terms_needed(c, r, epsilon) that takes as input a constant $c$, an $r$ satisfying $|r|<1$, and an $\epsilon>0$. The function should returns the smallest number of terms $N$ such that the remainder is less than $\epsilon$ in absolute value.

In [ ]:
# Tests
c = 1
r = 0.5
epsilon = 0.0001
N = geometric_terms_needed(c,r,epsilon)
remainder = 2 - geometric_partial_sum(c, r, N)
assert abs(remainder) < epsilon, "Remainder is not small enough."
remainder = 2 - geometric_partial_sum(c, r, N-1)
assert abs(remainder) >= epsilon, "Previous remainder is not big enough."

3. $\cosh x$

The function $cosh(x) = \frac{e^x + e^{-x}}{2}$ is analytic with Taylor series $$\cosh x = 1 + \frac{x^2}{2!} + \frac{x^4}{4!} + \ldots = \sum_{k=0}^\infty \frac{x^{2k}}{(2k)!}$$ What degree Taylor polynomial is needed to estimate $\cosh(1/2)$ with a truncation error less than $\epsilon>0$? Write a function cosh_half_degree(epsilon) that returns this degree. Write a second function cosh_half(epsilon) that returns the value of the degree $d$ Taylor polynomial evaluated at $\frac{1}{2}$.

Check that your answer achieves the desired accuracy with $\epsilon=10^{-8}$ by comparing your answer to the value of $cosh(\frac{1}{2})$ returned by the cosh function in the math library (i.e., m.cosh(0.5)).

(This is similar to Problem 5 from Chapter 2 of TAK.)

In [ ]:
# Some tests you should pass.
assert type(cosh_half_degree(0.1)) == int, "cosh_half_degree should return an integer."
assert type(cosh_half(0.5)) == float, "cosh_half should return a float."

4. Taylor Polynomials for the exponential function

Write a function exp_poly which takes as input a degree $d$ and returns a function that evaluates the Taylor polynomial of degree $d$ for the exponential function centered at zero.

Check your answers against the first few Taylor polnomials by hand.

(Hint: You are asked to do something similar to the sin_polynomial function in the Taylor Series notebook.)

In [ ]:
# Define $p$ to be the Taylor polynomial of degree 1:
p = exp_poly(1)
# Note that the Taylor polynomial should be p(x)=x+1
assert p(5.5)==6.5, "Unexpected result!"

5. Computing partial sums

For $n \geq 1$, let $a_n = \frac{2}{n(n+1)(n+2)}$. And consider the series $$\sum_{n=1}^\infty a_n = \sum_{n=1}^\infty \frac{2}{n(n+1)(n+2)}.$$

Write a function partial_sum_3 which takes as input an integer $N \geq 1$ and computes the partial sum consisting of the sum of the first $N$ terms in the series. Use Python floats for computing the series.

In [ ]:
# Some basic tests
assert abs(partial_sum_3(1) - 2/6) < 2**-50, \
    "The first partial sum is wrong."
assert abs(partial_sum_3(2) - (2/6+2/24)) < 2**-50, \
    "The second partial sum is wrong."

6. Mathematical limits

Print out the first 10 partial sums in the following format:
Partial sum 1 is 0.3333333333333333
Partial sum 2 is 0.41666666666666663
...

Conjecture the limit $L$ and store the value in the variable limit_prob_4. Recall that the $N$-th remainder is the difference $L-S_N$ where $S_N$ is the $N$-th partial sum $S_N$. Print out the first 10 values of the conjectured remainder in the following format:
Remainder 1 is ?.??????????????????
...

Conjecture a formula for the remainder. (Hint: You may also want to print the reciprocals of the remainders.) Give a rigorous inductive proof of this formula in a Markdown cell in the notebook. (Including a picture of your hand-written proof is okay.)

Explain why the proof of your remainder formula implies that your conjecture for the limit is correct.

7. Zeta(3)

For $s > 1$, the Reimann zeta function is given by the convergent series $$\zeta(s)=\sum_{n=1}^\infty \frac{1}{n^s}=\frac{1}{1^s}+\frac{1}{2^s}+\frac{1}{3^s}\ldots.$$ In the Series notebook, we calculated that $\zeta(2)=\frac{\pi^2}{6}$. Now we will approximate $\zeta(3)$.

Let $b_n=\frac{1}{n^3}$. Then $\zeta(3)=\sum_{n=1}^\infty b_n$. Define the partial sum and remainder $$S'_N=\sum_{n=1}^N b_n \quad \text{and} \quad R'_N=\sum_{n=N+1}^\infty b_n.$$

For $n$ large enough, find indices $l(n)$ and $u(n)$ such that $$a_{l(n)} \leq 2 b_n \leq a_{u(n)},$$ where $a_k = \frac{2}{k(k+1)(k+2)}$ as in the problems above. Use these inequalities to find lower and upper bounds for $R'_N$ in terms of the sequence of remainders $R_N$ from the series $\sum a_k$. Include a Markdown cell explaining what you have found.

Use your results to write a function zeta_3(epsilon) which takes as input a positive real number $\epsilon>0$ and returns a pair consisting of a lower bound $L$ and an upper bound $U$ for $\zeta(3)$ such that $U-L \leq \epsilon$. You may use floats and ignore the effects of round off error.

(Remark. To return a pair $(L,U)$ you can use the statement return (L,U) at the end of your function.)

In [ ]:
# Some simple tests
L,U = zeta_3(0.001)
assert type(L)==float, "The lower bound should be a float."
assert type(U)==float, "The upper bound should be a float."
assert L<=U, "The lower bound should be smaller than the upper bound."
assert U-L<=0.001, "The gap should be smaller than 0.001."