Processus stochastiques et équations aux dérivés partielles

République Algérienne Démocratique et Populaire

Ministère de l'Enseignement Supérieur et de la Recherche Scientifique

Université M'hamed BOUGARA Boumerdès

Faculté des Sciences Boumerdès

DEPARTEMENT DE MATHEMATIQUES

MEMOIRE DE MAGISTER
SPECIALITE : MATHEMATIQUES
OPTION : MODELES STOCHASTIQUES

Thème
PROCESSUS STOCHASTIQUES ET EQUATIONS AUX
DERIVEES PARTIELLES

Présenté par

HANECHE MOHAMED Soutenu publiquement le :29/06/2009

Devant le jury composé de:

Présidente :	Mme H. GUERBYENNE	Maitre de conférences A	USTHB Alger
Promoteur :	Mr K. KHALDI	Maitre de conférences A	UMBB Boumerdès
Co-promoteur :	Mr A. ABASSOV	Maitre de conférences A	UMBB Boumerdès
Examinateur:	Mr H. OSMANOV	Professeur	UMBB Boumerdès
Examinateur:	Mr S.MAKDECHE	Maître de conférences B	UMBB Boumerdès

Remerciement

J'exprime mes profonds remerciements à mon promoteur le Docteur KHALDI Khaled pour l'aide compétente qu'il m'a apportée, pour sa patience et son encouragement à finir ce travail.

Je tiens à remercier Monsieur ABASSOV Assim mon co-promoteur et Monsieur MAKDECHE Said pour leurs conseils et remarques pertinentes.

Ensuite, je tiens à exprimer mes remerciements aux membres du jury, Mme Hafida GUERBYENNE, Présidente ainsi que Messieurs Hamid OSMANOV et Said MAKDACHE, membres, qui ont accepté d'évaluer mon travail.

Certaines personnes ne peuvent être oubliées, mes professeurs: Monsieur OSMANOV Hamid, Monsieur DOUMAZE, et tous les enseignants du département de mathématiques, sans oublier tous les personnels administratif et de soutien du département de mathématiques

Je remercie enfin toute l'équipe de ma promotion et en particulier mon ami et collègue d'études TAZEROUTI Moussa.

J'ai également une pensée pour tous mes proches, famille et amis, qui ont fait preuve de beaucoup de patience et m'ont toujours encouragé. J'espère maintenant être plus disponible et leur rendre ce que je leur dois. Enfin, je souhaite remercier en particulier mes parents pour leurs efforts, encouragements, pour le temps qu'ils ont consacré pour moi et pour tout.

Abstract

The aim of this work is to show the relation between the partial differential equations of the second order and the stochastic processes of diffusion, and present some results obtained recently on the partial derivative equations by probabilistic methods.

These results provide a probabilistic method that we allow to avoid the complication of numerical methods and written the solution as expectation of functional of diffusion process. This work is presented in five chapters:

The chapter I, present the basic mathematics tools, the Brownian motion and the stochastic process solution of stochastic differential equation (SDE) well-known with noun of diffusion process i.e. that their future is not depending of any other state excepting the present state, is key notion of this study. We introduce a new character of integral, is the stochastic integral says Itô integral that allow to give a sense to the differential of Brownian motion, the important notion upon rest the SDE theories.

In the chapter II, we give the generality of partial differential equations (PDE) of second order and explain the method of finite difference method this method is used in case where the resolution by the analytic method is impossible.

In the chapter III, we exhibit the profound relation existed between the notion of partial differential equations and stochastic differential equations through the certain theory (Feynman-Kac), the generalization of this theory given a probabilistic interpretation of PDEs.

The chapter IV is an application which we help to comprehend the notions of the president chapter, we start by simulating the trajectory of Brownian motion, and next, we simulate the diffusion process and resolve a PDE by the probabilistic method.

The chapter V, it is an application in finance, where we applied the Black and Scholes formula by different methods.

????

ÉíäÇËáÇ ÉÌÑÏáÇ äã ÉíÆÒÌáÇ ÊÇÞÊÔãáÇ Ê?ÏÇÚã äíÈ ÉÏæÌæãáÇ ÉÞ?ÚáÇ ÑÇåÙÅ ìáÅ áãÚáÇ ÇÐå á?Î äã Ïåä

ÇåáÍ æ ÉíÆÒÌáÇ ÊÇÞÊÔãáÇ Ê?ÏÇÚãÈ ÉÞáÚÊãáÇ ÌÆÇÊäáÇ ÖÚÈ ìáÅ ÞÑØÊáÇ ìáÅ ÉÇÖ?ÇÈ

,

ÑÔäáá ÉíÆÇæÔÚáÇ ÊÇíáãÚáÇæ

. ÉíáÇãÊÍÇ ÞÑØÈ

: áæÕ ÉÓãÎ

ìáÅ

áãÚáÇ ÇÐå ãÓÞäí

,

ÇåáÂÇÔã

ÈäÌÊæ áÍáá ÉíãÞÑáÇ ÞÑØáÇ áÇãÚÊÓÇ ÈäÌÊ äã ÇääßãÊ ÌÆÇÊäáÇ åÐå

ÉÂÑÍ ÇåäíÈ äã æ , ÉÓÇÑÏáÇ åÐå í ÇåÌÇÊÍä íÊáÇ ÉíÖÇíÑáÇ ÊÇíÓÇÓ?Ç ÒÇÑÈÅ ìáÅ áæ?Ç áÕáÇ á?Î äã ÞÑØÊä

ÉíÆÇæÔÚáÇ ÊÇíáãÚáÇ ÇÖí ìãÓí íÐáÇæ,ÉíÆÇæ ÔÚáÇ ÉíáÖÇÊáÇ ÉáÏÇÚãáá áÍ æå íÐáÇ íÆÇæÔÚáÇ ÑíÓáÇ æ äÇíäæÑÈáÇ

. äÇíäæÑÈáÇ

ÉÂÑÍ

ìáÅ ÉÈÓäáÇÈ áÖÇÊáá ìäÚã ÁÇØÚÅ äã Çääßãí íÐáÇ æÊíÇ áãÇßÊ ìáÅ ÞÑØÊä Çã ,?????

ÍÑÔä ãË

, ÉíäÇËáÇ ÉÌÑÏáÇ äã ÉíÆÒÌáÇ ÊÇÞÊÔãáÇ Ê?ÏÇÚã áæÍ ÊÇíãæãÚ ìáÅ íäÇËáÇ áÕáÇ á?Î äã ÞÑØÊä

.Ê?ÏÇÚãáÇ åÐåá

?????? ??

ÇäíØÚÊ íÊáÇæ

,ÉíåÊäãáÇ ÞæÑáÇ ÉÞíÑØ íå æ Ê?ÏÇÚãáÇ åÐå áÍáÉíãÞÑáÇ ÞÑØáÇ ìÏÍÅ

ÉíÆÇæÔÚáÇ

ÉíáÖÇÊáÇ

Ê?ÏÇÚãáÇæ ÉíÆÒÌáÇ ÊÇÞÊÔãáÇ Ê?ÏÇÚã äíÈ ÉÏæÌæãáÇ ÉÞíãÚáÇ ÉÞ?ÚáÇ äíÈä

ËáÇËáÇ áÕáÇ í

?????

ÁÇØÚÅ äã

Çääßãí ÉíÑÙäáÇ åÐå ãíãÚÊ , ÓíÇ ????? ãÓÇ áãÍÊ íÊáÇ ßáÊ ÉÕÇÎ ÉÕÈæ

ÊÇíÑÙä á?Î äã

. íÈíÑÞÊ áÍ íáÇÊáÇÈæ ÉíÆÒÌáÇ ÊÇÞÊÔãáÇ

Ê?ÏÇÚãá ÉíáÇãÊÍÇ

äÇíäæÑÈáÇ ÉÂÑÍ ÑÇÓã áíËãÊÈ ÞíÈØÊáÇ ÇÐå ÄÏÈä æ, ÉÞÈÇÓáÇ áæÕáÇ í ÁÇÌ Çãá ÞíÈØÊ äÚ ÉÑÇÈÚ æå ÚÈÇÑáÇ áÕáÇ ,áãÚáÇ ÇÐåá ãå?Ç ÁÒÌáÇ ìáÅ áÞÊää ãË,ÑáíÇ ??íÑØ ÚÇÈÊÅÈ ßáÐæ ÉíÆÇæÔÚáÇ ÉíáÖÇÊáÇ ÉáÏÇÚãáá ÑÇÓã áíËãÊ ìáÅ áÞÊää ãË Çääßãí íÊáÇ ÉíÆÒÌáÇ ÊÇÞÊÔãáÇ ÉáÏÇÚã ÑÇÊÎä ÞíÈØÊáÇ áíåÓÊá , ÉíÆÒÌáÇ ÊÇÞÊÔãáÇ Ê?ÏÇÚãá ÉíáÇãÊÍ?Ç ÉãÌÑÊáÇ æåæ äã ÇåíáÚ áÕÍÊãáÇ ÌÆÇÊäáÇ äÑÇÞäæ ÉíåÊäãáÇ ÞæÑáÇ ÉÞíÑØÈ ÇåáÍä ãË,( ÆÇßã ÚØÞ ÉáÏÇÚ?) ÉíáíáÍÊáÇ ÞÑØáÇÈ ÇåáÍ

. ÉáãÚÊÓãáÇ ÞÑØ áÇ á?Î

áæÕáÇ í ÇåíáÚ

áÕÍÊãáÇ ÌÆÇÊäáÇ áÇãÚÊÓÇÈ

ßáÐæ íáÇãáÇ ãÇÙäáÇ í ÞíÈØÊ æå ÑíÎ?Ç æ ÓãÇÎáÇ áÕáÇ

ÇääßãÊ íÊáÇ ÞÑØáÇ

?????

ãíÏÞÊ á?Î äã,ÒáæÔ ß?È

ÌÐæãä í ÊÇÑÇíÎáÇ ãííÞÊ ìáÚ ÉÓÇÑÏáÇ ÒÂÑä ËíÍ ÉÞÈÇÓáÇ ,

ÊÇÞÊÔãáÇ ÉáÏÇÚãÈ

ÑæÑãáÇ

ÉÞíÑØ ÑíÎ?Ç íæ áÇíãæäíÈ ÉÞíÑØæ æáÑÇ íÊäæã ÉÞíÑØ í,ÊÇÑÇíÎáÇ ÉãíÞ ÏÇÌíÅ äã

. ÒáæÔ ß?Èá ÉíÆÒÌáÇ

Sommaire

INTRODUCTION : . 1

CHAPITRE I: Mouvement Brwonien et calcul stochastique

0. Introduction : 3

1. Processus aléatoires : 4

2. Processus gaussien : 5

3. Mouvement brownien : 6

3.1. Construction du mouvement brownien : 6

3.2. Régularité des trajectoires du mouvement brownien : 8

3.3. Mouvement brownien standard : 9

3.4. Transformations du mouvement brownien standard : 10

3.5. Semi-groupe du mouvement brownien : 10

3.6. Mouvement brownien multidimensionnel : 13

4. Martingales et temps d'arrêt : 13

4.1 Martingales : 13

4.2. Comportement d'une martingale à l'infini : 14

4.3 Temps d'arrêt : 16

4.4. Tribu du passé d'un temps d'arrêt : 17

5. Intégrales stochastiques : 19

5.1. Variation quadratique : 19

5.2. Intégrales stochastiques : 21

5.2.1 Propriétés de l'intégrale stochastique : 23

5.2.2. Extension de l'intégrale stochastique : 24

5.3. Formule d'Itô : 25

6. Equations différentielles stochastiques : 28

6.1. Introduction : 28

6.2. Solutions d'une équation différentielle stochastique : 29

CHAPITRE II : Généralité sur les EDPs et la méthode des différences finies

0. Introduction : 34

1. Généralité sur les équations aux dérivées partielles : 35

1.1. Définitions : 35

1.2. Conditions initiales et conditions aux limites : 36

1.3. Classification des EDPs du second ordre : 37

1.4. EDPs du second d'ordre à plusieurs variables indépendantes: 39

1.5. Besoins en termes de conditions initiales et aux limites : 39

2. Méthode des différences finies pour les EDPs paraboliques : 40

2.1. Problèmes du premier ordre en temps: équation de la chaleur : 41

2.2.1. Schémas numériques de discrétisation par différences finies : 42

2.2.2. Erreur de troncature, consistance et ordre d'un schéma : 45

2.2.3. Stabilité des schémas numériques : 46

2.2.4. Convergence des schémas : 50

CHAPITRE III : Interprétation probabiliste des EDPs

0. Introduction : 52

1. Générateur infinitésimal : 53

2. Interprétation probabiliste des EDP: 54

2.1 Equation de la chaleur : 54

2.2. Formule de Feynman-Kac multidimensionnelle: 56

2.3.Problème de Cauchy pour des opérateurs généraux : 57

2.4. Interprétation des EDP de type elliptique : 58

a)Problème elliptique sans frontière : 59

b) Problème elliptique avec condition de Dirichlet : 59

2.5. Interprétation des EDP de type parabolique : 60

a. Problème parabolique sans frontière : 61

b. Problème parabolique avec condition de Dirichlet : 62

3. Discrétisation de problèmes : 63

3.1. Le schéma d'Euler : 63

3.2. La méthode de Monte Carlo : 63

CHAPITRE IV : Application

0. Introduction : 67

1. Discrétisation du mouvement brownien : 68

2. Discrétisation d'un processus de diffusion: 71

2.1. Schéma d'Euler- Maruyama : 71

2.2. Exemple de discrétisation d'un processus de diffusion: 72

3. Application à l'interprétation probabiliste des EDPs : 75

3.1. Problème N°=1 : 75

3.1.1. Illustration numérique du problème par la méthode déterministe(La méthode des différences finies) : 76

3.1.2. Illustration numérique du problème par la méthode probabiliste : 80

3.2. Problème N° 2 : 84

3.2.1. Illustration numérique du problème par la méthode des différences finies : 85

3.2.2. Illustration numérique du problème par la méthode probabiliste : 88

CHAPITRE V : Application en finance

0. Introduction : 91

1. Modèle du prix de l'actif: 91

2. Formule de Black et Scholes : 93

a. La méthode de Monte Carlo : 95

b. La méthode Binomiale : 96

c. EDP de Black et Scholes : 98

i. Mèthode des différences finies pour l'EDP de Black et Scholes : 99

ii. Méthode probabiliste pour l'EDP de Black et Scholes : 100

CONCLUSION : 103

ANNEXE... ..104

BIBLIOGRAPHIE... 111

INTRODUCTION:

Le but de ce travail est de montrer le lien entre les équations aux dérivées partielles du second ordre et les processus stochastiques de diffusion ainsi que de présenter quelques résultats obtenus récemment sur les équations aux dérivées partielles par des méthodes probabilistes.

Ces résultats fournissent une méthode probabiliste qui nous permet d'éviter les complications des méthodes numériques et écrire la solution comme l'espérance d'une fonctionnelle d'un processus de diffusion. Ce travail est présenté en cinq chapitres :

Le chapitre I, présente les outils mathématiques de base ([5], [9], [11], [12], [16], [17], [19], [26], [29], [31], [37], [38], [40], [42], [45], [48], [51]) entre autres le mouvement Brownien et les processus stochastiques solution d'équations différentielles stochastiques (EDS) connu sous le nom de processus de diffusion qui sont markoviens i.e. que leur état futur ne dépend que de leur état présent, notion clé de cette étude. Nous introduisons un nouveau type d'intégrale qui est l'intégrale stochastique dite d'Itô qui permet de donner un sens à la différentielle d'un mouvement brownien, notion importante sur laquelle repose la théorie des EDS.

Dans le chapitre II, on donne des généralités sur les équations aux dérivées partielles (EDP) de second ordre ; ensuite nous présentons une méthode

numérique pour la résolution des EDP: la méthode des différences finies ([3], [4], [10], [18], [23], [32], [54]), méthode utilisée dans les cas où la résolution par les méthodes analytique est impossible. Dans ce cas la solution trouvée est une solution approchée.

Dans le chapitre III, Nous exhibons le lien profond existant entre la notion des équations aux dérivées partielles (EDP) et celle des EDS à travers des théorèmes, en particulier, ceux dit de Feynman-Kac ([6], [9], [15], [16], [23] , [24], [35], [40], [49], [52]). La généralisation des ces théorèmes nous permet de donner une interprétation probabiliste des EDPs, et par conséquent, une solution approchée sous forme d'une espérance d'une fonctionnelle.

Le chapitre IV, est une application des notions introduites dans les chapitres I, II et III. Nous commençons par simuler les trajectoires du mouvement Brownien, ensuite, la simulation de la diffusion par la méthode d'Euler ([2], [4], [28]) et enfin à la plus important dans ce travail, qui est l'interprétation probabiliste des EDPs. Dans cette partie, on choisit une EDP (parabolique) que nous pouvons résoudre analytiquement, pour ensuite la résoudre par la méthode des différences finies, et on compare la solution donnée par les méthodes utilisées ([1], [13], [25], [26], [30], [34], [35], [50]). Dans le chapitre V, nous abordons l'application en finance. On s'intéresse particulièrement à l'évaluation d'options dans le modèle de Black et Scholes ([29], [41], [43], [44], [46], [53]). On illustre les différentes méthodes pour obtenir la valeur de l'option, c'est à dire, la méthode de Monte Carlo, Binomiale et le passage par l'EDP de Black et scholes (parabolique).

CHAPITRE I :

Mouvement brownien et

calcul stochastique

0. Introduction:

Le mouvement brownien est associé à l'analyse de mouvements dont l'évolution au cours du temps est si désordonnée qu'il semble difficile de la prévoir, même pour un temps très court, tel le mouvement d'une particule microscopique en suspension dans un liquide et soumise à l'agitation thermique.

Celui-ci joue un rôle central dans la théorie des processus aléatoires, d'une part parce que, dans de nombreux problèmes appliqués, le mouvement brownien sert à modéliser les erreurs ou les perturbations aléatoires, et d'autre part parce que le mouvement brownien ou les processus de diffusion qui en découlent permettent de construire des modèles simples sur lesquels des calculs peuvent être faits.

Le calcul stochastique, ou calcul d'Itô, du nom d'un des pionniers en ce domaine, est en fait un calcul d'intégrale par rapport au mouvement brownien. Ce dernier étant une fonction qui n'est pas à variation finie, cette notion d'intégrale n'est pas usuelle et sa définition en est probabiliste. Elle permet en particulier de définir la notion d'équation différentielle stochastique qui est une équation obtenue par la perturbation aléatoire d'une équation différentielle ordinaire. Les solutions de ces équations définissent de nouveaux processus, appelés processus de diffusion, et qui sont à la base du calcul probabiliste moderne. Ces processus sont souvent markoviens, au sens où leur comportement futur, conditionnellement au passé, ne dépend en fait que de l'état présent. Cette propriété, dite de Markov, est souvent vérifiée dans la réalité, en particulier, en physique, dans les réseaux de télécommunication, ou en mathématiques financières. Ainsi, les processus de diffusion sont précieux dans la modélisation de nombreux phénomènes aléatoires. On verra par ailleurs qu'il existe des liens importants entre leur loi et certaines équations aux dérivées partielles. Ces liens sont à la base de beaucoup de développements récents liant des résultats d'analyse et des résultats probabilistes.

1. Processus aléatoires : Définition 1:

On appelle processus aléatoire X = (Xt)t"une famille de variables aléatoires indexée par IR₊, toutes ces variables étant définies sur le même espace de probabilité (fl, Jl, P). On suppose ici que chaque Xt est à valeurs réelles.

On peut également voir le processus comme une variable aléatoire X définie sur (fl, Jl, P) et à valeurs dans l'ensemble des fonctions de t 1-0 Xt de IR₊ dans IR .

La théorie moderne des probabilités repose sur les résultats fondamentaux de Kolmogorov qui permettent, en particulier, de construire sur cet espace de fonctions une tribu qui rend l'application X mesurable (et permet donc de parler de variable aléatoire). Kolmogorov montre également que la loi de cette variable est caractérisée

par ses lois marginales de dimension finie, définies comme étant les lois des k-uplets (Xt₁, ...,Xt_k) , pour tous temps (t1, ..., tk).

Définition 2:

Les variables aléatoires Xt -- X_s , t > s ~ 0, sont appelées accroissements du processus (Xt).

Définition 3:

On dit qu'un Processus Xt est à accroissements indépendants et stationnaires si : (Indépendance): Xt -- X_s 1 F_s^x = o^-(X_s, 0 s t), V t > s ~ 0

(Stationnarité) : Xt -- X_s-- Xt_+s -- X0, V t > s ~ 0

Pour de tels processus, donner la loi de Xt_+s -- X0,Vt > 0, ainsi que celle de X0 suffit à caractériser entièrement le processus.

2. Processus gaussien : Définition 4:

Le processus X est un processus gaussien si chaque famille finie (Xt1, ...,Xt_k) est un vecteur aléatoire gaussien.

Notons par m(t) = E(X_t) la « fonction » moyenne du processus, et par C(s, t) = cov(Xt,X_s) la fonction de covariance, les deux fonctions m et C déterminent complètement la loi de toute famille finie(Xt₁, ...,Xt_k), et donc aussi la loi du processus X, On a de manière évidente C(s, t) = C(t, s), réciproquement, on a le théorème suivant.

Théorème 1:

Soit une fonction m de 11₊ à valeurs réelles, et une fonction C de 11₊ x 11 à valeurs réelles, on suppose de plus que C est de type positif, c'est-à-dire que pour tous réels positifs (t1, ..., tk):

IC(tk,t1)tkt1 > -0

k,1

Alors, il existe un processus gaussien (Xt)t unique (en loi), tel que :

Vt E IR₊, m(t) = E(X_t)

V t, s E IR₊, C(s, t) = E ((Xt -- m(t))(X_s -- m(s)))

Ce théorème est une conséquence du théorème fondamental de Kolmogorov sur la construction des processus aléatoires.

Supposons maintenant que X soit un processus gaussien à trajectoires continues, au sens où, pour tout w, l'application t 1-0 Xt(w) est continue. On peut alors définir

t

l'intégrale : Yt = f o X_t(w) ds

et obtenir ainsi un nouveau processus Y. Puisque la limite de toute suite de lois normales est encore une loi normale (cela se prouve en utilisant la fonction caractéristique), on vérifie, en utilisant une approximation de Riemann de l'intégrale ci-dessus, que Y est encore un processus gaussien. Par un calcul on peut montrer que sa fonction moyenne m' et sa fonction de covariance C' sont données par :

t

me(t) = f m(s)ds

o

^sC'(t, s) = f t

du f C(u, v)dv

o o

3. Mouvement brownien :

3.1. Construction du mouvement brownien :

Le mouvement brownien, ou processus de Wiener, joue un rôle fondamental dans de nombreux domaines. Il fut introduit pour la première fois en 1827 par le botaniste Robert Brown en observant des mouvements de grains de pollen dans un liquide. Ensuit par Louis Bachelier en 1900 pour des applications à la finance et a de nouveau,

à l'heure actuelle, un rôle important en mathématiques financières. Il fut redécouvert peu après Bachelier par Einstein, et est devenu depuis un des outils majeurs de la modélisation en physique. On le note : W = (Wt)t,0 ou B = (Bt)t,0 , Il peut être construit de différentes manières. Les définitions les plus usuelles du mouvement brownien sont les suivantes :

Définition 5 : un processus gaussien

Le mouvement brownien est un processus gaussien centré (E(Wt) = 0) pour tout t, de covariance C(s, t) = min (s, t) .

Définition 6 : un processus à accroissements indépendants stationnaires

Le mouvement brownien est un processus à accroissements indépendants et stationnaires. Plus précisément, pour tous s, t 0, la variable Wt_+s -- Wt est indépendante des variables (W_r: r t), Wo = 0, et de plus la loi de l'accroissement Wt_+s -- Wt est la loi normale N (0, s). Elle ne dépend donc que de s.

Commentons les équivalences entre ces deux définitions. Si W est un processus gaussien centré de covariance C(s, t) = min (s, t), ses accroissements ont de manière évidente la distribution souhaitée. Puisque, pour r < t: E((Wt+s -- Wt)W_r) = E(Wt+sWr) -- E((W_tW_r) = 0, la variable Wt_+s -- Wt est indépendante des variables (W_r: r t), par les propriétés des vecteurs gaussiens. Le processus W est donc un processus gaussien à accroissements indépendants et stationnaires. Si, réciproquement, le processus W satisfait les propriétés de la définition 2, il est presque immédiat de vérifier que c'est un processus gaussien possédant les caractéristiques voulues.

Il existe une troisième approche classique du mouvement brownien, qui consiste à l'obtenir comme limite de marches aléatoires. Ce résultat fondamental est donné par le théorème de Donsker et est, en fait, au niveau des processus, une version du théorème usuel de la limite centrale.

Théorème 2 :

Soit une suite (Xn)n>0 de variables aléatoires réelles indépendantes, identiquement distribuées, avec E(X_n) = 0 et E(X_n²) = 1 . Soit S_n = E i<i.<n Xi. avec So = 0 .

Les processus des sommes normalisées Yt n = vin S[nt] (où [nt] désigne la partie

entière de nt ) convergent en loi, en tant que processus, vers le mouvement brownien.

3.2. Régularité des trajectoires du mouvement brownien :

Ainsi défini, le mouvement brownien n'a pas spécialement de régularité, en tant que fonction du temps à valeurs réelles. On va voir que, dans certains cas, le modèle mathématique permet de construire un processus « équivalent » au processus initial, et qui possède une propriété de continuité presque sure.

Définition 7 :

Deux processus X et X', définis sur le même espace de probabilité, sont dits des modifications l'un de l'autre si pour tout t , X = X' presque surement.

Le critère suivant, appelé critère de Kolmogorov, montre que sous certaines conditions, un processus peut avoir une modification qui est presque sûrement continue, au sens où, pour presque tout w, la trajectoire t 1-0 Xt(w) est continue. Théorème 3:

Si un processus X est tel qu'il existe trois constantes strictement positives a, /3, C avec, pour tous t et h :

E(|Xt+h. -- X_t|^a) <Ch¹⁺¹⁶ (1.1)

alors X admet une modification presque surement à trajectoires continues.

Dans le cas du mouvement brownien, la variable aléatoire Wt+h. -- Wt est gaussienne centrée de variance h, donc :

E(|Wt+h. -- W_t|⁴) < 3h2

Le critère de Kolmogorov s'applique et on obtient :

Théorème 4:

Il existe un processus W = (Wt)t,0 presque surement à trajectoires continues, à accroissements indépendants et stationnaires, et tel que, pour tout t, la variable aléatoire Wt est gaussienne centrée et de variance t.

Cette propriété peut également se comprendre grâce au théorème 1 de Donsker. Choisissons des variables X_n, qui prennent seulement deux valeurs -1 et 1, avec probabilité 1/2 pour chacune d'entre elles. Dans ce cas, les processus discontinus Yⁿ ont des sauts de taille #177;1/Vn, qui convergent, uniformément sur tout intervalle de temps borné, vers 0. On peut donc s'attendre à ce que, à la limite, le processus W soit continu.

Dans la suite, nous considérons toujours des mouvements browniens ayant des trajectoires continues et qui sont nuls au temps 0. On appelle souvent un tel processus un mouvement brownien standard.

En affinant le critère de Kolmogorov, on peut également montrer que les trajectoires du mouvement brownien sont höldériennes d'ordre a, pour a < 1/2, c'est-à-dire que, presque sûrement |Wt+h -- Wt| C|h|^a, pour une constante C. En dehors de ces résultats de continuité, les propriétés de régularité du mouvement brownien sont très mauvaises. On montre par exemple que :

Proposition 1 :

Les trajectoires t 1-0 Wt sont presque sûrement nulle part dérivables, et qu'elles sont même à variation infinie sur chaque intervalle de longueur positive. Cela signifie donc qu'il n'y a pas de mesure it,,(dx) telle que Wt(w) = it,,([0, t[). Pour cette raison, on ne pourra pas définir l'intégrale f f(s)dWt(w) qui serait l'intégrale de Stieltjes de la fonction f par rapport à la mesure (non existante) it,,(ds).

3.3. Mouvement brownien standard :

Définition 8:

Un mouvement brownien standard (abrégé M.B.S.) est un processus aléatoire à temps continu(Wt, t E IR₊) tel que :

i) Wo = 0 p. s,

ii) (Wt) est à accroissements indépendants et stationnaires

iii) (Wt)--N(0, t), Vt > 0.

Cette définition permet de démontrer la propriété que (W_t) est à trajectoires continues.

Remarque :

De cette définition, il suit que pour t ~ s ~ 0, Wt -- W_s--Wt__s--N(0, t -- s) c'est à dire :

E(Wt -- W_s) = 0 et E((Wt -- W_s)²) = t -- s.

En appliquant la loi des grands nombres, on trouve encore que /qtt --> 0 p. s, lorsque t --> 00. De plus, on a _/q,/: --N(0,1), pour tout t > 0.

3.4. Transformations du mouvement brownien standard :

Soit W un mouvement brownien standard. Alors les cinq processus ci-dessous sont également des mouvements browniens standards :

1) W_t⁽¹⁾ = --Wt (propriété de symétrie du mouvement brownien).

2) soit T E IL₊ fixé : W_t⁽²⁾ = Wt₊T -- WT, t E IL₊ (stationnarité).

3) soit T E IL₊ fixé : W_t⁽³⁾ = WT -- WT_t, t E [0, T] (renversement du temps).

4) soit a > 0 fixé : Wt (4) = 1 ,/a Wat, t E IL₊ (loi d'échelle).

5) Wt^(5) = tWi t

, t > 0 et Wt^(5): = 0 (inversion du temps).

3.5. Semi-groupe du mouvement brownien :

Définissons tout d'abord la notion de filtration, qui modélise l'évolution de l'information au cours du temps.

Soit un processus X = (Xt)t_>0 sur un espace de probabilité (12, Jl, P). Nous notons Ft la tribu engendrée par les variables aléatoires X_s pour s t, qui est la plus petite
tribu rendant toutes ces variables mesurables. Nous avons Ft c Jl, et également F_s c F_t si s t.

Définition 9 :

1) La famille croissante F = (Ft)t_>0 est appelée la filtration engendrée par le processus X, et est aussi notée F^" = (Ft")t>0.

2) Plus généralement, une filtration est une famille F = (Ft)t,0 de tribus de A satisfaisant F_s c Ft si s t.

Propriété de Markov :

Considérons un mouvement brownien W sur (12,A,P), et la filtration F = (Ft)t,0 qu'il engendre. Puisqu'il est à accroissements indépendants, la variable Y := Wt_+s -- Wt est indépendante de la tribu Ft, on a pour chaque fonction f borélienne bornée sur IR :

E(f(Wt+s)\Ft) = E(f(W_t + Y)\Ft)

= f IR 1/ ,2 its e'2/ ^2s f (W_t + x)dx (1.2)

Cette formule montre que conditionnellement à Ft, la loi de Wt_+s ne dépend pas de tout le passé (c'est-à-dire de toutes les variables W. pour r t), mais seulement de la valeur « présente » Wt du processus. On dira que le mouvement brownien est un processus de Markov.

Définition 10:

Un processus (Xt)t est un processus de Markov si, étant donné la filtration (F_t')_t engendrée par le processus, celui-ci vérifie la propriété de Markov, à savoir que pour tous s, t 0 et pour toute fonction f borélienne bornée sur IR :

E(f(X_t+s)\F_t') = E(f(Xt+s)\Xt) (1.3)

Dans le cas du mouvement brownien, les variables W. pour r t d'une part, et les variables W. pour r t d'autre part, sont indépendantes, conditionnellement à la valeur de Wt. De plus, la loi de Wt_+s sachant Ft dépend bien sûr de s, mais pas de t. On dit que le mouvement brownien est un processus de Markov homogène en temps.

Définition 11 :

Un processus (Xt)t est un processus de Markov homogène en temps si pour tous s, t 0 et pour toute fonction f borélienne bornée sur IR, il existe une fonction borélienne h_s telle que :

E(f(W_t+s)\F_t') = h_s (Xt) (1.4)

Si X = (Xt)t_>0 est un processus de Markov, son evolution est decrite par une famille (Ps,t(x, dy))o.s.tde probabilites de transition de Ir vers Ir, (ou plus

generalement de (E,E) vers (E,E) si X prend ses valeurs dans l'espace mesurable (E,E) . La probabilite P_s,t(x,
·) represente la loi conditionnelle de Xt quand Xt = x .

Grâce à la propriete de Markov et à la relation E( E(X\3)\~) = E(X\~), on montre la propriete de semi-groupe. Pour 0 r s t :

Pr,t(x^, A) = f Ir Pr,s(x,, ^dy) Ps,t6 ⁷ , A) (1.5)

qui s'ecrit plus rapidement P_r,t = P_r,sP_s,t .Observons que P_s,s(x, dy) est la mesure de Dirac 8_x(dy) au point x.

Dans le cas ou X est un processus de Markov homogène, la probabilite de transition

P_s,t depend seulement de la difference t - s. En ecrivant Pt = P~,t = P_s,s+t pour tout t , on definit alors le semi-groupe de transition.

Définition 12 :

Soit X un processus de Markov homogène. On appelle semi-groupe de transition de X la famille (Pt)t>0 d'operateurs positifs lineaires :

P_t: 4) e L'(Ir^t) 1-0 Pt4) e L'(Ir^t)

Pt4)(x) = E(4)(~t)\X0 = x) = .1^" 4)(y)Pt(x, dy)

Ir^d

qui satisfait Pt1 = 1 et la propriete de semi-groupe :

Po = Id, P_t+s = P_t 0 P_s , Vs,t > 0 (1.6)

Dans le cas du mouvement brownien, qui est un processus de Markov homogène, le semi-groupe (Pt)t_>0 est donne par :

Pt(x, dy) = ¹

V 27rt exp (-- (y_ 2:)2) dy (1.7)

comme le montre immediatement la formule (1.2).

En utilisant les relations (1.2) et (1.7), on obtient facilement les proprietes :

E(W_t\F_s) = W_s ; E(W_t²\~_s) = W_s² + t -- s (1.8)

3.6. Mouvement brownien multidimensionnel :

Définition 13:

Un mouvement brownien d-dimensionnel est une collection W = (W^t)1,t,d de dmouvements browniens à valeurs réelles W^t = (Wt^t)t~o, qui sont indépendants entre eux.

Ce processus est encore un processus de Markov homogène (et même un processus à accroissements indépendants). Son semi-groupe vaut alors :

Pt(x, clY) =

1 exi 11Y -- x11^2Vcl (V2mt) P 2t rY

où x et Y appartiennent à R^d, 11
·11désigne la norme euclidienne sur R^d, et clY la mesure de Lebesgue sur R^d.

4. Martingales et temps d'arrêt :

4.1 Martingales :

Nous allons maintenant définir et étudier en détails une classe fondamentale de processus, qui vérifient la propriété de la relation (1.3).

Définition 14:

Soit(Sl, c 4, P) un espace de probabilité, muni d'une filtration1F = (Ft)t~o.

· Un processus à valeurs réelles M = (Mt)t,o est une 1F-martingale si :

a) il est adapté à la filtration(Ft)t,o, ce qui veut dire que pour tout t, Mt est Ftmesurable ;

b) chaque variable Mt est intégrable, et :

s < t M_s = E(M_t\F_~) (1.9)

· On dit que M est une 1F-surmartingale (resp. une 1F-sousmartingale) si l'égalité cidessus est remplacée par :

M_s > E(M_t\F_~) (resp. M_s < E(M_t\F_~) ) (1.10)

En particulier, l'espérance E(Mt) d'une martingale, (resp. d'une surmartingale,
sousmartingale), est une fonction constante du temps (resp. décroissante,
croissante). De manière évidente, une martingale est un processus qui est à la fois

une surmartingale et une sous-martingale et si M est une surmartingale, alors --M est une sous-martingale.

Remarque 1 :

Le mot « martingale » vient du monde des jeux, et a été donné initialement à la suite d'une mauvaise interprétation des « martingales » soi-disant trouvées par les joueurs pour gagner à coup sûr.

En fait, le sens mathématique de ce terme, appliqué à un jeu, est le suivant. Supposons que la variable Mt, égale au gain du joueur s'il arrête de jouer au temps t, soit une martingale. L'espérance « conditionnelle » du gain du joueur, si celui-ci arrête de jouer en un temps t strictement supérieur à s, et sachant le passé jusqu'au temps s, est alors égale au gain M_s qu'il obtiendrait s'il arrêtait de jouer au temps s. Ainsi, le jeu est équitable en espérance ou, en d'autres mots, le joueur maximisera la moyenne de ses gains en ne jouant pas du tout !

Exemples :

Soit W = (1/10t,0 un mouvement brownien, et F = (Ft)t,0 la filtration qu'il engendre.

Exemple 1 : le processus W est une F -martingale. En effet, Wt_+s--Wrest centré et indépendant de Ft , et donc E(Wt_+s--Wt\Ft) = 0.

Exemple 2 : le processus Mt = Wt² -- t est une F -martingale. En effet :

Mt+s^-- Mt = (Wt+s^--Wt)² + ²(Wt+s^--Wt) -- s

~ ~ ~~~ ~~~~ ~ ~~

Exemple 3 : soit ~ ~ ~. Le processus ~~ t) est une F- martingale.

2

En effet :

~ ~ ~~~ ~~~~~~~~~~~ ~ ~~

~~~~

~ /~~ 2 ~~

A² 0^-2

le fait que si U est une variable N(0, a²), alors E(exp(A.U)) = exp ₂ .

4.2. Comportement d'une martingale à l'infini :

Soit M une surmartingale ou une sous-martingale pour la filtration F. On peut montrer qu'il existe une « version » de M pour laquelle toutes les trajectoires sont

continues à droite et ont des limites à gauche. Nous ne considérerons dorénavant que de telles versions. En fait, presque toutes les surmartingales et sousmartingales que nous rencontrerons dans ce travail auront des trajectoires continues.

Nous allons énoncer une propriété remarquable des martingales.

Proposition :

Si M est une martingale satisfaisant :

suptER+ E(|Mt|) < +00 (1.11)

alors t 1-> Mt(w) admet, pour presque tout w, une limite (qui peut être infinie) quand t 1-> +00 , et cette limite est notéet 1-> M,(w).

Observons que cette propriété ne s'applique pas à l'exemple 1 ci-dessus, pour

lequel E(|Wt|) = .12t h r, ni pour l'exemple 2, pour lequel E(|Mt|) = ct pour une
constante c > 0, tandis qu'il s'applique pour l'exemple 3, puisque E(IM_t^~~) =

E(M_t^~) = 1.

Une importante question est alors de savoir quand cette égalité (1.9) reste vraie pour s = 00 . La réponse (non triviale) est la suivante.

Théorème 5:

Mt = E(M_c\F_t) <=> Mt - M_oo dans L¹ (1.12)

Ces propriétés sont, comme on l'a vu grâce à l'exemple 3, strictement plus fortes que la relation (1.11). Elles sont en fait équivalentes au fait que M est une martingale uniformément intégrable, au sens suivant.

Définition 15:

La famille des variables aléatoires (Mt)t>0 est uniformément intégrable si :

lim sup t E(|Mt|¹|m~|>a) = 0 (1.13)

a-+oo t_ER₊

Inégalité de Doob :

L'inégalité suivante est fondamentale et spécifique aux martingales. Elle donne un contrôle en norme L^P du supremum d'une martingale, sur un intervalle de temps fini, en fonction de la valeur absolue de sa valeur terminale.

Théorème 6:

Soit M une martingale nulle en 0. On note. Mt^~ = sup .s^..t

IMsI.

Alors, pour tout entier p > 1,

E(( _MnP) ( )

^{P 13}

1 - p E(PMt|P)

4.3 Temps d'arrêt :

Cette notion joue un rôle très important en théorie des probabilités.

Définition 16:

Soit F une filtration. Une application T : 12 1-> [0, oo] est un temps d'arrêt si (T t} E Ftpour tout t 0.

Un temps d'arrêt est donc un temps aléatoire, tel que sur chaque ensemble

(w: T(w) t}, l'application w 1-> T(w) dépend seulement de ce qui s'est passé avant le temps t.

Un joueur honnête, qui ne peut pas anticiper sur les événements futurs, peut décider d'arrêter le jeu au temps aléatoire T uniquement si T est un temps d'arrêt. Un exemple trivial de temps d'arrêt est donné par T(w) = t pour tout w.

En dehors des temps constants, l'exemple fondamental de temps d'arrêt est le temps d'atteinte d'un ensemble borélien A par un processus x à trajectoires continues à droite et adapté à la filtration F. On définit plus précisément :

T = inf(t 0 ; x_t E A)

(avec la convention que l'infimum de l'ensemble vide vaut+oo).

Le fait que le temps d'atteinte de l'ensemble A est un temps d'arrêt est difficile à prouver pour un borélien A arbitraire, mais facile à prouver quand A est un ensemble ouvert. En effet, dans ce cas :

(T t} = (x_t E A} U _(urEQn[0,t] (x_r E A}) ce qui permet de conclure.

Considérons X et A comme ci-dessus et S un temps d'arrêt. Un nouvel exemple de temps d'arrêt est T = inf(t: t > S ; Xt E A). Ce temps est le temps d'atteinte de A par X, après le temps S.

En revanche, le dernier temps avant un temps fixé s où un processus adapté X visite un ensemble borélien A, défini par T = sup(t: t > s ; Xt E A) (où le supremum de l'ensemble vide est égal à 0) n'est pas un temps d'arrêt. En effet, T est mesurable, mais la valeur de T(w) dépend de tout ce qui s'est passé pour le processus X entre les temps 0 et s.

4.4. Tribu du passé d'un temps d'arrêt :

Définition 17:

Soit T un temps d'arrêt. On appelle tribu du passé T, la tribu notée FT égale à l'ensemble de tous les événements A E A tels que A n (T t) E Ft pour tout t.

Il est facile de voir que FT est une tribu, et cette notation est cohérente avec la notation Ft car si T est identiquement égal à t, il est clair que FT = Ft. La tribu FT est appelée tribu du passé de T, car elle peut être interprétée comme suit : A E FT si, quand on s'arrête au temps T, on sait si A est réalisé ou non.

Voici quelques propriétés simples. Ci-dessous, S, T, (R)n désignent des temps d'arrêt, et la filtration F est supposée vérifier que Ft =n_s>t F_s.

(U < t) E Ft , Vt = U est un temps d'arrêt.

S < T =FscFT

(S < T) E Fs n FT, (S T) E Fs n FT

R = inf_n(R_n) est un temps d'arrêt et FR =n_n FR_n

R = sup_n(R_n) est un temps d'arrêt.

Théorème d'arrêt :

Soit M une martingale. La propriété (1.9) peut facilement être étendue aux temps d'arrêt bornés.

Théorème 7:

Si S et T sont deux temps d'arrêt et si a E 11., alors :

Ms = E(MT\Fs) sur l'ensemble (S T a) (1.14)

En particulier, si T est un temps d'arrêt qui est borné, on a :

E(MT) = E(Mo) (1.15)

Quand Mt désigne de nouveau le gain d'un joueur au temps t , la propriété (1.15) peut être interprété comme suit. Quelle que soit la stratégie non anticipante que le joueur choisit pour arrêter le jeu, et s'il doit finir de jouer avant un temps déterministe donné (aussi grand que soit ce temps), alors la valeur espérée de son gain est constante et égale à son capital initial.

Observons que la relation (1.14) est, en général, fausse sur l'ensemble(S T}, et de même (1.15) est fausse si T n'est pas borné.

Par exemple si M = W est un mouvement brownien et si T = inf(t: Mt = 1), alors E(M0) = 0 < E(MT) = 1 Dans ce cas, le temps aléatoire T est presque sûrement fini, mais n'est pas borné et a même une espérance infinie.

En revanche, dans le cas d'une martingale uniformément intégrable, tout se passe bien.

Théorème 8:

Si M est une martingale uniformément intégrable, alors la relation (1.14) est satisfaite sur (S T}, et (1.15) a lieu pour tout temps d'arrêt T.

Remarquons que l'on peut considérer des temps d'arrêt qui peuvent prendre la valeur infinie, pourvu que l'on pose Mt = M₀ sur l'ensemble( T = 0}.

Une autre propriété importante, en lien avec la propriété (1.14), concerne le processus arrêté au temps T, défini par M_t^T = Mini,i(t,T).

On peut facilement déduire de (1.14) que si M est une martingale et T un temps d'arrêt, alors :

M^Test une martingale.

Définition 18:

Si M est une martingale et T un temps d'arrêt par rapport à la même filtration, on appelle martingale arrêtée au temps T la martingale MT.

5. Intégrales stochastiques :

On a vu que le mouvement brownien est à variation infinie, et que l'on ne peut donc pas definir une integrale de Stieltjes qui lui serait associee. On va toutefois voir qu'il est possible de definir une integrale d'une autre nature, definie dans un sens quadratique.

5.1. Variation quadratique :

Soit X un processus à valeurs reelles. On appelle « variation quadratique approximee » de X au niveau n le processus suivant :

V(X, n)_t = E__i (Xi/n -- X(i_1)/n)² (1.16)

où [nt]est la partie entière de nt.

Si X est un processus continu et à variation finie, au sens où X est la difference de deux processus Y et Z dont les trajectoires sont croissantes et de valeurs finies, il est facile de verifier que V(X,n)t tend vers 0 quand n tend vers l'infini.

Supposons maintenant que X = W soit un mouvement brownien. Chaque accroissement W/n -- W(i_1)/n a une loi normale .N(0, 1 /n ), et donc [nt]

n

E(V(W,n)_t) = . On en deduit que E(V(W,n)_t) tend vers t quand n tend vers

l'infini. Il est alors naturel et facile de prouver que V(W,n)t converge vers t dans L¹. Au vu de ce resultat, nous disons que la variation quadratique du mouvement brownien est (W, W)t = t. Cela montre en particulier que les trajectoires du mouvement brownien ont une variation infinie sur les intervalles finis.

On considère maintenant plus generalement une martingale continue M. Définition 19 :

La martingale continue M est dite de carre integrable si pour chaque t , Mt E L². On peut alors, grâce au resultat suivant, definir la variation quadratique de M. Proposition 20:

Soit M une martingale continue de carre integrable. La variable aleatoire

V(M,n)_t converge dans L¹, quand n tend vers l'infini, vers une variable

(M,M)t notée ( (M,M)t)tn. Le processus est croissant, continu, et est appelé la variation quadratique de M. De plus, il vérifie que :

M_t² -- (M,M)t est une martingale (1.17)

En fait, (M,M)t est l'unique processus croissant, continu et adapté, nul en 0, tel que l'on ait (1.17).

Remarque 2 :

On a donc :

(W, W)t = t ; Wt² -- t est une martingale.

En fait ceci est caractéristique du mouvement brownien d'après le théorème de Paul Lévy suivant.

Théorème 9 :

Toute martingale M continue de carré intégrable, et telle que M_t² -- t soit encore une martingale, est un mouvement brownien.

L'application M 1-0 (M,M)t se comporte comme une forme quadratique, ce qui explique la notation. Si a, b sont des réels, par la caractérisation ci-dessus, on obtient facilement que (aM, bM) = ab(M, M). On peut lui associer une sorte de « produit scalaire », en posant, pour chaque paire M,N de martingales de carré intégrable et continues :

(M, N) = 4i ((M + N, M + N) -- (M -- N, M -- N)) (1.18)

formule qui pourrait être comparée à la formule :

(x, 17) =

4 (11x + 1711² -- 11x -- 1711²)

1

pour le produit scalaire usuel et la norme euclidienne dans IRY.

Observons que MN -- (M, N) est une martingale, et (M, N) est l'unique processus Y continu, adapté, nul en 0, et à variation finie, tel que MN -- Y soit une martingale.

En particulier, si W = (W^t)i,t,d est un mouvement brownien d-dimensionnel, on a :

(W^t, Wl)_t = S^tlt (1.19)

où Sⁱⁱ = 1 si i = j et 8' = 0 sinon. C'est evident si i = j, et pour i # j le fait que (Wⁱ, W¹) = 0 vient du fait que le produit WⁱWⁱ est une martingale, comme produit de deux martingales independantes.

En fait, cette propriete caracterise le mouvement brownien. En effet, on peut enoncer le théorème de Paul Lévy suivant.

Théorème 10:

Toute martingale M d-dimensionnelle, continue, de carre integrable, telle que :

(M^t, Mⁱ)t = a^tit,

est un mouvement brownien d-dimensionnel.

5.2. Intégrales stochastiques :

Nous souhaitons maintenant donner un sens à l'integrale f_o^t HAM. quand W est un

mouvement brownien et H = (Ht)t,0 un processus, dont nous preciserons ces proprietes dans ce paragraphe.

Nous nous limiterons essentiellement à des integrants H qui sont continus en la

variable t. L'idee est alors d'obtenir C HAW_s comme limite de sommes de Riemann: I(H,n)_t = ? ~~~~~~ /~

~~~~ ~~~ ~ ~ ~~~~~~

~ ~ ~1.20~

~~~

La variable W/n -- 1/17(i_1)/_n a une taille d'ordre 1/vn, car elle est centree et de variance 1/n, et donc H(i_i) /n (Wi/n -- ¹4i_1/n) est egalement d'ordre 1/vn . La

taille de I(H,n)t devrait donc être vn. Mais, dans ce cas, les variables I(H,n)t ne pourront en general pas converger ; cela est coherent avec le fait que les trajectoires de W sont à variation infinie.

Pourtant, une sorte de « miracle » a lieu, quand on suppose de plus que le processus H est adapte à la filtration du mouvement brownien. Par simplicite, on supposera aussi que ce processus H est borne par une constante, mais cette hypothèse peut être allegee.

Dans ce cas, la variable

Y(n, i) = ¹(i-1) /n (Wi/n -- W(i-1)/n)

satisfait E(Y(n, i) /F(i-1) /n ) = 0

et E(Y(n, i)²/F(i-1) /n ) = H(i-1) /n , n c2 /n , de telle sorte que

[nt]

I(H,n)_t = 1 Y(n, i)

i=1

est centrée, de variance :

[nt] 1 E(Y(n, i)^2) i=1

< c^2t

Il n'est alors pas totalement déraisonnable de penser que la suite I(H,n)t converge et, effectivement, on peut montrer le théorème suivant.

Théorème 11:

Soit H un processus borné, continu et adapté à la filtration du mouvement brownien. Alors la suite :

[nt]

I(H,n)t = IH(i-i) /n (Wi/n -- W4i-1)/n)

i=1

converge dans L² , quand n tend vers l'infini, vers une limite, notée :

f

^t HAW, o et appelée l'intégrale stochastique de H par rapport à W sur l'intervalle [0, t].

Remarque 3 :

1. Si H n'est pas adapté, les sommes de Riemann ne convergent pas en général.

2. On observe que la relation (1.20) est une forme particulière de somme de Riemann.

Si on remplace

H(i-i) /n (Wi/n -- W(i^-1)/n)

par W

Ht(n,i)( i/n -- W(i-1)/n) ,

avec (i -- 1) /n t(n, i) i/n comme il est possible de le faire pour les

approximations par les sommes de Riemann pour les intégrales usuelles, la suite

associée I(H, n)t ne converge pas nécessairement et, si elle converge, la limite peut être différente def ^ot HsdWs.

3. La terminologie « intégrale stochastique » permet d'insister sur le fait que cette intégrale n'est pas une intégrale de Stieltjes usuelle, prise séparément pour chaque valeur de co, mais une limite dans L² .

4. Il est possible de définir l'intégrale f ot HsdWs pour des intégrants H bornés qui ne

sont pas continus en temps. Mais, dans ce cas, l'adaptation de H à la filtration n'est pas suffisante et il faut supposer plus de mesurabilité. (La propriété requise s'appelle la mesurabilité progressive).

5. Il est aussi possible de définir l'intégrale f _o^tH_sdW_s. pour des intégrants continus non bornés H, ou même pour des intégrants non bornés satisfaisants la mesurabilité progressive mais, dans ce cas, on doit supposer que l'intégrale f _o^tH_s²ds est finie pour tout t.

5.2.1 Propriétés de l'intégrale stochastique :

Les propriétés d'adaptation sont ici implicitement définies par rapport à la filtration du mouvement brownien W.

Théorème 12:

Soit H et K des processus bornés, continus et adaptés.

1. Pour tous réels a,13,

f o

t t t

(aH_s + 13K_s)dW_s = a f _HsdWs + 13 f KsdWs

o o

2. Le processus M défini par M_t = f₍^t₎ HsdWs est une martingale continue de carré intégrable, nulle en 0.

3. Si, de plus, Nt = f: K_sdW_s, on a :

(M,N)_t = f₍^t₎ H_sK_sds (1.21)

En particulier, on a l'isométrie fondamentale, donnée par la formule suivante :

²

E ((f: H_sdW_s)) = E(f o t H_s²ds ) (1.22)

5.2.2. Extension de l'intégrale stochastique :

Si l'on souhaite definir l'integrale stochastique d'un processus continu adapte H par

rapport à une martingale M de la forme Mt = f₍^t₎ Ks.dWs., il suffit de poser :

f : Hs.dMs. = f _o^t(Hs. Ks.)dWs. (1.23)

Plus generalement, on peut definir l'integrale stochastique d'un processus continu

adapte H par rapport à une martingale continue M, dès que (f (t) Hs.²d ^(M, N)s.) < 00
Cette integrale, encore notee par f :Hs.dMs., est obtenue comme limite des sommes
de Riemann (1.20) où M remplace W. Si Mt = f₍^t₎ Ks.dWs., on retrouve le processus

defini par la relation (1.23).

Nous avons les proprietes suivantes :

Proposition 2 :

Soit M, M' des martingales continues, H et K des processus continus, adaptes à F^m,

et H' continu et adapte à F^m' .

t t

1. H 1-> f o HdM_s.et M 1-> f o Hs.dMs. sont lineaires.

2. Si Nt = f _:Hs.dMs., alorsf o t Ks.dNs. = f : Hs.k^.dMs. .

3. Si N_t = f : Hs.dMs. et N' = f :H;dM;, on a :

(N, N')_t = f ^ot H_s.H;d (M, M')s. (1.24)

En particulier, la formule ci-dessous generalise (1.22) :

2

E ((f : Hs.dMs.) ) = E(f:H_s.²d (M, M)s.) (1.25)

4. Si H est borne, le processus N defini par Nt = f ot Hs. dMs. est une martingale

continue nulle en 0. Si H n'est pas borne, le manque possible d'integrabilite en fait ce qu'on appelle une « martingale locale », mais nous ne developperons pas ce point dans cette partie.

5.3. Formule d'Itô :

Quand t 1-0 x(t) est une fonction réelle, continue et à variation finie, la formule d'intégration par parties implique que pour toute fonction f continûment dérivable, on a :

f(x(t)) = f(x(0)) + f ot f' (x(s))dx(s) (1.26)

et il existe également une version multidimensionnelle de cette formule.

Celle-ci devient fausse quand la fonction x est remplacée par un mouvement brownien W, ou par une martingale continue M.

En effet, si la relation (1.26) était vraie, en prenant f(x) = x² , on obtiendrait

Wt² = 2 f : Ws.dWs. et, puisque le processus défini par l'intégrale stochastique est

encore une martingale, on pourrait en déduire que Wt² est une martingale. Or, on a vu au paragraphe précédent que la variation quadratique de W au point t est égale à t, et donc le processus W_t² -- t est aussi une martingale. Par différence, le « processus » t devrait alors être une martingale, ce qui est évidemment faux.

Ainsi la relation (1.26) est fausse pour le mouvement brownien. Pour obtenir une formule juste, on doit ajouter un terme de plus (appelé quelque fois terme de correction d'Itô) et supposer plus de régularité sur f.

Théorème 13:

Considérons une martingale continue X, de variation quadratique ( X, X). Soit f est une fonction de classe C², alors on a :

f(X_t) = f(Xo) + f o t f^'(Xs.)dXs. + ₂¹ f ot f" (Xs.) ^d(X, X)s. (1.27)

où f"est la dérivée seconde de f.

Cette formule est connue sous le nom de « formule d'Itô » et est extrêmement utile.

La formule d'Itô appliquée au mouvement brownien X = W et à donne

Wt² = t + 2 f : Ws.dWs. et on retrouve bien la propriété que Wt² -- t est une martingale.

Cette formule n'est pas suffisante pour les applications. On a souvent besoin de considérer des processus de la forme plus générale suivante.

Définition 21 :

Un processus X de la forme X = M + A , où M est une martingale continue et A un processus adapté, continu et localement à variation finie (à variation finie sur tout compact de IR₊ ), est appelé une semi-martingale continue.

En regroupant l'équation (1.26) pour A et (1.27) pour M, on obtient pour f de classe C^Z une version de la formule d'Itô pour les semi-martingales.

Théorème 14:

Soit X = M + A une semi-martingale, et f une fonction de classe C². On a:

~ ~~

~~~~~ ~ ~~~~~ ~ ~ ~~~~~~ ~~~ ~ ~ ~~~~~~ ~~~ ~ 2 1 ~ ~~~~~~~ ~~~, ~~~ ~1.28~

~ ~ ~

Donnons égalementt la forme « multidimensionnelle » de la formule d'Itô.. Théorème15::

Soit X = (X^t)i,t,d des semi-martingales de décompositionss X^L = M^t + A^t , oùu les M^t sont des martingales continues et les A^t des processus adaptés,, continus et localement àa variation finie. Soit f une fonction de classe C² sur Ir^t, dont les~~ dérivées partielles des premier et second ordres sont notées ~~ ~ et ~~~ b respectivement.

On a:_~

~~~~~ ~ ~~~~~ ~ ? ~~ ~~

~ ~

~~~~~ ~~~~ ~ ~ ~~ ~~~~~ ~~~~~

~~~ ~ ~ 0~1.29~

~ 2 1 ? ~ ~~~

~ ~ ~~~~~~ ~~~~, ~~~~

~,~~~ ~ s

Remarque 4 :

Par abus de notation, on écrit souvent cette équationn sous une forme différentielle, bien que ce soit purement formel :

^d

<Ant = fx^(t)(dA^lt + dmo i,j=1

1 ~2 ft7(Xt)d(1',, Mi)t (1.30)

Observons que ces formules impliquent en particulier que le processusf(X),, image
d'unee semi-martingale X par une fonction f de classe C², est encore une

semi-martingale, et que sa décomposition en somme d'une martingale et d'un processus localement à variation finie, les deux continus et adaptés, est donnée directement par (1.29). L'intégrale stochastique par rapport à M est la partie martingale, et le processus à variation finie est formé de la somme des autres termes. Exemple 1 :

Considérons deux semi-martingales continues X et Y avec :

X_t = M_t + At; Yt = Nt + Ct

M et N étant les parties martingales. En prenant f(x, y) = xy, on obtient immédiatement :

t t

X_tY_t = X0Y0 + 1 Y_s dX_s + 1 _XsdY~ + (M,N)_t

o o

Exemple 2 :

et ~~~, ~~ ~ ~~~ ~~~ ~ ~~

Appliquons (1.28) pour ~ ~ 2, et ~~ ~ ~ ~~, ~~ ~ ~ ~ ~ y). On

~ ~ ~ ~~

a ~~ ~ ~ ~~, ~~ ~ ~ 0, ~~ ~ ~ 0, ~~ ~ ~ ~, et ~~ ~ ~ ~~ , ~~~

~~ ~ ~~~ , ~~ ~ f On

obtient ainsi :

~

~~ ~ ~ ~~~ ~~~~ ~ ~~ 2 ~~ ~ 1 ~ ~ ~~~ ~~~~ ~ ~~ 2 s) dW_s

~

~

et on retrouve le fait que ~~ est une martingale, comme on l'avait déjà. Exemple 3:

Plus généralement, si H est un processus borné, continu, adapté, on définit :

~ ~

~~ ~ ~~~ ~~ ~~~~~ ~ 2 1 f H_sds)

~ ~

et on montre par la formule d'Itô que Z est solution de l'équation :

Z_t = 1 + f_o^t Z_sH_sdW_s; t 0

Cette équation est un exemple d'équation différentielle stochastique, équation que l'on étudiera au paragraphe des EDS.

Exemple 4 :

Un deuxième exemple d'équation différentielle stochastique, appelé le modèle non homogène de Black et Scholes, utilisé en finance, est donné par l'équation formelle:

dSt = S_t11(t)dt + Sta(t)dWt , SO = 1 (1.31)

Les fonctions 11 et a _?sont à valeurs réelles, avec a(t) > 0, W est un mouvement brownien, et S est le processus « inconnu ». C'est une équation linéaire, qui peut être résolue explicitement. La solution est de la forme St = exp(Mt + At), où M est une martingale continue et A un processus adapté, continu et localement à variation finie. Pour trouver M et A, on applique la formule d'Itô (1.30) à la fonction f(x, y) = e^x+Y , ce qui donne :

dSt = StdMt + StdAt +

1

2 Std(M, M)t

En identifiant avec (1.31), on obtient dMt = a (t)dWt, d'où Mt = f : a (t)dWt. On en déduit que (M, M)_t = f_o^t a(s)²ds, d'où dA_t + ₂¹ a(t)² dt = 11(t)dt et

At = f ot (11(s)-- ₂¹ a(s) ² ) ds Finalement une solution de (1.31) est dormée par :

t

1

S_t = exp (1 a (s)dW_s. + f (11(s)-- t

₂ a(s)²) ds)

~ ~

6. Equations différentielles stochastiques :

6.1. Introduction :

De manière informelle, on appelle équation différentielle stochastique une équation différentielle ordinaire perturbée par un terme stochastique. Plus précisément, c'est une équation du type suivant :

dXt = b(t, Xt)dt + a(t, Xt)dWt, X0 = x0 , (1.32)

Dans cette équation, dWt est la « différentielle » d'un mouvement brownien W, et b, a sont les coefficients de l'équation (ce sont des fonctions de 11₊ x 11 dans 11), et x0 E 11 est la valeur initiale. Tous ces termes sont donnés. La notation (1.32) est la plus usuelle, mais elle n'a pas de sens (car la différentielle dWt n'a pas de sens). Le sens mathématique sera donné sous forme d'une équation intégrale.

Définition 22 :

Rechercher une « solution » de l'équation (1.32) consistera à rechercher un processus X = (Xt)t_>0 satisfaisant l'équation intégrale :

t t

X_t = f ₀ b(s,X_s)ds + f ₀ a(s,X_s)dW_s (1.33)où la seconde intégrale est une intégrale stochastique.

L'équation (1.32) ou l'équation (1.33) étaient jusqu'à présent unidimensionnelles. On peut également définir une équation d-dimensionnelle de la manière suivante. Le processus inconnu X = (X^t)i,t,d est une famille de processus à valeurs réelles X^t = (Xt^t)t,0, la condition initiale x0 = (x^t0)i,t,d appartient à R^d , le mouvement brownien W = (W^t)i<t<ti est q-dimensionnel, et les coefficients ont les dimensions appropriées, soit b = (b^t)i,t,d et a = (a^t,i)i,t,d,i,i<ti, où les coefficients b^t et at,i sont des fonctions de R₊ x R^d dans R. On écrit encore l'équation sous les formes (1.32) ou (1.33), mais cela signifie maintenant que l'on a :

Xt^t = x^t0 + f :b^t(s,X_s)ds + Eci=i f :a^t,i (s,X_s)dW_s^t , i = 1, ...., d (1.34)

Définition 23 :

Quand les coefficients b et a ne dépendent pas du temps et sont seulement des fonctions définies sur R^d , on dit que l'équation est homogène.

Le coefficient b est appelé le coefficient de dérive, tandis que a_? est le coefficient de diffusion. Un processus qui résout l'équation (1.32), ou de manière équivalente (1.34), est appelé processus de diffusion ou, plus simplement, une diffusion. On notera par EDS le terme « équation différentielle stochastique ».

6.2. Solutions d'une équation différentielle stochastique :

Notre première tâche est de donner un sens précis au mot « solution » de l'équation différentielle stochastique (1.33).

Définition 24 :

Soit un mouvement brownien W = (W^t)i<t<ti q-dimensionnel sur un espace de probabilité (12,A,P) , et F = (Ft)t,0 la filtration qu'il engendre. Les coefficients b et a, de même que la condition initiale x0 E R^d , sont donnés. Nous appelons solution forte de (1.34) tout processus d-dimensionnel X = (X^t)i,t,d, continu, adapté à la filtration, et tel que (1.34) a lieu.

Observons que l'adaptation de X est nécessaire pour que l'intégrale stochastique dans (1.34) ait un sens.

Remarque 5 :

Il existe une autre notion de solution, qui en un sens est plus naturelle, bien que moins intuitive. Plus précisément, nous appelons solution faible de l'EDS (1.34) la loi d'un processus X qui résout l'équation (dans un certain espace de probabilité muni d'un mouvement brownien). On peut montrer que si l'on a unicité de la solution forte, alors, on a unicité de la solution faible.

Le premier problème à résoudre concerne l'existence et l'unicité d'une solution. Considérons tout d'abord le cas d'une équation différentielle ordinaire, de la forme:

dXt = b(t, Xt)dt , X0 = x0 (1.35)
où le coefficient b et la condition initiale x0 sont donnés.

Bien sûr, un tel problème est purement déterministe. Un résultat classique énonce que (1.35) admet une et une seule solution dès que :

11b(t,x) -- b(t,Y)11 C11x -- Y11 , 11b(t, 0)11 C (1.36)

où C est une constante, et 11
·11est la norme euclidienne dans Wⁱ.

La première condition dans (1.36) est appelée condition de Lipschitz. Revenons maintenant à l'EDS (1.34). On peut prouver essentiellement le même résultat. Théorème 16 :

~ (1.37)

Sous les hypothèses de la définition 24 et si, de plus, les coefficients b et a^- satisfont, pour tous t, x, Y :

11b(t, x) -- b(t, ^Y)11 + 110(t, x) -- a(t, Y)11 C11x -- Y11

11b(t,0)11 + 110(t, 0)11 C

(où 11.11 désigne la norme euclidienne dans Wⁱ ou W^ig et C une constante positive), alors pour tout T > 0, il existe une et une seule solution forte X dans l'espace :

t.T

Li = [Xcontinus, adaptés; E (sup11X_t11²) < 00}

Remarque 6 :

1) L'unicité est comprise au sens presque sûr : si X et X' sont deux solutions fortes, alors, en dehors d'un ensemble de P-mesure nulle, on a Xt = X; pour tout t.

2) Le fait que l'on obtienne une solution dans un espace de type L ² n'est pas surprenant : c'est en effet dans ce cadre que l'on a développé le calcul stochastique. Nous avons déjà rencontré, l'EDS non homogène linéaire suivante :

dXt = pt(t)Xtdt + a(t)XtdWt, X0 = 1 (1.38)

où pt_? et a_? sont des fonctions continues bornées sur IR₊ . Elle a une solution explicite donnée par :

t

1

X_t = exp (1 a₍s)dW_s. + 1 (pt(_s)-- t

2 a(s)²) ds)

~ ~

on a ici d = q = 1 et les coefficients b(t,Xt) = pt(t)x et a(t,x) = a(t)x

satisfont à (1.37). Le théorème implique alors que la solution explicite est l'unique solution forte de cette équation.

Quand pt(t) = pt et a(t) = a sont constants, la solution est appelée «mouvement brownien géométrique avec drift ».

Comme pour les équations différentielles ordinaires, les hypothèses (1.37) de lipschitzianité et de bornitude sont suffisantes, mais pas nécessaires, pour obtenir l'existence et l'unicité de la solution. Par exemple, on a le résultat plus fort suivant, spécifique à la dimension 1.

Théorème 17 :

Dans le cas d'une EDS unidimensionnelle, si on a :

|b(t,x) -- b(t,Y)| + |a(t,x) -- a(t,Y)|^2 C|x -- Y| |b(t, 0)| + |a(t,0)| C

~ (1.39)

alors il existe une unique solution forte.

Ici, a n'est plus lipschitzienne, mais seulement höldérienne de rapport 1/2. Remarque 7 :

Très souvent, on considère des EDS homogènes. Pourtant, en pratique (spécialement en finance), les EDS non homogènes sont utiles. Dans ce cas, l'hypothèse (1.37) est un peu restrictive, car souvent, les coefficients « explosent » quand le temps tend vers l'infini. Il est alors utile de savoir que l'existence et l'unicité sont préservées, quand (1.37) est remplacée par :

Ilb(t, x) -- b(t, y)Il + Ila(t, x) -- a(t, y)Il CtIlx -- yIl

I

Ilb(t, 0)Il + Ila(t, 0)Il Ct

où t 1-0 Ct est une fonction croissante à valeurs finies. Processus d'Ornstein-Uhlenbeck :

Le mouvement brownien a ete construit pour modeliser le deplacement d'une particule microscopique, en suspension dans un liquide, soumise à l'agitation thermique. Une critique importante faite à cette modelisation est que les accroissements sont independants et ne dependent pas de la vitesse de la particule au debut de chaque periode. Un modèle plus sophistique, tenant compte de la vitesse de la particule, conduit alors au processus suivant.

Définition 25 :

On appelle processus d'Ornstein-Uhlenbeck le processus X à valeurs reelles, solution de l'EDS :

dXt = --bXtdt + cdw_t , Xo = xo (1.40)

où xo E R , b E R et c > 0.

Les coefficients de cette equation satisfont (1.37), donc il y a une solution et une seule. On peut montrer à l'aide du calcul stochastique que X a la forme explicite :

Xt = _xoe-bt + c f :e^-b(t-s)dw_s (1.41)

Nous observons que Xt est egal à l'integrale stochastique d'une fonction deterministe par rapport à un mouvement brownien. Nous en deduisons que Xt est une variable aleatoire gaussienne, comme limite (à une constante additive près) des approximations de Riemann de l'integrale stochastique, qui sont gaussiennes. On en deduit même que le processus X est gaussien. Il est facile de calculer sa fonction moyenne m(t) et sa fonction de covariance C(s, t) qui, pour s t , valent :

m(t) = xoe^-bt , C(s, t) = 2^C2b (e^-b(t-s) -- e-b(t+s))

(1.42)

Modèle de Vasicek :

Ce modèle, introduit par Vasicek en 1977, sert à modéliser un taux de placement à court terme. Ce taux fluctue suivant l'EDS décrite ci-dessous, qui généralise la précédente. Ici, x0 E I1 , b, b' E I1 et c > 0 :

dXt = (--bXt + b)dt + cdWt , X0 = x0 (1.43)

Les coefficients de cette équation satisfont (1.37), donc il y a une solution et une seule. De nouveau, on peut trouver une forme explicite pour la solution :

t

X_t = b' + (x -- b') e^-bt + c 1 e^-b(t_-^s)dW_s (1.44)

⁰ b 0

et, par le même argument que pour le processus de Ornstein-Uhlenbeck, nous voyons que cette unique solution forte X est un processus gaussien.

Remarque 8 :

Une caractéristique commune aux exemples précédents est que les coefficients dépendent linéairement de X, ou n'en dépendent pas du tout. C'est pourquoi on peut trouver une solution explicite. Dans le cas de coefficients plus généraux, il n'y a en général pas de solution explicite, comme pour une équation différentielle ordinaire.

CHAPITRE II :

Généralités sur les EDPs et la

méthode des différences finie

0. Introduction:

Le métier d'ingénieur est toujours besoin d'utilisation des logiciels de modélisation. Ces logiciels résolvent des équations telles que les équations aux dérivées partielles, mais la résolution ici est faite par une méthode discrète. Donc ces équations différentielles ne peuvent en général pas être résolues de façon exacte. Elles sont résolues de façon approchée, à l'aide des méthodes numériques.

Les méthodes numériques ne donnent pas la solution véritable du problème que l'on cherche à résoudre. Des méthodes numériques mal employées peuvent conduire à des résultats totalement faux, Il est indispensable pour un ingénieur de posséder des notions de base sur les méthodes numériques utilisé pour résoudre son problème.

Dans ce chapitre on présente une méthode numérique pour la résolution des EDPs,

La méthode des différences finies.

1. Généralité sur les équations aux dérivées partielles :

1.1. Définitions :

· Equation aux dérivées partielles :

Une équation aux dérivées partielles fait intervenir plusieurs variables indépen- dantes (temps et espace pour les équations de l'ingénieur), ainsi que les dérivées partielles de la variable dépendante ( c.-à-d. la solution recherchée) par rapport à ces variables indépendantes.

Exemple : l'équation de convection (parfois appelée advection)

=

ax 0 (2.1)

ac

+
at +u

ac

est une EDP, la variable dépendent est c , les variables indépendantes sont le temps t et l'espace x. La grandeur u (homogène à une vitesse) peut (ou non)être fonction de t ,x et c.

· Ordre d'une EDP :

L'ordre d'une EDP est défini exactement de la même façon que pour une EDO : c'est l'ordre le plus élevé parmi toutes les dérivées partielles de l'EDP.

Exemple : l'EDP (2.1) est une EDP d'ordre 1 (par rapport à t ou x ). En revanche, l'EDP suivant :

ac

+
at +u

ac

+ D

ax

ax2 = 0 (2.2)

a²c

est une EDP d'ordre 2 car sa dérivée partielle d'ordre le plus élevé est une dérivée seconde(en l'occurrence par rapport à x).

· EDP linéaires, quasi-linéaires et non- linéaires :

EDP linéaire : est une EDP qui ne fait intervenir que des combinaisons linéaires des dérivées partiales de la variable dépendante.

Exemple : l'équation suivante est linéaire :

a^2u A+2B ax^2

a^2u ~~~ ~~~~ ~ ~ ~2.3~ ~~~ ~ ~ ~~, ~, ~~ ~~ , ~~ ay) = 0

où les coefficients A, B et C sont dépend de X et Y (connues), F est une fonction lineaire et U est la fonction inconnue (à déterminer).

EDP quasi-linéaire : est une équation linéaire par rapport aux dérivées partielles d'ordres supérieurs.

Exemple : l'équation précédente (2.3) est quasi-linéaire si les coefficients A, B et C sont en fonction de U et/ou de X, Y et F peut être non linéaire.

EDP non-linéaire : est une EDP où l'une des dérivées partielles intervient comme argument d'une fonction non-linéaire.

Exemple : l'équation suivante:

Co

~~~~ ~~~~ ~ ~~ ~~~~ ~2.4~ ~~~~ ~ ~~ ~~~ ~~~ ~ ~~ ~~~ ~~~ ~ 0

est non-linéaire si au moins des fonctions Co(x) et C1(x) n'est pas linéaire par rapport à x.

1.2. Conditions initiales et conditions aux limites :

Contrairement aux EDOs, les conditions initiales ne suffisent pas à assurer l'unicité de solution. Il faut également fournir des conditions aux limites. Pour certains types d'équations (ex. EDPs elliptiques), le concept de condition initiale n'a pas de sens.

Les conditions initiales et les conditions aux limites se distinguent de la manière suivante :

· Une condition initiale s'applique pour une valeur donnée (et unique) d'une variable indépendante. A partir de cette condition initiale, il est possible de déduire la solution pour toutes les autres valeurs de la variable indépendante.

· Une condition aux limites est appliquée en tout point de la frontière du domaine sur lequel on souhaite résoudre l'équation (et non en un point unique). Il n'est pas possible de déterminer la solution en partant simplement d'un seul point de la limite de domaine et en progressant à l'intérieur de celui-ci, car la solution est également conditionnée par sa valeur en tous les autres points de la frontière.

1.3. Classification des EDPs du second ordre :

Les EDPs du second ordre représente une classe importante des EDPs du monde de l'ingénierie, on traite dans quel suit des EDPs second ordre quasi-linéaires de la forme générale :

a²u A + B
aX²

a²u a²u

+

aY2 D

au

+ E

aX

aY + F = 0 (2.5)

au

où A ,B,C,D,F sont des fonctions de X,Y et u. X et Y sont les variables
indépendantes (ce pourrait être t,x,y etc) de l'EDP et u la variable dépendante.

Selon la valeur des coefficients A ,B et C, l'EDP est dite hyperbolique, parabolique ou elliptique.

EDPs hyperbolique :

Une EDP du type (2.5) est dite hyperbolique si son discriminant A= B² -- 4AC est strictement positif.

Exemple : l'équation suivante est du type hyperbolique :

a²u

at²

A. ax2

² = 0 (2.6)

En effet, l'équation (2.5) peut être mise sous la forme (2.6) en posant X = t, Y = x, A = 1, C = --A.², et B = D = E = F = 0. Il est facile de vérifier que A= B² -- 4AC > 0. A noter que A. est homogène à une vitesse.

EDPs paraboliques :

Une EDP du type (2.5) est dite paraboliques si son discriminant A= B² -- 4AC est nul.

Exemple : l'équation suivante est du type paraboliques :

ax2 = 0, v > 0 (2.7)

au

-- v

at

a²u

En effet, l'équation (2.5) peut être mise sous la forme (2.7) en posant X = t, Y = x , D = 1, C = v, et A = B = E = F = 0. On vérifiera que A= B² -- 4AC = 0.

EDPs elliptique :

Une EDP du type (2.5) est dite elliptique si son discriminant A= B² -- 4AC est strictement négatif.

Exemple : l'équation suivante est du type elliptique :

0y2 = Q (2.8)

a²u ax2 +

a²u

En effet, l'équation (2.5) peut être mise sous la forme (2.8) en posant X = t, Y = x , A = C = 1, et B = D = E = F = 0, on peut vérifier que A= B² -- 4AC < 0.

· Une aide mnémotechnique :

L'<< astuce >> suivante peut etre utilisé pour déterminer la nature d'une EDP du seconde ordre : il suffit, dans l'équation originale, de remplacer les dérivées a^Pu / a XP (p = 1,2) par X^P et a^Pu / a Y^P (p = 1,2) par Y^Pet le second membre par une constante. Ainsi (2.5) devient :

AX²+BXY + Y² + DX + EX + F = Cst (2.9)

L'équation (2.9) est l'équation d'une courbe conique dans le plan(X, Y). Si cette courbe est une ellipse, l'EDP est elliptique ; si la courbe est une parabole, l'EDP est parabolique ; enfin, si la courbe est une hyperbole l'EDP est hyperbolique.

Exemple : en appliquant la méthode ci-dessus à l'équation(2.6), on obtient :

T² -- A.²X² = Cst (2.10)
qui est l'équation d'une courbe hyperbolique dans le plan(X, T).

1.4. EDPs du second d'ordre à plusieurs variables indépendantes:

Le principe de classification des EDPs reste le même quand on se trouve en présence de 3(ou 4) variable indépendantes. Par exemple, l'EDP (2.6) peut être généralisée à deux dimensions d'espace :

a²U a²u a²u

A.2 A² = 0 (2.11)

at2 x ax2 Y ay2

où A._x² et 4 sont les vitesses de propagation des ondes les directions x et y (elles sont égales dans un milieu isotrope). En utilisant la transformation exposée dans le paragraphe précédent. On obtient :

T² -- A._x²X²--A._Y²Y² = Cst (2.12)

1.5. Besoins en termes de conditions initiales et aux limites :

· EDPs hyperboliques :

Les EDPs hyperbolique abordées dans ce chapitre comprendrons en général une variable de temps et une (ou deux) d'espace. Elles seront de la forme (2.6). On cherche une solution de (2.6) en tout point d'un domaine de calcul SI = [x1,x2] et pour des dates t 0. Pour pouvoir déterminer la solution de (2.6) de façon unique, il faut connaître :

· la valeur de U en tout point du domaine de calcul à t = 0 (condition initiale) ;

· la valeur de U à tout date t 0 aux limites du domaine (condition aux limites) ;

Comme il y a des EDPs hyperbolique n'ont besoin que d'une condition à la limite.

· EDPs paraboliques :

La plupart des EDPs paraboliques de l'ingénieur sont d'ordre 1 par rapport au temps et d'ordre 2 par rapport à une (ou plusieurs) dimension(s) d'espace. C'est le cas de l'équation (2.7). Les conditions aux limites nécessaires à l'existence et l'unicité de la solution sont les mêmes pour que les EDPs hyperboliques.

· EDPs elliptique :

Les EDPs elliptiques étudiées dans ce qui suit impliquent 2 dimensions de l'espace et aucune de temps. Cela signifie que la solution ne dépend pas du temps. Dans ce cas, seules des conditions aux limites sont nécessaires.

2. Méthode des différences finies pour les EDP paraboliques :

Dans ce paragraphe nous avons aborderons la résolution numérique de problèmes d'équations aux dérivées partielles d'évolution par la méthode des différences finies. Nous étudierons un type de problèmes d'évolution du premier ordre en temps, dénommés également problèmes paraboliques.

Les équations intervenant dans ces problèmes sont constituées pour partie d'une combinaison de dérivées partielles par rapport à la variable temporelle dont nous détaillerons le traitement numérique et pour partie d'une combinaison de dérivées partielles par rapport à la variable spatiale ; le problème pouvant être posé dans un domaine , monodimensionnel, bidimensionnel ou tridimensionnel; pour simplifier l'exposé nous considérerons que le domaine est le segment [0, 1], le cas bi et tridimensionnel ne présentant pas de difficultés majeures.

2.1. Problèmes du premier ordre en temps : équation de la chaleur :

· Position du problème :

On considère une barre métallique de longueur unité; on suppose que cette barre est soumise à un apport de chaleur f (x, t ) par unité de longueur et de temps et

que, de plus, la température u (x, t) de la barre est maintenue à zéro à chacune de ses extrémités. On désigne par C la capacité thermique et par K le coefficient de diffusion de chaleur ; rappelons que le coefficient de diffusion de chaleur exprime qu'en tous les points de la barre les températures ont tendance à s'uniformiser. En supposant que la dimension transversale de la barre est négligeable par rapport à sa dimension longitudinale, la modélisation de ce problème conduit à déterminer u (x, t) en chaque point x E[0,1] et pour tout instant t E[0,T ], T représentant l'horizon d'intégration, solution de l'équation de la chaleur :

~ ~ ~~~~, ~~

~~ ~ ~ ~~~~~, ~~

~~~ ~ ~~~, ~~, ~ ~ ~0,1~ ; ~ ~ ~0, ~~

~ ~~0, ~~ ~ ~~1, ~~ ~ 0, ~ ~ ~0, ~~

~

u(x, 0) = uo(x), x E [0,1]

avec uo(x) température initiale de la barre.

Remarque :

Les conditions aux limites de type Dirichlet homogènes u (0, t) = u (1, t) = 0 correspondent, évidemment, au fait que les extrémités de la barre sont maintenues à une température nulle.

Remarque :

Dans le cas d'un matériel non homogène, le problème précédent est remplacé par :

~ ~~~~ ~~~~, ~~ ~~ ~~~~~ ~~~~, ~~

~~ ~ ~ ~~~, ~~, ~ ~ ~0,1~ ; ~ ~ ~0, ~~

~~ ~ ~

~

~ ~~0, ~~ ~ ~~1, ~~ ~ 0, ~ ~ ~0, ~~

~

~ ~ ~ ~0,1~
~~~, 0 ~ ~~~~~, Dans la suite, on supposera pour simplifier que C = 1,K = a (constante positive) ;

on ne considérera de plus que les conditions aux limites de Dirichlet homogènes, ce qui conduira à la formulation du problème modèle suivant :

~~~~~, ~~

~~ ~ ~ ~~~~~, ~~

~~~ ~ ~~~, ~~, ~ ~ ~0,1~, ~ ~ ~0, ~~

~ ~~0, ~~ ~ ~~1, ~~ ~ 0, ~ ~ ~0, ~~

~

u(x, 0) = uo(x), x E [0,1]

2.2.1. Schémas numériques de discrétisation par différences finies :

La discrétisation du problème précédent consiste à remplacer par une technique appropriée le problème continu par un système linéaire algébrique. L'approximation s'effectue par la manière suivante:

~

on divise l'intervalle [0,1] en (n + 1) intervalles de longueur h = .

n+1

on divise l'intervalle [0, T ] en M intervalles de temps k, tels que T = Mk. On pose xi = ih, pour i E f 1 , 2 , ..., n}, t_m = mk pour m E (0, 1, . . . , M} ; on écrit le problème modèle au point (xi, t_m ) avec 1 i n et m 0:

(au ~~ ~ ~ ~~~

~~~~ ~~~, ~~~ ~ f(xi,t_i)

~~~~~,~~~ ~~~~~~,~~~

On remplace et par des quotients différentiels faisant intervenir

~~ ~~~

les valeurs de u aux points xi, 0 i n + 1, et aux instants t_m, m > 0 . Exemple : on écrira :

a^2u(xi, tm)		~~~~~~, ~~~ ~ 2, ~~~ ~ ~~~~~~, ~~~ ~ + O(h^2) h^2	ax^e

	au(xi, tm)		u(xi, tm+i) -- u(xi, tin) = + O(k) k	et on obtient ainsi :
at

u(xi, tm+i) -- u(xi, tm)

~

~

~~ ~~~~~~, ~~~ ~ 2, ~~~ ~ ~~~~~~, ~~~

~ ~ O, ~~ ~ ~~~~, ~~~

~~

~ 1 ~ ~ ~ ~ ~~ ~ ~ 0

u(xo, tm) = u(xn+1, tm) = 0 u(xi, 0) = uo(xi)

~

~

Cette écriture permet de définir les schémas numériques définissant um r-r^, u(xi, t_m), en traitant de manière particulière le terme temporel. La façon la plus générale de définir les schémas aux différences finies est de considérer une combinaison convexe de l'équation de la chaleur aux instants t_m et t_m+1 ; plus

précisément, soit 0 E [0,1] un nombre réel ; écrivons l'approximation de (1 - 0) fois le problème modèle considéré à l'instant t_m et de 0 le même problème à l'instant t_m+1, ce qui permet de définir le 0 -- schema suivant :

uro_. -- 2ui^m + ur 1r _Er -- 24¹⁺¹ + uTV

ui um-- a ((1 -- 0) + 0u

k h² h²

= (1 -- 0)fi^m + 0fi^m+1 , 1 < i < n et m > 0

~ ug = u_n^p ₊₁ = 0, p = m et p = m + 1

~ ui^o = uo(ih)

avec f_i^p = f(xi, t_p), pour p = m ou p = m + 1

Pour 0 = 0, le schéma précédent définit un schéma explicite, c'est-à-dire que les valeurs um étant connues (par la condition initiale lorsque m = 0, et pour m > 0 par le calcul de l'itéré en temps précédent), les valeurs sont ur^l déterminées par la récurrence :

1

ur + ¹ = k f iⁿ' + ur + a (ur^. 1-- 24 ' + 4¹: 1) , 1 < i < n et m > 0 um = u^m = 0 o n+1
ui^o = uo(ih)

avec a =

ak h^2.

Si l'on définit les vecteurs de dimension n par :

F^m = (kf~ ^m, ... , kfnm)t U^m = (u^m, ,un^m)^t et la matrice B de dimension n par :

2 --1

--1 2

B =
·..

)

--1 2 --1

--1 2

le schéma explicite s'écrit matriciellement :

U^m+1 = F^m + (I -- aB)U^m

Pour 0 # 0, le schéma précédent définit un schéma implicite, c'est-à-dire qu'à chaque pas de temps, en supposant connu le vecteur U^m (toujours par la condition

ainsi que le schéma rétrograde implicite, également à trois niveaux, suivant :

~~ ~~~ ~ ~~ ~~~ ~~~ ~ ~~ ~~~ ~ ~~~~

~

~~~~

~ ~ ~~

~ ~ ~

~ 2 ~~

1

~ ~~ ~ ~ ~~~~

~~~ ~ 0

~ ~ ~~ ~ ~ ~~~~~~

~ ~~ ~

~ ~ ~ ~ ~~ ~ ~ 0

initiale lorsque m = 0 et pour m > 0 par le calcul de l'itéré en temps précédent), le vecteur U^m+1 est obtenu en résolvant le système linéaire suivant:

(I -- a0B)U^m+1 = (1 -- 6)F^m + OF' + (I -- a(1 -- 6)B)U^m

Notons immédiatement que, pour 0 E0, 1], la matrice (I -- a0B) est inversible ; dans le cas du problème monodimensionnel en espace, la résolution de ce système s'effectue par une simple adaptation de la méthode de Gauss, (pour l'adaptation de cet algorithme au cas des matrices tridiagonales). Pour des valeurs particulières du paramètre 0, on retrouve des schémas numériques particuliers ; ainsi pour 0 = 1, on obtient le schéma purement implicite qui s'écrit matriciellement :

(I -- aB)U^m+1 = F^m+1 + U^m

alors que pour 0 = 0,5, on obtient le schéma implicite de Crank-Nicolson:

1 a

(I --^a₂ 13) U^m+1 = 2 (F^m + F^m+1) + -- 2 /3) U^m

Remarque :

Il existe d'autres types de schémas numériques pour résoudre le problème modèle, par exemple le schéma saute-mouton défini par :

~ ~

~~ ~~~ ~ ~~

~ ~ ~~ ~

2 ~ ~ ~~~~

~ ~ 2 ~ ~ ~~~~

1 ~ ~ ~ ~ ~~ ~ ~ 0

~~

~ ~ ~~~~

~ ~ 0

~ ~~

~

u° = uo(ih)

C'est un schéma à trois niveaux explicite qui, bien que plus précis que les schémas précédents, n'est pas intéressant sur le plan numérique car il est instable, c'est-à-dire sensible à la propagation des erreurs systématiques d'approximation, de chute ou de troncature. Citons également le schéma explicite à trois niveaux proposé par Du Fort et Frankel:

~~~

~ 1 ~

~3 ~~~~ ~ ~ ~~~~

~~~ ~ 2

2 ~~ ~ ~ ~~~~

2 ~~ ~~~ ~ 2 ~ ~ 1 ~ ~~

~ ~ ~~

1 ~ ~ ~ ~ ~~ ~ ~ 0

~ ~~ ~~~ ~ ~~~~

= 1.10(ih)

~~~ ~ 0

~ ~ ~~ ~

2.2.2. Erreur de troncature, consistance et ordre d'un schéma :

~

Pour chacun des schémas précédents, on remplace dans le schéma le terme ~~ par

u(xi, t_m) et l'on définit l'erreur de troncature, notée, comme la différence entre le premier et le second membre de ces quantités. Ainsi pour le 0 -- schéma considéré au paragraphe précédent, on a l'expression suivante de l'erreur de troncature :

~~~~, ~~~~~ ~ ~~~~, ~~~ ~ ~ ~~1 ~ ~~ ~~~~~~, ~~~ ~ 2, ~~~ ~ ~~~~~~, ~~~

~~ ~ ~ ~

~~

~~~~~~, ~~~~~ ~ 2, ~~~~~ ~ ~~~~~~, ~~~~~

~ ~ ~ ~ ~1 ~ ~~~~~~, ~~~ ~ ~~~~~, ~~~

~~

Remarque :

On définit l'erreur de troncature associée au schéma saute-mouton, au schéma rétrograde et au schéma de Du Fort et Frankel de manière analogue.

On adapte à la situation des problèmes d'évolution les notions d'erreur de troncature, de consistance et d'ordre de schéma comme suit :

Définition 1 : On appelle erreur de troncature E, la quantité définie par : E = Max(E1 ⁿ, 1 < i < n,1 < m < M}

Définition 2 :

Un schéma numérique est consistant si E --> 0 lorsque h --> 0 et k --> 0. Définition 3 :

Un schéma est précis à l'ordre p en espace et à l'ordre q en temps, s'il existe une constante C', indépendante de h et k telle que |E| Cⁱ(h^P + k ^q)

On a le résultat suivant:

Théorème 1 :

On suppose que la solution du problème modèle est suffisamment régulière ; alors pour 0 E [0, 1], 0 # 0,5, le 0-schéma est d'ordre 2 en espace et d'ordre 1 en

temps. Pour 0 = 0,5, le schéma de Crank-Nicolson est d'ordre 2 en espace et en temps.

Remarque : Il résulte du résultat précédent que le schéma explicite et le schéma purement implicite sont d'ordre 2 en espace et d'ordre 1 en temps.

Remarque :

On vérifie que le schéma saute-mouton et le schéma rétrograde sont d'ordre 2 en espace et en temps. L'étude du schéma de Du Fort et Frankel est plus délicate; en

effet, ce dernier n'est consistant que si le rapport kh --> 0 lorsque h --> 0 et k --> 0.
De plus, si l'on veut que le schéma de Du Fort et Frankel soit d'ordre 2 en espace, il faut prendre k = 0(h²).

2.2.3. Stabilité des schémas numériques :

La stabilité d'un schéma numérique pour l'approximation de la solution d'une équation aux dérivées partielles est une question importante et, apparemment, nouvelle ; lorsque le pas de discrétisation tendait vers zéro, le système linéaire pouvait devenir mal conditionné, c'est-à-dire que la solution calculée du système discret pouvait être dénaturée ; en effet, la question qui se pose peut être formulée de la façon suivante : compte tenu d'une part des erreurs systématiques introduites lors de la construction du schéma d'approximation et d'autre part compte tenu également de la mauvaise représentation des nombres en machine qui induisent des erreurs d'arrondi ou de troncature, est-ce que l'erreur va rester bornée et tendre vers zéro lorsque les paramètres de discrétisation h et k vont tendre vers zéro ?

Il existe de nombreuses méthodes pour étudier la stabilité d'un schéma numérique ; nous en étudierons deux, connues sous le nom de :

méthode matricielle ;

méthode de Fourier ou méthode du coefficient d'amplification.

La méthode matricielle : consiste à considérer qu'un schéma est stable si, lorsqu'une erreur est introduite sur les conditions initiales (par exemple, erreur d'arrondi ou de troncature), l'erreur sur la solution calculée à un instant fixé t n'est pas amplifiée ; l'idéal est, bien entendu, qu'elle ne soit pas amplifiée du tout et que

supi,m(Iui^m - u(xi, t_m)1) reste borné quels que soient h et k. L'étude de la stabilité d'un schéma revient donc en fait à celle de l'amplification d'une perturbation sur la condition initiale U^°; on va vérifier que l'erreur sur la solution cherchée satisfait le schéma numérique homogène associé, ce qui va permettre de dégager des critères de stabilité. A ce stade, on peut donc donner la définition intuitive suivante de la stabilité :

Definition :

Un schéma numérique est stable si la solution du schéma homogène associé est bornée en tout point xi et t_m quels que soient h et k.

Pour illustrer cette approche intuitive, considérons le 0 - schéma dont la formulation est la suivante :

(I - a0B)U^m+1 = (1 - 0)F^m + 0F^m+1 + (I - a(1 - 0)B)U^m

et soit U^m+1 la solution exacte du schéma numérique et V^m+~ celle du schéma perturbé : on peut traduire l'imprécision par la relation suivante :

V^m = U^m + E^m, Vm > 0

avec E^m erreur au pas m.

En remplaçant U^m+1 par V^m+~ dans le schéma précédent, on obtient :

(I - a0B)E^m+1 = (I - a(1 - 0)B)E^m

et on vérifie ainsi que le terme d'erreur est solution du schéma numérique homogène associé. Cette relation peut encore s'écrire :

E^m+1 = (I - a0B)^-1(I - a(1 - 0)B)E^m = PE^m =
· = Pm+1E° avec P = (I - a0B)^-1(I - a(1 - 0)B).

Le schéma numérique sera donc stable si la norme de la matrice P^mreste bornée lorsque m -> 00. La matrice B étant symétrique, on vérifie aisément que la matrice P est normale (c'est-à-dire que P^tP = PP^t ) ; on peut alors exprimer la relation précédente dans la base des vecteurs propres de la matrice P ; soit A la matrice des valeurs propres de la matrice P. La relation E^m+1 = P^m+1E^°, s'écrit alors E~^m+~ = A^m+1E^°, où E~^m+~ est l'écriture de E^m+1 dans la base des vecteurs propres, et une condition suffisante de stabilité est donc :

|uk(P) | < 1,V k E (1,...,n )

avec uk(P) valeurs propres de la matrice P, c'est-à-dire les coefficients diagonaux de la matrice .4.

On vérifie facilement que :

1 -- a(1 -- 0)A.k(B)

uk(P) 1 + a0A.k(B) '

= V k E (1,...,n )

avec A.k(B) valeurs propres de la matrice B; on vérifie alors que la

condition - 1 < uk(P) < 1, V k E (1, ... , n ), se traduit par la condition suivante :

ak a = < h^2

1 , 0 E [0,1] 2(1 -- 20)

qui représente donc la condition de stabilité. On peut résumer l'étude précédente par le théorème suivant :

Théorème 2 :

Le 0-schéma est inconditionnellement stable si 0 > 0,5 ; si 0 < 0,5, ce schéma est stable si :

ak

a = <

h²

1

2(1 -- 20)

Corollaire 1 :

Le schéma purement implicite (0 = 1) et le schéma de Crank-Nicolson (0 = 0,5) sont inconditionnellement stables. Le schéma purement explicite (0 = 0) est stable à condition que :

1

2

ak

a = <

h²

Remarque :

Le schéma saute-mouton est instable ; le schéma de Du Fort et Frankel et le schéma rétrograde sont inconditionnellement stables.

La méthode matricielle présentée ci-dessus est applicable dans le cas où la matrice B est symétrique. Or il existe des situations où cette matrice ne vérifie pas cette condition ; c'est le cas par exemple du problème de convection-diffusion d'évolution. Dans ce cas, le calcul du coefficient d'amplification permet de s'affranchir de cette

contrainte et fournit une méthode d'analyse plus générale de la stabilité. La définition de la stabilité étant inchangée, considérons toujours le 0-schéma homogène :

ur^l -- -- a[(1 -- ~~~~~~~

~ -- 2ur + ~~~~

~ ~

~~~~~~~

~~~ -- 24¹⁺¹ + ~~~~

~~~~~ = 0,1 < i < n et m > 0

que l'on peut encore formellement écrire en repassant à la variable d'espace continue x et en mettant dans le premier membre les termes à l'instant discret (m + 1) et dans le second membre ceux à l'instant discret m :

u^m+1(x) -- a0(u^m+1(x + h) -- 2u^m+1(x) + u^m+1(x -- h))

= u^m(x) + al1 -- 0) ~~~~~

~~~~~ + h) -- 24¹⁺¹(x) + ~~~~

~~~~~ -- h)) , m > 0

Considérons la transformée de Fourier du schéma ainsi écrit ; on rappelle que cette transformation est définie par :

~~

1

~~~~~ ~ _v2~ ~ ~~~~~~~~~~~~~ ~~, ~~ ~ ~1

~~

On vérifie aisément par un calcul direct que :

~~

1

v2~ ~ ~~~ ~ ~~~~~~~~~~~ ~~ ~ ~~~~~~~~~~~~~~ _00

et le schéma discret en temps s'écrit alors :

It^m+1(0 = c(a)u^m(

où c() est le coefficient d'amplification défini par :

1 + 2a(1 -- 0)(cos(~h) -- 1)

~~~~ ~1 -- 2a0(cos(~h) -- 1)

~

1 -- 4a(1 -- 0)sin2 ~~~ 2 ~

1 + 4a0sin² (^Q1₂ )

=1

asin2 ~~~

2 ~

1 + a0sin2 ~~~

2 ~

Comme précédemment, une condition suffisante de stabilité est donnée par --1 c(0 1 ; or on a évidemment c() < 1; la vérification de la seconde
inégalité conduit à la condition suivante :

ak a= < h^2

1 , ~ ~ ~0,1~ 21 ~ 2

qui correspond à la condition obtenue par la méthode matricielle ; on retrouve donc les résultats obtenus et énoncés dans le théorème 2 et le corollaire 1.

2.2.4. Convergence des schémas :

Soit e^m , le vecteur erreur au temps discret m et défini par :

e^m = (er, ... , ei ^m, ... , ein^t

avec :

~~ ~ ~~ ~ -- u(xi,t_m), 1 i n,m 0

Définition :

Un schéma numérique est convergent si :

|~T -- u(xi, t_m)| --) 0 lorsque h --) 0 et k --) 0 ,Vi e (1,...,n),Vm 0.

Considérons le 0- schéma ; en considérant d'une part l'écriture du schéma et d'autre part celle de l'erreur de troncature et en soustrayant membre à membre, on vérifie aisément que les composantes de l'erreur satisfont le schéma numérique, avec un second membre égal à --Elⁿ, soit :

~~~ ~ ~~ ~

~ ~~ ~ ~ 2 ~ ~ ~~~~

~

~ ~~1 ~ ~~ ~~~~

~ ~~

~ ~~~

~~~~

~~~ ~ 2 ~~~ ~ ~~~~

~~ ~ ~~~ ~, 1 ~ ~ ~ ~ ~~ ~ ~ 0

~ ~~

~ ~ ~ ~ ~ 1

~ ~~ ~ ~~

~ ~~~~ ~ 0, ~ ~ ~

~

soit matriciellement :

(I -- a0B)e^m+1 = E^m + (I -- a(1 -- 0)B)e^m

où les composantes du vecteur E^msont égales à l'erreur de troncature, au signe près ; en multipliant scalairement par e^m+1, il vient :

((I -- a0B)e^m+1, e^m+1) = (E^m, e^m+1) + ((I -- a(1 -- 0)B)e^m, eⁿ⁺¹),Vm 0 Or, la matrice(I -- a0B) étant symétrique définie positive, on peut minorer le second membre; si, de plus, on note par p , la norme matricielle de (I -- a(1 -- 6)B) , on obtient finalement :

eP ~ = 0

Ile^m+1112 ^filEm12 _+pllem112,vm 0

avec /3 inverse de la plus petite valeur propre de la matrice (I -- a0B).

En posantT^m = Ile^m112,Vm 0, on obtient finalement :

T^m+1 /311E^m112 + pT^m, Vm 0

Supposons, de plus, que le schéma vérifie la condition de stabilité, ce qui a pour conséquence que la quantité T^m = Ile^m112, Vm 0, est borné. Pour pouvoir obtenir une majoration de la norme de l'erreur, on utilise le lemme suivant dont le résultat s'obtient par récurrence :

Lemme 1 : Si l'on a une suite qui vérifie la relation de récurrence :

T^m+1 pT^m + a,Vⁿ 0,Vm 0, p # 1

alors :

p^m -- 1

T^m < pT^° + (5 IlE^m112,Vm 0

p -- 1

En appliquant ce résultat à notre situation, compte tenu du fait que E^° = 0 , on obtient :

~~~ ~ ~~~~~ ~ ~ ~~ ~ 1

~ ~ 1 , ~~ ~ 0

Si, de plus, le schémaa numériquee est consistant, l'inégalitée précédentee nous permet de vérifierr qu'ill est convergent. On peut donc énoncerr le résultatt suivant ^:Théorèmee 3 : Si le 0 - schémaa est stable et consistant alors ilt est convergent.

Le résultatt précédentt exprime que la stabilitée et la consistance constituent une condition suffisante pour la convergence. En fait le théorèmee de Lax( le théorèmee de Lax prévoitt que dans un problèmee bien posé,, et avec un schémaa numériquee consistant, la stabilitée est une condition nécessairee et suffisante pour la convergence), nous donne le résultatt suivant :

Théorèmee 4 :

La stabilitée est une condition nécessairee et suffisante pour qu'unn schémaa numériquee soit convergent vers la solution du problème,, si ce dernier est discrétisée par un schémaa consistant.

CHAPITRE III :

Interprétation probabiliste des

EDPs

0. Introduction:

Le but de ce chapitre est de montrer les liens qui peuvent exister entre la théorie des processus stochastiques et les équations aux dérivées partielles (EDP). Les processus stochastiques utilisés sont des processus possédant la propriété de Markov. L'idée principale est de montrer que l'espérance mathématique de fonctionnelles de ces processus fournit une représentation probabiliste de solutions de certaines équations. Nous donnons les outils nécessaires pour cette présentation: les processus de diffusions obtenues comme solutions d'équations différentielles stochastiques à partir du processus de Wiener. Ces diffusions nous permettent, de représenter les solutions d'EDP du second ordre, Nous exprimons ici les solutions d'EDP du second ordre paraboliques avec divers types de conditions aux bords, comme espérances de fonctionnelles de processus de Markov correctement choisis. Nous étudions d'abord l'exemple le plus simple des dérivées partielles linéaires, l'équation de la chaleur, ensuite les techniques de martingales, de martingales exponentielles qui conduisent à la généralisation connue sous le nom d'équations de Feynman-Kac.

1. Générateur infinitésimal :

Il est fondamental pour beaucoup d'applications que nous pouvons associer un

operateur différentiel de second ordre A à une diffusion d'Itô Xt ; la relation de base entre A et Xt est que A est le générateur de la diffusion d'Itô Xt.

Définition 1:

Soit (Xt} une diffusion d'Itô homogène dans W^L, le générateur infinitésimal A de Xt est définie par :

Af(x) = lim

~~~

Ex[fVD] -- f(x)

_t , x E

L'ensemble des fonctions ( f: 1l --> 11 ) telle que la limite existe en x est noté par DA(x), ou la notation : l'ensemble des fonctions pour qui la limite existe pour tout x E Rⁿ.

Pour montre la relation entre A et les coefficients a et b on a le théorème suivant : Théorème 1:

Soit Xt la diffusion d'Itô homogène:

dXt = b(Xt)dt + Xt)dWt

Si f e cg(Ir) alors f e DA et

2 ~~~~~~~,~~~~ ~~~ 1

(x) (3.1)

~~~~~~

(x) +

~

~,~

Proposition 1:

Soit Xt la diffusion d'Itô non homogène:

dXt = b(t,Xt)dt + t, Xt)dWt

a² f

+ ₂ I(aa^T)i_, j(t, x) 1

(t, x) (3.2)

~~~~~~

~,~

Si f E c_o^l,20R⁺ x IV) alors f E DA_t et

~~~~~, ~~ ~ ~ ~~~~, ~~ ~~ ~~, ~~

oxi

~

Exemples 1:

1. Le mouvement brownien de dimension n et bien sûr la solution de l'EDS :

dXt = dWt

(3.3)

² Ox2

au

1 a2u

at .

c'est à dire nous avons b = 0 et a = I_n la matrice identité; donc le générateur de Wrest:

~~ ~ ₂ ~ ~~~

1 ; ~ ~ ~~~~, ... ~~~ ~ ~~ ~~~~~

~~~ ~

~

1

i.e. A = 2 ?, ou ? est le laplacien.

2. Soit W le Mouvement Brownien à une dimension et soit X = ~~~

~~~ la solution de

l'EDS :

~

~~~ ~ ~~ ; ~~~0~ ~ ~~ ~~~ = dWt ; X2(0) = xo i.e. dX = bdt + cidWt ; X(0) = ~~~

~~~.

Avec b = (₀¹) and a = CD, le générateur A de X est donné par :

1

⁰²f

2

ax^e

^Af = at^of

+

Remarque : tous les espaces cités dans ce chapitre sont définis dans l'annexe A. 2. Interprétation probabiliste des EDP:

2.1 Equation de la chaleur :

Nous allons présenter quelques rapports entre les diffusions et les équations aux dérivées partielles. Grâce à la formule d'Itô, il est possible de donner une interprétation probabiliste à certaines équations aux dérivées partielles, ce qui permet ainsi de prouver l'existence de solutions.

L'exemple le plus simple est celui de l'équation de la chaleur en dimension 1. Considérons une barre métallique infinie, assimilée à l'axe réel. Cette barre est chauffée et on note f(x) sa température à l'instant t = 0 et pour la position x sur la barre. Soit u(t, x), la température de la barre au temps t et à la position x. Avec un choix approprié d'unités, donc la fonction u est solution de l'équation de la chaleur en dimension 1 :

avec la condition initiale, u(0,x) = f(x), x E IR, que l'on suppose continue. On a alors le théorème suivant :

Théorème2 :

1) Si u est une fonction continue sur [0, oo[de classe C^1,2 sur [0, oo[ x IR, et solution de (3.3), alors :

u(t, x) = E_x(f(Wt))

(3.4)

, tf(t,x) E [0,oo[ x IR^,+co

= f-co f _37)p(t; x, y)dy

où y 1-0 p(t; x, y) est la densité du mouvement brownien issu de x au temps t. Cela entraîne donc l'unicité d'une telle solution.

2) Supposons qu'il existe a > 0 tel que la fonction f vérifie f _IR,e^-ax2|f (x)| dx < 8, alors, pour tout 0 < t < 21 , la fonction u définie par (3.4) est dérivable à tous
ordres et est solution de (3.3).

Cette fonction u est donc l'unique solution de l'équation de la chaleur (3.3), qui soit une fonction continue et de classe C^1,2 sur [0, oo[ x IR,.

Preuve :

La preuve de la deuxième assertion utilise la formule explicite de la densité p et consiste juste en une vérification des dérivations sous le signe somme.

Le plus intéressant est la preuve de la première assertion, qui va pouvoir s'étendre à des cas plus généraux. L'interprétation probabiliste introduit un retournement du temps. Soit T un réel fixé. Pour montrer l'unicité de la solution u, on applique la formule d'Itô à la fonction (t,x) 1-0 u(T -- t,x) et au mouvement brownien W issu

de x.

On obtient alors :

u(T -- t, Wt) = u(T, x) + M_t

où M est une martingale nulle en 0.

Le terme à variation finie s'annule, du fait que u est solution de (3.3). En prenant l'espérance dans chaque terme de cette égalité, au temps t = T, on obtient finalement que :

u(T, x) = E_x(u(0, WT)) = E_x(f(WT))

De ce theorème, on deduit egalement une solution de l'equation de la chaleur retrograde, pour T > 0 fixe, donnee par :

at 1 a2v

+ = 0 (3.5)

at ² ax2

avec la condition terminale, v(T, x) = f(x), x E IR

Si nous sommes dans les conditions precedentes, l'unique solution de classe C^1,2 de (3.5) est donnee par :

v(t, x) = u(T -- t, x) = Ex(f(WT_t))

2.2. Formule de Feynman-Kac multidimensionnelle:

On va generaliser l'approche precedente à d'autres equations paraboliques. Plus precisement, nous allons tout d'abord donner une representation probabiliste de la solution d'une equation retrograde, apparaissant classiquement en mathematiques financières. Pour T > 0 fixe, introduisons l'equation :

v 2

-- a + kv = 1 Av + g , (t, x) E [0, T) x IR^a

at

(3.6)

v(T, x) = f(x) ,x E I^a

R

pour des fonctions k: IR^~ --> [0, +00), g: IR ^a x [0, +00) --> IR et f: IR ^a --> IR , les fonctions k, g et f etant supposees continues et bornees. Une solution v est dite solution du problème de Cauchy pour l'equation retrograde (3.6) avec potentiel k et lagrangien g. Elle admet la representation probabiliste suivante.

Théorème 3 :

Supposons que v soit une fonction continue sur [0, T) x IR ^d, de classe C1,2 sur [0, T) x IR ^a, solution de (3.6). Alors v admet la representation probabiliste : v(t, x) = E_x [f (WT_t)exp f-- f _o^{T t} k (W_s.)d.s}

+ f ₀^T_t g(t + u,W_u)exp{-- t_o^u k (W_s.)d.s} du ] ;0 tT,x E IRa (3.7) où P_x est la loi du mouvement brownien W issu de x, et E_x designe l'esperance sous P_x.

Ce theorème donne en particulier un resultat d'unicite pour les solutions de (3.6), continues sur [0,T] x IR ^a et de classe C^1,2 sur [0,T) x IR a.

Remarque :

Si on n'est pas sous des hypothèses qui assurent l'existence d'une solution de classe ~^1,2 au problème de Cauchy, alors la fonction définie par le terme de droite de (3.7) peut ne pas être de classe C^1,2.

Les hypothèses sur les fonctions k,g et f peuvent être allégées. Il faut qu'elles assurent l'existence des termes stochastiques qu'on introduit dans la démonstration.

Un ensemble de conditions suffisantes pour que la solution v de classe C^1,2 existe et soit à croissance polynomiale est donné dans Friedman [16]. Si on a ellipticité uniforme de la matrice de diffusion, bornitude de a = ww^*, b, k, höldérianité de a, b, k et g, et croissance polynomiale de f et g, alors on a existence d'une telle fonction v, et ainsi existence, unicité, et une forme explicite de la solution v par le théorème précédent.

On déduit alors du théorème 2 une représentation d'une solution de classe C^1,2 de l'équation suivante :

at 2

+ kv = ¹Au + g , (t, x) E (0, oo) x ^d (3.8)

at

vérifiant la condition initiale :

u(0,x) = f(x) ,x E ~

Corollaire 1 :

Supposons que f: ^a --> , k: a --> [0, +oo) et g: ^a x [0, +oo) --> sont

continues bornées et que la fonction u, continue sur [0, +oo) x ^d et de classe ~1,2 sur [0, oo) x ^d, est solution de (3.8). Alors u admet la représentation stochastique :
u(t,x)= E_x [f (Wt)exPf--f o k (Ws)ds}

(3.9)

+ f _o^t g (t -- 0, We)exp f-- f _o ^e k (W_s)ds} 01;0 < t < oo,x E 2.3. Problème de Cauchy pour des opérateurs généraux :

On regarde ici la situation la plus générale que l'on peut obtenir par cette approche. On considère un temps arbitraire 7' > 0 , des fonctions k(t,x): d --> [0, +oo), g(t,x): ^a x [0, +oo) --> et f(x): a --> continues et telles qu'il existe L > 0 et ~ > 0 avec :

If(x)I L(1 + 11x11^2A) ou bien f(x) 0 ; x E: Rd

Ig(t, x)I L(1 + 11x11^2A) ou bien g(t, x) 0 ; x E: IV; 0 t T (3.10)

On note At l'operateur defini par la relation (3.2), b et a satisfaisant les hypothèses de lipschitzianite et de croissance lineaire, et soit X le processus solution de l'equation differentielle stochastique (1.33) associee. On note Et_,x l'esperance sous laquelle Xt = x. On a alors une forme generale du theorème de Feynman-Kac.

Théorème 4 :

Sous les hypothèses precedentes, soit v une fonction continue sur [0. T] x R^d et de classe C^1,2 sur [0, T) x R ^d, satisfaisant au problème de Cauchy :

v(T,

x) =

f(x)

,x E

R d

et à la condition de puissance polynomiale :

1 -- at + kv = ¹ A_tv + g , (t, x) E [0,T) x R d at (3.11)

max Iv(t, x)I M(1 + 11x11²) ;x E: R^d (3.12)

0.T

pour M > 0 et it 0. Alors v admet la representation probabiliste :

(3.13)

v(t, x) = Et, ^[f _(XT)exp f-- f r k (s^, Xs)ds}

+ f r g(u, X_u)exp{-- f tu k (s, X _s)ds} du ] ; 0 t T, x E IV^I Comme dans le cas precedent, la preuve repose sur la formule d'Itô appliquee à la semi-martingale :

s

v(t, x)expf-- 1 k (u, X_u)dul ; s E [t, T]

t

2.4. Interprétation des EDP de type elliptique :

Soit b de R^ddans R^det a de R^d dans les matrices d x d mesurables et bornees telles qu'il existe K > 0 tel que :

11b(x) -- b(Y)11 + 110(x) -- a(Y)11 k11x -- Y11,vx,Y E Rd

Les normes etant les mêmes qu'au theorème 16 (chap. I). On suppose que la matrice a = aa^T verifie :

3a> 0 / (a(x)Y,Y) a11Y11²

et que a est à derivees bornees.

W étant un P-Ye- processus de Wiener, on dénote par X^x la diffusion solution de l'EDS :

t t

X^x = x + 1 b(X_s.^x)ds + 1 _o-(Xs.x)dWs. , t 0

o o

dont le générateur L s'écrit :

L =

avec

d d

;114

1=1

)^.1

Pour une fonction c de R^d dans R continue, bornée et négative, on définit :

~

L = L + cl

a) Problème elliptique sans frontière :

Supposons de plus que c(x) -- /3 < 0 avec /3 assez grand pour qu'il existe 0 > 0 tel que :

(--Lu,u) 0iiuii² (3.14)

pour tout u dans H¹(R^d), L~ étant considéré comme un opérateur de H^l(R^d) dansH^-1(R^d).

Pour f continue bornée sur R^det de carré intégrable, on considère l'équation : Lu(x) = --f(x), x E R^d (3.15)
La solution de ce problème admet une représentation probabiliste suivante :

~ ~

u(x) = E 1 (Xexp {1

n c (X_s.^x)ds} dt , Vx E R^d (3.16)

o o

b) Problème elliptique avec condition de Dirichlet :

Soit l) un ouvert borné de R^d, à frontière al) assez régulière (C² par exemple). La condition de coercivité (1.14) devient :

(--Lu,u) eiiuii2 _a(D) (3.17)

pour tout u de Ha(D) (nulle au bord). L'inégalité (3.17) n'étant demandée que sur D, on peut encore affaiblir les hypothèses sur c en incluant notamment le cas L = A et c = 0. On dénote par T_x le temps d'atteinte du fermé D^C:

T_x = infft > 0: X^x D} (3.18)

Pour f continue sur D

= D U oD, on considère le problème : f Lu(x) = --f (x), x E D tu(x) = 0, x E oD (3.19)

dont l'unique solution u dans Ha(D) admet la représentation probabiliste suivante, à l'aide de la diffusion X^x:

~~ ~

u(x) = E 1 (Xexp {1

n c (X_s^x)ds}dt ,Vx E D~ (3.20)

~ ~

Cette formule se démontre comme (3.16) en remplaçant T par le temps d'arrêt T_x dont on montre que son espérance est finie. Si, au lieu de la condition de (3.19), nous considérons la condition inhomogène u(x) = (P(x) pour x E oD, il suffit d'ajouter le terme :

~

E ^[(P(X~x)exP {1 c (X_s^x)ds}1

~

à (3.20), pour obtenir la solution de (3.19) avec cette nouvelle condition au bord.

2.5. Interprétation des EDP de type parabolique :

On se replace ici dans le cadre des diffusions à coefficients dépendant du temps, c'est-à-dire que l'on se donne deux applications, b de I_E x R^V dans I^V et a de I_E x R^V dans les matrices d x d, que nous supposerons bornées donc vérifiée (avec a = aa^T); il existe K > 0 tel que pour tout t E It_Eet tout x, Y E I^V:

f IIb(t,x)II + IIa(t,x)II² K²(1 + IIxII²)

(3.21)

I IIb(t,x) -- b(t,Y)II + IIa(t,x) -- a(t,x)II KIIx -- YII

nous supposerons de plus que les dérivées spatiales de a sont bornées. A tout (t, x ) de IL_E x I^V, on associe le processus [X_s^x,t,s > t} solution de :

{

= _x

dX_s^t,x = b(s, X_s^t,x~ds + a(s, X_s^t,x~dW_s, s > t X,x

Pour T > 0 fixé, on se donne c, u~ et f tels que :

c E Cb([0, T] x 1^d) , c(t, x) 0

u⁷ E L²(I^d) n Cb(R^d)

f _E L²([0, T] x I^d) n Cb([0, T] x I^d) Le générateur de (X^x,t) est donné par :

d d d

1 L(t) = ₂11 _a^a_xi(ctii(t,x) ^a ) + cti(t, x) ^a

ax. axi

.1

1=1 j=1 i=1

d ,

octt -

Avec cti(t,x) = bi(t, x) -- 1 ax ^-I (t_' x)

-

i=1

on définit L_(t) par :

L_(t) = L(t) + c(t, x)I

a. Problème parabolique sans frontière :

On considère l'équation parabolique rétrograde (avec dormée finale en T ) sur IER^d :

u(T,x) =

11

(x),

x E Ir^l Ce problème a une unique solution :

l

au _

+ L(t)u = --f, 0 t T (3.22)

u dans L²([0,T]: H¹(^d)) n C([0,T]: L²(^d)) qui se représente de la manière suivante :

u(t, x) = E f r f (s, Xst,x)exp{f _t^s c (o, X,^t9^,x)d0} ds

+E [17(X:^,x)exp {f:^. c (s, X_s^t,x)ds}1

(3.23)

La démonstration se fait encore en appliquant la formule d'Itô à :

s

1p_s = u(t,X_s^t,x)exp ff co,X,^t9^,x)del

entre s = t et s = T, en prenant l'espérance et en utilisant (3.22).

Dans le cas homogène, où les coefficients (b, a, c, f ) ne dépendent pas de t, on peut poser v (t,x ) = u (T - t,x ) ; il résulte de (3.23) avec t = 0 que :

~

~~~, ~~ ~ ~ ~ ~~~, ~~ ~,~~~~~~~ ~

~ ~ ~~, ~~ ~,~~~~~ ~~

~

~3.24~

~~ ~~~~~~ ~,~~~~~ ~~ ~

~ ~ ~~, ~~ ~,~~~~~~

qui donne la solution au point (T, x) du problème :

{^ay = L(t)v + f , t > 0

(3.25)

v(0, x) = u(x)

Dans le cas f = 0, cette formule est appelée formule de Feynman-Kac énoncé précédemment ; dans le cas du mouvement brownien, en tombe dans l'exemple de l'équation de la chaleur.

~3.26~

Ce problème a une unique solution u dans :

L²([0, T]: Ha(D)) n c([0, T^.]: L²(D))

En utilisant le processus (X^x,t,s Odéfini précédemment, on introduit le temps

d'arrêt Tt_,x défini par :

b. Problème parabolique avec condition de Dirichlet :

Soit D un ouvert borné de R^d à frontière al) de classe C² .Avec les données du paragraphe précédent, on s'intéresse au problème :

~

~~ ~~ ~ ~~~~~~ ~ ~~, ~~~~ ~0, ~~ ~ ~
~~~, ~~ ~ 0, ~~~ ~0, ~~ ~ ~~
~~~, ~~ ~ ~~~~~, ~~~~ ~

T_x = inf{s > 0: X_s^t,x D} u(t, x) admet alors la représentation :

u(t, x) = E ~ ~~~, ~~

~~,~~~ ~,~~~~~~~ ~

~ ~ ~~, ~~ ~,~~~~~ ~~

~

~3.27~

~~ ~~~~~,~~~~~~~~~ ~,~~~~~ ~~ ~

~ ~ ~~, ~~ ~,~~~~~~

pour tout (t, x) E [0, T ] x D.

Cette formule se démontre comme (3.16) en appliquant ce coup-ci la formule d'Itô, de s = t à s = Tt_,x A T, à iY's pour une suite (u_n) régulière convergente vers u.

On peut ainsi considérer le problème non homogène où la deuxième condition de (3.26) est remplacée par u (t, x ) = (/) (t, x ) sur [0, T ] x aD ; il faut alors ajouter à (3.27) le terme :

_rTtx

E _{[II{Ttx<70(rt,x,X,t)exp} fic (0, X L^,ide

t ll

A nouveau, si (b, a, c, f , (/) ) ne dépendent pas de t, on peut poser v (t, x ) = u (T- t, x ) pour obtenir une représentation de la solution de :

av

= L(t)v + f sur [0, T ] x D

at

avec les conditions v = (/) sur ]0, T] x aD et v(0, x) = u(x) dans D.

3. Discrétisation de problème :

On a vu dans la partie précédente de ce chapitre que l'interprétation des EDPs conduit à écrire la solution sous forme E[f(X~)] avec la quantité f(Xt(w)) dépend de temps t et du hasard w. Le processus X est une diffusion, il n'est alors pas possible de calculer exactement E[f(X_~)] et l'on a naturellement recours à la simulation numérique. Tout d'abord une discrétisation temporelle de la dynamique de X permit de générer une variable aléatoire Xt dont la loi est proche de celle de Xt. On applique ensuite la méthode de Monte Carlo : la moyenne arithmétique de N copies indépendantes de la variable f(X^~_t) converge vers E[f(X^~_t)].

3.1. Le schéma d'Euler :

La méthode d'Euler qui représente la première étape dans la simulation, permet de donner une variable aléatoire X qui est proche en loi de X, où X = (X^x) est le processus de Markov solution de l'équation différentielle stochastique.

Soit T un horizon de temps fini, on peut sans perte de généralité prendre T = 1; on cherche à estimer la loi de Xf. = Xf qui en général n'est pas connue. Pour ce faire on approche X^xpar son schéma d'Euler d'ordre n > 1, disons Xn^,x défini de la façon suivante: on considère la subdivision régulière

y_n = (0 = ~_~ⁿ <
·
·
· < t_nⁿ _i < t_nⁿ = 1) de l'intervalle [0,1] c.-à-d. t_nⁿ = k /n et on pose X,^x = x, et pour tout k = (0,1 ...,n -- 1) et t E [tIⁿ, q₊1]

_

g n,x

_t -- gⁿr^,x + b (g _tk r^,x) (t -- 4¹) + a (g _tr^,x) (W_t -- Wtkn)

tk

Le schéma d'Euler est simulable. C'est une petite perturbation X^xque l'on peut expliquer de la façon suivante :

· On part à la date 0 de la valeur vraie x.

· Sur l'intervalle[0, tⁿ], on gèle les coefficients de l'EDS en leurs valeurs exactes b(x) et a (x) à gauche de cet intervalle et on calcule la valeur gⁿe en el^l de la

tl

solution de cette nouvelle petite EDS.

· Sur l'intervalle[t~ⁿ, tⁿ], on gèle les coefficients de l'EDS en leurs valeurs exactes b(Xⁿe) et CY aⁿe ) à gauche de cet intervalle et on calcule la valeur re en q de

tl tl t,

la solution de cette nouvelle petite EDS.

· Et ainsi de suite...

De plus on a les variable aléatoire W_til^._ct_+l -- W_tI_ct sont mutuellement indépendantes et

de même loi gaussienne centrée de variance q₊₁ -- ~~ⁿ.

3.2. La méthode de Monte Carlo :

On voit que les EDPs peuvent être interprétés à l'aide de processus de Markov bien choisis : on interprète u à l'aide du générateur infinitésimal du semi-groupe de transition d'un processus de Markov (Xt, t 0). Les motivations de cette démarche peuvent être d'ordre théorique et/ou numérique. En effet, en particulier lorsque (Xt, t 0) est solution d'une équation différentielle stochastique, le calcul stochastique permet parfois d'obtenir des résultats d'existence, d'unicité ou de régularité de la solution de l'EDP.

D'autre part, dès que l'on peut écrire la solution de l'EDP sous la forme d'une espérance du type u(t) = EF(X) avec F fonctionnelle sur l'espace des trajectoires de X entre 0 et t, on peut chercher à développer une méthode de Monte-Carlo pour approcher u(t) même si on ne sait pas simuler des trajectoires exactes de X : il suffit de construire un processus proche (en loi) de X, comme on a vu dans la paragraphe

précédente (schéma d'Euler), en simuler un grand nombre de trajectoires entre 0 et t, évaluer la fonctionnelle F le long de chaque trajectoire simulée et enfin moyenner toutes les valeurs obtenues.

Donnons un exemple élémentaire. Considérons l'équation de la chaleur :

at

(t, x) = vAu(t,x), t f(t,x) E ]0,7] x Ir^l (3.28)
at avec pour condition initiale u(0,
·) = uo(
·) une fonction mesurable bornée. Le

paramètre v est strictement positif et est appelé « paramètre de viscosité » en mécanique des fluides ou « volatilité » en finance.

On vérifie facilement d'après qui précédent que la fonction :

V(t,x) E ]0,7] x Ir^l, u(t,x) = Euo(x + V2vWt)

où (We) est un mouvement brownien standard à valeurs dans Ir^l, satisfait (3.28) ainsi que ¹imt,o u(t,x) = uo(x) en tout point de continuité de uo. Par application de la loi des grands nombres, on peut donc approcher u(t, x) par :

N

)

_N1 1

uo(x + V2vtgj(w)

j=1

où les (gj(w)} forment une famille de variables aléatoires gaussiennes indépendantes, à valeurs dans Ir^l, centrées et de matrice de covariance Id _Ird . Cet algorithme est très simple à mettre en oeuvre : on sait effectuer des tirages gaussiens indépendants à l'aide d'appels à un générateur de nombres pseudo-aléatoires uniformément répartis. La vitesse de convergence est décrite par des théorèmes- limite tels que le théorème de limite centrale, la loi du logarithme itéré, l'inégalité de Berry-Essen : la convergence est d'ordre 1/VN, elle est donc lente. Toutefois, le coût de l'algorithme croît seulement linéairement avec la dimension d de l'espace puisqu'on simule Nd trajectoires d'un mouvement brownien unidimensionnel standard, et ce coût est indépendant du paramètre v.

Typiquement, les méthodes de Monte-Carlo pour des équations aux dérivées partielles elliptiques ou paraboliques peuvent permettre de traiter des problèmes

extrêmes, en très grande dimension ou avec de très faibles viscosités, lorsqu'il serait difficile, ou extrêmement coûteux, d'utiliser des algorithmes classiques.

A l'aide de schéma d'Euler on peut écrire :

N

)

u(t,x)- 1

N Luo(~~r,x(wi)

i=1

CHAPITRE IV :

Application

0. Introduction:

On applique dans ce chapitre quelques formules et méthodes citées dans les trois chapitres précédents. Cette application repose sur la simulation numérique. Pour commencer nous abordons la simulation des trajectoires de mouvement brownien qui est la base des autres simulations stochastiques de notre travail. Il faut discrétisé le problème, cette discrétisation se faisant au niveau des variables, par exemple sur la variable de temps pour le mouvement Brownien.

Après la simulation du mouvement Brownien, on procède à la simulation de processus de diffusion, que l'on utilise dans notre problème principal concernent l'interprétation probabiliste des EDPs.

1. Discrétisation du mouvement brownien :

Le mouvement Brownien standard est une variable aléatoire W(t) qui dépend continument de temps (t E [0, T]) et satisfait les quartes conditions connues(Chap1).

Pour la computation proposée il est utile de considérer le mouvement Brownien discrétisé, où W(t) est spécifie aux valeurs discrètes de t; ainsi nous donnons At = T /N pour un certain entier positive N et soit W~ une notation de W(ti) avec ti = jLt.

La première condition est W_o = 0 avec une probabilité 1, la 2^eme et la 3^eme sont données par :

Wi = Wi + AWi , j = 1,2, ...,N,

où OW~ est une variable aléatoire indépendante de -VAt.7V(0,1).

Le programme (PROG1) en MATLAB nous donne la simulation d'une trajectoire du Mouvement Brownien standard sur l'intervalle [0,1] avec N = 1000 ; le résultat de la simulation est illustré dans la figure (Fig.01) qui représente 1000 points (ti ,Bi) en les joignant par interpolation linéaire:

Remarque : tous les algorithmes des programmes (en MATLAB) cités dans ce chapitre sont donnés dans l'annexe B.

W(t)

-0.2

-0.4

0.8

0.6

0.4

0.2

1.4

1.2

0

1

simulation d'une trajectoire d'un mouvement brownien

0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

Fig 01: t

-2.5 -2 -1.5 -1 -0.5 0 0.5 1 1.5

Fig 02:

1

0.8

0.6

0.4

0.2

0

-0.2

-0.4

-0.6

-0.8

simulation d'une trajectoire d'un mouvement brownien dans 2D

Dans la figure (Fig.2) on donne la représentation d'une trajectoire de mouvement Brownien dans deux dimensions. On voit ici que le point de trajectoire marche aléatoirement dans tous les sens du plan.

Avec les mêmes démarches suivies dans la simulation d'une trajectoire d'un mouvement Brownien on peut simuler M trajectoires, le programme (PROG2) sur MATLAB, nous a donné la simulation de M = 1000 trajectoires, le résultat est illustré dans la figure(Fig.03) :

-2

-4

-5

-1

4

0

3

5

2

1

-3

W(t)

0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

Fig.03: 1000 trajectoires de mouvement Brownien

t

2. Discrétisation d'un processus de diffusion:

2.1. Schéma d'Euler- Maruyama :

Une des simples méthodes de discrétisation de processus de diffusion est l'approximation d'Euler, appelée parfois Euler- Maruyama. On considère un

processus de diffusion X = (Xt, to t T} satisfait l'EDS :

dX_t = a( t, Xt)dt + b( t, X_t)dW_t (4.1)

dans [ to, T] avec la valeur initiale Xt_o = Xo.

Pour une discrétisation to = To
·
·
· TN = T de l'intervalle de temps [ to, T],

l'approximation d'Euler est un processus stochastique continu Y = (Y(t), to t T}
qui satisfait le schéma itératif :

Y_n+i = Y_n + c(r_n, Y_n)(r_n+i--T_n) + b(r_n, Y_n)(W_n+i--W_n) (4.2)

pour n = 0,1,2 ..., N -- 1 avec la valeur initiale Yo = X0, où nous avons écrit :

Yn = Y(rn)

Pour la valeur de l'approximation de temps de la discrétisation x_n, nous écrirons aussi : A_n= ;₁+1 -- x_n pour le n^tenie pas de temps et appelé S = max_n A_n le pas de temps maximale. Pour ce paragraphe on considère un temps de discrétisation équidistant x_n = to + nS avec S = A_nE A= (T -- to) /N pour N entier assez grand pour que DE [0,1] .

Quand le coefficient de diffusion est identiquement égale à zéro, l'itératif de schéma (4.2) réduit au schéma d'Euler déterministe pour l'équation différentielle ordinaire (EDO) :

dx

dt = a(t, x) (4.3)

La séquence (Y_n, n = 1,2., ..., N} des valeurs de l'approximation d'Euler (4.2) aux instants de temps discrétisé (x)8 = (r_n, n = 1,2, ..., N} peut avoir le même schéma

simulé dans le cas déterministe, la différence principale est que nous avons besoin à générer le pas aléatoire :

^AWn = Wrn-Ft--Wrn (4.4)

pour n = 0,1, 2 ..., N -- 1 de mouvement brownien W , on sait que ces pas sont des variables aléatoires gaussiens indépendants avec un moyenne E(AW_n) = 0 et de variance E(AW_n²) = A_n, donc on peut utiliser des nombres pseudo-aléatoires pour les incréments de mouvement brownien.

_n+i = Y_n + c/..A_n + b.AW_n

r

(4.5)

Y0 = X0

Pour simplifier nous écrivons : f = f(r_n,Y_n) pour une fonction définie sur IR⁺ x IR^l et n = 0, 1, 2 ..., N -- 1, alors nous pouvons écrire le schéma d'Euler (4.2) dans la forme abrévié :

La structure récurrente de schéma d'Euler qui évalue les valeurs approximatives du processus de diffusion aux instants de discrétisation est toujours la clé du Succès d'implémentation par ordinateur.

2.2. Exemple de discrétisation d'un processus de diffusion:

Nous illustrons l'aspect de simulation pour approximer le processus de diffusion, en examinant un simple exemple avec simple coefficient, on considère le processus X = (Xt, t 0) qui satisfait l'EDS linéaire :

dXt = aXtdt + bXtdWt (4.6)

Pour t E [0, 7] avec X0 E IR la valeur initiale, c'est un processus de diffusion avec drift a(t,x) = ax et un coefficient de diffusion b(t,x) = bx .

On peut avoir la solution exacte de cette EDS explicitement :

X_t = XoexpKa --¹1₂b²)t + bWt1 (4.7)

Pour t E [0, T] et W = (We, t 0) un Mouvement Brownien, connaître la solution de (4.6) explicitement donne la possibilité de comparer l'approximation d'Euler avec celle de la solution exacte, et calculer l'erreur.

Pour simuler la trajectoire de l'approximation d'Euler pour un temps de discrétisation, nous simplement donnons la valeur initiale Y0 = X0 et procéder par récurrence la prochaine valeur :

Y_n+i = Y_n + aY_nA_n + bY_nAW_n

pour n = 1,2, ... , accorder au schéma d'Euler avec drift et diffusion. Ici AW_n est la distribution gaussienne N(0, A_n) du mouvement brownien W sur le subintervalle x_n < t < T_n+1.

Pour la comparaison, nous pouvons utiliser (4.7) pour déterminer la valeur de la solution exacte pour un même échantillon de mouvement brownien, obtenu :

~

X,,, = Xoexp Ka -- 1 1₂ b²) r_n + b 1 Wi-il

~

Une application numérique de cet exemple en prenant l'intervalle de temps [0,1] , le pas A= 10^-2 et pour le processus de diffusion a = 0.5, b = 0.5 et X0 = 1. A travers le programme (PROG3) nous avons le résultat de la simulation de la solution exacte et la solution approximée par la méthode d'Euler de l'EDS (4.6). Le résultat est illustré dans la figure (Fig.04), où on peut voir que la trajectoire de solution approximée approche celle de la solution exacte.

solution exacte solution approximée

1.15

X(t)

1.1

1.05

1

0.95

0.9

0.85

0.8

0.75

0.7

0.65

0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

t

Fig.04:solution exacte et approximée d'un processus de diffusion

3. Application à l'interprétation probabiliste des EDPs :

Dans cette partie on applique la méthode probabiliste sur quelques exemples des EDPs paraboliques, Pour voir l'efficacité de la méthode on applique en parallèle à ces EDPs une méthode déterministe (différences finies), on veut dire par efficacité, estce-que la méthode nous a donné une approximation acceptable de la solution ? On commence par poser le premier problème :

3.1. Problème N°=1 (équation de la chaleur):

On prend comme premier problème à traiter l'équation de la chaleur. Nous choisissons un problème que nous savons résoudre analytiquement et par différences finies pour essayer par une comparaison de voir l'efficacité et les avantages de la méthode probabiliste.

L'équation s'écrit comme suit :

lau(x, t) a ²u(x, 0

= a dans D = [0,T] x [0, xf]

at ax² (4.8)

u(x, 0) =

f(x)

avec :

1 -- domaine D = [0,0.1] x [0,1] et a = 1.

conditions aux limites

u(1,)

= 0

condition initiale f (x) = sinn-x . _fu(0, 0 = 0 (4.9)

On sait résoudre cette équation analytiquement et la solution (exacte) est :

u(x, t) = sin(n-x)exp(--n-²t)

La figure (Fig.05) présente la solution analytique u(x, t) en fonction de t et x:

Chapitre IV

Application

0.6

0.08

0.04

t

0

0

1

0.8

0.6

0.4

0.2

0

1

0.8

u(xt)

0.1

0.06

0.02

x

0.4

0.2

Fig.05:solution analytique en fonction de x et de t

3.1.1. Illustration numérique du problème par méthode déterministe (la méthode des différences finies) :

Cherchons l'approximation numérique de la solution u(x, t) qui satisfait l'équation (4.8) avec la condition initiale et les conditions aux limites (4.9). La méthode des différences finies nous permet de donner cette approximation ; on commence par discrétiser le problème. On divise le domaine de l'espace [0,xf] à M subintervalles, chacun de longueur Ax = xf/M, et on divise le domaine de temps [0,T] à N subintervalles, chacun de durée At = T/M, en utilisant la méthode d'Euler explicite (chap. II) on obtient :

k

ui₊1 - 214 +14_1

Ax²

k+1 k

ui -ui

= (4.10)

At

a

avecui^k

r-r^, u(xi,tk), la méthode explicite nous permet de résoudre le problème par

itération :

k+1 At

u_i = r(u+₁ + ut_₁) + (1 - 2r)u^k-, avec r = a

ex²

On a d'après le chapitre II que ce schéma est consistent, donc il suffit de choisir ?t et ?x pour que la condition de stabilité soit vérifiée, de cette façon on obtient un

schéma convergent. La condition de stabilité dans ce cas est que r ~ ~ .

~

Le programme (PROG.4) fournit une fonction qui calcule la solution à tous les points (tk,xt) explicitement, avec condition initiale et condition aux limites quelconques noté dans le programme par f(x) , bx0 et bxf :

Le programme (PROG.5) définit les fonctions : condition initiale et conditions aux limites (4.9), le pas de discrétisation ?t et ?x . Il suffit de donner par exemple M = 20, N = 300, donc r = 0.133 , et on trouve la solution:

La figure (Fig.06) donne l'illustration graphique de la solution approximée en fonction de t et x, par une première comparaison on voit la grande ressemblance entre la figure (Fig.05) et (Fig.06) :

On continue la comparaison par le calcul de l'erreur de troncature :

E = MaxtE~^k,1 < i < M,1 < k < N1

E = 1.4797 x 10^-4

0.06

0.04

0.2

0.4

x

t

0.02

0

0

Fig.06:solution approximée par la méthode des difference finies

1

0.8

u(x,t)

0.6

0.4

0.2

0

1

0.8

0.6

0.1

0.08

Puisque on connaît la solution exacte de l'équation (4.8), alors on peut calculer le vecteur erreur défini par :

I 1 k 1

et = k_.ei ,
·
·
· , _et ,
·
·
· , et

t )

avec :

e^k = lu~^k -- u(xi,tk)1, 1 < i < M, 1 < k < N

où u(xi, tk) est la solution analytique (exacte) au point (xi, tk) et ui^kest la solution approximée par la méthode des différences finies au même point.

La figure (Fig.07) illustre la solution exacte et approximée au long de l'axe t et pour un point d'espace fixé x = 0.5. On remarque que les deux graphes sont presque confondus ce qui veut signifier que l'erreur est acceptable; ce que montre la figure (Fig.08) qui illustre le vecteur d'erreur e(t, x) pour x = 0.5.

430

0.9

0.8

0.7

0.6

0.5

0.4

1

solution approximée solution exacte

0 0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.08 0.09 0.1

t

Fig.07:solution exacte et approximée par la méthode des différences finies pour x=0.5

x 10^-4

e(t)

0.5

1.5

0

1

0 0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.08 0.09 0.1

t

Fig.08: L'erreur de la méthode pour x=0.5

On calcule l'erreur quadratique pour xi = 1 / 2 par :

i1 VN

_e quad _ (e c)² (4.11)

_

N Li k=1 ^E

dans ce cas squad = 1.1699 x 10^-4

La Table.1 donne l'erreur maximum (troncature) et le temps (7'E) d'exécution de programme pour des résolutions différentes de l'équation (4.8) par la méthode des différences finies en changeant les valeurs de N et M (toujours pour x = 0.5):

N	90	600	1000	2000	8000	1.8x10^4
M	20	20	50	100	200	300
At	0.0013	1.66x10^-4	10^-4	5x10^-5	1.25x10^-5	5.5x10^-6
Ax	0.05	0.05	0.02	0.01	0.005	0.0033
r	0.44	0.066	0.25	0.5	0.5	0.5
ema x	0.0013	4.54x10^-4	6.05x10^-5	6.05x10^-5	1.51x10^-5	6.72x10^-6
7'E(s)	0.031	0.25	0.656	4.82	144.7	1076.53

Table.1

On constate de cette table que : chaque fois qu'on diminue la valeur de Ax il faut diminuer la valeur de At pour obtenir la stabilité de solution, cette situation provoque un système linéaire de grand dimension qui nécessite un temps d'exécution plus grand, mais l'erreur devient plus petite.

3.1.2. Illustration numérique du problème par la méthode probabiliste :

On passe par suite à l'application de la deuxième méthode sur le même problème, donc on s'intéresse à l'interprétation probabiliste de la solution de l'équation (4.8) pour t E [0, 7'], ce problème est traité dans le chapitre III, et par l'application de la formule de Feynman-Kac on peut écrire la solution sous forme :

u(x, t) = E _x[f(X_~)] (4.12)

où le processus stochastique sous-jacent (Xt) est la solution de :

X_t = x + f _o^t V2 dW_s ,0 t 7' (4.13)

La solution numérique est obtenue par la simulation de Monte-Carlo; l'application de schéma d'Euler donne une solution approximée de l'équation (4.13). On simule N trajectoires de cette solution, et par la fonction de condition initiale on calcule les valeurs f(X~), en moyennant ces valeurs on trouve :

~

u(x,t) ~ 1 ~

~~~

L'illustration graphique de cette solution est donnée par la figure (Fig.09) :

0

0

0.08

t

1.2

0.8

0.6

0.4

0.2

0

-0.2

1

x 0.8

0.6

0.4

0.2

0.1

0.06

0.04

0.02

1

u(x,t)

Fig.09: solution approximée par la méthode probabiliste

Le programme (PROG.6) donne la solution approximé par la méthode probabiliste, où le processus (X1 ~) est discrétisé avec même pas de temps utilisé à la première méthode (D.F), c.-à-d. on prend n = 300 et donc ?t = T /n = 1 / 300.

La figure (Fig.10) illustre la solution probabiliste et la solution exacte au long de l'axe de temps et pour un point d'espace fixé x = 0.5 , on peut voir que la solution probabiliste approche de la solution exacte avec des perturbations.

U(tpX)

0.9

0.8

0.7

0.6

0.5

0.4

1

solution probabiliste solution exacte

0 0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.08 0.09 0.1

t

Fig.10:solution exacte et probabiliste pour x=0.5

On calcule l'erreur de troncature :

E = MaxtE ^k,1 k N)

E = 0.0227

L'erreur de troncature de cette méthode est plus grande que celle de la première. Ceci provient de l'ordre de convergence de la méthode de Monte Carlo.

Le vecteur d'erreur est donné par la formule :

e (xi, tk) = lUex(xi, tk) -- Uprob(xi, tk)I

où Uex est la solution exacte et Uprob est la solution par la méthode probabiliste.

La figure (Fig.11) illustre le vecteur e(x, t) pour une valeur de x fixé comme dans la première méthode x = 0.5.

e(tx)

0.025

0.015

0.005

0.02

0.01

0

0 0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.08 0.09 0.1

Fig.11:l'erreur de la méthode probabiliste pour x=0.5 t

On calcule l'erreur quadratique : s qUad = 0.0098

On remarque que l'erreur quadratique est supérieure que celle de la première méthode. Ce qui signifie que la solution approximée par méthode des différences finies est meilleure que celle de la méthode probabiliste.

On peut affaiblir la valeur d'erreur par l'augmentation de la valeur de N (le nombre des trajectoires). La Table.2 donne l'erreur maximale des résultats de simulation de solution pour des valeurs différentes de N (toujours pour x = 0.5 ) :

N	500	1000	2000	8000	10000
emax	0.0110	0.0131	0.0098	0.0089	0.0065
TE (s)	0.141	0.734	3.313	215.21	330.82

Table.2

D'après la Table.2 on constate que la convergence de cette méthode est très lente. Mais le temps d'exécution de programme est inferieur à celle de la première méthode, car dans la méthode probabiliste on peut trouver la solution pour un x fixé directement, contrairement à la première méthode il faut trouver toute la solution (c'est à dire la matrice u(t,x)) et en suite extraire la solution pour un x fixé.

Interprétation des résultats :

Malgré que l'approximation de la solution par la méthode des différences finies soit meilleure que celle de la méthode probabiliste, on peut dire que la méthode probabiliste donne une approximation acceptable de la solution.

Puisque la méthode est basée sur la simulation de Monte Carlo, la convergence de la méthode est lente, et l'approximation sera faible de celles des méthodes déterministes.

L'utilisation de cette méthode est avantagée dans le cas de résolution des EDPs où la résolution par les méthodes déterministes implique une résolution des systèmes linéaires à grandes dimensions.

3.2. Problème N° 2 :

On prend comme deuxième problème l'équation aux dérivées partielles suivante:

{ au 1 0²u 1 au 1

x² (x, t)--x (x, t) + u(x, t) = 0 , dans D = [0,T] x IR,

at (x, t) + 2 ax² 2 ax ' 2 ' (4.14)

u(x, 0) = f(x)

{-- domaine D = [0,1] x [0,1] .

condition initiale f (x) = x² .

(4.15)

)u(0,t = 0

conditions aux limites {u(1, u(1, t) = exp (-- t / 2 )

On peut résoudre cette équation analytiquement et la solution (exacte) est :

u(x, t) = x²exp(--t/2)

La figure (Fig.12) illustre la solution en fonction de x et t :

Chapitre IV

Application

0.4

0.2

0.2

0.8

t

0

0

1

0.8

0.6

0.4

0.2

0

1

x

0.8

0.6

1

u(x,t)

0.4

0.6

Fig.12:solution analytique(exacte) en fonction de x et t

3.2.1. Illustration numérique de problème par la méthode des différences finies :

Cherchons l'approximation de la solution par la méthode des différences finies, en utilisant la méthode explicite (chap. II). On obtient le schéma suivant :

1

ur = 0-1B i -- r2Ai)u^li,__i + (1 -- 2 At -- riBi + r2Ai) ut -- r2Atu^li,₊1

avec r1 = At/Ax , r2 = AtiAx², Bt_. = -- ¹₂ x i , At = ¹₂ xi² .

On donne dans le programme (PROG.7) une fonction qui fait cette discrétisation et donne la solution approximée, la figure (Fig.13) représente cette solution en fonction de x et t.

Pour étudier l'erreur de la méthode, on fixe un point d'espace, comme dans le premier problème pour x = 0.5, la figure (Fig.14) illustre la solution exacte et la solution approximée par la méthode des différences finies, et la figure (Fig.15) illustre l'erreur commise par la méthode, d'après cette erreur on déduit que l'approximation est bonne, avec une erreur maximale emax = 1.5169 X 10~~

87

0

0

Chapitre IV

Application

0.6

0.4

0.2

0.4

0.2

0.8

t

0.6

x

1

0.8

0.6

0.4

0.2

0

1

0.8

u(x,t)

Fig.13:solution approximée par la méthode des différences finies

u(x,t)

0.19

0.18

0.17

0.16

0.15

0.14

0.13

0.12

0.11

0.2

0.1

solution approximée solution exacte

0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 t 1

Fig.14:solution exacte et approximée par la méthode des différences finies pour x=0.5

1

e(x,t)

0.8

0.6

0.4

0.2

X 10^-4

1.6

1.4

1.2

0

1

Fig.15:le vecteur erreur pour x=0.5

0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 t 1

3.2.2. Illustration numérique de problème par la méthode probabiliste :

La formule de Feynman-Kac permet de représenter la solution de problème :

1

u(x, t) = E_x [exp (- ₂ t) f(X_~)i , t E [0, T] (16)

avec f(x) = x² et le processus sous-jacent:

Xt = X0 -

1 i^t t 2J dt + f Xs dWs , X0 = x o o

Le programme (PROG.8) donne l'approximation de la solution en passant par la méthode de Monte Carlo :

La figure (Fig.16) représente la solution exacte et la solution approximée par la méthode probabiliste et la figure (Fig.17) représente l'erreur commise par cette méthode, avec N = 100 trajectoires.

L'erreur maximale emax = 0.0094, mais pour la plupart des points, l'erreur est inferieure à 0.005, avec une erreur absolue maximale égale à 5.9%, qui laisse une bonne approximation.

Interprétation des résultats :

Dans le deuxième problème l'importance de la méthode probabiliste est bien apparue, car l'application de la méthode des différences finies sur cet exemple est difficile à cause des coefficients non constant, c'est-à-dire, dépendent de la variable de l'espace, cette dépendance rend la solution non stable parfois, la méthode probabiliste grâce à la formule Feynman-Kac nous donne une bonne représentation de la solution, et plus on augmente le nombre des trajectoires dans la méthode de Monte Carlo on obtient une solution plus précise avec une erreur faible.

u(x,t)

0.26

0.24

0.22

0.18

0.16

0.2

solution exacte solution apprximée

0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

t

0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

Fig.17: le vecteur erreur de la méthode

t

0.01

e(x,t)

0.009

0.008

0.007

0.006

0.005

0.004

0.003

0.002

0.001

0

Fig.16:solution exacte et solution approximée par la méthode probabiliste por x=0.5

CHAPITRE V :

Application en finance

0. Introduction:

Une option est un produit financier qui donne le droit à son acheteur d'effectuer une opération portant sur un sous-jacent (qui peut être une action, une obligation, une cargaison de pétrole) avant la date d'échéance.

1. Modèle du prix de l'actif:

La théorie des modèles de Black et scholes pour un actif est un processus stochastique (EDS):

dS = jiS_~dt + aS~dW~

La figure (Fig1) illustre ces échantillons d'un actif, avec it = 0.1 et 0- = 0.3.

A la date d'échéance t = T, le prix de l'actif est une variable aléatoire donné par la densité :

exp (--(~o9(x /S0 ) -- (¹ -- 0-²/ 2 )T)²)

20-²T

f (x) = ' , pour x > 0

x0-V 2n-T

avec f(x) = 0 pour x 0, pour confirmer on donne la figure (Fig.02), qui illustre l'histogramme où on prend la valeur finale de l'actif S(T) pour 10⁴ échantillons, et la courbe présente la densité f(x).

at)

2.4

2.2

0.8

0.6

0.4

1.8

1.6

1.4

1.2

2

1

0 10 20 30 40 50 60 70 80 90 t 100

Fig.01:100 trajectoire d'un actif

0. Formule de Black et Scholes :

La formule de Black-Scholes permet de calculer la valeur théorique d'une option à partir des cinq données suivantes :

· St: la valeur actuelle de l'action sous-jacente.

· t : le temps qui reste à l'option avant son échéance T.

· K: le prix d'exercice fixé par l'option.

· r : le taux d'intérêt sans risque.

· a : la volatilité du prix de l'action.

Le prix théorique d'une option d'achat (call), qui donne le droit mais pas l'obligation d'acheter l'actif S à la valeur K à la date T, est caractérisé par son payoff :

(ST -- K)⁺ = max(ST -- K; 0)

Le prix de l'option est donné par l'espérance sous probabilité risque neutre du payoff terminal actualisé :

C = E( payo f f X e^-rT)

De plus, pour employer le modèle du prix de l'actif, Black et scholes imposent un nombre de supposition de simplification aux options de marché, alors ils utilisent la formule suivante pour les valeurs de l'option européenne au temps t et le prix d'actif S :

C(S, t) = S.Ar(di) -- Ke^-r(T-~).Ar(d2)

où

1

log(S/K)+ (r + ₂ o^-2)(T -- t)

di =

et .Ar(.) est une fonction de la distribution de loi normale .Ar(0.1) :

.Ar(x) =

1 r s^2 e-^2 ds 112n-

Le paramètre r dans la formule est le taux d'intérêt. Si le prix d'actif aujourd'hui (temps zéro) est So la valeur de black et scholes call option est C(So, 0). Le programme (PROG2) fourni une fonction qui donne la valeur de l'option à la date t et le prix de l'actif S.

Il donne la valeur de call C = 1.1447 pour t = 0, quand on change les valeurs de t on obtient autres valeurs de call par exemple :

t = 1 --> C = 1.0972 ; t = 1.5 --> C = 1.0731 ..... a. La méthode de Monte Carlo :

Une propriété plus utilisée pour obtenir la valeur de l'option de Black et Scholes est d'estimer la moyenne de payoff, sous la condition de risque it = r, autrement dit, nous pouvons reproduire la valeur de l'option par fixer it = r dans le modèle de l'actif et computer la moyenne de payoff de tous les trajectoires de l'actif, dans la pratique, ceci peut être faire par la simulation de Monte Carlo (moyenner le payoff sur un grand nombre des trajectoires de l'actif). Pour le call européen on a seulement besoin de connaitre le prix de l'actif à l'échéance, ainsi on peut prendre At = T à chaque échantillon, et on donne le programme (PROG3) :

Ici, Pi est le payoff du i^eme trajectoire de l'actif. La moyenne totale P_moy est l'estimation de Monte Carlo de la valeur de l'option. La variance P~ar est utilisée pour donner l'intervalle de confiance à 95% suivant :

[P_moy -- 1.96 X P_var/M ; P_moy + 1.96 X P_var/M] Après l'exécution on trouve :

P_moy = 1.1453 avec l'intervalle de confiance égal [1.1435 ; 1.1471]

On rappelle qu'on a trouvé avec ces paramètres dans la première méthode que C = 1.1447 est approché à la valeur trouvée dans la deuxième méthode. La figure (Fig03) montre comment l'approximation de la méthode de Monte Carlo varie avec le nombre d'échantillon M, ici on prend S = 10, K = 9,r = 0.06, a = 0.1 et T = 1, les croix dans la figure représentent l'approximation de Monte Carlo et les lignes horizontales représentent l'intervalle de confiance de cette approximation, la valeur de Black et Scholes est représenté par la ligne verticale discontinu.

b. La méthode Binomiale :

On commence par simplifier le modèle du prix de l'actif, on discrétise l'intervalle de temps de façon équidistante 0 = to < ti <
·
·
· < tM = T avec ti_. = iAt. Donnons le prix de l'actif au temps zéro to, on pose que le prix de l'actif au temps ti

fait un mouvement vers le bas avec dS0 ou un mouvement vers le haut avec uSo, où d < 1 et u > 1, alors au temps t2 les mêmes mouvements bas/haut sont faits avec les trois possibilités d²S0, duS0 ou u²S0, et on continue cette restriction, au temps ti = iLt on obtient i + 1 possibilités donné par : S_nⁱ = d^i-nuⁿSo, 0 < n < M.

A la date d'échéance ti = tm = T, il y a M + 1 possibilités du prix de l'actif fS_n^m}_n^m _o. Soit f C_ni_n^m_ 0 correspond aux payoffs à la date d'échéance pour une option d'un call européen, nous savons que:

C_n^m = max(S_n^m -- K, 0) , 0 < n < M

La méthode binomiale procède par une marche arrière sur l'intervalle de temps. La valeur de l'option Ct, correspond au prix de l'actif Sn à la date ti est déterminée comme la moyenne de deux prix d'actif S_nⁱ⁺¹ et SI+i au temps ti₊1. La formule est :

C_nⁱ = e^-rAt(p.q⁺+¹1 + (1 -- p)S_nⁱ⁺¹ ) , 0 < n < i , 0 < i < M -- 1

Le paramètre p est considéré comme la probabilité d'un mouvement vers le haut du prix de l'actif, cette formule permet de marcher en arrière vers le temps zéro et déterminer la valeur de l'option C8. Il faut choisit les paramètres de la méthode telle que le modèle binomial de l'actif s'accorde avec la version de Black et Scholes quand At --> 0. Une fois At soit fixé, on aboutit à deux équation pour les trois paramètres restants, et par conséquence plusieurs possibilités des solutions, le choix par exemple est :

_erAt --d

d = A -- 11A² -- 1 , u = A + 1124² -- 1 , p = _{u --d} ,

où A = 2 ¹ ( e-rAt -- e(r+0-2)At). k

Le programme (PROG4) donne la valeur de l'option par la méthode binomiale, avec les mêmes paramètres utilisés dans les méthodes précédentes, la valeur approximative de l'option W = 1.1448 s'accorde avec la valeur de Black et Scholes C = 1.1447 .

c. EDP de Black et Scholes :

La formule de Black et Scholes pour la valeur d'une option européenne (call) intervient comme la solution d'une équation aux dérivées partielles (EDP). L'EDP est d'une forme parabolique avec les conditions de Dirichlet, soit x = T -- t la marche temps vers la date d'échéance, alors l'EDP prend la forme suivante:

²aC 1 a C as a-²S² rrS + C = 0 (5.1)

ax 2 aS² as

avec condition initiale :

C(S, 0) = max(S(0) -- K, 0)

et conditions aux limites :

C(0,x) = 0 , C(S, x ) r-r^, S -- Ke' , pour S grand

dans l'intervalle S > 0 et 0 < x < T.

On suit la même procédure du chapitre IV, pour résoudre cette EDP par les deux méthodes (différences finies et probabiliste) :

i. Différences finies pour l'EDP de Black et Scholes :

On subdivise l'intervalle de S , 0 S L, et on utilise le maillage de différences

finies {jh, ik} avec les pas h = L/N_S et k = T7N_t, la méthode explicite nous permet d'écrire l'équation (5.1) par :

V.i+1 -- Vi J J k

1 (V^i 1 --2VJ^i + V ) V^i -- Vi a2(jh)2 J+^-1-- r(jh) J+1 J1 + rV^i = 0 2 h^2 h J

Avec VJⁱ r-r^, C(jh,ik)

où par la formule matricielle :

_yi+1 = Fyi + pi, pour 0 i Nt -- 1,

où F E IR(Ns-1)x(Ns-1) est tridiagonale et pⁱ E 1^(Ns-1) est déterminé par les conditions aux limites.

Le programme (PROG5) donne les valeurs de call V à la date t et au prix S, avec les paramètres T = 1, L = 10, K = 4,r = 0.03, a = 0.5.

Application en finance

Chapitre V

10

0.6

6

0.8

4

0.4

t

8

S

La figure (Fig.4) illustre le résultat de ce programme :

10

8

6

4

2

0

0

0.2

XS t)

1 2

0

Fig.04:valeurs d'une option européenne(call) par la méthode des différences finies

ii. Méthode probabiliste l'EDP de Black et Scholes :

L'interprétation probabiliste de l'EDP de Black et Scholes peut donner la valeur de l'option avec les mêmes paramètres précédents, le programme (PROG6) fait cette interprétation. On donne premièrement le processus sous-jacent :

t t

S(t) = So + f rS(s)ds + f o^-S(s)dW_s

o o

et la solution sera :

C(S, t) = E((S(t) -- K, 0)⁺exp (r(7' -- t))) Les valeurs sont données dans la figure (Fig.05):

Application en finance

Chapitre V

0

0

10

8

6

4

2

0

1

C(S,t)

0.8

10

0.6

t

S

8

0.4

6

0.2

4

2

Fig.05:valeurs d'une option européenne (call) par la méthode probabiliste

conclusion

CONCLUSION

Dans ce travail nous présentons une des utilités du mouvement Brownien, qui est l'interprétation probabiliste des EDPs. A partir des propriétés qui sont concrétisées soit dans la notion de mouvement brownien soit par la généralisation de cette notion où on trouve la notion de processus stochastique de diffusion, la formule de Feynman-Kac nous permet de donner la représentation de la solution cherchée.

Pour bien présenter l'utilité de la méthode probabiliste on fait une application sur un exemple d'EDP parabolique qui nous savons résoudre analytiquement et par la méthode des différences finies, la comparaison fait entre les solutions approximatives de chaque méthode et la solution exacte montre qu'on peut adopter la méthode probabiliste comme une méthode de résolution, mais les avantages de cette méthode apparus dans les différents cas où on veut résoudre un problème en dimension élevée (par exemple supérieure à 4) les méthode classique conduisent à l'invention de systèmes linéaires d'une taille telle qu'elles deviennent impraticable, et la méthode probabiliste est souvent utilisée.

De même, la méthode probabiliste est souvent préférable lorsque l'on cherche les valeurs de la solution en certains points du domaine de calcul seulement : le cas des calculs de prix d'option en finance est typique puisque l'on ne s'intéresse qu'à une ou à quelques valeurs des prix.

ANNEXE

A. Principaux espaces fonctionnels :

1. Espaces de fonctions « régulières » (au moins continues)

Soit 12 un ouvert de IV , K un compact de Ir. Les fonctions sont ici à valeurs réelles ou complexes.

C(12) Espace des fonctions continues sur 12. Muni de la famille de

semi- normes : pK(u) = sup_xEK|u(X)| pour tout
compact K c 12 , c'est un espace de Fréchet.

Cb(12) Espace des fonctions continues bornées sur 12. C'est un

espace de Banach pour la norme : llull = sup_xEK|u(X)|.

Co(12) Espace des fonctions continues sur 12 et tendant vers zéro

au bord de 12. C'est un espace de Banach pour la norme précédente.

C_c(12) Espace des fonctions continues sur 12 à support compact

dans 12.

C^k(12), k E N Espace des fonctions de classe C^k, ou encore espace des

fonctions dont toutes les dérivées d'ordre < k existent et sont continues. C'est un espace de Fréchet pour la famille de semi-normes :

sup

xEK

llullk,K = sup

|p|k

C"(12) espace des fonctions infiniment dérivables ; espace de

Fréchet, avec les semi-normes :

llflla,K = sup

|D^a f(X)| , a E N^n, K c 12

xEK

Co"(12) = D(12) espace des fonctions infiniment dérivables à support compact. Une suite (fk) de D(12) tend vers zéro dans D(12) si U supp fk = K est borné et :

sup |D^a fk(X)| --> 0, tea E Nn

xEK

2. Espaces de fonctions intégrables :

Soit p un réel,p 0 , et encore 12 ouvert de Ir.

Les fonctions considérées étant à valeurs réelles ou complexes, on note :

L^P(12) l'espace des (classes de) fonctions mesurables sur 12 telles que la

fonctionx E 12 --> |f(x)|^P soit intégrable sur 12 . C'est un espace de Banach pour la norme :

₁¹/_P

IlfIlP = [f |f(x)|^Pdx

~

3. Espaces de Sobolev : Soit m E N, p E 41, +0(4, s E ~.on note :

H^m(12) l'ensemble des fonctions f telles que fEL^P(12) avec

D^a f E L^P(12) , Va E Nⁿ vérifiant |a| = (a1 +
·
·
· + a_n m),

avec, bien sûr, D^a fdéfini au sens des distributions.

C'est un espace de Hilbert pour le produit scalaire :

(f, 9) = 1 f D^a f(x)
· Da ~~~~~~~~~~ 9(x) dx

|ct|m ⁿ

H~^m(12) l'adhérence D(12) de dans H^m(12).

H^-m(12) l'espace dual de H^m(12). C'est un espace de Hilbert pour la

norme duale.

Algorithme 01 : simulation d'un mouvement brownien Variables {Déclaration des donné}

[to, T] : tableau des {réels L'intervalle de temps} N : entier {la taille de trajectoire}.

Début

dt = T/N {Discrétisé l'intervalle de temps} Wo = 0 {Initialisation de trajectoire}

Pour i = 1 jusqu'à N {La boucle de simulation}

Simuler dW(i) par la loi N(0,1) W(i) = W(i -- 1) + -Vdt dW(i)

Finpour

Fin.

B. Algorithmes :

Algorithme 02: simulation de M trajectoire d'un mouvement brownien Variable {Déclaration des donné}

[to, T] : tableau des réels {l'intervalle de temps}

N : entier {la taille de trajectoire}

M : entier {nombre des trajectoires}

Début

dt = T/N {Discrétisé l'intervalle de temps}

Wo(i) = 0 pour 1 i M {Initialisation des trajectoires}

{ La boucle de simulation} Pour j = 2 jusqu'à N

Pour i = 1 jusqu'à M

Simuler dW(i, j) par la loi N(0,1)

W(i,j) = W(i,j -- 1) + -Idt dW(i, j)

Finpour

Finpour

Fin.

Algorithme 03: simulation d'une diffusion Variable

[t0, T] : tableau des réels {l'intervalle de temps}

N : entier {la taille de trajectoire}

a, b , X0 : réels {Constantes}

Début

dt = T/N {Discrétisé l'intervalle de temps}

{La boucle de simulation de la solution exacte X_ex} X_ex(1) = X0 {Initialisation de trajectoire}

Pour j = 1 jusqu'à N

Simuler dW(j) par la loi N(0,1)

W(j) = W(j -- 1) + Vdt dW(j)

Xex(j) = XoexP [(a -- ₂¹ b) + (bW(j))1

Finpour

{La boucle de simulation de la solution approximée par la méthode

d'Euler-Maruyama Xapp}

Xapp(¹) = X0 {Initialisation de trajectoire}

Pour j = 1 jusqu'à N

Xapp(j) = Xapp(j -- 1) + a Xapp(j -- 1)dt + bXapp(j -- 1)W(j)

Finpour

Fin.

Algorithme 4 : résolution d'une EDP parabolique (équation de la chaleur)
par la méthode des différences finies

Variable {Déclaration des donné}

[to, T] : tableau des réels {l'intervalle de temps} [xo,xf] : tableau des réels {l'intervalle de l'espace} M , N : entier {la taille de trajectoire}

f : fonction {Condition initiale; fonction de }.
b_x,b_xf : fonction {Conditions aux limites ; fonction de t } a : réel {constant}

Début

dt = T/N {Discrétisé l'intervalle de temps}

dx = xf/M {Discrétisé l'intervalle de temps par}

dt ₂1

si ( r = a dx2 < ) alors {Condition de stabilité}

{La boucle de la solution aux points (to, x(i)) avec (i) E [xo,xf] } Pour i = 1 jusqu'à M + 1

u(to^,x(i)) = f(x(i))

Finpour.

{La boucle de la solution aux points(t(k), xo) et (t(k),xf) } Pour k = 1 jusqu'à N + 1

u(t(k), xo) = bx(t(k)) u(t(k), xf) = b_xf(t(k))

Fin.

{La boucle de la solution aux points(t(k), x(i)) }

Pour k = 1 jusqu'à N

Pour i = 2 jusqu'à M

u(t(k), x(i)) = r [u(t(k), x(i + 1)) + u(t(k), x(i -- 1))]

+(1 -- 2r)u(t(k), x(i))

Finpour Finpour

Finsi

Fin.

Algorithme 5 : résolution d'une EDP parabolique (équation de la chaleur)

par la méthode probabiliste.

Variable

[t0, T] : tableau des réels {l'intervalle de temps} n : entier {la taille de trajectoire}

x0 : réel {La valeur initiale}

N : entier { La taille de trajectoire}.

Début

dt = T/N {Discrétisé l'intervalle de temps}

algorithme 1 {Simulation de mouvement Brownien}

{La boucle de simulation de la solution au point (t(k), x0) } Pour k = 1 jusqu'à N

Pour i = 1 jusqu'à n

X(t(k), i) = sin(n^- (x0 + Ai2 W(k, i)) Finpour

Finpour

Pour k = 1 jusqu'à N

Fin

u(t(k), x0) = lsi E ls_ 1 X (t (k), i) {L'approximation de Monte Carlo} Finpour

BIBLIOGRAPHIE:

[1] Søren Asmussen ,Peter W. Glynn. Stochastic Simulation: Algorithms and Analysis. Springer Science+Business Media, LLC 2007.

[2] V. Bally, D. Talay. The law of the Euler scheme for stochastic Differential Equations (I): Convergence rate of the distribution function. Probability Theory and Related Fields, 104 (1), 1996.

[3] V. Bally, D. Talay. The law of the Euler scheme for stochastic Differential Equations (II) : convergence rate of the density. Monte Carlo Metods and Applications, 2:93-128, 1996.

[4] Timothy J. Barth, Michael Griebel, David E.Keyes, Risto M. Nieminen, Dirk Roose, Tamar Schlick. Numerical Solution of Partial Differential Equations on Parallel Computers. SpringerVerlag Berlin Heidelberg, Printed in The Netherlands. 2006.

[5] Fred Espen Benth , Giulia Di Nunno ,Tom Lindstrøm, Bernt Øksendal,Tusheng Zhang. Stochastic Analysis and Applications. Springer-Verlag Berlin Heidelberg 2007.

[6] N. Bouleau. Processus stochastiques et Applications. Hermann, 2000.

[7] Z. Brzezniak, T. Zastawniak. Basic Stochastic Processes. A course through exircises. Springer-Verlag London, printed in great Britain. 1999.

[8] D. Dacunha-Castelle, M. Du. o. Probabilités et statistiques, tome 1, Problèmes à temps fixe. Masson, 1982.

[9] R. Dautray et al. Méthodes probabilistes pour les équations de la physique. Eyrolles, 1989.

[10] DEBEAUMARCHÉ (G.). - Introduction aux équations aux dérivées partielles linéaires. Techniques de l'Ingénieur, traité Sciences fondamentales AF 162 p. 1-21 (1999).

[11] J.F .Delmas, B. Lapeyre. Processus Aléatoires. CERMICS, 2000.

[12] R. M. DUDLEY. Real analysis and probability. Cambridge University Press 2004.

[13] Vincenzo Capasso, David Bakstein. An Introduction to Continuous-Time Stochastic Processes, Theory, Models, and Applications to Finance, Biology, and Medicine. Birkhauser Boston.2005.

[14] YuX Egorov, M. A. Shubin (Eds.). Partial Differential Equations VI, Elliptic and Parabolic Operators. Springer-Verlag Berlin Heidelberg.

[15] Jean-Pierre FOUQUE. Relations entre probabilités et équations aux dérivées partielles. Techniques de l'Ingénieur, traité Sciences fondamentales AF 565.

[16] A, Friedman. Stochastic differential equations and application, vol1 and vol2. Academic press.1975,1976.

[17] J. Gani, C.C. Heyde, P. Jagers, T.G. Kurtz. Probability and Its Applications. SpringerVerlag London Limited 2008.

[18] Gary.M. Lieberman. second order parabolic differential equations. By world scientific publishing 1996.

[19] I. I. Gihman, A. V. Skorohod. Stochastic differential equations. Springer-Verlag, 1972.

[20] Vincent Guinot, Bernard Cappelaere. Méthodes Numériques Appliquées, Résolution numérique des équations différentielles de l'ingénieur. Polytech'Montpellier STE 2. 2005/2006.

[21] Desmond. j. Higham, Nicholas. j. Higham. MATLAB Guide. second edition. university of Strachclyde, Elasgow, Scotland. University of Manchester, Manchester, England. 2005.

[22] Franck Jedrzejewski. introduction aux méthodes numériques, deuxième édition. springer-verlag france 2001 pour la 1ère édition.

[23] Yu. Kabanov, R. Lipster, J. Stoyanov. From Stochastic Calculus to Mathematical Finance, The Shiryaev Festschrift. Springer-Verlag Berlin Heidelberg 2006.

[24] I. Karatzas, S. Shreve. Brownian Motion and Stochastic Calculus. Springer Verlag, New-York, 1988.

[25] Jaan Kiusalaas. numerical methods in engineering with MATLAB. The Pennsylvania State University. cambridge university press. 2005.

[26] Andrew Knight. Basics of Matlab and Beyond. Boca Raton London New York Washington D. C. Chapman & hall/crc. 2000.

[27] Fima c Klebaner. Introduction to stochastic calculus with applications, second edition. Monash University, Australia. Imperial college press,2005.

[28] P. E. Kloeden, E. Platen. Numerical Solution of Stochastic Differential Equations. Springer Verlag, 1999.

[29] A. G. Malliaris, W. A. Brock. Stochastic methods in economics and finance. Elsiever Science 1982.

[30] John H.Mathews, Kurtis D.Fink. Numerical methods using MATLAB. Prentice Hall.1999.

[31] Sylvie MÉLÉARD. Mouvement brownien et calcul stochastique. Techniques de l'Ingénieur, traité Sciences fondamentales AF 566

[32] K. W. Morton, D.F. Mayers. Numerical Solution of Partial Differential Equations, An Introduction. Cambridge University Press 2005.

[33] D. Lamberton, B. Lapeyre. Probabiltés numériques et théorie des options. In N. Bouleau and D. Talay, editors, Probabilités numériques, Collection Didactique, 180-186, INRIA,1986.105.

[34] David P. Landau, Kurt Binder. A Guide to Monte Carlo Simulations in Statistical. Second Edition. Cambridge university press.2005.

[35] B. Lapeyre, E. Pardoux, R.Sentis. Méthodes de Monte-Carlo pour les équations de transport et de diffusion. Mathématiques et applications 29, Springer-Verlag, 1997.

[36] Mario Lefebvre. Applied Stochastic Processes. springer 2007.

[37] Don s. Lemons. An Introduction to Stochastic Processes in Physics. the johns hopkins university press baltimore and london.2002.

[38] O. Lévêque. Cours de Probabilités et Calcul Stochastique. EPFL, 2004.

[39] M. Loéve . Porbability theory I,4^th edition. springer-verlag.1977.

[40] B. Øksendal. Stochastic differential equations, an introduction with application. Springer- Verlag, 2000.

[41] P. Protter. A partial introduction to financial asset pricing theory. Stochastic Processes and their Applications 91 , 2001, 169-203.

[42] L.C.G. Rogers, David Williams. Diffusion, Markov processes and martingales , volume 2.

₂ed

Ito calculus edition. Cambridge university press.2000.

[43] Steven E. Shreve. Stochastic calculus for finance II, continuous-time models. Springer Science+Business Media,Inc 2004.

[44] Steven E. Shreve. Stochastic calculus for finance I, The Binomial asset pricing model. Springer-verlag new York.LLC 2004.

[45] Galen R. Shorack. Probability for Statisticians. springer-verlag New York.2000.

[46] Dieter Sondermann. Introduction to Stochastic Calculus for Finance, A New Didactic Approach. Springer-Verlag Berlin Heidelberg 2006.

[47] Pierre SPITERI. Méthode des différences finies pour les EDP stationnaires. Techniques de l'Ingénieur, traité Sciences fondamentales AF 500.1999.

[48] David Stirzaker. Stochastic Processes and Models. Oxford University Press, 2005.

[49] D. Talay. Probabilistic Numerical Methods for Partial Differential Equations: Elements of Analysis. Probabilistic models for nonlinear partial differential equations, Lecture notes in mathematics 1927, 1996.

[50] D. Talay, E. Peirano. Domain decomposition by stochastic methods. Fourteenth Internatinal Conference on Domain Decomposition methods.2003.

[51] Albert Tarantola. Probability And Measurements. Université de Paris, Institut de Physique du Globe 4, place Jussieu; 75005 Paris; France. December 3, 2001.

[52] Henk C. Tijms. A First Course in Stochastic Models. Vrije Universiteit, Amsterdam, the Netherlands. John Wiley & Sons Ltd. 2003.

[53] Won Young Yang. Wenwu Cao.Tae-Sang Chung. Applied numerical methods using Matlab. John Morris Chung-Ang University, Korea. Pennsylvania State University. The University of Auckland, New Zealand.2005.

[54] B. Ycart. Introduction aux équations différentielles stochastiques. LMC-IMAG, 1998.