Lecture 5 : Solving nonlinear systems¶

May, 2022 - François HU

Master of Science - EPITA

This lecture is available here: https://curiousml.github.io/

Table of contents¶

Objective
Bisection method
Fixed-point iteration method
Newton's method
Exercices
- Exercice 1: bisection algorithm
- Exercice 2: fixed_point algorithm
- Exercice 3: newton algorithm

Solving $A\times X = B$
- (optional) Exercice 4: (linear algebra) LU decomposition
- (optional) Exercice 5: (linear algebra) Matrix inversion

1.1. Nonlinear systems¶

In linear systems: $$ \begin{bmatrix} a_{0,0} & a_{0,1} & \cdots & a_{0,n-1} \\ a_{1,0} & a_{1,1} & \cdots & a_{1,n-1} \\ & \vdots & \\ a_{n-1,0} & a_{n-1,1} & \cdots & a_{n-1,n-1} \end{bmatrix} \times \begin{bmatrix} x_0\\ x_1\\ \vdots\\ x_{n-1} \end{bmatrix} = \begin{bmatrix} b_0\\ b_1\\ \vdots\\ b_{n-1} \end{bmatrix}\iff Ax - b = 0 $$

In nonlinear systems: $$ f(x) = 0 $$

We wish to find a solution to the non linear equation $f(x) = 0$ with $x = \begin{bmatrix} x_0\\ x_1\\ \vdots\\ x_{n-1} \end{bmatrix}$ the unkwown variable in the form of a "floating point" number. An equation of this type may have none, one or more solutions. We will assume that it has at least one which we will try to determine.

1.2. Condition number¶

We quantity the condition number that measures how sensitive the output of a function is on its input:

$$ \text{change in output} = \text{condition number}\times \text{change in input} $$

for a linear system $Ax = b$: a condition number gives a bound on how inaccurate the solution $x$ will be after approximation. If $A$ is invertible, the associated condition number is: $$\kappa(A) = \lvert\lvert A^{-1}\lvert\lvert_2\cdot\lvert\lvert A\lvert\lvert_2$$

for a nonlinear function $f$: a condition number measures how sensitive a $f$ is to changes (or errors) in the input. The associated condition number is:
- for one variable: $$ \kappa(x) = \left\lvert \dfrac{xf'(x)}{f(x)} \right\lvert $$
- for several variables: $$ \kappa(x) = \dfrac{\lvert\lvert x\lvert\lvert_2 \cdot \lvert\lvert J(x)\lvert\lvert_2 }{\lvert\lvert f(x)\lvert\lvert_2} $$ with where $J ( x )$ is the Jacobian matrix

1.3. Order of convergence¶

a solution $x^*$ is called simple if $f(x^*) = 0$ and $f'(x^*) \neq 0$

order/rate of convergence: The rate or order of convergence for $x_k$ obtained at the $k$-th iteration, is characterized by $q$ (order) or $r$ (rate) such that

$$ \lim\limits_{k\to\infty} \dfrac{ x_{k+1} - x^* }{ \lvert\lvert x_{k} - x^* \lvert\lvert^q } = r $$

for each iteration we are $q$ times more precise.
if $q=1$ the convergence is called a linear convergence,
$q=2$ is called quadratic convergence,
$q=3$ is called cubic convergence.

2. Bisection method ¶

We want to solve $f(x) = 0$, knowing that the root $r$ belongs to $[a, b]$ and that $f(a)\times f(b) < 0$ (in other words, $f(a)$ and $f(b)$ have opposite signs). We assume that we have an unique solution in this interval.

We consider the middle $m = \dfrac{b+a}{2}$ of the interval $[a, b]$

If $f(m)$ is of the sign of $f(a)$ then replace $a$ by $m$, otherwise $b$ by $m$

The operation is repeated until the interval size corresponds to the desired accuracy

Pseudo-code¶

bisection(f, a, b, eps)
    Input : function f, extremities a and b, precision eps
    Output : m such that |b-a| < eps

    while b - a >= eps do:
        m = (a+b)/2
        if f(a) x f(m) > 0 do:
            a = m
        else:
            b = m

Remarks:

At each iteration the interval decreases by half

The convergence is linear ($q=1$)

3. Fixed-point iteration method ¶

Given a function $g$, we call a fixed point, a value $x$ such that: $g(x) = x$

Schema for finding a fixed point, with iteration: $$ x_{k+1} = g(x_k) $$

Solving the equation $f(x) = 0$ can be done via a fixed point iterative approach. For one equation, there can be several approaches. But not all of them converge.

Example: for $f(x) = x^2 -x -5$, we set $g(x) = x^2 - 5$ or $g(x) = \sqrt{x+5}$.