• Floating Point Arithmetic
  • Double Precision Floating Point
    • Special Numbers
    • Subnormal numbers
    • Error
  • Arithmetic Operations
    • Multiplication
    • Addition
    • Subtraction
• Condition of a problem
• Stability of an Algorithm
• Finding Roots
  • Intermediate Value Theorem
  • Secant Method
    • Convergence
  • Newton's Method
    • Convergence
  • Bisection Method
    • Convergence
  • Fixed point iteration
    • Uniqueness
    • Convergence
  • Implementation
  • Condition Number
  • Matrices
• Finding Minimums
  • Implementation
    • Golden Section Search
    • Parabolic Interpolation
    • Solving for minimum of parabola
• Interpolation
  • Power basis
  • General matrix form
  • Lagrange basis
• Newton Form Interpolating Polynomial
  • Divided difference
  • Interpolation with Derivatives
  • Error
  • Horner's Method
• Chebyshev Interpolating points.
• Piecewise Interpolation
  • Hermite Interpolation
• Numerical Integration
• Midpoint Method
  • Error
    • MVT for Integrals
• Trapezoidal Method
  • Error
• Simpson's Method
  • Error
• Degree of Precision
• Composite Midpoint Rule
• Composite Trapezoidal Rule
• Composite Simpson's Rule
• Gaussian Quadrature
• Approximating Functions
  • Properties of Norms
  • Least Squares
  • Chebyshev Polynomial
    • Alternative Representation
• todo

• ascher and greif, a first course in numerical methods

Floating Point Arithmetic

Double Precision Floating Point

binary64 - specified by IEEE 754

• 1 sign bit
• 11 exponent bits

• 52 fraction bits

  64 bits total

  \(fl(x) = (-1)^{d_s}(1 + d_1 d_2 ... d_52) \cdot 2^{exp - 1023}\)

  This shifted exponent is known as a biased exponent and gives us a range of exponents \([-1023, 1024]\). However, we reserve -1023 and 1024 to represent NaN and Inf, respectively.

  \(realmin = (1.0000...0) \cdot 2^{-1022}\) \(realmax = (1.1111...1) \cdot 2^{1023}\)

Special Numbers

• \(0(0...0)_2 \cdot 2^{1...1} = +Inf\)
• \(1(0...0)_2 \cdot 2^{1...1} = -Inf\)
• \(0(\neq 0)_2 \cdot 2^{1...1} = +NaN\)
• \(1(\neq 0)_2 \cdot 2^{1...1} = -NaN\)
• \(0(0...0)_2 \cdot 2^{0...0} = +0\)
• \(1(0...0)_2 \cdot 2^{0...0} = -0\)

#+end_examples

Subnormal numbers

\(fl(x) = (-1)^{d_s} (0 + d_1 d_2 d_3 ... d_n)_2 \cdot 2^{exp}\)

• \(0(\neq = 0)_2 \cdot 2^{0...0} = subnormal\)

Smallest exp - -1074 Largest exp - 1023

Error

In Octave, exp() gives the distance between a number and the next smallest floating point number.

Arithmetic Operations

Let \(x,y \in \mathbb{R}\)

Multiplication

\(f(x)\) can be expressed as \(x(1 + \delta)\) where \(|\delta| \leq \eta\)

\(f(x \cdot y) = (x(1 + \delta_x) \cdot y(1 + \delta_y)){(1 + \delta_m)}{rounding after multiplication}\) = xy(1 + δ_x + δ_y + δ_m) + {O(η²)}{≤ c ⋅ η²}

\(\frac{|x \cdot y - f(x \cdot y)|}{|x \cdot y|} = |\delta_x + \delta_y + \delta_m| + O(\eta^2) \leq 3 \eta + O(\eta^2)\)

Addition

\(f(x + y) = \[x(1 + \delta_x) + y(1 + \delta_y)\]\underbrace{(1 + \delta_a)}_{\text{error after rounding addition result}} &= x(1 + \delta_x)(1 + \delta_a) + y(1 + \delta_y)(1 + \delta_a) = x(1 + \delta_x + \delta_a) + y(1 + \delta_y + \delta_a) + O(\eta^2) \\ &= (x + y) + x(\delta_x + \delta_a) + y(\delta_y + \delta_a) + O(\eta^2)\)

\(|\frac{f(x + y) - (x + y)}{x + y}| \leq |\frac{x( \delta_x + \delta_a)}{x + y}| + |\frac{x( \delta_x + \delta_a)}{x + y}|\)

If \(x\) and \(y\) have same sign, \(|\frac{x}{x+y}| \leq 1\) and \(|\frac{y}{x+y}| \leq 1\), then $ ≤|δ_x + δ_a| + |δ_y + δ_a| ≤ 4 η$

But if \(x\) and \(y\) have different sign, then \(|x - y|\), can be very small and potentialy causes the result to blow up.

We conclude that the result of any floating point arithmetic operation must be equal to the result using infinite precision then rounding to \(t\) binary digits.

Subtraction

\(100.0 - 99.99 \Rightarrow 1.000 \cdot 10^2 - 0.9999 \cdot 10^2 = 0.001 \cdot 10^2\)

\(fl(x-y) = .001 \cdot 10^2 = 0.1\) \(\text{relative error} = \frac{0.1 - 0.01}{0.01} = 9\)

with guard digit

\(100.0 - 99.99 \Rightarrow 1.0000 \cdot 10^2 - 0.9999 \cdot 10^2 = 0.01 \cdot 10^2\)

Condition of a problem

Stability of an Algorithm

So far, we've defined two types of error. These errors require that we have the exact answer to the problem, \(y\).

For some problems, we don't have access to the exact answer, so we instead compute a different kind of error.

\(y_n = \int_0^1 \frac{x^n}{x + 10} dx\) - the thing were trying to calculate

\(y_n + 10y_{n-1} = \frac{1}{n}, n \geq 1\) - the problem \(y_0 = \ln(11) - \ln(10)\)

when we calculate \(y_n\) using floating point arithmetic, this is the algorithm

Finding Roots

\(y_{i+1} = -f(y_i)/f'(y_i)\)

Intermediate Value Theorem

Intermediate value theorem

Secant Method

Convergence

Newton's Method

Convergence

Bisection Method

Convergence

Fixed point iteration

This can be interpreted \(x_{k+1} = g(x_k)\) When the solution \(x^*\) is found, \(x^* = g(x^*)\)

There are many functions \(g\) such that \(x^* = g(x^*)\), like \(g(x) = x + \alpha f(x)\)

Assume \(g(a) - a > 0\) and \(g(b) - b < 0\)

From IVT, there exists some \(x^* \in [a,b]\) such that \(g(x^*) - x^* = 0\)

Uniqueness

Convergence

Implementation

Root finding scripts combine multiple methods to achieve speed and reliability in finding roots of \(f\)

1. check \(f\) at an initial guess \(x_1\) and finds two points \(a,b\) around \(x_1\) such that the signs of \(f(a)\) and \(f(b)\) are different
2. iterate newtons method a few times until error stops decreasing significantly
   • (for example, iterate until \(\text{curr. err} \leq \frac{\text{prev. err.}}{2} = |x_{k+1} - x_k | \leq \frac{|x_k - x_{k-1}|}{2}\)

Condition Number

We can't use the relative condition number, \(k_r(x) = \frac{|x| \cdot |f'(x)|}{|f(x)|}\), since \(k_r \rightarrow \infty\) as \(f(x) \rightarrow 0\).

Instead we use the absolute condition number to quantify the problem.

\[\log_{10} | \thing{\frac{y - \hat{y}}{y}}{10^{-p}} | = \log_{10} k_r + \log_{10} | \thing{\frac{x - \hat{x}}{x}}{10^{-16}} |\]

\(p = 16 - \log_{10} k_r\)

When finding roots, we find \(x\) given \(y\), so \(x = f^{-1}(y) \rightarrow k_r = \frac{|y| \cdot |f^{-1}'(y)|}{ |f^{-1}(y)| }\).

\(f^{-1}'(y) = \frac{d}{dy} f^{-1} (y) = \frac{dx}{dy} = \frac{1}{dy/dx} = \frac{1}{f'(x)}\)

\(k_{r,f^{-1}} = \frac{|f(x)|}{|x| \cdot |f'(x)|}\)

\(\hat{y} - y = (\hat{x} - x) f'(\psi)\) \(|\hat{y} - y| = |(\hat{x} - x)| |f'(\psi)|\) where \(\psi \in [x,\hat{x}]\)

Approximate \(f'(\psi)\) with \(f'(x)\)

\(|\hat{y} - y| \approx |f'(x)| \cdot |\hat{x} - x|\)

\(|\hat{x} - x| = \frac{1}{|f'(x)|} |\hat{y} - y|\)

\(k_a = \frac{1}{|f'(x)|}\)

Matrices

\[Ax = b\]

\(x_1 + x_2 = 3\) \(x_1 - x_2 = 1\)

\[\begin{matrix} 1 & 1 \\ 1 - 1 \end{matrix} \cdot \begin{matrix} x_1 \\ x_2 \end{matrix} = \begin{matrix} 3 \\ 1\]

Relative condition: A(x + δ x) = b + δ b

\[\frac{|| \delta x ||}{|| x ||} = k_r \cdot \frac{|| \delta b ||}{|| b ||}\]

\(k_r = || A || \cdot || A^{-1} ||\)

\[|| A ||_1 = max_{1 \leq j \leq n} \sum_{i=1}^n |a_{ij}|\] - max column sum norm \[|| A ||_\infty = max_{1 \leq i \leq n} \sum_{j=1}^n |a_{ij}|\] - max row sum norm

Finding Minimums

f = @(x) x.^2 - x - 2; # parabola with minimum at x = .5

x = fminbnd(f, -2, -3)

Assume we are trying to minimize a function \(\phi (x)\), where \(x^*\) is the location of the minimum (minimizer).

Assume \(\phi (x)\) is twice differentiable.

Writing \(\phi\) as a taylor series: \[\phi (x) = \phi (x^*) + (x - x^*) \phi '(x^*) + \frac{(x - x^*)^2}{2} \phi ''(x^*) + O((x - x^*)^3)\]

Since \(x^*\) is a minimizer, \(\phi '(x^*) = 0\), so \[\phi (x) = \phi (x^*) + \frac{(x - x^*)^2}{2} \phi ''(x^*) + O((x - x^*)^3)\]

If \(\phi ''(x^*) > 0\), then \(x^*\) is a minimizer, and if \(\phi ''(x^*) < 0\), then \(x^*\) is a maximizer. If \(\phi ''(x^*) = 0\), there is not enough information to know.

Implementation

A common minimum finder employs two methods:

1. Golden section search
2. Parabolic interpolation

f = @(x) x.^2 + 4 * cos(x);
fplot(f, [0 3]);
options = optimset('Display', 'iter');
[xmin, fmin] = fminbnd(f, 0,3, options);

octave:4>     [xmin, fmin] = fminbnd(f, 0,3, options);

Func-count     x          f(x)    Procedure
    1        1.85410    2.319570 initial
    2        2.29180    2.611785 golden
    3        1.83004    2.323648 parabolic
    4        1.88882    2.316881 parabolic
    5        1.89619    2.316809 parabolic
    6        1.89554    2.316808 parabolic
    7        1.89549    2.316808 parabolic
    8        1.89549    2.316808 parabolic
    9        1.89549    2.316808 parabolic
   10        1.89549    2.316808 parabolic
   11        1.89549    2.316808 golden

Optimization terminated:
the current x satisfies the termination criteria using OPTIONS.TolX of 1.000000e-08

Golden Section Search

Assume \(f\) is a unimodal function (only one minimum).

\(\tau = \frac{\sqrt{5} - 1}{2} \approx 0.618\)

1. Begin with two points, \(x_1 = a + (1 - \tau)(b - a)\), \(x_2 = a + \tau (b - a)\)

   case 1: if \(f(x_1) < f(x_2)\), then there must be a minimum in the range \([a, x_2]\)

   case 2: if \(f(x_1) > f(x_2)\), then there must be a minimum in the range \([x_1, b]\)

2. Find two new points \(x_1, x_2\) in the new range, noticing that the unused existing point is already in the right place, ie:

   case 1: \(x_2' = x_1\)

   case 2: \(x_1' = x_2\)

The new length of the interval after each iteration is \(\tau\) of its previous length.

Parabolic Interpolation

1. Start with 3 points, \(u, v, w\), such that v is the minimizer: $f(v) < f(u), f(w)$$

2. Draw a parabola \(p(x)\) that goes through the points \((u,\ f(u))\ ,(v,\ f(v))\ ,(w,\ f(w)))\)

3. Find the minimum of the parabola, \(m\)

4. case 1: if \(m < v\), then \(u'=u,\ v'=m\ ,w' = v\)

   case 2: if \(m > v\), then $u

Solving for minimum of parabola

\(L(x) = \frac{f(u)(x - v)(x - w)}{(u-v)(u-w)} + \frac{f(v)(x - u)(x - w)}{v - i)(v - w)} + \frac{f(w)(x - u)(x - v)}{w - u)(w - v)}\)

\(L'(x) = 0\)

\(p = (v - u)^2(f(w) - f(v) - (v - w)^2(f(u) - f(v))\)

Interpolation

Interpolation allows us to:

• data fitting - if we only have a limited number of sample points, we can find a function that fits the data
• approximate functions that are difficult to compute - we can compute a function at a handful of points then interpolate to fill in the gaps

Power basis

Given a set of coordinates \((x_i, y_i), \cdots, (x_n, y_n)\), find a function that interpolates the values. We limit the the function to a polynomial of \(n\) degrees to keep the function simple.

Let \(P_n(x) = a_1 x^n + \cdots a_n x + a_{n+1}\) be an nth degree polynomial.

\[\underbrace{\left[ \begin{matrix} x_1^n & \cdots & x_1 & 1 \\ x_2^n & \cdots & x_2 & 1 \\ \vdots & \vdots & \vdots & \vdots \\ x_{n+1}^n & \cdots & x_{n+1} & 1 \end{matrix} \right]}_{V} \left[ \begin{matrix} a_1 \\ a_2 \\ a_3 \\ a_3 \end{matrix} \right] = \left[ \begin{matrix} y_1 \\ y_2 \\ \vdots \\ y_{n+1}\end{matrix} \right]\]

If \(x_i \neq x_ x_j\), then \(Va = y\) has a unique solution. \(a\) contains the coefficients of the interpolation polynomial.

The basis here is \({x^n, x^{n-1}, \cdots, x, 1}\) and is known as a power basis.

#+begin_examples

1. (1,16), (3,21), (5,15), (15,12). Fit to a 3rd degree polynomial: \(P_3(x) = a_1 x^3 + a_2 x^2 + a_3 x + a_4\)

   \[\left[ \begin{matrix} 1 & 1 & 1 & 1 \\ 27 & 9 & 3 & 1 \\ 125 & 25 & 5 & 1 \\ 216 & 36 & 6 & 1 \end{matrix} \right] \left[ \begin{matrix} a_1 \\ a_2 \\ a_3 \\ a_3 \end{matrix} \right] = \left[ \begin{matrix} 16 \\ 21 \\ 15 \\ 12\end{matrix} \right]\]

   \(V\) is the Vander Monde matrix, and if \(x_i \neq x_j\), then \(V\) is non-singular

#+end_exapmles

However Vander Monde matrices are known for being ill-conditioned, so this is not a good method for computing interpolating polynomials.

General matrix form

Let the basis functions be \(\phi_1 (x), \phi_2 (x), \cdots, \phi_{n+1} (x)\)

\(y = c_1 \phi_1 (x) \cdot + \phi_{n+1} (x)\)

\[\left[ \begin{matrix} \phi_1 (x_1) & \cdots & \phi_n (x_1) & \phi_{n+1} (x_1) \\ \phi_1 (x_2) & \cdots & \phi_n (x_2) & \phi_{n+1} (x_2) \\ \vdots & \vdots & \vdots & \vdots \\ \phi_1 (x_{n+1}) & \cdots & \phi_n (x_{n+1) & \phi_{n+1} (x_{n+1}) \end{matrix} \right] \left[ \begin{matrix} a_1 \\ a_2 \\ a_3 \\ a_3 \end{matrix} \right] = \left[ \begin{matrix} y_1 \\ y_2 \\ \vdots \\ y_{n+1}\end{matrix} \right]\]

Lagrange basis

\[L_i(x) = \frac{(x - x_1)(x - x_2) \cdot (x - x_{i-1}) (x - x_{i+1}) \cdot (x - x_{n+1})}{(x_i - x_1) (x_i - x_2) \cdots (x_i - x_{i-1}) (x_i - x_{i+1}) \cdots (x_1 - x_{n+1})}\]

\(L_i\) has the property that \[L_i(x_j) = \begin{cases} 1 & i = j \\ 0 & i \neq j \end{cases}\]

\(L(x) = y_1 L_1 (x) + y_2 L_2 (x) + \cdots + y_{n+1} L_{n+1} (x)\)

Each \(L_i (x)\) is an nth degree polynomial.

\[L(x) = \sum_1^{n+1} y_i L_i (x)\]

Where \[L_i(x) = \frac{\prod_1^{n+1} (x - x_i)}{\prod_1^{n+1} (x - x_j)}\]

Let \(\psi(x) = \prod_{j = 1}^{n + 1} \frac{f(x_i)}{ x - x_j} w_i\)

\[w_i = \frac{1}{\prod_{j = 1,\ j \neq i}^{n+1} (x_i - x_j)}\]

$Li(x) =

\[L = \sum_i=1^{n+1} y_i \frac{\phi (x)}{x - x} \cdot w_ifdA \]L9x) ϕ(x) \

\(L(x) = \frac{\sum_{i=1}^{n+1} \frac{z_i x - x_i}}\)

Newton Form Interpolating Polynomial

The divided difference approaches the derivative if all \(x_i\) begin to approach the same point: \[\lim_{(x_0, \cdots, x_n) \rightarrow (x, \cdots, x)} f[x_0, \cdots, x_n] = f^{(n)}(x)\]

So the Newton form polynomial is equivalent to a Taylor polynomial if all \(x_i\) begin to approach the same point: \[\lim_{(x_0, \cdots, x_n) \rightarrow (x, \cdots, x)} f[x_1] + f[x_1, x_2](x - x_1) + \dots + f[x_1, x_2, \dots, x_n](x - x_1) \dots (x - x_n) = f(x) + (\varepsilon - x) f'(x) + (\varepsilon - x)^2 \frac{f(x)}{2!} + \cdots + (\varepsilon - x)^n \frac{f(x)^(n)}{(n!)}\]

Divided difference

Interpolation with Derivatives

\[f[x_1, x_2] = \frac{f[x_2] - f[x_1]}{x_2 - x_1}\]

\[\lim_{x_2 \rightarrow x_1} f'(x_2)\]

so define \(f[x_1, x_1] = f'(x_1)\),\(\underbrace{f[x_1, x_1,\dots , x_1]}_{\text{k+1 times}} = \frac{f^{(k)}}(x_1){k!}\)

Error

\[e_n(x) = |f(x) - P_n(x)| = \frac{f^{(n+1)}(\zeta) }{(n+1)!} \phi(x)\], where \(\phi(x) = (x - x_1) \cdots (x - x_{n+1})\) and \(\frac{f^{(n+1)}(\zeta)}{(n+1)!} = f[x_1, x_2, \cdots, x_{n+1}, x]\)

for \(n=1\), \(P_1(x) = f[x_1] + f[x_1, x_2](x - x_1)\).

Error in the polynomial is driven by two terms:

• the derivative of f
• placement of interpolated points - equally spaced points usually drive this term high and lead to a large error. See chebyshev polynomial for information about ideal point placement.

Horner's Method

\(P = c_{n+1} (x - x_n)\)

\(P = \left[ P + c_n \right] (x - x_{n-1})\)

$P = [P+c_n − 1] (x - x_n-2)

jth iteration

\(P = \left[ P + c_j \right] (x - x_{j-1})\)

Pseudo code

Chebyshev Interpolating points.

Piecewise Interpolation

Hermite Interpolation

Assume we are interpolating a function \(f(x)\). We can perform nth degree Hermite interpolation at \((x_1, x_2, \cdots, x_n)\) if we know \(f(x_i), f'(x_i), \cdots, f^{(n)}(x_i)\).

Numerical Integration

Suppose we want to integrate

\[I = \int_a^b f(x) dx\]

We can replace \(f(x)\) by an interpolating polynomial and solve numerically.

\[Q = \int_a^b P_n(x) dx\]

where the error is \[E = \int_a^b f(x) - P_n(x) dx\]

Midpoint Method

Let \([a, b]\) be a small interval.

Approximate \(f\) over \([a,b]\) with \(f(\frac{a+b}{2})\).

\(Q = \int_a^b f\left( \frac{a+b}{2} \right) dx = f \left(\frac{a+b}{2} \right) (b - a)\)

This is equivalent to using an interpolating polynomial of degree 0.

midpoint method drawing

Error

\(f(x) = f(m) + f'(\varepsilon) \cdot (x - m) + \frac{f''(\varepsilon)}{2!}(x - m)^2\), where \(\varepsilon \in [a,b]\)

The first part of the error term is

\(\int_a^b f'(m) (x - m) dx = 0\)

The second part of the error term is

\(\text{error} = \frac{1}{2} \int_a^b f''(\varepsilon) (x - m)^2 dx\)

Since \((x-m)^2\) doesn't change sign, we can use the mean value theorem for integrals to conclude,

\(\text{error} = f''(\eta) \frac{1}{2} \int_a^b (x-m)^2 dx = \frac{f''(\eta) (b-a)^2}{24}\), for \(\eta \in [a,b]\)

MVT for Integrals

Assume \(g(x) \in C[a,b]\). \(\phi(x)\) is integrable that is nonnegative or nonpositive on \([a,b]\). Then there is a point \(\eta \in [a,b]\) such that \(\int_a^b g(x) \phi(x) dx = g(\eta) \inta^b \phi(x) dx\)

Trapezoidal Method

Let \([a, b]\) be a small interval.

Approximate \(f\) over \([a,b]\) with \(f(a) + f[a,b](x - a)\).

\[\begin{aligned} Q &= \int_a^b \left[ f(a) + f[a,b] (x - a) \right] dx \\ &= f(a)(b - a) + f[a,b] \frac{(x - a)^2}{2}_\eval{a,b} \\ &= \frac{b-a}{2} \cdot \left[ f(a) + f(b) \right] \end{aligned}\]

trapezoidal method drawing

Error

\(f(x) = f(a) + f[a,b](x - a) + \frac{f''(\varepsilon)}{2} (x-a)(x -b)\)

\(\text{error} = \int_a^b \frac{f''(\varepsilon)}{2} \underbrace{(x - a)}_{\geq 0}\underbrace{((x - b)}_{\leq 0}\), where \((x - a)(x - b)\)

Since \((x - a)(x - b)\) is nonpositive,

\[\text{error} = \frac{f''(\eta)}{2} \int_a^b (x - a)(x - b) dx = -f''(\eta)}\frac{(b - a)^3}{12}\]

Simpson's Method

Replace \(f(x)\) by a 2nd order polynomial in Lagrange form.

\[Q = \int_a^b P_2(x) dx = \frac{h}{6} \left[ f(a) + 4 f(m) + f(b) \right]\]

Error

\(\text{error} = \frac{f^{(4)}(\eta) (b - a)^5}{90 \cdot 2^5}\)

Degree of Precision

The degree of precision

Composite Midpoint Rule

Composite Trapezoidal Rule

Composite Simpson's Rule

Gaussian Quadrature

Assume we are computing \(I = \int_{-1}^1 f(x) dx\). With Gaussian quadrature method, we use carefully chosen points \([x_1, \cdots, x_{n+1}] \in [-1, 1]\) to make the integral precision as high as possible.

We begin by defining a family of orthogonal polynomials known as Legendre polynomials

We use the roots of \(P_k(x)\) as the interpolation points for quadratic interpolation.

Approximating Functions

\[min ||f - v||_2^2 = \int_a^b \left( f(x - \sum_{j=0}^{n+1} c_j \phi_j(x) \right)^2 dx\]

Differentiate with respect to \(c_k\)

\[\frac{d}{dk} \int_a^b \left( f(x - \sum_{j=0}^{n+1} c_j \phi_j(x) \right)^2 dx = 2 \int_a^b \left( f(x - \sum_{j=0}^{n} c_j \phi_j(x) \right) ( - \phi_k(x) ) dx = 0\]

\[\int_a^b f(x)(- \phi_k(x)) dx = \int_a^b \sum_{j=0}^{n} \ c_j phi_j(x) \phi_k(x) dx = \sum_{j=0}^n c_j \int_a^b \phi_j(x) \phi_k(x) dx\]

so

\[\underbrace{\sum_{j=0}^n c_j \int_a^b \phi_j(x) \phi_k(x) dx}_{\sum_{j=0}^n c_j B_{kj}} = \underbrace{\int_a^b f(x)(- \phi_k(x)) dx}_{b_k}\]

If we choose the monomial basis,

\(c_i = \frac{b_i}{B_{ii}}\)

where \(b_i = \int_a^b w(x) f (x) \phi_i(x) dx\)

Properties of Norms

1. \(|| g || \geq 0\), \(||g|| = 0 \Rightarrow g(x) = 0\) on \([a,b]\)

2. \(||\alpha g|| = |\alpha|||g||\)

3. \(||f + g||_2^2 \leq \left( ||f||_2 + ||g||_2 \right)^2\)

Least Squares

\[Va \approx y\]

We can get the residual error by $r = y - Va\[ Since $Va^*$ is a projection of $y$ to a lower dimension, $r$ is orthogonal to $Va^*$ and therefore $V^T$. $V^T(y - Va^*) = 0 \Longrightarrow V^Ty = V^TVa$ ** Orthogonal Polynomials $\phi_{-1} = 0$ $\phi_0 = 1$ \]ϕ_j(x) - x ϕ_j-1(x) = b ϕ_j-1 - c ϕ_j-1} + d ∑_k=0^j-3\[ three term recurrence relation: $\phi_j(x) = (x - b_j) \phi_{j-1}(x) - c_j \phi_{j-2}(x)$ $(\phi_j, \phi_{j-1}) = 0 = (x \phi_{j-1}, \phi_{j-1}) - \b_j (\phi_{j-1}, \phi_{j-1}) - c_j (\phi_{j-2}, \phi_{j-1})$ Multiply by $\phi_{j-1}$, integrate from a to b with weight for $w$ to find $b_j$ \]b_j = = \[ Repeat with $\phi_{j-2}$ to find $c_j$ \]c_j = \[ ** Legendre Polynomial $w(x) = $ $\phi_1(x) = (x - b_1) \phi_0 - c_1 \phi_{-1} = (x - b_1) 1 - c_1 0 = x - b_1$ From our earlier results, \]b_j = = \[ \]c₂ = = = = 1/3\[ \]ϕ₂ = x² - 1/3\[ *** Orthogonality \]P_n(x) = (1 - x²)ⁿ\[ $P_0(x) = 0$, $P_1(x) = 1$ We will show $\int_{-1}^1 P_j(x) P_k(x) dx = 0$ where $j < k$ \]∫_-1¹ (1 - x²)^j ⋅ (1 - x²)^k dx

Chebyshev Polynomial

Family of orthogonal polynomials on the interval \([-1, 1]\) with weights \(w(x) = \frac{1}{\sqrt{1 - x^2}}\).

i.e. the weights approach infinity when \(x \rightarrow \pm 1\), which puts more weight at the ends of the interval.

Alternative Representation

\(T_j = \cos(j \arccos(x)) = \cos(j \theta)\), where \(\theta = \arccos(x)\) ⇒ \(\cos(\theta) = x\)

which has roots \(j \theta = (2k - 1) \frac{\pi}{2} \Rightarrow \theta_k = \frac{2k - 1}{j} \frac{\pi}{2}\) for \(k = 1,2, \cdots, j\)

The function has extreme points at \(\cos(j \theta) = \pm1 \Rightarrow \theta = \frac{k \pi}{j}\)

\(T_j = \cos(j \theta)\)

\(T_{j+1} = \cos((j + 1) \theta) = \cos(j\theta + \theta) = \cos(j \theta) \cos( \theta) - sin(\j \theta) \sin(\theta)\)

\(T_{j-1} = \cos((j - 1) \theta = \cos(j \theta - \theta) = \cos(j \theta) \cos( \theta) + sin(j \theta)\sin(\theta)\)

\(T_{j+1} + T_{j - 1} = 2 \cos(j \theta) \cos(\theta) = 2 T_j (\theta) \cdot x = 2x t_j(x) - T_{j-1}(x)\)

\(T_0 = 1\)

\(T_1 = x\)

\(T_2 = 2x^2 - 1\)

\(T_3 = 2x(2x^2 - 1) - x = 4x^3 - 3x\)

\(\phi_j = \prod_{k=1}^j (x - \cos(\theta_k))\) for \(\theta_k = \frac{2k - 1}{j} \frac{\pi}{2}\)

\(\phi_j = \frac{1}{2^{j-1}} T_j(x)\)

todo

• algo vs problem
• problem conditioning
• backwards stability
• find k of integration prob
• types of convergence
  • linear
  • quadratic
  • forbinius
  • 3.5
• hessian matrix
• read chapter 10
• review parabolic interpolation

• piecewise polynomial interpolation

• chebyshev polynomial ch12

• everything up to piecewise polynomial interpolation

• read chebyshev polynomial 10.6