Numerical Algorithms

Introduction

In today’s lab, we are going to use recurrence relations to analyze a numerical algorithm. Numerical algorithms are used for solving problems in computational mathematics. For example, your computer doesn’t naturally “just know” how to compute things like \(\sin{x}\), \(\sqrt{x}\), or \(\log x\) – people need to write algorithms for computing these quantities. Making these algorithms reliable and fast is one of the priorities of the field called numerical analysis, which sits at the intersection of mathematics and computer science.

Warmup (~50 mins)

Problem 1

Introduction

Consider the following algorithm. This algorithm is used to compute a familiar mathematical quantity. It does so by repeatedly computing the sequence defined by the recurrence relation

\[ \begin{aligned} x_{n+1} = \frac{1}{2}\left(x_n + \frac{a}{x_n}\right) \end{aligned} \]

for some input number \(a\). We do this until the difference \(|x_{n+1} - x_n|\) is smaller than some \(\epsilon\). The algorithm outputs the final value of \(x\) after this threshold is reached. Here’s a more formal description of this algorithm in pseudocode:

Algorithm 1

Inputs: nonnegative real number \(a\), small positive real number \(\epsilon > 0\).

Initialize \(x = a + 1\), \(\Delta = 2\epsilon\).
While \(\Delta > \epsilon\):
- \(x_\mathrm{new} \gets \frac{1}{2} \left(x + \frac{a}{x}\right)\).
- \(\Delta \gets |x_\mathrm{new} - x|\).
- \(x \gets x_\mathrm{new}\).
Return \(x\).

Part A

Write a function in Python or Java that implements Algorithm 1. The user should be able to specify the values of \(a\) and \(\epsilon\).

Notes:

It is possible but not required (or easier) to do this with a recursive function.
Set the default value of \(\epsilon\) at \(10^{-16}\). In Python this can be done by writing epsilon = 1e-16 in the appropriate location.
It’s fun, though not required, to print the value of \(x\) at each iteration.

Part B

What familiar mathematical quantity does this function approximately compute? You don’t have to reason about this function mathematically – just try a few values of a and make an informed guess.

Note: proving the answer to this question, and describing how well the function approximates its target, is a classical question in numerical analysis.