Master 2 Mathématiques Fondamentales Mars-Juin 2021

Techniques d’intelligence artificielle pour

les équations différentielles autonomes

Stage de Master 2 - Rapport

Stage encadré par Philippe CHARTIER, Mohammed LEMOU et Florian MEHATS

Maxime BOUCHEREAU 1 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

Remerciements
Je remercie Philippe CHARTIER, Mohammed LEMOU et Florian MEHATS pour l’encadrement,
l’aide et les conseils apportés tout au long du stage, mais aussi pour l’intérêt qu’ils m’ont
fait porter pour ce domaine.

Maxime BOUCHEREAU 2 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

Introduction
L’intelligence artificielle, outre l’application à des problèmes complexes tels que le traite-
ment d’image, du signal ou encore la classification des données [1],[2], peut s’appliquer à
l’étude de problèmes d’évolution. En effet, il est possible d’utiliser des techniques d’IA
sur des équations différentielles - voire aux dérivées partielles lorsqu’elles sont discrétisées
en espace - autonomes ẏ = f (y) afin de retrouver le champ de vecteurs f via l’étude de
trajectoires [3],[4] y compris lorsque ces dernières sont perturbées par un bruit [5], afin
d’imiter la trajectoire (discrétisée en temps) d’une solution d’EDO [6] ou lorsque l’on veut
simplifier un modèle physique afin de réduire le coût de calcul [7]

L’objectif de ce stage consiste en l’application de méthodes d’intelligence artificielle

à des équations différentielles, c’est-à-dire que, étant donné une EDO ẏ = f (y), il s’agit
de retrouver le champ de vecteurs f de manière approchée à partir de données sur la
trajectoire de la solution partant d’une donnée initiale y0 , et en connaissant une donnée
”finale” y1 . En réalité, il est nécessaire de connaı̂tre plus d’une trajectoire (donc on se
(k) (k)
donne plusieurs conditions initiales y0 et plusieurs données y1 ) afin de pouvoir calculer
f . C’est là que les méthodes d’IA entrent en jeu, puisque les données sont souvent nom-
breuses. Par ailleurs, f est retrouvé de manière approchée, c’est-à-dire que des méthodes
numériques pour les EDO entrent en jeu [8],[9]. Une part du travail est consacrée à l’étude
de la précision avec laquelle le champ de vecteurs f est retrouvé et le lien avec l’ordre de
la méthode numérique, en particulier lorsque l’on rajoute des points intermédiaires entre
y0 et y1 , permettant d’apporter plus de précision.

La première partie traite de l’introduction des méthodes d’intelligence artificielle pour

l’étude d’EDO. Il y est introduit la définition de fonction de perte Loss, qui mesure l’écart
entre des points de la trajectoire réelle et de la trajectoire calculée par une IA, la notion de
réseaux de neurones, notamment le PMC [2], ainsi que le gradient stochastique, technique
d’optimisation largement employée dans le cadre du stage afin de minimiser la fonction
Loss [1],[10].

La seconde partie se concentre sur le cas linéaire et aborde quatre méthodes numériques
classiques des EDO de petit ordre [8],[9] permettant de retrouver le champ de vecteurs.
Cette partie accorde une place importante à la convexité de la fonction Loss, critère très
apprécié en optimisation.

La troisième partie propose toujours d’étudier l’approximation de f dans le cas linéaire,

en fonction de certaines méthodes numériques mieux adaptées aux algorithmes d’optimisa-
tion que l’on peut trouver dans la littérature [3],[8],[9]

Enfin, la quatrième partie aborde le cas où f n’est plus linéaire, et met en évidence une
nouvelle structure de réseau de neurones [1][7], adaptée à l’étude de systèmes dynamiques
complexes utilisé dans d’autres disciplines comme la physique ou la biologie [3],[7]. Les
méthodes numériques employées reprennent celles utilisées dans le cas linéaire [3],[8],[9].

Maxime BOUCHEREAU 3 Université Rennes 1

Table des matières

1 Intoduction et outils 5
1.1 Position du problème . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5
1.2 Réseaux de neurones . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6
1.2.1 Couche de neurones . . . . . . . . . . . . . . . . . . . . . . . . . . . 6
1.2.2 Perceptron Multi-Couche (PMC) . . . . . . . . . . . . . . . . . . . 7
1.3 Algorithme du gradient stochastique (SGD) . . . . . . . . . . . . . . . . . 8

2 Cas linéaire: Première formulation 9

2.1 Position du problème . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9
2.2 Méthode d’Euler Explicite . . . . . . . . . . . . . . . . . . . . . . . . . . . 10
2.3 Méthode d’Euler Implicite . . . . . . . . . . . . . . . . . . . . . . . . . . . 13
2.4 Méthode du Point Milieu . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15
2.5 Méthode de Rung-Kutta 2 . . . . . . . . . . . . . . . . . . . . . . . . . . . 18

3 Cas linéaire: Nouvelle formulation et méthodes supplémentaires 24

3.1 Position du problème . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24
3.2 Convexité des fonctions Loss . . . . . . . . . . . . . . . . . . . . . . . . . . 25
3.3 Méthodes à deux points . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26
3.4 Méthodes Multi-pas (explicites) à pas constant . . . . . . . . . . . . . . . . 30
3.4.1 Ordre de convergence . . . . . . . . . . . . . . . . . . . . . . . . . . 30
3.4.2 Exemples . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 33
3.5 Méthodes Multi-pas à pas non constant . . . . . . . . . . . . . . . . . . . . 37
3.5.1 Ordre de convergence . . . . . . . . . . . . . . . . . . . . . . . . . . 37
3.5.2 Exemples . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39
3.6 Généralisation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41
3.6.1 Ordre de convergence . . . . . . . . . . . . . . . . . . . . . . . . . . 41
3.6.2 Exemples . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 43

4 Cas non linéaire 46

4.1 Structure du réseau de neurones . . . . . . . . . . . . . . . . . . . . . . . . 46
4.2 Ordre de convergence . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 47
4.2.1 Première formulation . . . . . . . . . . . . . . . . . . . . . . . . . . 48
4.2.2 Nouvelle formulation . . . . . . . . . . . . . . . . . . . . . . . . . . 51

4
Partie 1

Intoduction et outils

1.1 Position du problème

Soit f : Rd → Rd un champ de vecteurs continu et localement lipschitzien. Soient h >
0, N ∈ N∗ et t0 = 0 < t1 < · · · < tN = h. On considère le problème de Cauchy suivant:

ẏ = f (y)
y(0) = y0

et on suppose que, pour tout j ∈ [[0, N ]], y tj = φtj (y0 ) est connu, où φ est le flot associé
h
à l’équation différentielle ẏ = f (y).

L’idée du travail réalisé est de retrouver le champ de vecteurs f à partir des données
yj de la solution à différents temps. Concrètement, il s’agit de retrouver le champ de
vecteurs f à partir d’observations sur la trajectoire partant de y0 .
(k)
Bien entendu, plusieurs observations sont menées, à partir de K données initiales y0 ,
(k) (k)
pour k ∈ [[0, K − 1]], donnant ainsi, pour tout j ∈ [[1, N ]], y tj = φtj (y0 ).
h
La stratégie est d’approcher le champ de vecteurs f aussi précisément que possible
afin de pouvoir reconstituer toute la trajectoire. C’est là que les techniques d’IA entrent
en jeu:

Figure 1.1: Illustration de la technique d’IA

- Les données y tj sont entrées dans un réseau de neurones (NN), muni de paramètres.
h

5
Master 2 Mathématiques Fondamentales Mars-Juin 2021

- Il en ressort un champ de vecteurs approché fapp (calculé en fonction des paramètres de

N N ), ainsi qu’une valeur en h prédite de la solution, notée yˆ1 , faisant intervenir fapp au
moyen d’un schéma numérique donné. Par exemple, si on opte pour le schéma d’Euler
explicite, on aura yˆ1 = y0 + hfapp (y0 ).

- On évalue la différence yˆ1 − y1 à l’aide le la fonction de perte (Loss). C’est cette

fonction que l’on va chercher à minimiser en les paramètres du réseau de neurones,
faisant intervenir des méthodes d’optimisation

On dispose de K données, mais on ne les utilise pas toutes afin d’optimiser la fonction
de perte. Une partie seulement le sera (K0 < K données). C’est la phase d’apprentissage
(ou entraı̂nement), qui consiste à minimiser cette fonction:

K0 −1
1 X 1 (k)
2
LossT raining = yˆ1 (k) − y1
K0 k=0 h2

Le reste des données servira à vérifier que les paramètres du réseau de neurones ap-
ˆ(k) (k)
prochent correctement f , et que les autres différences y1 − y1 ne sont pas trop grandes
non plus. C’est la phase de test, qui consiste à vérifier que la fonction suivante n’est pas
trop grande par rapport à la fonction LossT raining :

K−1
1 X 1 (k)
2
LossT est = 2
yˆ1 (k) − y1
K − K0 k=K h
0

Remarque. - Lorsque la fonction LossT est devient trop grande par rapport à LossT raining ,
on dit qu’il y a surapprentissage.

- En notant p ∈]0, 1[ la proportion de données qui servent à l’entraı̂nement, on a K0 =

bpK0 c

1.2 Réseaux de neurones

Dans cette section, nous allons détailler le fonctionnement d’un réseau de neurones, tel
qu’il est utilisé dans le cadre du stage (ce qui ne correspond donc pas au cas général).

1.2.1 Couche de neurones

Une couche de neurones peut être vue comme une application Rd → Rζ qui envoie un
vecteur vers un autre vecteur, et dont les paramètres sont les suivants:

- Une matrice W ∈ Mζ,d (R), appelée matrice poids

- Un vecteur b ∈ Rζ appelé poids du biais (le biais a la valeur 1)

- Une fonction σ : R → R appelée fonction d’activation

On dit que la couche comporte ζ neurones, et voici son fonctionnement:

Maxime BOUCHEREAU 6 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

- On prend un vecteur x ∈ Rd , appelé entrée

- On calcule le vecteur y = W x + b ∈ Rζ

- On calcule le vecteur x0 = Σ(y) ∈ Rζ en appliquant la fonction d’activation σ à chacune

des composantes de y: Σ((y1 , . . . , yζ )T ) = (σ(y1 ), . . . , σ(yζ ))T

- On obtient un vecteur x0 ∈ Rζ , appelée sortie

Figure 1.2: Illustration du fonctionnement de la couche de neurones

Remarque. Pour tout j ∈ [[1, ζ]], le vecteur (wj,1 , . . . , wj,d ) ∈ Rd est appelé poids du
neurone j.

1.2.2 Perceptron Multi-Couche (PMC)

Un PMC est un type de réseau de neurones utilisé dans le cadre du stage. Il se compose
d’une succession de couches dans lesquelles la sortie de la j-ème couche correspond à
l’entrée de la (j + 1)-ème couche:

Figure 1.3: Illustration du fonctionnement d’un PMC comportant J couches de neurones,

l’entrée, la sortie et le poids du biais associés à la couche j sont respectivement notés x(j) ,
x0(j) et b(j)

Dans le cadre du stage, en notant W (j) la matrice poids associée à la couche j, nous
obtenons une fonction Loss qui dépendra des paramètres W (1) ,...,W (J) ,b(1) ,...,b(J) . La
stratégie est d’optimiser en ces paramètres.

Maxime BOUCHEREAU 7 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

1.3 Algorithme du gradient stochastique (SGD)

On cherche à minimiser la fonction en les paramètres du réseau de neurones, que l’on note
W:

K0 −1
1 X
LossT raining (W ) = L(W, k)
K0 k=0

2
1 ˆ(k) (k)
où L(W, k) = h2
y1 − y1 .

On se donne un réel α > 0, appelé taux d’apprentissage (learning rate en anglais).

L’algorithme du gradient stochastique est donné par la boucle suivante, et s’effectue en
un nombre d’itérations que l’on note Niter :

(W )0 est choisi au hasard

Pour tout t ∈ [[0, Niter − 1]]:

- it ∼ U([[0, K0 − 1]])
- (W )t+1 = (W )t − α∇L((W )t , it )

Cet algorithme présente l’avantage de ne pas avoir à calculer le gradient de LossT raining
à chaque itération, ce qui serait coûteux, en particulier lorsque K0 devient grand, ce qui
est le cas lorsque nous avons beaucoup de données, tout en assurant de descendre dans
suffisamment de directions afin de faire converger l’algorithme.

Remarques. - Dans le cadre du stage, il a été observé que l’algorithme converge si les
fonctions L(·, k) sont (strictement) convexes sur le domaine sur lequel se trouvent les
(W )t .

- Il existe des variantes de l’algorithme SGD. Par exemple, il est possible de calculer la
moyenne de plusieurs gradients de fonctions L(·, k) à chaque itération (petits lots), ou
bien de faire décroı̂tre le taux d’apprentissage au cours des itérations.

Maxime BOUCHEREAU 8 Université Rennes 1

Partie 2

Cas linéaire: Première formulation

2.1 Position du problème

Dans cette partie, on va se restreindre au cas linéaire, i.e. de la forme:

ẏ = Ay
y(0) = y0
où A ∈ Md (R)

Par ailleurs, on considèrera que N = 1, c’est-à-dire que l’on connaı̂t la solution du

(k) (k)
problème de Cauchy qu’en t = 0 et t = h, cette dernière étant donnée par: y1 = ehA y0 .

Ainsi, on utilise un réseau neuronal comportant une couche de d neurones, sans biais,
ce qui se traduit par une matrice poids W ∈ Md (R). Ainsi, nous avons:

ˆ(k) (k0)
y1 = Ah y0
où Ah est la matrice d’itération, qui dépend du schéma numérique employé:
- Pour la méthode d’Euler Explicite, on a:

Ah = Id + hW

donnant ainsi:

1
(k) 2
L(W, k) = Id + hW − ehA y0
h2
- Pour la méthode d’Euler Implicite, on a:

Ah = (Id − hW )−1

donnant ainsi:

1
(k) 2
L(W, k) = (Id − hW )−1 − ehA y0
h2

9
Master 2 Mathématiques Fondamentales Mars-Juin 2021

- Pour la méthode du Point Milieu, on a:

 −1  
h h
Ah = Id − W Id + W
2 2

donnant ainsi:

 −1   ! 2
1 h h hA (k)
L(W, k) = 2 Id − W Id + W −e y0
h 2 2

- Pour la méthode de Runge-Kutta 2, on a:

h2 2
Ah = Id + hW + W
2
donnant ainsi:

2
h2 2
 
1 hA (k)
L(W, k) = 2 Id + hW + W − e y0
h 2

Remarque. Si on écrit une fonction Loss de la façon suivante:

1 (k)
2
L(W, k) = 2
M y0
h
alors la condition pour que LossT raining (W ) = 0 ⇔ M = 0 (minimisant ainsi la Loss)
(k)
est d’avoir d données y0 linéairement indépendantes. Ainsi, peu de données suffisent
dans le cas linéaire.

2.2 Méthode d’Euler Explicite

Proposition (Fonction Loss - Méthode d’Euler Explicite)
Soit k ∈ [[0, K − 1]]. La fonction
1
(k) 2
L(·, k) : W 7→2
Id + hW − ehA y0
h
est convexe sur Md (R) et est presque partout strictement convexe sur Md (R).

Démonstration. L(·, k) est une forme quadratique en W , donc est une application régulière.
Tout d’abord, calculons la différentielle de L(·, k). Soient W, H ∈ Md (R):

1
(k) 2
L(W + H, k) = Id + h(W + H) − ehA y0
h2
1
(k) 2 2
(k) (k) (k)
2
= 2 Id + hW − ehA y0 + h Id + hW − ehA y0 |Hy0 i + Hy0
h h

Maxime BOUCHEREAU 10 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

Donc, on retient le terme linéaire en H pour obtenir la différentielle (le terme restant
est un terme quadratique):
2
(k) (k)
dL(W, k) : H 7→ h Id + hW − ehA y0 |Hy0 i
h
- Tout d’abord, montrons la convexité. Soient W1 , W2 ∈ Md (R):

1
(k) 2
L(W1 , k) + dL(W1 , k) · (W2 − W1 ) = Id + hW1 − ehA y0
h2
2
(k) (k)
+ h Id + hW1 − ehA y0 |(W2 − W1 )y0 i
h
1
(k) 2
6 2 Id + hW1 − ehA y0
h
2
(k) (k)
+ h Id + hW1 − ehA y0 |(W2 − W1 )y0 i
h
2
(k)
+ (W2 − W1 )y0
1
(k) (k)
2
6 2
Id + hW1 − ehA y0 + h(W2 − W1 )y0
h
1
(k) 2
6 2 Id + hW2 − ehA y0
h
6 L(W2 , k)

Donc L(·, k) est convexe sur Md (R).

- Montrons la stricte convexité presque partout sur Md (R).

Cette propriété se montre en reprenant le raisonnement précédent, et en utilisant le fait

(k)
que, pour presque tous W1 , W2 ∈ Md (R), y0 ∈ Rd :

2
(k)
(W2 − W1 )y0 > 0

puisque les matrices non inversibles sont de mesure nulle dans Md (R) (pour la mesure
(k)
de Lebesgue), et que y0 6= 0 presque partout sur Rd . Ainsi, on a, pour presque tous
(k)
W1 , W2 ∈ Md (R), y0 ∈ Rd :

L(W1 , k) + dL(W1 , k) · (W2 − W1 ) < L(W2 , k)

ce qui montre la stricte convexité presque partout sur Md (R)

La stricte convexité presque partout fournit l’existence d’au plus un minimiseur WM in

qui annule la fonction L(·, k), et à fortiori LossT raining , cette matrice étant la matrice
recherchée pour approximer A.

Maxime BOUCHEREAU 11 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

Proposition (Méthode d’Euler Explicite - Ordre de convergence)

L’apprentissage pour la méthode d’Euler Explicite est d’ordre 1:

WM in = A + O(h)
h→0

Démonstration. On rappelle que:

1
(k) 2
L(W, k) = Id + hW − ehA y0
h2
On a:

ehA − Id
WM in =
h
On fait un développement limité à l’ordre 2:

1
Id + hA + O(h2 ) − Id


WM in =
h→0 h
WM in = A + O(h)
h→0

Figure 2.1: Courbe de convergence en échelle log-log pour la méthode d’Euler Explicite,
obtenue avec une matrice A prise au hasard dans Md ([−10, 10]) dans le cas d = 2,
Niter = 500 et α = 10−3

Maxime BOUCHEREAU 12 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

2.3 Méthode d’Euler Implicite

Proposition (Fonction Loss - Méthode d’Euler Implicite en dimension 1)
(k)
Soit R > 0. Soient k ∈ [[0, K − 1]] et y0 =6 0. La fonction
 2
1 1 hA (k)2
L(·, k) : W 7→ 2 −e y0
h 1 − hW
est strictement convexe sur [−R, R] dès que h est assez petit.

Démonstration. On a, pour tout W ∈ R\ h1



 2
1 1 hA (k)2
L(W, k) = 2 −e y0
h 1 − hW
Cette fonction est deux fois différentiable sur R\ h1 et on a:


 
∂L 2 1 hA (k)2
(W, k) = 2
−e y0
∂W h(1 − hW ) 1 − hW
2
 
∂ L 2 2 hA (k)2
(W, k) = − e y 0
∂W 2 (1 − hW )3 1 − hW
Prenons h < R1 . Soit W ∈ [−R, R]. On a W 6 R < h1 , donc 1
1−hW
> 0. Donc L(·, k)
est strictement convexe si et seulement si:

2
− ehA > 0
1 − hW
(stricte positivité de la dérivée seconde). Or, nous avons:

2
− ehA > 0 ⇔ 1 − hW < 2e−hA
1 − hW
⇔ 2e−hA + hW − 1 > 0
Or, un développement limité en h assure que:

2e−hA + hW − 1 = 2(1 − hA) + hW − 1 + O(h2 )

h→0
= 1 + h(W − 2A) + O(h2 )
h→0

Donc, si h est assez petit, on a:

2e−hA + hW − 1 > 0
Donc L(·, k) est strictement convexe sur [−R, R].

Étudions maintenant une matrice WM in telle que LossT raining (WM in ) = 0:

Maxime BOUCHEREAU 13 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

Proposition (Méthode d’Euler Implicite - Ordre de convergence)

L’apprentissage pour la méthode d’Euler Implicite est d’ordre 1:

WM in = A + O(h)
h→0

Démonstration. On rappelle que:

1
(k) 2
L(W, k) = (Id − hW )−1 − ehA y0
h2

On a comme unique possibilité pour WM in :

Id − e−hA
WM in =
h

On fait un développement limité à l’ordre 2:

1
Id − Id + hA + O(h2 )


WM in =
h→0 h
WM in = A + O(h)
h→0

1
Remarque. En dimension 1, le point critique de L(·, k) vérifie 1−hW
− ehA = 0. donc, au
voisinage de ce point critique:

(k)2
∂ 2L 2y0 1
2
(W, k) = 3
·
∂W (1 − hW ) 1 − hW
(k)2
2y0
=
(1 − hW )4
> 0

donc on y a de la convexité stricte, la convergence vers ce point critique, minimiseur

via la méthode du gradient y est garantie. De façon plus générale, en partant d’un point
quelconque, la méthode de gradient va converger vers le point critique, et ce grâce à la
convexité sur tout segment lorsque h est assez petit.

Maxime BOUCHEREAU 14 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

Figure 2.2: Courbe de convergence en échelle log-log pour la méthode d’Euler Implicite,
obtenue avec une matrice A prise au hasard dans Md ([−10, 10]) dans le cas d = 2,
Niter = 500 et α = 10−3

2.4 Méthode du Point Milieu

Proposition (Fonction Loss - Méthode du Point Milieu en dimension 1)
(k)
Soit R > 0. Soient k ∈ [[0, K − 1]] et y0 6= 0. La fonction
!2
1 1 + h2 W hA
 2
(k)
L(·, k) : W 7→ 2 h
− e y0
h 1 − 2W
est strictement convexe sur [−R, R] dès que h est assez petit.
2
Démonstration. On a, pour tout W ∈ R\ h
:

!2
1 1 + h2 W 
(k)
2
L(W, k) = 2 − ehA y0
h 1 − h2 W
2
Cette fonction est deux fois différentiable sur R\ h
et on a:

!
∂L 2 1 + h2 W hA
 2
(k)
(W, k) =
2 h
−e y0
∂W h
h 1 − 2W 1 − 2W
!
∂ 2L 2  2
(k) 2 1 + h2 W hA
 2
(k)
2
(W, k) =
2 y0 +
3 h
−e y0
∂W h
1 − 2W h
1 − 2W 1 − 2W

Pour tout W 6= h2 :

Maxime BOUCHEREAU 15 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

2 !
∂ 2L
  
2 h h h  2
(k)
(W, k) =
4 1 − W + 1 + W − 1 − W ehA y0
∂W 2 1 − h2 W 2 2 2

2
Si h < R
, alors W 6 R < h2 , donc 1
h 4 > 0. Donc L(·, k) est strictement convexe
( 2 W)
1−
si et seulement si:

 2  
h h h
1− W + 1 + W − 1 − W ehA > 0
2 2 2

(stricte positivité de la dérivée seconde).

On fait un développement limité au voisinage de h = 0:

 2  
h h h h
1 − W + 1 + W − 1 − W ehA = 1 − hW + 1 + W + O(h2 )
2 2 2 h→0 2
 
h
1 − W 1 + hA + O(h2 )


−
2
 2  
h h h h h
1− W +1+ W − 1− W ehA = 2 − W − 1 + W − hA + O(h2 )
2 2 2 h→0 2 2
 2  
h h h
1− W +1+ W − 1− W ehA = 1 − hA + O(h2 )
2 2 2 h→0

Donc, si h est assez petit:

 2  
h h h
1− W + 1 + W − 1 − W ehA > 0
2 2 2

Donc L(·, k) est strictement convexe sur [−R, R].

On s’intéresse à une matrice WM in qui minimise LossT raining en l’annulant.

Proposition (Méthode du Point Milieu - Ordre de convergence)

L’apprentissage pour la méthode du Point Milieu est d’ordre 2:

WM in = A + O(h2 )
h→0

Démonstration. On rappelle que:

 −1   ! 2
1 h h (k)
L(W, k) = 2 Id − W Id + W − ehA y0
h 2 2

Maxime BOUCHEREAU 16 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

On a ainsi:

 −1  
h h
Id − WM in Id + WM in = ehA
2 2

On applique le logarithme matriciel:

( −1  )
h h
hA = log Id − WM in Id + WM in
2 2

On utilise les propriétés du logarithme (valables pour les matrices):

    
1 h h
A = log Id + WM in − log Id − WM in
h 2 2

On fait un développement limité à l’ordre 3:

(  2  2 !)
1 h 1 h h 1 h
A = WM in − WM in + O(h3 ) − − WM in − − WM in + O(h3 )
h→0 h 2 2 2 2 2 2
1
= hWM in + O(h3 )
h→0 h
= WM in + O(h2 )
h→0

Donc on a:

WM in = A + O(h2 )
h→0

1+ h W
Remarque. En dimension 1, le point critique de L(·, k) vérifie 2
1− h W
− ehA = 0. donc, au
2
voisinage de ce point critique:

(k)2
∂ 2L 2y0
(W, k) =
2
∂W 2 1 − h2 W
> 0

Donc on y a de la stricte convexité, ce qui garantit la convergence de la méthode de

gradient vers ce point critique. De façon plus générale, en partant d’un point quelconque,
la méthode de gradient va converger vers le point critique, et ce grâce à la convexité sur
tout segment lorsque h est assez petit.

Maxime BOUCHEREAU 17 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

Figure 2.3: Courbe de convergence en échelle log-log pour la méthode du Point Milieu,
obtenue avec une matrice A prise au hasard dans Md ([−10, 10]) dans le cas d = 2,
Niter = 500 et α = 10−3

2.5 Méthode de Rung-Kutta 2

Proposition (Fonction Loss - Méthode de Runge-Kutta 2 en dimension 1)
(k)
Soient k ∈ [[0, K − 1]] et y0 6= 0. La fonction
2
h2 2

1 hA (k)2
L(·, k) : W →
7 1 + hW + W − e y0
h2 2
1
est strictement convexe sur R si et seulement si A < 0 et h > |A|
log(2)

Démonstration. On a, pour tout W ∈ R:

2  
h2 2

1 hA (k)
2
L(W, k) = 2 1 + hW + W − e y0
h 2

Cette fonction est deux fois différentiable sur R et on a:

h2 2
  
∂L 2 hA (k)
2
(W, k) = (1 + hW ) 1 + hW + W − e y0
∂W h 2
2
 
∂ L hA 3 2 2  (k) 2
(W, k) = 2 2 − e + 3hW + h W y0
∂W 2 2
  
3 2 1 hA (k)
2
= 2 (1 + hW ) + − e y0
2 2

Maxime BOUCHEREAU 18 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

Une condition nécéssaire et suffisante de stricte convexité est que, pour tout W ∈ R,
∂2L
∂W 2
(W, k) > 0, i.e. ∆ < 0, où ∆ est le discriminant associé à la fonction polynomiale de
∂2L
degré 2 ∂W 2 (·, k). ∆ est donné par:

(k)4
∆ = 12h2 y0 2ehA − 1

Donc:

∂ 2L
(·, k) > 0 sur R ⇔ ∆ < 0
∂W 2
⇔ 2ehA − 1 < 0
⇔ hA < − log(2)

- Si A > 0, c’est impossible

log(2)
- Si A < 0, h > |A|

Remarque. Ce résultat nécéssite de ne pas prendre h trop petit et on n’a pas de convexité
globale quand A > 0

Proposition (Fonction Loss - Méthode de Runge-Kutta 2 en dimension 1)

(k)
Soit R > 0. Soient k ∈ [[0, K − 1]] et y0 6= 0. Si on a:

    
1 1 1 5
h < min 1− √ , log
R 2 |A| 4

alors la fonction:
2  
h2 2

1 hA (k)
2
L(·, k) : W →
7 1 + hW + W − e y0
h2 2
est strictement convexe sur [−R, +∞[.

1
log 45 donc ehA 6 eh|A| < 54 , soit 21 − ehA > − 34


Démonstration. - h < |A|
 
1 1
- De plus, si h 6 R 1 − 2 , alors, pour tout W ∈ [−R, 0[:
√

   
1 1 1 1
h6 1− √ = √ −1
|W | 2 W 2

 
1 1
h 6 1− √
|W | 2
 
1 1
h 6 √ −1
W 2

Maxime BOUCHEREAU 19 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

Donc hW > √1 − 1 i.e. 1 + hW > √1 . Ainsi, on a:

2 2

1
(1 + hW )2 >
2

Donc, pour tout W ∈ [−R, 0[:

 

∂ 2L 3
2 1 
hA 
(W, k) = 2 (1 + hW ) + − e >0

∂W 2  2 | {z } |2 {z }

1
>2
>− 34

Donc L(·, k) est strictement convexe sur [−R, 0[.

1
log 54 6 |A|
1 1


- Enfin, h < |A| log (2) donc hA 6 h|A| < log(2) i.e. ehA < 2 soit 2
− ehA >
− 23 .
Donc, si W > 0:
 

∂ 2L 3
2 1 
hA 
(W, k) = 2 (1 + hW ) + − e >0

∂W 2  2 | {z } |2 {z }

>1
>− 32

Donc L(·, k) est strictement convexe sur [0, +∞[.

Donc L(·, k) est strictement convexe sur [−R, +∞[.

Figure 2.4: Lorsque h n’est pas assez petit, l’algorithme du gradient stochastique donne
la convergence vers le mauvais minimiseur de la Loss

Maxime BOUCHEREAU 20 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

Figure 2.5: Lorsque h est suffisamment petit, l’algorithme du gradient stochastique donne
la convergence vers le bon minimiseur de la Loss, et ce grâce à la convexité.

Proposition (Fonction Loss - Méthode de Runge-Kutta 2 en dimension supérieure)

(k)
Soit R > 0. Soient k ∈ [[0, K − 1]] et y0 6= 0. La fonction
2
1 h2 2 (k)2
L(·, k) : W 7→ 2
Id + hW + W − ehA y0
h 2
est convexe sur B||·|| (0, R) où || · || est une norme matricielle quelconque.

Démonstration. On a, pour tout W ∈ Md (R):

T 
h2 2 h2 2
 
1 (k)T hA hA (k)
L(W, k) = y Id + hW + W − e Id + hW + W − e y0
h2 0 2 2
1 (k)T 2 T

2

(k)
= y h(W − A) + O(h ) h(W − A) + O(h ) y0
h→0 h2 0
2
(k)
= (W − A)y0 + O(h)
h→0

Donc, si h est assez petit, L(·, k) ressemble à une forme quadratique en W → 0d , donc
on a de la convexité autour de cette matrice.

On s’intéresse à une matrice WM in qui minimise LossT raining en l’annulant.

Proposition (Méthode de Runge-Kutta 2 - Ordre de convergence)

L’apprentissage pour la méthode de Runge-Kutta 2 est d’ordre 2:

WM in = A + O(h2 )
h→0

Maxime BOUCHEREAU 21 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

Démonstration. On rappelle que:

2
h2 2
 
1 hA (k)
L(W, k) = 2 Id + hW + W − e y0
h 2

On a ainsi:

h2 2
Id + hWM in + WM in = ehA
2
On applique le logarithme matriciel:

h2 2
 
hA = log Id + hWM in + WM in
2
h2 2
 
1
A = log Id + hWM in + WM in
h 2

On fait un développement limité à l’ordre 3:

( 2 )
2 2

1 h 2 1 h 2
A = hWM in + WM in − hWM in + WM in + O(h3 )
h→0 h 2 2 2
h2 2 h2 2 h3 3 h4 4
 
1 3
= hWM in + WM in − WM in − WM in − WM in + O(h )
h→0 h 2 2 2 8
1
= hWM in + O(h3 )
h→0 h
= WM in + O(h2 )
h→0

Donc on a:

WM in = A + O(h2 )
h→0

Remarque. En dimension 1, le point critique de L(·, k) vérifie 1 + hW + 21 h2 W 2 − ehA = 0

i.e. 21 [1 + (1 + hW )2 ] = ehA donc, au voisinage de ce point critique:

∂ 2L (k)2
(W, k) = 2y 0 (1 + hW )2 > 0
∂W 2
Donc on y a de la stricte convexité, ce qui garantit la convergence de la méthode de
gradient vers ce point critique. De façon plus générale, en partant d’un point quelconque,
la méthode de gradient va converger vers le point critique, et ce grâce à la convexité sur
tout segment lorsque h est assez petit.

Maxime BOUCHEREAU 22 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

Figure 2.6: Courbe de convergence en échelle log-log pour la méthode de Runge-Kutta

2, obtenue avec une matrice A prise au hasard dans Md ([−10, 10]) dans le cas d = 2,
Niter = 500 et α = 10−3

Maxime BOUCHEREAU 23 Université Rennes 1

Partie 3

Cas linéaire: Nouvelle formulation et

méthodes supplémentaires

3.1 Position du problème

Dans cette partie, on va toujours se restreindre au cas linéaire:

ẏ = Ay
y(0) = y0
où A ∈ Md (R)

En revanche, on considèrera que N > 1 quelconque, c’est-à-dire que l’on connaı̂t la

solution du problème de Cauchy en des temps 0 = t0 < t1 < . . . < tN = h, i.e. on connaı̂t
(k) (k)
la solution en des points intermédiaires, cette dernière étant donnée par: y tj = etj A y0 .
h

Ainsi, on utilise toujours un réseau neuronal comportant une couche de d neurones,

sans biais, ce qui se traduit par une matrice poids W ∈ Md (R). Cependant, nous allons
considérer ces formes pour le vecteur yˆ1 (k) :

N
X N
X
(k) (k) (k)
yˆ1 = αj y t j + W βj hj y tj (3.1)
h h
j=0 j=0

avec h0 + · · · + hn = h. Ainsi, la fonction L(·, k) est donnée, pour tout W ∈ Md (R) par:

N N 2
1 X (k) (k)
X (k)
L(W, k) = 2 αj y tj − y1 + W βj hj y tj
h j=0 h
j=0 h

- Pour la méthode d’Euler Explicite, on a:

(k) (k)
yˆ1 (k) = y0 + hW y0

donnant ainsi:

1
(k) 2
L(W, k) = Id + hW − ehA y0
h2

24
Master 2 Mathématiques Fondamentales Mars-Juin 2021

- Pour la méthode d’Euler Implicite, on a:

(k) (k)
yˆ1 (k) = y0 + hW y1

donnant ainsi:

1
(k) 2
L(W, k) = (Id − hW )ehA − Id y0
h2

- Pour la méthode du Point Milieu, on a:

(k) h  (k) (k)


yˆ1 (k) = y0 + W y0 + y1
2

donnant ainsi:

  2
1 hA h hA (k)
L(W, k) = 2 e − Id − W (Id + e ) y0
h 2

3.2 Convexité des fonctions Loss

(k)
L’écriture de y1 sous la forme (3.1) permet d’obtenir une fonction L(·, k) sous la forme
d’une forme quadratique en W :

1 (k) 2
L(W, k) = u + W v (k)
h2

Proposition (Fonction Loss - Forme quadratique en W )

Soit k ∈ [[0, K − 1]]. La fonction
1 (k) 2
L(·, k) : W 7→
2
u + W v (k)
h
est convexe sur Md (R) et est presque partout strictement convexe sur Md (R).

Démonstration. Cette propriété se montre de la même manière que pour la fonction Loss
associée à la méthode d’Euler Explicite dans la partie précédente, en utilisant la différen-
tielle de la fonction.

Remarque. Dans la plupart des cas rencontrés, la fonction Loss s’écrit sous la forme:

1 (k) (k)
2
L(W, k) = 2
By0 + W Cy0
h

Maxime BOUCHEREAU 25 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

où B, C ∈ Md (R). Ainsi, la propriété précédente de stricte convexité donne l’existence

d’au plus un minimiseur, que l’on peut expliciter, dans le cas où la matrice C est inversible,
par:

WM in = −BC −1

car ce minimiseur annule la fonction L(·, k).

3.3 Méthodes à deux points

(k) (k)
Dans cette sous-partie, on suppose que seules y0 et y1 sont données (N = 1). On peut
en déduire les trois méthodes principales: Euler Explicite, Euler Implicite et Point Milieu.

Proposition (Méthodes d’Euler et du Point Milieu - Ordres de convergence)

- L’apprentissage pour la méthode d’Euler Explicite est d’ordre 1:

WM in = A + O(h)
h→0

- L’apprentissage pour la méthode d’Euler Implicite est d’ordre 1:

WM in = A + O(h)
h→0

- L’apprentissage pour la méthode du Point Milieu est d’ordre 2:

WM in = A + O(h2 )
h→0

Démonstration. - Pour la méthode d’Euler Explicite, on raisonne comme dans la propo-

sition de la partie précédente.

- Pour la méthode d’Euler Implicite, même si la Loss est différente, on a toujours:

Id − ehA
WM in =
h

et un développement limité à l’ordre 2 donne le résultat souhaité, comme dans la partie

précédente

- Pour la méthode du Point Milieu, on a:

2
WM in = (Id + ehA )−1 (ehA − Id )
h

Maxime BOUCHEREAU 26 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

qui annule L(·, k). On fait un développement limité à l’ordre 3:

−1 
h2 2 h2 2
 
2 3 3
Wmin = 2Id + hA + A + O(h ) hA + A + O(h )
h→0 h 2 2
−1
h2 2 h2 2
   
1 h 3 3
= Id + A + A + O(h ) hA + A + O(h )
h→0 h 2 4 2
" #
2
h2 2 h2 2 h2 2
  
1 h h 3 3
= Id − A − A + A+ A + O(h ) hA + A + O(h )
h→0 h 2 4 2 4 2
h2
  
1 h
= Id − A + O(h3 ) hA + A2 + O(h3 )
h→0 h 2 2
2 2
 
1 h 2 h 2 3
= hA − A + A + O(h )
h→0 h 2 2
1
hA + O(h3 )


=
h→0 h
= A + O(h2 )
h→0

Figure 3.1: Courbe de convergence en échelle log-log pour la méthode d’Euler Implicite
(Loss convexe), obtenue avec une matrice A prise au hasard dans Md ([−10, 10]) dans le
cas d = 2, Niter = 500 et α = 10−3

Maxime BOUCHEREAU 27 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

Figure 3.2: Courbe de convergence en échelle log-log pour la méthode du Point Milieu
(Loss convexe), obtenue avec une matrice A prise au hasard dans Md ([−10, 10]) dans le
cas d = 2, Niter = 500 et α = 10−3

La proposition suivante donne l’odre maximal que l’on peut atteindre avec la forme
(3.1) pour yˆ1 (k) :

Proposition (Ordre de convergence maximal avec deux points)

Si on prend yˆ1 (k) de la forme:

 
(k) (k) (k) (k)
yˆ1 (k) = α0 y0 + α1 y1 + hW β0 y0 + β1 y1

Alors, l’apprentissage est, au mieux, d’ordre 2.

Démonstration. On a, pour tout k ∈ [[0, K − 1]]:

1 
 (k) 2
L(W, k) = 2
α0 Id + (α1 − 1)ehA + hW β0 Id + β1 ehA y0
h
On a ainsi:

1
−1
β0 Id + β1 ehA α0 Id + (α1 − 1)ehA


WM in = −
h
1
∼ − [(β0 + β1 )Id + β1 hA]−1 [(α0 + α1 − 1)Id + (α1 − 1)hA] (3.2)
h→0 h
On doit avoir α0 + α1 = 1 pour que le troisième facteur de (3.2) soit en O(h) et
”compense” le premier facteur en h1 . on a ainsi:
 
α0 β1 2 α0
WM in ∼ Id − hA + O(h ) A ∼ A
h→0 β0 + β1 β0 + β1 h→0 β0 + β1

Maxime BOUCHEREAU 28 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

On doit avoir α0 = β0 + β1 . L’objectif est de monter en ordre. Comme on a:


α0 + α1 = 1
β0 + β1 = α0
On obtient:

1 −1 
β0 Id + β1 ehA α0 Id − α0 ehA

WM in = −
h
On fait un développement limité à l’ordre 4:

−1
h2 2 h3 3

4
WM in = (β0 + β1 )Id + β1 hA + β1 A + β1 A + O(h )
h→0 2 6
2 3
 
h 2 h 3 h 4 4
· α0 A + A + A + A + O(h )
2 6 24
α0
En posant b = β0 +β1
, on a:

−1
h2 2 h3 3
 
α0 4
WM in = Id + b hA + A + A + O(h )
h→0 β0 + β1 2 6
2 3
 
h 2 h 3 h 4 4
· A + β1 A + A + A + O(h )
2 6 24
−1
h2 2 h3 3
 
4
WM in = Id + b hA + A + A + O(h )
h→0 2 6
2 3
 
h 2 h 3 h 4 4
· A + A + A + A + O(h )
2 6 24
En écrivant WM in sous forme d’une série de Taylor, on a, à l’ordre 4:

   
1 2 2 1
Wmin = A+ − b hA + b − b + h2 A3
h→0 2 6
 
3 3 2 7 1
+ −b + b − b + h3 A4 + O(h4 )
2 12 24
Si on veut de l’ordre 2, on doit avoir b = 21 , et les conditions pour avoir de l’ordre 2
sont:
 
 α0 + α1 = 1  α0 = 2β0
β0 + β1 = α 0 ⇔ α0 + β1 = α0
β0 = β1 β1 = β0
 

Ce qui correspond à la méthode du Point Milieu.

En plus de ces conditions, si on veut avoir de l’ordre 3, on doit avoir:

1
b2 − b + = 0
6
ce qui est impossible avec b = 21 , donc on ne peut pas monter en ordre.

Maxime BOUCHEREAU 29 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021


Remarques. - si α0 + α1 − 1 = α1 − 1 = 0, alors (α0 , α1 ) = (0, 1):

 
(k) (k) (k)
yˆ1 (k) = y1 + hW β0 y0 + β1 y1

1 
(k) (k)
 2
L(W, k) = hW β y
0 0 + β y
1 1
h2
Donc WM in = 0, ce qui n’approche pas forcément A
- Si β0 + β1 = 0 et α0 + α1 = 1, alors, dans l’expression (3.2), on a:

1 − α1
WM in = Id + O(1) si β1 6= 0
h→0 hβ1
α0
= Id + O(1)
h→0 hβ1

Cela impose α0 = 0, donc α1 = 1, et on se ramène au cas précédent, et si β1 = 0, alors

β0 = 0 et yˆ1 (k) ne dépend pas de W .

3.4 Méthodes Multi-pas (explicites) à pas constant

3.4.1 Ordre de convergence
jh
On suppose que, pour tout j ∈ [[0, N ]], tj = N
, et que la solution y est connue en les
points tj , i.e. pour tout k ∈ [[0, K − 1]]:

(k) j (k)
yj = e N A y0
N

Posons alors:

N −1 N −1
ˆ(k) X (k) h X (k)
y1 = αi y1− i+1 + W βi y1− i+1
i=0
N N i=0
N

On obtient ainsi:

N −1 N −1 2
1 X (k) (k) h X (k)
L(W, k) = 2 αi y1− i+1 − y1 + W βi y1− i+1
h i=0 N N i=0
N

et, comme minimiseur:

N −1
!−1 N −1
!
N X N −i−1
hA
X N −i−1
hA
WM in = − βi e N αi e N − ehA
h i=0 i=0

Maxime BOUCHEREAU 30 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

Proposition (Méthode Multi-pas à pas constant - Ordre de convergence)

Si on a:

-
N
X −1
αi = 1
i=0

- Pour tout q ∈ [[1, p]]:

N
X −1
iq αi − qiq−1 βi = (−1)q
i=0

Alors:

WM in = A + O(hp )
h→0

Démonstration. On a:

N −1
!−1 N −1
!
N X N −i−1
hA
X −i
hA hA
WM in = − βi e N αi e N −e
h i=0 i=0
N −1
!−1 (N −1 )!
N N −1
hA
X
− Ni hA N −1
hA
X
− Ni hA h
A
WM in = − e N βi e e N αi e −e N
h i=0 i=0

N −1
hA
Les matrices commutent (polynômes en A), donc on peut simplifier par e N :

N −1
!−1 N −1
!
N X N −i−1
hA
X −i
hA h
A
WM in = − βi e N αi e N −e N
h i=0 i=0

On fait un développement limité à l’ordre p du numérateur et du dénominateur:

N p 
−1 q q !−1
N X
X i h
WM in = − βi − A + O(hp+1 )
h→0 h i=0 q=0 N q!
−1 Xp 
N  q  q  !
X i 1
· αi − A − A hq + O(hq+1 )
i=0 q=0
N N
p   q "N −1
# !−1
N X1 A X
= − (−i)q βi hq + O(hp+1 )
h→0 h q=0 q! N i=0
p
"N −1 # !
X1 A  q X
q q p+1
· (−i) αi − 1 h + O(h )
q=0
q! N i=0

Maxime BOUCHEREAU 31 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

Dans le deuxième facteur, on sépare la somme en deux parties: un premier terme pour
q = 0 et un second terme pour q allant de 1 à p:

p  q "N −1
# !−1
X 1 A X
WM in = −N (−i)q βi hq + O(hp+1 )
h→0
q=0
q! N i=0
−1 p   q "N −1
N
# !
1X 1 X1 A X
· αi − + (−i)q αi − 1 hq−1 + O(hp )
h i=0 h q=1 q! N i=0
p  q N " −1
# !−1 N −1
!
N X1 A X X
= − (−i)q βi hq + O(hp+1 ) · αi − 1
h→0 h q=0 q! N i=0 i=0
p  q "N −1
# !−1
X 1 A X
− N (−i)q βi hq + O(hp+1 )
q=0
q! N i=0
p
" −1 # !
 q−1 N
A X1 A X
· (−i)q αi − 1 hq−1 + O(hp )
N q=1 q! N i=0

à la dernière ligne, on fait le changement d’indice q 7→ q − 1 dans la somme:

p  q "N −1
# !−1 N −1
!
N X 1 A X X
WM in = − (−i)q βi hq + O(hp+1 ) · αi − 1
h→0 h q=0
q! N i=0 i=0
p−1  q " N −1
# !−1
X 1 A X
+ A − (−i)q βi hq + O(hp )
q=0
q! N i=0
p−1
"N −1 # !
X1 A   q
1 X
q+1 q p
· (−i) αi − 1 h + O(h )
q=0
q! N q + 1 i=0

Donc, si on a ces conditions:

-
N
X −1
αi − 1 = 0
i=0

- Pour tout q ∈ [[0, p − 1]]:

−1 −1
N N
!
X 1 X
− (−i)q βi = (−i)q+1 αi − 1
i=0
q+1 i=0

Alors on a:

WM in = A (Id + O(hp ))
h→0
WM in = A + O(hp )
h→0

Maxime BOUCHEREAU 32 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

en particulier, la deuxième condition est équivalente à:

−1 −1
N N
!
X 1 X
∀q ∈ [[1, q]], − (−i)q−1 βi = (−i)q αi − 1
i=0
q i=0

N
X −1
⇔ ∀q ∈ [[1, q]], (−1)q iq αi + q(−1)q−1 iq−1 βi = 1
i=0

N
X −1
⇔ ∀q ∈ [[1, q]], iq αi − qiq−1 βi = (−1)q
i=0


(k) (k) (k)
Remarque. Concrètement, pour un problème avec N + 1 points donnés y0 , y 1 , . . . , y1 ,
N
pour avoir l’ordre p, il faut que:

α0 + α1 + α2 + · · · + αN −1 = 1
α1 + 2α2 + (N − 1)αN −1 − (β0 + · · · + βN −1 ) = −1
.. ..
. .
α1 + 2 α2 + · · · + (N − 1) αN −1 − p β1 + 2 β2 + · · · + (N − 1) βN −1 = (−1)p
p p p−1 p−1

3.4.2 Exemples
- Méthode de Nyström: On se place dans le cas N = 2 (1 point intermédiaire donné), et
on prend les coefficients:

(α0 , α1 , β0 , β1 ) = (0, 1, 2, 0)

Ainsi, avec:

ˆ(k) (k) h  (k) 

y1 = y0 + W 2y 1
2 2

l’apprentissage est d’ordre 2:

WM in = A + O(h2 )
h→0

- Une méthode d’ordre 3: On se place dans le cas N = 2 (1 point intermédiaire donné),

et on prend les coefficients:

(α0 , α1 , β0 , β1 ) = (−4, 5, 4, 2)

Maxime BOUCHEREAU 33 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

Figure 3.3: Courbe de convergence en échelle log-log pour la méthode de Nyström (Loss
convexe), obtenue avec une matrice A prise au hasard dans Md ([−10, 10]) dans le cas
d = 2, Niter = 500 et α = 10−3

Ainsi, avec:

ˆ(k) (k) (k) h  (k) (k)


y1 = −4y 1 + 5y0 + W 4y 1 + 2y0
2 2 2

l’apprentissage est d’ordre 3:

WM in = A + O(h3 )
h→0

Maxime BOUCHEREAU 34 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

Figure 3.4: Courbe de convergence en échelle log-log pour la méthode Multipas à pas
constant d’ordre 3 (Loss convexe), obtenue avec une matrice A prise au hasard dans
Md ([−10, 10]) dans le cas d = 2, Niter = 500 et α = 10−3

- Méthode de Milne: On se place dans le cas N = 4 (3 points intermédiaire donné), et

on prend les coefficients:

 
8 4 8
(α0 , α1 , α2 , α3 , β0 , β1 , β2 , β3 ) = 0, 0, 0, 1, , − , , 0
3 3 3

Ainsi, avec:

 
ˆ(k) (k) h 8 (k) 4 (k) 8 (k)
y1 = y0 + W y3 − y1 + y1
4 3 4 3 2 3 4

l’apprentissage est d’ordre 4:

WM in = A + O(h4 )
h→0

Maxime BOUCHEREAU 35 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

Figure 3.5: Courbe de convergence en échelle log-log pour la méthode de Milne (Loss
convexe), obtenue avec une matrice A prise au hasard dans Md ([−10, 10]) dans le cas
d = 2, Niter = 500 et α = 10−3

- Une méthode d’ordre 5: On se place dans le cas N = 3 (2 points intermédiaire donné),

et on prend les coefficients:

(α0 , α1 , α2 , β0 , β1 , β2 ) = (−18, 9, 10, 9, 18, 3)

Ainsi, avec:

ˆ(k) (k) (k) (k) h  (k) (k) (k)


y1 = −18y 2 + 9y 1 + 10y0 + W 9y 2 + 18y 1 + 3y0
3 3 3 3 3

l’apprentissage est d’ordre 5:

WM in = A + O(h5 )
h→0

Maxime BOUCHEREAU 36 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

Figure 3.6: Courbe de convergence en échelle log-log pour la méthode Multipas à pas
constant d’ordre 5 (Loss convexe), obtenue avec une matrice A prise au hasard dans
Md ([−10, 10]) dans le cas d = 2, Niter = 1000 et α = 10−4

3.5 Méthodes Multi-pas à pas non constant

On suppose la solution connue aux temps 0 = t0 < t1 < · · · < tN = h, où cette foix on ne
suppose plus forcément que tj = jh
N
:

(k) (k)
y tj = etj A y0
h

3.5.1 Ordre de convergence

On considère les polynômes interpolateurs de Lagrange:
Y t − ti
∀j ∈ [[0, N ]], Lj,N (t) = ∈ RN [t]
tj − ti
06i6N
i 6= j

Ainsi, pour tout i ∈ [[0, N ]], Lj,N (ti ) = δi,j

Soit ainsi:

N
X
(k) (k)
PN (t) = Lj,N (t)W y tj ∈ R[t]
h
j=0

(k) (k)
t 7→ PN (t) doit approcher t 7→ W etA y0 sur [0, h] par interpolation aux points t0 , t1 , · · · , tN .
Posons, pour tout k ∈ [[0, K − 1]]:

Maxime BOUCHEREAU 37 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

Z h
(k) (k) (k)
yˆ1 = y0 + PN (t)dt
0
N Z
X h
(k) (k)
= y0 +W Lj,N (t)dty tj
j=0 0 h

On a ainsi:

N Z h  2
1 hA (k)
X (k)
L(W, k) = 2 (Id − e )y0 + W Lj,N (t)dt etj A y0
h j=0 0

Proposition (Méthode Multi-pas à pas non constant - Ordre de convergence)

Si A ∈ GLd (R) (ce qui est le cas presque partout dans Md (R)), l’apprentissage pour
cette méthode est d’ordre N + 1:

WM in = A + O(hN +1 )
h→0

Démonstration. Tout d’abord, notons que:

N Z
!−1
X h
Lj,N (t)dtetj A ehA − Id


WM in =
j=0 0

Soit la matrice de polynômes:

N
X
QN = etj A Lj,N (X) ∈ RN [X]d×d
i=0

Soit t ∈ [0, h]. Par le théorème de Taylor, il existe ξ ∈]t, h[ tel que:

1
etA − QN (t) = AN +1 eξA πN (ξ)
(N + 1)!
avec:

N
Y
πN (ξ) = (ξ − tj )
j=0

Or, on a:

Z h −1
ehA − Id


WM in = QN (t)dt
0
Z h Z h −1
tA tA
ehA − Id


= e dt − e − QN (t)dt
0 0

Maxime BOUCHEREAU 38 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

Donc, pour tout t ∈ [0, h]:

1
etA − QN (t) 6 ||A||N +1 eh||A|| hN +1
(N + 1)!

Or, on a:

Z h
etA dt = A−1 ehA − Id


0

Notons que toutes les matrices utilisées commutent entre elles, puisque ce sont des polynômes
en A, donc on a:

 Z h −1
−1 hA tA
ehA − Id



WM in = A e − Id − e − QN (t)dt
0
 Z h −1
hA

−1 hA

−1 tA
ehA − Id


= A e − Id Id − A e − Id e − QN (t)dt
0
 Z h −1
hA

−1 tA
= A Id − A e − Id e − QN (t)dt
0

Or, on a:

−1 1
A ehA − Id ∼ Id
h→0 h
Z h
etA − QN (t)dt = O(hN +2 )
0 h→0

Donc on a:

−1
WM in = A Id + O(hN +1 ) = A + O(hN +1 )
h→0 h→0

3.5.2 Exemples
- Une méthode d’ordre 3
On prend N = 2 (1 point intermédiaire), et on connaı̂t les solutions en 0, h3 et h.
L’apprentissage est d’ordre 3:

WM in = A + O(h3 )
h→0

Maxime BOUCHEREAU 39 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

Figure 3.7: Courbe de convergence en échelle log-log pour la méthode Multipas à pas
variable d’ordre 3 (Loss convexe), obtenue avec une matrice A prise au hasard dans
Md ([−10, 10]) dans le cas d = 2, Niter = 500 et α = 10−3

- Une méthode d’ordre 4

On prend N = 3 (2 points intermédiaires), et on connaı̂t les solutions en 0, h4 , h2 et h.
L’apprentissage est d’ordre 4:

WM in = A + O(h4 )
h→0

Figure 3.8: Courbe de convergence en échelle log-log pour la méthode Multipas à pas
variable d’ordre 4 (Loss convexe), obtenue avec une matrice A prise au hasard dans
Md ([−10, 10]) dans le cas d = 2, Niter = 500 et α = 10−3

Maxime BOUCHEREAU 40 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

3.6 Généralisation
3.6.1 Ordre de convergence
Dans cette sous-partie, on veut déterminer l’ordre d’apprentissage directement à partir
de l’ordre de la méthode numérique, tant qu’elle est sous la forme (3.1). Par exemple,
vu que les méthodes d’Euler sont d’ordre 1, l’apprentissage est d’ordre 1, pour le Point
Milieu, qui est une méthode d’ordre 2, l’apprentissage est lui aussi d’ordre 2.

Théorème (Ordre de convergence - Cas linéaire)

Si on prend yˆ1 de la forme (3.1):

N
X N
X
yˆ1 = αj y t j + W βj hj y tj
h h
j=0 j=0

où la méthode est d’ordre p, c’est-à-dire que, si on pose:

N N
∼ X X
y1 = αj y t j + A βj hj y tj
h h
j=0 j=0

on a:

∼
y1 − y1 = O(hp+1 )
h→0

et si on suppose que:

N
X hj
βj 6= 0
j=0
h

Alors, l’apprentissage est d’ordre p:

WM in = A + O(hp )
h→0

Démonstration. On a cette expression pour la fonction L(·):

N N 2
1 X X
L(W ) = 2 αj y tj − y1 + W βj hj y tj
h j=0 h
j=0
h

Remarque. Par souci d’allègement des notations, l’indice (k) de la donnée a été supprimé.
(k)
Dans le cas d’une vraie donné y0 , on écrirait L(W, k).

Ainsi, on a:

Maxime BOUCHEREAU 41 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

N
! N
X X
WM in βj hj y tj = y1 − αj y tj
h h
j=0 j=0
N
! N N
!
X X X
(WM in − A) βj hj y tj = y1 − αj y tj − A βj hj y tj
h h h
j=0 j=0 j=0
∼
= y1 − y1

WM in peut-être construite à partir de A et de etj A par somme, produit ou inversion,

donc à partir de polynômes en A, c’est pour cette raison que l’on peut faire commuter les
facteurs des produits. Ainsi, on a:

N N
X ∼ X
βj hj WM in etj A y0 = y1 − y1 + βj hj Aetj A y0
j=0 j=0
N
! N
!
X X ∼
WM in βj hj etj A y0 = A βj hj etj A y0 + y1 − y1
j=0 j=0
N
!−1
X  ∼

WM in y0 = Ay0 + βj hj etj A y1 − y1
j=0
N
!−1
X hj 1 ∼

WM in y0 = Ay0 + βj etj A y1 − y1
j=0
h h

PN hj
Comme j=0 βj h
6= 0, on a:

N N
!
X hj X hj
βj etj A −→ βj Id
j=0
h h→0
j=0
h

qui est une matrice inversible (perturbation d’un multiple non nul de l’identité). Par
∼
ailleurs, comme on a y1 − y1 = O(hp+1 ), on obtient ainsi:
h→0

WM in = A + O(hp )
h→0

Maxime BOUCHEREAU 42 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

3.6.2 Exemples
Une méthode d’ordre 2
Si on prend cette méthode numérique:

h
yˆ1 = y0 + W (y0 + 2y 1 + y1 )
4 2

∼ h
y1 = y0 + A(y0 + 2y 1 + y1 )
4 2

Alors on a:

∼ h  h

Id − ehA y0 + A Id + 2e 2 A + ehA y0


y1 − y1 =
4
h  h
2
Id − ehA y0 + A Id + e 2 A y0


=
4
On a donc:

∼
y1 − y1 = O(h3 )
h→0

qui est de l’ordre 3 localement, ainsi, l’apprentissage est d’ordre 2:

Wmin = A + O(h2 )
h→0

Figure 3.9: Courbe de convergence en échelle log-log pour l’exemple, obtenue avec une
matrice A prise au hasard dans Md ([−10, 10]) dans le cas d = 2, Niter = 1000 et α = 10−3 .

Maxime BOUCHEREAU 43 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

Une méthode d’ordre 4

Si on prend cette méthode numérique:

h
yˆ1 = y0 + W (y0 + 3y 1 + 3y 2 + y1 )
8 3 3

∼ h
y1 = y0 + A(y0 + 3y 1 + 3y 2 + y1 )
8 3 3

Alors on a:

∼ h  h 2h

Id − ehA y0 + A Id + 3e 3 A + 3e 3 A + ehA y0


y1 − y1 =
8
h  h
3
Id − ehA y0 + A Id + e 3 A y0


=
8

On a donc:

∼
y1 − y1 = O(h5 )
h→0

qui est de l’ordre 5 localement, ainsi, l’apprentissage est d’ordre 4:

Wmin = A + O(h4 )
h→0

Maxime BOUCHEREAU 44 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

Figure 3.10: Courbe de convergence en échelle log-log pour l’exemple, obtenue avec une
matrice A prise au hasard dans Md ([−10, 10]) dans le cas d = 2, Niter = 1000 et α = 10−3 .

Maxime BOUCHEREAU 45 Université Rennes 1

Partie 4

Cas non linéaire

Dans cette partie, nous allons nous intéresser au cas général:


ẏ = f (y)
y(0) = y0
(k) (k) (k)
et on suppose que l’on connaı̂t K données y0 (conditions initiales) et y1 = φh (y0 )
(solution au temps h, où h est le pas de la méthode numérique choisie).

4.1 Structure du réseau de neurones

Afin d’avoir une idée du réseau de neurones le plus simple à créer, nous allons utiliser le
théorème d’universalité:

Théorème (Théorème d’universalité)

Soit g ∈ C(Rd , R), soit σ : R −→ R une fonction continue non polynomiale. Alors,
pour tout ε > 0, il existe gapp ∈ Mσ telle que, pour tout Ω ⊂ Rd compact:

||g − gapp ||L∞ (Ω) 6 ε

avec:

Mσ = V ect x 7→ σ(w · x + b), w ∈ Rd , b ∈ R



Remarques. - σ est appelée fonction d’activation. On peut par exemple choisir une sig-
moı̈de:

1
σ(x) =
1 + e−x

mais on retiendra, dans le cadre du stage, la tangente hyperbolique:

σ(x) = tanh(x)

46
Master 2 Mathématiques Fondamentales Mars-Juin 2021

- Vu que notre champ de vecteurs est à valeurs dans Rd et que le théorème s’applique à
une fonction à valeurs dans R, on doit appliquer ce théorème à chaque composante de
f.
Ainsi, avec ce théorème, on peut construire un réseau de neurones qui approchera
correctement notre champ de vecteurs f par un champ approché fapp . Voici la structure
du réseau de neurones, qui comporte deux couches:

- Une première couche de ζ neurones, et donc les paramètres sont W (1) ∈ Mζ,d (R) (ma-
(1)
trice poids des ζ neurones), ainsi qu’un biais de poids w0 ∈ Rζ

- Une seconde couche de d neurones, et donc les paramètres sont W (2) ∈ Md,ζ (R) (matrice
poids des ζ neurones), ne comportant pas de biais.

Figure 4.1: Illustration du fonctionnement du réseau de neurones permettant d’approcher

fapp pour le cas général. La fonction Σ : Rζ → Rζ applique la fonction σ à chacune des
composantes du vecteur qu’elle prend en entrée

Ainsi, le champ de vecteur approché fapp s’écrit:

 
(2) (1) (1)
fapp (y0 ) = W Σ W y0 + w0

Par exemple, pour la méthode d’Euler Explicite, on a cette valeur prédite en h:

ˆ(k)
 
(k) (k)
y1 = y0 + hfapp y0

ainsi que cette fonction loss:


(1) (1)
 1 (k) (k)
 
(k)
2
L W , w0 , W (2) , k = 2 y0 − y1 + hfapp y0
h

4.2 Ordre de convergence

Dans cette sous-partie, on considèrera une certaine méthode numérique multi-pas. On
considère également que le champ de vecteur f recherché peut être approché de manière
aussi précise que l’on souhaite par un champ fapp (en d’autres termes, on peut utiliser

Maxime BOUCHEREAU 47 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

autant de neurones que l’on veut dans le théorème d’universalité). on considère alors un
vecteur yˆ1 donnant la valeur prédite en h du problème de Cauchy:


ẏ = f (y)
y(0) = y0

à l’aide du champ de vecteurs fapp . Ce vecteur yˆ1 utilise une certaine méthode
numérique, en faisant intervenir fapp , ainsi que y tj , via notre réseau de neurones.
h

On introduit la fonction Loss:

1
L(fapp ) = |yˆ1 − y1 |2
h2

∗ ∗
Soit fapp ∈ C 0 (Rd , Rd ) tel que L(fapp ) = 0 (on suppose que l’on peut avoir une précision
aussi fine que possible, en pratique, il est compliqué d’avoir une Loss nulle). L’objectif
de cette sous-partie est d’établir un lien entre l’ordre de la méthode numérique utilisé et
∗
l’ordre d’apprentissage, c’est-à-dire étudier la différence entre les champs de vecteurs fapp
(champ de vecteurs appris) et f (champ de vecteurs réel).

4.2.1 Première formulation

Dans cette section, on souhaite déterminer l’ordre d’apprentissage directement à partir de

l’ordre de la méthode numérique, tant que cette dernière est de la forme:

N N
X X hj  
yˆ1 = αj y t j + h βj fapp y tj (4.1)
j=0
h
j=0
h h

Maxime BOUCHEREAU 48 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

Théorème (Ordre de convergence dans le cas général - Première formulation)

Si on prend yˆ1 de la forme (4.1):

N N
X X hj  
yˆ1 = αj y t j + h βj fapp y tj
j=0
h
j=0
h h

où la méthode est d’ordre p, c’est-à-dire que, si on pose:

N N
∼ X X hj  
y1 = αj y tj +h βj f y tj
j=0
h
j=0
h h

on a:

∼
y1 − y1 = O(hp+1 )
h→0

et si on suppose que:

N
X hj
βj 6= 0
j=0
h

Alors, l’apprentissage est d’ordre p, autrement dit, pour tout y0 ∈ Rd , et pour h assez
petit:

∗
fapp (y0 ) = f (y0 ) + O(hp )
h→0

Démonstration. On a:

N N
X hj
X  
yˆ1 = αj y t j +h βj fapp y tj
j=0
h
j=0
h h

hj
avec h
constant (la position des tj ne change pas par rapport à h). Soit la fonction
Loss:

1
L(fapp ) = |yˆ1 − y1 |2
h2
N N 2
1 X X hj  
= 2 αj y tj − y1 + h βj fapp y tj
h j=0 h
j=0
h h

On a ainsi:

N N
X hj ∗   X
h βj fapp y tj = y1 − αj y tj
j=0
h h
j=0
h

Maxime BOUCHEREAU 49 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

Donc:

N N N
X hj  ∗   X X hj  
h βj f − f y tj = y1 − αj y t j − h βj f y tj
j=0
h app h
j=0
h
j=0
h h

∼
= y1 − y1
Ainsi, on a:

N
X hj  ∗   1 ∼

βj fapp − f y tj = y1 − y1
j=0
h h h
∗

1 ∼

ϕh fapp − f (y0 ) = y1 − y1
h
où l’application (linéaire) ϕh est donnée par:

ϕh : C 0 (Rd , Rd ) −→ C 0 (Rd , Rd )
PN hj

g 7−→ j=0 βj h g φtj (y0 )
On a ainsi:

N
!
hj X
ϕh −→ Id βj
h→0
j=0
h
hj
qui est bien une application inversible puisque N
P
j=0 βj h =6 0. On a ainsi:

1
∗
fapp (y0 ) = f (y0 ) + ϕ−1 (φh − Fh ) (y0 )
h h
où Fh est le champ de vecteurs donné par:

N
! N
X X hj ∗

Fh = αj φtj Id + h βj fapp ◦ φtj
j=0 j=0
h
Or, on a:

1
ϕ−1
h (Fh − φh ) (y0 ) ∼ PN hj
(φh − Fh ) (y0 )
h→0
j=0 βj h
Donc on a:

∗ 1
fapp (y0 ) = f (y0 ) + PN hj
(φh − Fh ) (y0 )
h→0 h βj
j=0 h
1  ∼

= f (y0 ) + PN hj
y1 − y1
h→0 h βj
j=0 h
D’où:

∗
fapp (y0 ) = f (y0 ) + O(hp )
h→0


Maxime BOUCHEREAU 50 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

4.2.2 Nouvelle formulation

Dans cette section, on souhaite déterminer l’ordre d’apprentissage directement à partir de
l’ordre de la méthode numérique, tant que cette dernière est de la forme:

N N
!
X X hj
yˆ1 = αj y tj + hfapp βj y tj (4.2)
j=0
h
j=0
h h

Théorème (Ordre de convergence dans le cas général - Nouvelle formulation)

Si on prend yˆ1 de la forme (4.2):

N N
!
X X hj
yˆ1 = αj y tj + hfapp βj y tj
j=0
h
j=0
h h

où la méthode est d’ordre p, c’est-à-dire que, si on pose:

N N
!
∼ X X hj
y1 = αj y tj + hf βj y tj
j=0
h
j=0
h h

on a:

∼
y1 − y1 = O(hp+1 )
h→0

et si on suppose que:

N
X hj
βj 6= 0
j=0
h

Alors, l’apprentissage est d’ordre p, autrement dit, pour tout y0 ∈ Rd , et pour h assez
petit:

∗
fapp (y0 ) = f (y0 ) + O(hp )
h→0

Démonstration. On a:

N N
!
X X hj
yˆ1 = αj y tj + hfapp βj y tj
j=0
h
j=0
h h

hj
avec h
constant (la position des tj ne change pas par rapport à h). Soit la fonction
Loss:

Maxime BOUCHEREAU 51 Université Rennes 1

Master 2 Mathématiques Fondamentales Mars-Juin 2021

1
L(fapp ) = |yˆ1 − y1 |2
h2
N N
! 2
1 X X hj
= 2 αj y tj − y1 + hfapp βj y tj
h j=0 h
j=0
h h

On a ainsi:

N
! N
∗
X hj X
hfapp βj y tj = y1 − αj y tj
j=0
h h j=0
h

Donc:

N
! N N
!
 ∗  X hj X X hj ∼
h fapp − f βj y tj = y1 − αj y tj − hf βj y tj = y1 − y1
j=0
h h j=0
h
j=0
h h

Ainsi, on a:
N
!
 ∗  X hj 1 ∼

fapp − f βj y tj = y1 − y1 = O(hp )
j=0
h h h h→0

De plus, puisque l’on a:

N
X hj
βj 6= 0
j=0
h

PN hj
alors on peut changer j=0 βj h
y tj en y0 . En effet:
h

N N
!
X hj X hj
βj y tj = βj φtj (y0 )
j=0
h h j=0
h

et:

N N
!
X hj X hj
βj φtj −→ βj Id
j=0
h h→0
j=0
h

PN hj
qui est inversible, donc j=0 βj h φtj est inversible quand h est assez petit (c’est une
∼
perturbation d’un multiple non nul de l’identité). De plus, le terme y1 − y1 reste d’ordre
p + 1.

Maxime BOUCHEREAU 52 Université Rennes 1

Bibliographie

[1] S. Mallat, Notes des Cours de Stéphane Mallat Chaire ”Sciences des Données” du
Collège de France, Cours 2019: L’apprentissage par réseaux de neurones profonds.
Notes sur tout le cours 2019 par J-E. Campagne

[2] Site OpenClassrooms d’initiation au Deep Learning,

[Link]

[3] M. Raissi, P. Perdikaris, G.E. Karniadakis, Multistep neural networks for data-driven
discovery of nonlinear dynamical systems, arXiv preprint, arXiv:1801.01236.

[4] Steven L. Brunton, Joshua L. Proctor, and J. Nathan Kutz. “Discovering governing
equations from data by sparse identification of nonlinear dynamical systems”. In: Pro-
ceedings of the National Academy of Sciences 113.15 (2016), pp. 3932–3937. ISSN:
0027-8424

[5] Duong Nguyen et al. EM-like Learning Chaotic Dynamics from Noisy and Partial
Observations. 2019. arXiv: 1903.10335

[6] Ricky T. Q. Chen et al. “Neural Ordinary Differential Equations”. In: Advances in
Neural Information Processing Systems. Ed. by S. Bengio et al. Vol. 31. Curran As-
sociates, Inc., 2018.

[7] F. Regazzoni, L. Dedè, and A. Quarteroni. “Machine learning for fast and reliable so-
lution of time-dependent differential equations”. In: Journal of Computational Physics
397 (2019), p. 108852. ISSN: 0021-9991.

[8] J-P. Demailly, Analyse numérique et équations différentielles, coll. Grenoble Sciences,
1999.

[9] E. Hairer, S. P. Nørsett and G. Wanner. Solving Ordinary Differential Equations I,

2008, ISBN 978-3-540-56670-0.

[10] Pierre Gillot, Akka Zemmari1, Jenny Benois-Pineau et Yurii Nesterov. Algo-
rithmes de Descente de Gradient Stochastique avec le filtrage des paramètres pour
l’entraı̂nement des réseaux à convolution profonds.

53