Multiple subjects explained in french 67

MÉTHODE DE DICHOTOMIE

**Hypothèses :**

- f continue sur [a, b]

- f(a) · f(b) < 0 (signes opposés → zéro garanti par le TVI)

**Algorithme :** on coupe l’intervalle en deux, on garde le sous-intervalle qui contient le zéro, on répète n fois.

```python

# -*- coding: utf-8 -*-

import numpy as np

import matplotlib.pyplot as plt

def f(x):

return x\*\*2 - 2

fig, ax = plt.subplots()

X = np.linspace(1, 2)

ax.plot(X, f(X))

ax.spines[“left”].set_position(“zero”)

ax.spines[“bottom”].set_position(“zero”)

ax.spines[“top”].set_visible(False)

ax.spines[“right”].set_visible(False)

# Dichotomie brute

a, b = 1, 2

n = 20

for i in range(n):

m = (a + b) / 2

if f(m) \* f(b) > 0:

    b = m

else:

    a = m

print(f"étape {i} : x{i} = {m}, f(x{i}) = {f(m)}")

ax.plot(m, 0, "r.")

# Dichotomie sous forme de fonction

def dicho(a, b, f, n=20):

"""

Trouve un zéro de f dans \[a, b\] par dichotomie.

Paramètres

----------

a, b : float      Bornes (f(a)\*f(b) < 0)

f    : callable   Fonction

n    : int        Nombre d'itérations (défaut 20)

Retourne

--------

float : approximation du zéro

"""

for i in range(n):

    m = (a + b) / 2

    if f(m) \* f(b) > 0:

        b = m

    else:

        a = m

return m

x1 = dicho(1, 2, f)

x2 = dicho(-2, -1, f)

print(x1, x2)

```

**Explication ligne par ligne :**

- `m = (a+b)/2` → on calcule le milieu

- `if f(m)*f(b) > 0` → m et b sont du même signe → le zéro est entre a et m → on rétrécit par `b = m`

- sinon → le zéro est entre m et b → `a = m`

- `ax.plot(m, 0, “r.”)` → on visualise chaque approximation sur l’axe x

-–

## MÉTHODE DE NEWTON

### PREUVE 1 — Dérivée numérique

On part du développement de Taylor :

```

f(x+h) = f(x) + h·f’(x) + (h²/2)·f’'(x) + …

f(x-h) = f(x) - h·f’(x) + (h²/2)·f’'(x) + …

```

En soustrayant les deux :

```

f(x+h) - f(x-h) = 2h·f’(x) + O(h³)

```

Donc :

```

f’(x) ≈ [ f(x+h) - f(x-h) ] / (2h)

```

C’est la formule des différences finies centrées. L’erreur est en O(h²), bien meilleure que la formule simple f’(x) ≈ [f(x+h)-f(x)]/h (erreur en O(h)).

### PREUVE 2 — Formule de récurrence x_{n+1}

On cherche x tel que f(x) = 0. En x_n, on trace la tangente à la courbe :

```

tangente : y = f(x_n) + f’(x_n) · (x - x_n)

```

On cherche où cette tangente coupe l’axe x (y = 0) :

```

0 = f(x_n) + f’(x_n) · (x_{n+1} - x_n)

f’(x_n) · (x_{n+1} - x_n) = -f(x_n)

x_{n+1} = x_n - f(x_n) / f’(x_n)

```

**Hypothèses pour Newton :**

- f dérivable sur l’intervalle

- f’(x_n) ≠ 0 à chaque étape

- x0 suffisamment proche du zéro (sinon diverge)

```python

# -*- coding: utf-8 -*-

import numpy as np

import matplotlib.pyplot as plt

def f(x):

return x\*\*2 - 2

def fprime(x, h=0.01):

"""Dérivée numérique par différences finies centrées."""

return (f(x + h) - f(x - h)) / (2 \* h)

def newton(x0, f, n=20):

"""

Trouve un zéro de f par la méthode de Newton.

Paramètres

----------

x0 : float     Point de départ

f  : callable  Fonction

n  : int       Nombre d'itérations

Retourne

--------

float : approximation du zéro

"""

for i in range(n):

    print(f"étape {i} : x{i} = {x0}, f(x{i}) = {f(x0)}")

    x0 = x0 - f(x0) / fprime(x0)

return x0

fig, ax = plt.subplots()

X = np.linspace(0, 5)

ax.plot(X, f(X))

ax.spines[“left”].set_position(“zero”)

ax.spines[“bottom”].set_position(“zero”)

ax.spines[“top”].set_visible(False)

ax.spines[“right”].set_visible(False)

x1 = newton(1, f)

x2 = newton(-1, f)

print(x1)

```

**Explication :**

- `fprime(x)` calcule f’(x) numériquement avec la formule prouvée ci-dessus

- `x0 = x0 - f(x0)/fprime(x0)` applique la récurrence prouvée ci-dessus

- Converge en quelques itérations seulement (beaucoup plus rapide que dichotomie)

-–

## APPLICATION — Factorisation de polynômes

```python

# -*- coding: utf-8 -*-

import numpy as np

import matplotlib.pyplot as plt

def dicho(a, b, f, n=20):

for i in range(n):

    m = (a + b) / 2

    if f(m) \* f(b) > 0:

        b = m

    else:

        a = m

return m

def p(x):

return x\*\*3 - 5\*x\*\*2 + 3\*x + 8

def q(x):

return x\*\*3 + 2\*x\*\*2 + 10\*x - 20

fig, ax = plt.subplots()

X = np.linspace(-2, 4)

ax.plot(X, p(X), label=“p”)

ax.plot(X, q(X), label=“q”)

ax.spines[“left”].set_position(“zero”)

ax.spines[“bottom”].set_position(“zero”)

ax.spines[“top”].set_visible(False)

ax.spines[“right”].set_visible(False)

ax.legend()

# Zéros de p

x1 = dicho(-1, 0, p)

x2 = dicho(2, 3, p)

x3 = dicho(3, 4, p)

print("Zéros de p : ", x1, x2, x3)

# Zéros de q

x4 = dicho(1, 2, q)

print("Zéros de q : ", x4)

print(f"p(x) = (x - {x1:.4f})(x - {x2:.4f})(x - {x3:.4f})")

print(f"q(x) = (x - {x4:.4f})(…)")

```

**Explication :** on lit sur le graphe les intervalles où la courbe change de signe, puis on applique dicho sur chacun. Si x1, x2, x3 sont les trois racines d’un polynôme de degré 3 : p(x) = (x−x1)(x−x2)(x−x3).

-–

## RÉGRESSIONS

### PREUVE — Régression linéaire (y = Ax)

On minimise S(A) = Σ(y_i − A·x_i)²

```

dS/dA = 0

-2 · Σ x_i·(y_i − A·x_i) = 0

Σ x_i·y_i = A · Σ x_i²

```

Donc :

```

A = Σ(x_i · y_i) / Σ(x_i²)

```

-–

### PREUVE — Régression affine (y = Ax + B)

On minimise S(A, B) = Σ(y_i − A·x_i − B)²

En dérivant par rapport à A et B séparément et en égalisant à 0 :

```

∂S/∂A = 0 → Σ(x_i·y_i) = A·Σ(x_i²) + B·Σ(x_i)

∂S/∂B = 0 → Σ(y_i) = A·Σ(x_i) + n·B

```

On note x̄ = Σx_i/n et ȳ = Σy_i/n. En résolvant le système :

```

A = [ Σ(x_i·y_i) − n·x̄·ȳ ] / [ Σ(x_i²) − n·x̄² ]

B = ȳ − A·x̄

```

-–

### PREUVE — Régression exponentielle (y = λ·Aˣ)

Changement de variable : on pose Y’ = ln(y)

```

ln(y) = ln(λ) + x·ln(A)

Y’ = b + a·x avec a = ln(A), b = ln(λ)

```

-> C’est une régression affine sur les couples (x_i, ln(y_i)).

Une fois a et b trouvés, on retransforme :

```

A = exp(a), λ = exp(b)

```

-–

### PREUVE — Régression monomiale (y = c·xᵏ)

Changement de variable : on pose X’ = ln(x) et Y’ = ln(y)

```

ln(y) = ln(c) + k·ln(x)

Y’ = ln(c) + k·X’

```

-> C’est une régression affine sur les couples (ln(x_i), ln(y_i)).

Une fois la pente k et l’ordonnée b trouvées :

```

k = pente A, c = exp(B)

```

-–

### Code du professeur — Régression exponentielle

```python

# -*- coding: utf-8 -*-

import numpy as np

import matplotlib.pyplot as plt

def r(X, Y):

"""Coefficient de corrélation entre X et Y."""

n = len(X)

num  = sum(X \* Y) - n \* sum(X)/n \* sum(Y)/n

den1 = sum(X\*\*2)  - n \* (sum(X)/n)\*\*2

den2 = sum(Y\*\*2)  - n \* (sum(Y)/n)\*\*2

return num / np.sqrt(den1 \* den2)

# —Données—

X = np.array([0, 2, 4, 6, 8, 10, 12, 14, 16, 18])

Y = np.array([1406, 855, 618, 422, 282, 189, 138, 84, 56, 37])

n = len(X)

# —Coefficients de corrélation—

print(“lin-lin”, r(X, Y))

print(“lin-log”, r(X, np.log(Y)))

# print(“log-log”, r(np.log(X), np.log(Y)))

# —Calcul de la régression y = A * B**x—

# Changement de variable : Y’ = ln(Y)

Ynew = np.log(Y)

a = (sum(X * Ynew) - n * sum(X)/n * sum(Ynew)/n) / (sum(X**2) - n * (sum(X)/n)**2)

b = sum(Ynew)/n - a * sum(X)/n

A, B = np.exp(b), np.exp(a)

# —Erreur quadratique moyenne—

err = sum((Y - A * B**X)**2) / n

print(“Erreur :”, err)

# —Représentation graphique—

fig, ax = plt.subplots()

ax.plot(X, Y, “ro”, label=“Données”)

T = np.linspace(0, 20, 100)

ax.plot(T, A * B**T, label=“Régression”)

ax.set_title(f"y = {A:.2f} · {B:.4f}^x")

ax.legend()

```

**Explication :**

- `r(X, np.log(Y))` teste si lin-log est le meilleur modèle (r proche de ±1 = bon)

- `Ynew = np.log(Y)` applique le changement de variable prouvé ci-dessus

- `a` et `b` sont calculés avec les formules de régression affine prouvées ci-dessus

- `A, B = np.exp(b), np.exp(a)` retransforme vers le modèle y = A·Bˣ

-–

## LOTKA-VOLTERRA

**Système :**

```

dx/dt = a·x − b·x·y − e·x

dy/dt = c·x·y − d·y − e·y

```

**Avec limite environnement (K) et pêche (e) :**

```

dx/dt = a·x·(1 − x/K1) − b·x·y − e1·x

dy/dt = c·x·y − d·y·(1 − y/K2) − e2·y

```

**Discrétisation (schéma d’Euler) :**

```

x_{t+1} = x_t + h · (a·x_t·(1−x_t/K1) − b·x_t·y_t − e1·x_t)

y_{t+1} = y_t + h · (c·x_t·y_t − d·y_t·(1−y_t/K2) − e2·y_t)

```

```python

# -*- coding: utf-8 -*-

import numpy as np

import matplotlib.pyplot as plt

# SANS PÊCHE

x0, y0, n, h = 2, 8, 50000, 0.001

a, b, c, d, e = 1, 0.05, 0.04, 0.2, 0

X, Y = [x0], [y0]

for t in range(n):

x = X\[t\] + h \* (a\*X\[t\] - b\*X\[t\]\*Y\[t\] - e\*X\[t\])

y = Y\[t\] + h \* (c\*X\[t\]\*Y\[t\] - d\*Y\[t\] - e\*Y\[t\])

X.append(x if x > 0 else 0)

Y.append(y if y > 0 else 0)

fig, ax = plt.subplots()

ax.plot(X, Y, label=“sans pêche”)

ax.plot(x0, y0, “r*”, label=f"initial ({x0}; {y0})")

ax.set_xlabel(“Proies”)

ax.set_ylabel(“Prédateurs”)

ax.legend()

# AVEC PÊCHE + LIMITE ENVIRONNEMENT

x0, y0, n, h = 2, 8, 50000, 0.001

a, b, c, d = 1, 0.05, 0.04, 0.2

e1, e2 = 0, 0.1 # taux de pêche

K1, K2 = 10, 5 # limites de l’environnement

X, Y = [x0], [y0]

for t in range(n):

x = X\[t\] + h \* (a\*X\[t\]\*(1 - X\[t\]/K1) - b\*X\[t\]\*Y\[t\] - e1\*X\[t\])

y = Y\[t\] + h \* (c\*X\[t\]\*Y\[t\] - d\*Y\[t\]\*(1 - Y\[t\]/K2) - e2\*Y\[t\])

X.append(x if x > 0 else 0)

Y.append(y if y > 0 else 0)

ax.plot(X, Y, label=“avec pêche”)

ax.legend()

```

**Explication :**

- `h=0.001` → pas de temps très petit = simulation précise

- `n=50000` itérations × h=0.001 = 50 unités de temps simulées

- `x if x>0 else 0` → les populations ne peuvent pas être négatives

- Le **diagramme de phase** (X vs Y) montre si le système converge vers un équilibre ou tourne en cycles

- `K1, K2` limitent la croissance : sans prédateurs, x ne peut pas dépasser K1

- `e1, e2` modélisent la pêche : augmenter e2 peut faire disparaître les prédateurs

-–