MAT321 Precept 1¶

Q1: Fast fixed-point¶

Suppose $x^* \in [ a, b]$ is a fixed point of $g$ and $h$, and $g$, $h \in C^2([a,b])$ (which implies that $g'', h''$ are bounded on $[a,b]$) and $g'(x^*)\neq h'(x^*)$.

  1. Show that there exists constant $\alpha$, $\beta$ such that the fixed point iteration:

$$ x_{k+1} = \alpha g(x_k) + \beta h( x_k) $$ converges quadratically for initial points $x_0$ close enough to $x^*$. That is, show that it converges to $x^*$ locally, and that the convergence is quadratic.

  1. Do the fixed points iterations $x_{k+1} = g(x_k)$ and $x_{k+1} = h(x_{k+1})$ need to converge for the iteration above to converge?

Task 1: Implement the Bisection Method¶

Problem Description¶

The bisection method is a root-finding method that repeatedly bisects an interval and then selects a subinterval in which a root must lie for further processing. This is done by verifying where the sign change occurs in the function.

Instructions:¶

  • Implement the bisection method as a Python function my_bisection(f, a, b, tol, n).
  • The inputs are:
    • f: A Python function where the root is being sought.
    • a, b: Interval endpoints such that the function changes sign between a and b.
    • tol: A tolerance for stopping the method.
    • n: Maximum number of iterations.
  • The method should return after (b - a < tol) or after a maximum number of iterations n.
In [ ]:
import numpy as np
def my_bisection(f, a, b, tol, n):
    return a

# Test case for bisection method
f = lambda x: np.sin(x - 2.3)
root = my_bisection(f, 1, 3, 1e-14, 100)
assert np.abs(root - 2.3) < 1e-14, f"Bisection method failed, got {root}"
print("Bisection method test passed.")

Task 2: Compute $\sqrt{a}$¶

Problem Description¶

Use your bisection method to compute $\sqrt{a}$ for a positive value of $a$.

In [ ]:
def compute_sqrt(a):
    return a

# Test case for sqrt computation using bisection
a = 2
sqrt_a = compute_sqrt(a)
assert np.abs(sqrt_a - np.sqrt(a)) < 1e-14, f"Sqrt computation failed, got {sqrt_a}"
print("Sqrt computation test passed.")

Task 3: Implement Newton's Method¶

Problem Description¶

Newton's method is a way to find (x) such that (f(x) = 0). It is defined as simple iteration on the function:

$$ x_{k+1} = x_k - \frac{f(x_k)}{f'(x_k)} $$

Instructions:¶

  • Implement Newton's method as a Python function newton_method(f, df, x0, tol, n).
  • Apply Newton's method to the function $ f(x) = x^2 - a $ to compute $\sqrt{a}$.
  • Plot the error on each iteration of Newton’s method until the error is around (10^{-16}).
In [ ]:
import matplotlib.pyplot as plt
# Test case for Newton's method
a = 4
f = lambda x: x**2 - a  # Root at x = 2
df = lambda x: 2 * x
# root_newton = newton_method(f, df, 3, 1e-6, 100)


def newton_method(f, df, x0, tol, n):
    return x0

root_newton = newton_method(f, df, 3, 1e-14, 100)
assert np.abs(root_newton - np.sqrt(a)) < 1e-14, f"Newton's method failed, got {root_newton}"
print(f"Computed sqrt({a}) = {root_newton}, Test passed.")

# Use this function to plot the error
# help(plt.semilogy)

Task 4: Chebyshev Polynomials¶

Problem Description¶

The Chebyshev polynomials of the first kind are defined by the recurrence relation:

$$ T_0(x) = 1, \quad T_1(x) = x, \quad T_{n+1}(x) = 2xT_n(x) - T_{n-1}(x) $$

Instructions:¶

  • Write a function mycheb(n, x) that computes the (n)-th Chebyshev polynomial at a value $x \in [-1, 1] $ using the above recurrence relation.
In [ ]:
from numpy.polynomial.chebyshev import Chebyshev 

def mycheb(n, x):
    Tn = 1
    return Tn

# Test case for Chebyshev polynomials
cheb_poly = mycheb(3, 0.5)
expected_value = Chebyshev.basis(3)(0.5)  # T_3(0.5) should be -0.125
assert abs(cheb_poly - expected_value) < 1e-14, f"Test failed, got {cheb_poly}"
print("Chebyshev polynomial test passed.")