INFO-F-305 - TP 9

Equations aux différences II

Jacopo De Stefani - jdestefa@ulb.ac.be

Github: https://github.com/jdestefani/ulb-infof305

Slides: https://jdestefani.github.io/ulb-infof305

Rappels Théoriques - Equation aux différences affines

Equation affine aux différences, linéaire, d’ordre 1 et à coefficients constants

$x(k+1) = a x(k) + b$

• $ a \in \mathbb{R}, a \neq 0$
• $ b \in \mathbb{R}$

Rappels Théoriques - Equation aux différences affines

Théorème 7.9

Pour chaque valeur $x^0 = x(0) $ l'equation a une solution unique

• $a=1$ $\Rightarrow$ $x(k) = x(0) + kb$

  • $b=0$ $\Rightarrow$ Aucun équilibre, trajectoires divergentes
  • $b \neq 0$ $\Rightarrow$ Toutes les valeurs de $x$ sont des valeurs d'équilibres, trajectoires constantes

• $a \neq 1$ $\Rightarrow$ $x(k) = (x(0) - \alpha)a^k + \alpha$

  • $\bar{x} = \alpha = \frac{b}{1-a}$ seul équilibre
  • $\vert a \vert < 1$ $\Rightarrow$ Toutes trajectoires convergent vers $\alpha$
  • $\vert a \vert > 1$ $\Rightarrow$ Toutes trajectoires divergentes sauf celle pour $x(0) = \alpha$
  • $a = -1$ $\Rightarrow$ Equilibre $\alpha = \frac{b}{2}$ et toutes solutions oscillantes, à part celle constante pour $x(0) = \alpha = \frac{b}{2}$

Rappels Théoriques - Equation aux différences non linéaires

Systèmes discrets non-linéraires d'ordre 1

$x(k+1) = f(x(k))$

• $ f \in C^{\infty}$
• $ x(k) \in \mathbb{I} \subset \mathbb{R}$
• $x(1) = f(x(0))$
• $x(2) = f(x(1)) = f(f(x(0))) = f^2(x(0))$
• $\cdots$
• $x(k) = f(x(k-1)) = \cdots = f^k(x(0))$

N.B

• $x(k) = f^k(x(0)) \neq (f(x(0)))^k \neq f^{(k)}(x(0))$

Rappels Théoriques - Equation aux différences non linéaires

$x(k+1) = f(x(k))$

Etat d'équilibre

$\bar{x}$ est un état d'équilibre du système ${f,I}$ si $\bar{x} = f(\bar{x})$

• En pratique - Graphique: Intersection entre $y=f(x)$ et $y=x$ sur le plan $x,y$

Théorème 7.10

Soit $I=[a,b]$ un intervalle clos et borné et $f: I \mapsto I$. Alors il existe toujours $\bar{x}$ et $\bar{x} = f(\bar{x})$

Rappels Théoriques - Cycles

Un cycle d'ordre $s$ est un ensemble de $s$ valeurs différents $\{\bar{x}_0,\cdots,\bar{x}_{s-1}\}$ un système discret $\{f,I\}$ telles que:

$\bar{x}_1 = f(\bar{x}_0)$ $\bar{x}_2 = f(\bar{x}_1)$ $\cdots$ $\bar{x}_0 = f(\bar{x}_{s-1})$

La quantité $s$ est le période de l'orbite.

Théorème 7.11

La paire $\{\bar{x}_0,\bar{x}_{1}\}$ du système $\{f,I\}$ est un cycle d'ordre $2$ si et seulement si $\bar{x}_0$ ,$\bar{x}_{1}$ sont des points d'équilibre du système $\{f^2,I\}$ mais pas de $\{f,I\}$.

Rappels Théoriques - Stabilité

Stabilité de l'équilibre

• Stable: Si pour chaque $\varepsilon$, il existe $\delta > 0$ telle quelle $\vert x(0) - \bar{x} \vert < \delta \Rightarrow \vert x(k) - \bar{x} \vert < \varepsilon$
• Globalement asymptotiquement stable/attractif: Si pour chaque $x(0) \in I$ $\lim\limits_{k \rightarrow \infty} x(k) = \bar{x}$
• Localement asymptotiquement stable/attractif: Si il existe $\eta$ telle quelle pour chaque $x(0) \in I \cap \{x-\eta,x+\eta\}$ $\lim\limits_{k \rightarrow \infty} x(k) = \bar{x}$

Théorème 7.12

Soit $\bar{x}$ un équilibre du système $\{f,I\}$ où $f \in C^1$:

• Instable: $\vert f'(\bar{x}) \vert > 1$.
• Aucune information: $\vert f'(\bar{x}) \vert = 1$.
• Localement asymptotiquement stable/attractif: $\vert f'(\bar{x}) \vert < 1$.

Rappels Théoriques - Stabilité ordre 1

Théorème 7.13

Soit $\bar{x}$ un équilibre du système $\{f,I\}$ où $f \in C^2$ et $f'(\bar{x})=1$:

• $f''(\bar{x}) > 0$ $\Rightarrow$ $\bar{x}$ inférieurement asymptotiquement stable et supérieurement instable (où répulsif)
• $f''(\bar{x}) < 0$ $\Rightarrow$ $\bar{x}$ supérieurement asymptotiquement stable et inférieurement instable (où répulsif)

Théorème 7.14

Soit $\bar{x}$ un équilibre du système $\{f,I\}$ où $f \in C^3$ et $f'(\bar{x})=1$ et $f''(\bar{x})=0$:

• $f'''(\bar{x}) > 0$ $\Rightarrow$ $\bar{x}$ instable
• $f'''(\bar{x}) < 0$ $\Rightarrow$ $\bar{x}$ localement asymptotiquement stable

Rappels Théoriques - Stabilité ordre 1

Théorème 7.15

Soit $\bar{x}$ un équilibre du système $\{f,I\}$ où $f \in C^3$ et $f'(\bar{x})=-1$:

• $2 f'''(\bar{x}) + 3(f''(\bar{x}))^2 > 0$ $\Rightarrow$ $\bar{x}$ localement asymptotiquement stable
• $2 f'''(\bar{x}) + 3(f''(\bar{x}))^2 < 0$ $\Rightarrow$ $\bar{x}$ instable

Exercice 1

Trouver les cycles du système suivant :

$x(k + 1) = x(k)^2 − 1$

Solution - Exercice 1

$f(x) = x^2 - 1$

$f^2(x)=f(f(x))= f(x^2 - 1) = (x^2 - 1)^2 - 1 = x^4 - 2x^2$

x = -2 : 0.01 : 2;
f = @(x) x.^2 - 1;

hold on;
plot(x,f(x))
plot(x,f(f(x)))
plot(x,x)
set(gca, "xaxislocation", "origin");
set(gca, "yaxislocation", "origin");
xlabel("x(k)");
ylabel("x(k+1)");

$\bar{x} = f(\bar{x})$ $\Rightarrow$ $\bar{x} = \bar{x}^2 - 1$ $\Rightarrow$ $\bar{x}^2 - \bar{x} - 1$

$\bar{x} = f^2(\bar{x})$ $\Rightarrow$ $\bar{x} = \bar{x}^4 - 2\bar{x}^2$ $\Rightarrow$ $\bar{x}^4 - 2\bar{x}^2 - \bar{x}$

c_f = [1 -1 -1]; # Coeffiecients de x=f(x)
c_f_2 = [1 0 -2 -1 0]; # Coeffiecients de x=f^2(x)
roots(c_f) # Racines de x=f(x)
roots(c_f_2) # Racines de x=f^2(x)
setdiff(single(roots(c_f_2)),single(roots(c_f))) # Racines de f^2 mais pas de f

ans =

  -0.61803
   1.61803

ans =

   1.61803
  -1.00000
  -0.61803
   0.00000

ans =

  -1
   0

plot_system_simulation(f,0:0.01:1,0.5,5)

plot(1:max_steps,x_k);
xlabel("k");
ylabel("x(k)");

error: 'max_steps' undefined near line 1 column 8

Exercice 2

$x(k + 1) = 6x(k)^2(1 − x(k))$

où $x(k) \in [0, 1]$.

1. Calculer analytiquement $f^2(x)$ et tracer les graphiques respectifs de $f(x)$ et de $f^2(x)$
2. Trouver les points d’équilibre
3. Etudier la stabilité des points d’équilibre
4. Pour chaque point d’équilibre, simuler graphiquement et numériquement une trajectoire de 5 étapes qui soit compatible avec l’analyse de stabilité
5. Trouver graphiquement s’il existe un cycle d’ordre 2.

1. Calculer analytiquement f^2(x) et tracer les graphiques respectifs de f(x) et de f^2(x)

$f(x) = 6x^2(1 − x)$

$f^2(x)=f(f(x))= f(6x^2(1 − x)) = 6(6x^2(1 − x))^2(1- (6x^2(1 − x)) = 216x^4(1-2x+x^2)(1 - 6x^2 + 6x^3)$

x = 0 : 0.01 : 1;
f = @(x) ((6*x.^2).*(1-x));

hold on;
plot(x,f(x))
plot(x,f(f(x)))
plot(x,x)
set(gca, "xaxislocation", "origin");
set(gca, "yaxislocation", "origin");
xlabel("x(k)");
ylabel("x(k+1)");

1. Trouver les points d’équilibre

$\bar{x} = f(\bar{x})$ $\Rightarrow$ $\bar{x} = 6\bar{x}^2(1-\bar{x})$ $\Rightarrow$ $-6\bar{x}^3 + 6\bar{x}^2 - x$

• $x^{(1)} = \frac{-6-\sqrt{12}}{12}$
• $x^{(2)} = \frac{-6+\sqrt{12}}{12}$
• $x^{(3)} = 0$

c_f = [-6 6 -1 0]; # Coeffiecients de x=f(x)
x_eq = roots(c_f) # Racines de x=f(x)

x_eq =

   0.78868
   0.21132
   0.00000

1. Etudier la stabilité des points d’équilibre

$f(x) = 6x^2(1 − x)$ $\Rightarrow$ $f'(x) = 12x(1 − x) + -6x^2$

pkg load symbolic
syms x;
f_sym = 6*x^2*(1-x)
diff(f_sym)

f_sym = (sym)

     2        
  6*x *(1 - x)

ans = (sym)

       2               
  - 6*x  + 12*x*(1 - x)

diff_f = @(x) -6*x.^2+12*x.*(1-x)
x_eq
diff_f(x_eq)

diff_f =

@(x) -6 * x .^ 2 + 12 * x .* (1 - x)

x_eq =

   0.78868
   0.21132
   0.00000

ans =

  -1.73205
   1.73205
   0.00000

Théorème 7.12

• $\vert f'x^{(1)} \vert $ > 1 $\Rightarrow$ $x^{(1)}$ instable
• $\vert f'x^{(2)} \vert $ > 1 $\Rightarrow$ $x^{(2)}$ instable
• $\vert f'x^{(3)} \vert $ < 1 $\Rightarrow$ $x^{(3)}$ asymptotiquement stable

# 4. Pour chaque point d’équilibre, simuler graphiquement et numériquement une trajectoire de 5 étapes qui soit compatible avec l’analyse de stabilité.
# Pour x(0) = 0.05
plot_system_simulation(f,0:.001:.05,0.05,5)

# 4. Pour chaque point d’équilibre, simuler graphiquement et numériquement une trajectoire de 5 étapes qui soit compatible avec l’analyse de stabilité.
# Pour x(0) = 0.1
plot_system_simulation(f,0:.01:1,0.25,5)

# 4. Pour chaque point d’équilibre, simuler graphiquement et numériquement une trajectoire de 5 étapes qui soit compatible avec l’analyse de stabilité.
# Pour x(0) = 0.7
plot_system_simulation(f,0:.01:1,0.7,5)

# 5. Trouver analytiquement s’il existe un cycle d’ordre 2.
syms x;
f_sym = 6*x^2*(1-x)
f2_sym = 216*x.^4.*(1-2*x+x.^2).*(1 - 6*x.^2 + 6*x.^3)
[c1,t1] = coeffs(f_sym)
[c2,t2] = coeffs(f2_sym)

f_sym = (sym)

     2        
  6*x *(1 - x)

f2_sym = (sym)

       4 / 2          \ /   3      2    \
  216*x *\x  - 2*x + 1/*\6*x  - 6*x  + 1/

c1 = (sym) [-6  6]  (1x2 matrix)
t1 = (sym 1x2 matrix)

  [ 3   2]
  [x   x ]

c2 = (sym) [1296  -3888  3888  -1080  -432  216]  (1x6 matrix)
t2 = (sym 1x6 matrix)

  [ 9   8   7   6   5   4]
  [x   x   x   x   x   x ]

c_f = [double(c1) -1 0]; # Coeffiecients de x=f(x) => f(x) - x = 0
c_f_2 = [double(c2) 0 0 -1 0]; # Coeffiecients de x=f^2(x) => f^2(x) - x = 0
roots(c_f) # Racines de x=f(x)
roots(c_f_2) # Racines de x=f^2(x)
setdiff(single(roots(c_f_2)),single(roots(c_f))) # Racines de f^2 mais pas de f

ans =

   0.78868
   0.21132
   0.00000

ans =

   1.05393 + 0.00000i
   0.86603 + 0.00000i
   0.78868 + 0.00000i
   0.60287 + 0.00000i
  -0.35946 + 0.00000i
   0.21132 + 0.00000i
  -0.08169 + 0.12935i
  -0.08169 - 0.12935i
   0.00000 + 0.00000i

ans =

  -0.08169 - 0.12935i
  -0.08169 + 0.12935i
  -0.35946 + 0.00000i
   0.60287 + 0.00000i
   0.86603 + 0.00000i
   1.05393 + 0.00000i