FR 2734069 A1 19961115 FR 9506085 A 19950512 La présente invention concerne une méthode pour prédire, par une technique d'inversion de scénarios, l'évolution de la production d'un gisement souterrain et notamment d'un gisement recélant des hydrocarbures. ETAT DE LA TECHNIQUE: Les techniques d'inversion sont utilisées de façon extensive dans le domaine de l'ingénierie de gisement. On peut citer, à titre d'exemple, les domaines d'utilisation usuels suivants: - Au stade du laboratoire, on s'en sert pour déterminer différents paramètres représentatifs du comportement des roches vis-à-vis des fluides. Ces paramètres peuvent être, par exemple: la perméabilité absolue et relative, les courbes de capillarité etc. -Les techniques d'inversion sont également couramment utilisées pour interpréter des tests de puits. Les paramètres de l'inversion sont alors, par exemple, la perméabilité d'un ou plusieurs faciès, les limites géométriques d'une structure géologique, les coefficients affectant l'indice de productivité d'un puits, etc. -Pour l'étude d'un gisement, l'inversion est utilisée pour caler la réponse d'un simulateur numérique sur les mesures disponibles de production (ou "historique de production"). Les paramètres peuvent être, par exemple, la porosité des roches, leurs perméabilités absolue et relative, les indices de productivité des puits, etc. Technique de calage de l'historique de production d'un gisement Les techniques d'inversion constituent un domaine de recherche très actif en simulation de gisements. Ce processus de calage d'un historique de production par ajustement des paramètres d'un modèle de siumation est également désigné par les spécialistes par "History Matching". Il consiste essentiellement à trouver un ensemble de paramètres régissant les équations d'écoulement des fluides, qui, intégrées dans un simulateur numérique, permettent de retrouver les données observées ou relevées. Les données observées à chaque puits sont, à titre d'exemple, principalement la pression, la composition des fluides ou le débit des différentes phases. On peut utiliser aussi toute combinaison de ces variables et notamment la fraction d'eau ("water cut"), le G.O.R ("Gas Oil Ratio"), etc.Les solutions au problème inverse, permettant de reconstituer l'historique de production d'un gisement par un processus d'ajustement, peuvent, d'une façon générale, être multiples. Il importe donc d'incorporer les connaissances géologiques disponibles de façon à restreindre l'espace des solutions possibles. La prise en compte de cette connaissance géologique initiale rend plus réalistes et plus fiables les prévisions de production à partir d'une simulation. Pour construire le modèle initial, on y incorpore toute l'information disponible: données brutes ou interprétées, études géologiques, mesures sismiques, etc. La connaissance physique d'un gisement peut être intégrée en considérant par exemple: - la structure des unités sédimentologiques; - les limites de variations des valeurs pétrophysiques (porosité, perméabilités etc.) associées aux lithofaciès; et - l'information statistique portant sur les moyennes, les écarts-types, la corrélation spatiale, etc. Le procédé dit de calage d'historique comporte classiquement les étapes suivantes: a) Un modèle de simulation est construit en se basant sur la connaissance géologique initiale et en intégrant le maximul de données disponibles. b) Des paramètres significatifs du modèle sont sélectionnés pour le processus d'inversion en considérant: la connaissance du comportement du gisement en fonction de ces paramètres, leur influence qualitative sur la production, les incertitudes initiales associées à ces paramètres, etc. c) on ajuste avec les paramètres du modèle pour restituer l'historique de production ou les données observées. Un ensemble de valeurs de paramètres étant donné, une simulation directe permet de comparer les résultats prévus avec les observations. La méthode la plus couramment utilisée est un processus d'essais et 3 2734069 d'erreurs, dans lequel intervient l'expérience de l'ingénieur de réservoir: les valeurs des paramètres sont ajustées en fonction de la connaissance du gisement et de la compréhension de son comportement dynamique. Il est parfois possible d'accélérer cette étape d'ajustement par une procédure automatique permettant une modification itérative de la valeur des paramètres pour obtenir un meilleur accord entre le résultat du calcul et les données d'observation. L'inconvénient de cette méthode classique, en grande partie empirique, est de ne pas directement quantifier l'influence de la connaissance initiale du gisement dans le processus d'inversion. Il est en effet possible de caler un historique de production tout en s'éloignant des grandeurs physiques admissibles. De plus, les incertitudes sur les prévisions de production ne sont pas directement quantifiées durant le processus d'inversion. Des études de sensibilité sont parfois menées sur les valeurs des paramètres obtenus par l'ajustement, pour appréhender leur influence sur les prévisions, mais ces études ne font pas partie intégrante du processus d'inversion. Le processus de calage d'historique peut être formalisé de la manière suivante: On modélise les phénomènes physiques par un modèle mathématique désigné par g, correspondant ici au simulateur des écoulements de fluides dans le gisement. Ce modèle dépend de np paramètres (vecteur x), et du temps t. Les phénomènes physiques modélisés peuvent alors être écrits sous la forme d = g(.,t). Si.dm représente un ensemble de nm observations correspondant aux mesures de d à différents temps ti, on peut alors écrire la condition suivante pour que le modèle restitue les observations: dm = {g(x,ti), i= 1, nm} Pour définir un critère d'adaptation du modèle aux observations, on peut définir une fonctionnelle, par exemple sous la forme classique de type "moindres carrés", avec l'expression suivante: Fm(.) = 1/2 Zi= 1,nm {( di - g(x, ti))2} Les paramètres à retenir sont ceux qui minimisent la fonctionnelle Fm(x) en tenant compte des contraintes imposées aux paramètres. Il faut noter que la qualité des prévisions permises par cette procédure peut être parfois sujette à caution pour deux raisons principales. La première c'est que la connaissance initiale ou a priori sur le réservoir n'est pas quantifiée si bien qu'il est possible qu'elle disparaisse durant l'opération de minimisation. La deuxième est qu'il peut y avoir plusieurs solutions qui minimisent la fonctionnelle. Aussi, pour accroître la qualité des prévisions, il est important de prendre en compte de l'information disponible sur le modèle géologique ainsi que des statistiques éventuelles caractérisant ces paramètres, durant l'opération d'ajustement. Technique d'inversion de type Bayesien d'un modèle de simulation Pour prendre en compte cette information initiale, on peut utiliser un formalisme d'inversion connu des spécialistes sous le nom "d'inversion Bayesienne", et décrit par exemple par: Floris F.J.T. et al, in "Flow Constrained Reservoir Characterization using Bayesian Inversion", 4th ECMOR, Norway, 7-lO june 1994. Ce formalisme connu des spécialistes permet de prendre en compte des informations a priori sur les paramètres sous forme de modèles probabilistes, les erreurs faites sur les mesures avec leurs probabilistes associées, et aussi les erreurs dues à l'utilisation d'un simulateur numérique pour modéliser les écoulements. Avec les moyens de calcul couramment utilisés, le formalisme Bayesien est malheureusement très coûteux en temps de calcul pour être mis en oeuvre et il est peu pratique quand on s'intéresse à la totalité d'un champ, car il impose que le simulateur d'écoulements "balaye" l'ensembles de paramètres possibles pour évaluer la probabilité de chaque solution. Le principe de base de cette inversion Bayesienne appliquée à la simulation des écoulements dans un gisement, est de modifier le modèle géologique initial en utilisant des données de production dans le but d'améliorer la qualité de ce modèle initial. Les principales étapes en sont: a) La construction d'un modèle géologique initial. f s'agit d'un modèle "a priori" de type probabiliste. Le gisement est décrit sous la forme d'un jeu de paramètres avec leurs densités de probabilité associées p(x). Des valeurs moyennes pour ces paramètres a priori (notée Xa) et des écarts-types (notés ax) peuvent être utilisés pour modéliser cette fonction de densité probabiliste (ou "probability density function", notée pdf. Avec des incertitudes de type gaussiennes, modélisées par un opérateur de covariance noté Cx, la fonction pdf a priori a pour expression: p(x) = Cte exp {- (a - _x)T Cx-1 (xa - x)} Pour exprimer cette information a priori dans une fonction objectif, il faut rajouter au terme précédemment défini Fm(x) un nouveau terme, noté Fx(x), pouvant s'exprimer ainsi dans le cas d'incertitudes gaussiennes Fx(x)): 1/2 {-x(a - x)T Cx-1 (La - x) b) l'utilisation des données observées: Les données d'observation sont utilisées pour le calcul de la densité de probabilité a posteriori p(x/d=dm) en accord avec la règle de Bayes: p(/d=dnm) = (p(d=dm/Jx). p(x)) / p(d=dm), o p(_=_dm/x) est la fonction de probabilité ("likelihood function" en notation anglo-saxonnes) qui quantifie les différences entre les valeurs observées dm et les valeurs calculées d pour un jeu donné de paramètres x. Quand les incertitudes sur les mesures et sur le modèle sont gaussiennes, et modélisées par un opérateur de covariance Cm, cette fonction de probabilité a pour expression: p(d=4nm/x) = Cte exp{m-(dm-d)T Cm-1 (dmrn-d)} Pour obtenir une interprétation statistique en termes de densités de probabilité pdf, il faut dans la pratique, déterminer la fonction a posteriori pour la gamme complète de paramètres. Si on s'intéresse uniquement au maximum de probabilité a prosteriori, les paramètres correspondant à cet optimum sont ceux qui minimisent la fonction objectif F=Fm () + Fx(). En résumé, pour chaque jeu de paramètres, la probabilité associée p(x_) est accrue quand la réponse du modèle g(x,t) fournit des valeurs proches de celles dm observées, c'est-à-dire quand la fonction de probabilité p(=-nm/Xj) est grande. La méthode classique de calage d'historique et la technique d'inversion Bayésienne, rappelées et discutées ci-dessus, sont utilisées pour adapter un modèle de simulation aux données réelles mesurables. Les inconvénients majeurs de la première méthode sont la possibilité de perdre une partie de l'information initiale et une grande incertitude sur la fiabilité des prévisions. Le formalisme Bayésien permet de mieux caractériser les incertitudes sur les paramètres géologiques et prend en compte l'information géologique a priori. Cependant, l'inversion de type Bayésien est, en pratique, très difficilement applicable à des cas d'étude réels, le coût dû au nombre de simulations nécessaires étant prohibitif. De plus, ces méthodes n'ont pas pour objet de répondre directement à la question de la quantifiction incertitudes liées aux prévisions de la production. Les incertitudes sont classiquement estimées à posteriori, après la phase d'inversion. La prévision des incertitudes est cependant cruciale pour établir un programme d'exploitation d'un gisement, faisant appel à des critère économiques, basés sur la fiabilité des prévision. La méthode selon l'invention permet d'intégrer directement au processus d'inversion 1 aquantification des incertitudes de production dans une même procédure. Cette méthode répond donc pour objectif à la fois de permettre: - un calage optimal d'un modèle de simulation sur un historique de production, en prenant en compte toute l'information a priori; et - une quantification directe des incertitudes sur les prévisions de production et tests d'hypothèses sur l'évolution de ces prévisions de production. La méthode permet à un ingénieur de gisement de tester des "scénarios" prévisionnels de production et ainsi de mieux valider sont modèle par rapport à la connaissance initiale du gisement. Il est ainsi possible de vérifier si un scénario de production est compatible avec les données géologiques, ou de quantifier les incertitudes de production provenant des incertitudes géologiques. Des prévsions concernant l'évolution de la production d'une zone souterraine recelant des fluides telle qu'un gisement d'hydrocarbures, peuvent être obtenues par lé méthode selon l'invention en réalisant les étapes suivantes: a) Dans une première étape classique, déjà décrite pour les méthode d'inversion usuelles, on construit un modèle géologique initial défini par une ensemble de paramètres, en intégrant toutes les données disponibles. A ce stade, les méthodes d'inversion usuelles permettent de déterminer les valeurs de paramètres qui permettent d'ajuster la réponse du modèle sur l'historique de production; b) On définit un ou plusieurs scénarios possibles dévolutions de la production. Pour chacun de ces scénarios, on ajoute aux données observées (ou données "réelles") des données fictives, à des instants futurs donnés correspondant à ces hypothèses; c) Pour chaque scénario, on ajuste les paramètres dudit modèle pour restituer à la fois les mesures réelles provenant de l'historique de production et les données rajoutées à des temps futurs et l'on applique une procdure d'inversion pour obtenir, pour chaque scénario, un nouveau jeu de paramètres x caractérisant un modèle géologique modifié et correspondant aux hypothèses du scénario; d) Pour chacun de ces modèles géologiques modifiés, on fait alors des simulations directes pour établir de nouvelles prévisions de production. il est possible d'utiliser des procédures classiques dans la phase d'inversion. Il est cependant avantageux d'utiliser une procédure d'inversion automatique dans laquelle on minimise une fonction objectif en y incluant l'information initiale associée au modèle géologique, et de faire appel au formalisme Bayésien pour prendre en compte ces données initiales. Dans le but, la méthode peut comporter une combinaison d'un formalisme d'inversion de type bayésien classique et d'un algorithme d'optimisation basé sur l'utilisation de la méthode des gradients qui permet d'obtenir les dérivées des prévisions de production par rapport aux paramètres. La combinaison de cette méthode avec un formalisme de type bayésien, permet de réduire dans des proportions importantes le nombre d'opérations nécessaires pour obtenir l'optimum de probabilité, en alliant les avantages de l'inversion Bayésienne (prise en compte de l'information initiale) et l'efficacité des algorithmes d'inversion utilisant les gradients. La méthode selon l'invention permet donc principalement: - de traduire l'incertitude inhérente au modèle géologique sous la forme d'une quantification de l'incertitutde relative à la production. - de tester si des hypothèses de prévisions sont compatibles avec la connaissance initiale du gisement. Enfin, la combinaison d'un formalisme de type bayésien avec une méthode des gradients permet de traduire cette méthode sous forme d'algorithmes plus efficaces que ceux utilisés antérieurement. D'autres caractéristiques et avantages de la méthode selon l'invention, apparaîtront à la lecture de la description ci-après de modes de réalisation et de tests de validation décrits à titre d'exemples non limitatifs, en se référant aux dessins annexés o: - la Fig. 1 montre schématiquement un modèle géologique formé de de deux couches Cl, C2 en section verticales; - la Fig.2 montre la densité de probabilité a priori p(kl, k2) ("a priori pdf"), caractérisant le modèle géologique initial, en fonction des perméabilités K1, K2 des deux couches C1, C2 du modèle de la Fig. 1; - la Fig.3 montre comment varie la fonction de ressemblance ou de probabilité ("likelihood function") p(P=Pm/K1, K2) caractérisant la probabilité d'une simulation en fonction de ces même paramètres; - la Fig. 4 montre la variation de la densité de probabilité a posteriori p(kl, k2/P- Pm) ("a posteriori pdf") obtenue par inversion de type bayésien; - la Fig. 5 montre le processus d'ajustement automatique, basé sur la combinaison du formalisme de type bayésien et de la méthode des gradients, permettant d'obtenir l'optimum de la densité de probabilité a posteriori. Le processus est illustré en partant de quatre points de départ différents SP1, SP2, SP3, SP4; - la Fig. 6 montre, dans l'espace des paramètres, les modèles géologiques obtenus pour une prévision optimale, une prévision maximale et une prévision minimale; et - la Fig. 7 montre la courbe FPa représentative des prévisions de variation au cours du temps, de la pression au fond du puits producteur, encadrée par les courbes FPM et FPm correspondant respectivement aux prévisions maximales et au minimales. DESCRIPTION DE LA MÉTHODE La méthode d'ajustement de scénarios de production selon l'invention comporte les étapes suivantes 1) On construit tout d'abord un modèle géologique initial de la zone du gisement étudiée, en incluant l'information dont on dispose a priori. Une méthode d'inversion d'un type connu et notamment la procédure d'inversion de Bayésienne précédemment décrite peut être utilisée pour effectuer un premier calage d'historique de production. Le modèle ajusté peut être utilisé pour obtenir une prévision moyenne. 2) Ce modèle initial étant établi, on définit un ou plusieurs scénarios possibles d'évolution de la production. Ces scénarios correspondent à des hypothèses sur l'évolution de la production au cours du temps. Pour chacun de ces scénarios, on ajoute des données fictives, correspondant à ces hypothèses, à des instants futurs donnés. Ces données fictives peuvent traduire par exemple une possible venue d'eau à une période donnée ou encore des hypothèses sur la récupération d'huile. 3) On ajuste les paramètres du modèle initial pour caler à la fois l'historique de production (données observées) et les données additionnelles correspondant aux hypothèses de production. Pour chaque scénario, une procédure d'inversion classique est utilisée de façon à obtenir un nouveau modèle géologique, défini par un nouveau jeu de paramètres x, et correspondant aux hypothèses antérieures. 4) On vérifie directement la probabilité de chaque scénario créé en utilisant un formalisme d'inversion de type Bayesien par exemple, avec la valeur de la fonction pdf a posteriori. ) Pour chacun de ces scénarios, des simulations directes à partir des modèles obtenus permet d'établir de nouvelles prévisions de production. Plusieurs cas peuvent se produire. En considérant le choix des paramètres et les contraintes sur le modèle géologique, on peut trouver un jeu de valeurs de paramètres permettant de reproduire les données mesurées et ajoutées. Cette information peut être très utile pour l'ingénieur de gisement. Dans le cas, par exemple, o le scénario créé a pour objet de tester si une venue d'eau prématurée peut se produire à un puits déterminé, et si la procédure d'évaluation du scénario conduit à un jeu de valeurs des paramètres physiquement admissibles, alors l'ingénieur de gisement peut prendre des décisions relatives à l'exploitation du gisement pour prévenir ce risque, ou décider de faire de nouvelles mesures pour vérifier si le risque est réel. Autre cas envisageable, celui o la procédure d'ajustement de scénario n'est possible qu'au prix d'une distorsion importante du modèle géologique initial. L'ingénieur de gisement peut alors être conforté dans ses prévisions, et être assuré qu'un tel événement est peu probable. Dans la même hypothèse de la prévention d'une venue d'eau prématurée, il peut estimér que le risque est faible avant un certain laps de temps. Une autre application de la méthode est la quantification directe des incertitudes de production, par la recherche de scénarios extrêmes. Cette application fait l'objet de la description ci-après. Comparaison d'un scénario pessimiste et d'un scénario optimiste correspondant à une estimation d'incertitude. La fonction objectif précédente F introduite dans la description des techniques classiques d'inversion utilise l'information obtenue depuis la mise en production du gisement jusqu'au temps présent. Plusieurs jeux de valeurs de paramètres peuvent être des solutions qui minimisent la fonction objectif. Gérer les incertitudes dans le processus de prévision depuis l'instant présent (t=tt) de la période de production jusqu'à un instant futur (t=tf) peut revenir à choisir parmi plusieurs solutions possibles, le jeu de valeurs de paramètres qui donne le "meilleur" scénario, et celui qui donne le "pire". Les scénarios peuvent être définis en utilisant un critère il dépendant de la statégie de production dans la période prévisionnelle envisagée, par exemple entre tt et tf. Le critère peut être par exemple la quantité totale d'huile produite, le débit total maximal, l'instant de la venue d'eau, le G.O.R etc. La meilleure façon de résoudre ce problème, est d'appliquer la méthode selon l'invention à deux problèmes d'optimisation correspondant respectivement le premier, au scénario le"meilleur", et le deuxième, au "pire". Dans le premier cas, on trouve le jeu de paramètres qui donne à la fois un bon calage de l'historique de production du gisement (c'est-à-dire le jeu de paramètre qui minimise Fm + Fx) et une valeur maximale pour le critère: Min(Fm + Fx) + Max(Critère). Dans le deuxième, on trouve le jeu de paramètres qui combine un bon ajustement de l'historique avec une valeur minimale du critère: Min(Fm + Fx) + Min(Critère). Formalisme de la méthode: Un instant futur tf. est choisi et de nouvelles données sont ajoutées aux données de production existantes. Chacune des nouvelles valeurs correspond à un scénario prévisionnel de production. On détermine les scénarios prévisionnels de production min/max à cet instant futur tf. A cet effet, les valeurs des paramètres sont ajustées pour obtenir un calage optimal de l'historique de production tout en déformant la courbe de prévision de production de façon que le minimum ou le maximum soit atteint au temps futur t=tf. Pour contrôler la possibilité de chaque scénario, on impose une contrainte d'iso- probabilité. Le problème se réduit à un problème d'optimisation avec des contraintes non- linéaires. La fonction objectif utilisée est F= Fm+Fx+Ff Fm + Fx représente la fonction objectif précédente incluant les données a priori. Ff a pour expression: Ff = 1/2 ( df - g(x,tf))2; o df représente les données additionnelles au temps tf; et g(,tf), les prévisions de production au temps tf. Ff est une contrainte additionnelle utilisée pour obtenir les scénarios de min/max. Les données ajustée sont df pour satisfaire les deux contraintes suivantes: 1) Fm + Fx 2) g(.,tf) représente une prévision maximale ou minimale. La première contrainte est équivalente à une contrainte de probabilité relative sur une densité de probabilité a posteriori, en se référant au maximum de probabilité. Un rapport de probabilité a posteriori (ppr) étant donné: ppr = p(x/d=dm) / p(.oo/d=dm), o xOo est la position de l'optimum. (valeurs des paramètres correspondant au maximum de probabilité).Le critère de norme admissible de calage de l'historique devient alors: Fhmc = Fm(xoo) + Fx(oo) - log(ppr) Pour minimiser la fonction objectif (F = Fm + Fx + Ff),, on peut utiliser un algorithme d'optimisation à convergence quadratique, tel que l'algorithme de Gauss- Newton, bien connu des spécialistes. Cet algorithme peut être utilisé pour conduire une procédure automatique d'ajustement de l'historique, en recherchant simultanément le maximum ou le minimum des prévisions de production. Un simulateur de gisement multi-usage peut être avantageusement utilisé pour conduire cette procédure en faisant intervenir la méthode des gradients. On obtient en effet, par la méthode des gradients, directement la sensibilité du modèle aux paramètres de l'inversion. Cette méthode est bien connue des spécialistes et décrite par exemple par: Anterion F. et al., "Use of parameter gradients for reservoir history matching", SPE 18433, Houston, TX, February 6-8, 1989. En utilisant les valeurs de gradients, l'algorithme d'optimisation est d'une grande efficacité pour traiter le problème dans une méthode itérative globale, incluant la recherche de scénarios extrêmes et le calage optimal de l'historique de production. La méthode selon l'invention permet d'utiliser un algorithme d'inversion combinant les avantages du formalisme Bayesien et l'efficacité des algorithmes d'optimisation basés sur l'utilisation gradients. La méthode permettant de combiner le formalisme Bayésien et la méthode des gradients est décrite ci-après. Formalisme d'inversion de type bayésien modifié Pour conduire les procédures d'inversion, on peut utiliser avantageusement une méthode d'inversion particulière dérivée du formalisme de type bayésien, qui permet une recherche directe et beaucoup plus rapide que le formalisme classique, de la position de l'optimum de la densité de probabilité a posteriori. Ce but est atteint en définissant une fonctionnelle à optimiser et en combinant le formalisme du type bayésien. avec la méthode des gradients. L'information initiale associée au modèle géologique peut ici être incluse dans cette fonctionnelle à minimiser. De plus, moyennant certaines suppositions, il est possible d'obtenir une approcimation de la densité de probabilité a posteriori autour de cette solution optimale. Pour inclure cette connaissance initiale dans la fonctionnelle à minimiser (notée F), on introduit un nouveau terme Fx dépendant de l'écart entre une moyenne Xa donnée a priori et la valeur x des paramètres: F = (Fm + Fx)/2, avec Fm = (dm-d)T Cm-1 (dm-d), et Fx =(xa - x)T Cx-1 (xa-x). Des opérateurs de covariance sont introduits pour modéliser les incertitudes gaussiennes sur les mesures (matrice Cm) ainsi que la densité de probabilité a priori sur les paramètres (matrice Cx). La formulation proposée constitue une extension de la fonction objectif classique des moindres carrés.Il y a un lien direct entre la fonction objectif et la règle de Bayes. La fonction pdf a posteriori peut s'exprimer sous la forme: p(xJ/d=4nm) = Cte. exp (-F) = Cte [exp (-Fm). exp(-Fx)), avec p(d=dmJx) = Ctel {exp (-Fm)} p(x) = Cte2 {exp (-Fx)}. On peut considérer ainsi que les paramètres obtenus pour le maximum de la fonction pdf a posteriori, désignés par xoo, sont ceux qui minimisent la fonction objectif F. Cette correspondance est bien connue des spécialistes de l'inversion. La procédure modifiée se distingue par: - l'utilisation de la méthode des gradients pour minimiser la fonction objectif, qui, selon l'invention, permet d'obtenir des scénarios de production extrêmes avec un critère d'iso-probabilité; et par - l'utilisation de la méthode des gradients combinée au formalisme de type bayésien pour obtenir, par linéarisation, une approximation de la densité de probabilité a posteriori au voisinage de l'optimum. Ce procédé permet en effet de satisfaire demanière efficace les contraintes d'iso-probabilité a posteriori dans les algorithmes d'inversion. Quand la convergence est atteinte dans la procédure d'optimisation et que l'on a trouvé les valeurs optimales xoo en utilisant la méthode des gradients. En supposant que la densité de probabilité a posteriori est localement gaussienne, on la modélise par un opérateur de covariance a posteriori. Le calcul de cet opérateur de covariance donne une approximation de la fonction "densité de probabilité a posteriori". Cas expérimental On a testé la validité de la méthode selon l'invention sur un cas expérimental synthétique proposé ci-après qui démontre bien l'amélioration apportée en combinant le formalisme d'inversion Bayésien avec une méthode de gradient, et la possibilité de quantification des incertitudes sur les prévisions par la recherche de scénarios extrêmes. Après avoir défini une représentation géologique synthétique, on a réalisé une simulation numérique pour obtenir des données de production synthétiques, en utilisant un simulateur industriel de gisement d'un type usuel). On choisit une représentation géologique initiale qui constitue l'information dont on dispose a priori sur le gisement. Le modèle géologique est ici une section verticale comportant deux couches (Fig. 1l) C1, C2 de hauteur 7,5 m et d'épaisseur m, avec des perméabilités horizontales K1, K2 constantes de 200 mD et 100 mD respectivement. La grille de simulation est constituée de dix cellules pour chaque couche dans la direction horizontale, et la longueur totale s'élève à 500m. On injecte de l'eau, sur le côté gauche du réservoir, dans la coupe inférieure, et un débit total de liquide (eau + huile) est produit sur le côté droit de la couche supérieure. On impose le même débit volumique (200 m3/j ou 1250B1/j) à l'injecteur et au puits producteur afin de maintenir la pression moyenne et d'éviter que la pression de l'huile ne descende au-dessous du point de bulle. Les fluides initiaux qui se trouvent dans le réservoir forment un réservoir pétrolifèere comportant 12 % d'eau irréductible. Une simulation numérique est mise en oeuvre durant une période de production de 600 jours. Au cours d'une première période (de t = 0 à t = 300 jours), seule la phase huile est produite. L'eau n'arrive au puits producteur qu'après une durée d'environ 300 jours. Le volume d'eau produit augmente au cours de la dernière période: 10 % d'eau à t = 380 jours et 77,7 % d'eau à la fin de la simulation (t = 600 jours). On produit des mesures synthétiques en utilisant les résultats de cette simulation de référence. Seule la première période de production est utilisée de t = 0 à t = 300 jours): on admet que le temps o se produit une arrivée d'eau est inconnu pour conserver une grande part d'incertitude relative à la caractérisation initiale du réservoir. La mesure choisie est la pression de fond dans le puits producteur, avec cinq valeurs observées entre 25 et 300 jours. On choisit une gamme d'incertitude globale pour décrire à la fois les incertitudes de mesure et les erreurs de simulation numérique. Ces incertitudes sont modélisées au moyen de densités de probabilité gaussiennes centrées sur les résultats de la simulation de référence et avec am = 138 kPa (20 psi) pour les écarts-types. Ces écarts-types sont mis en oeuvre dans la fonction objective dans les termes diagonaux de l'opérateur de covariance Cm relatif aux mesures. Si l'on considère le premier terme de la fonction objectif Fm qui quantifie l'écart entre les observations dm (dans le cas présent, la pression de fond) et la réponse de simulation g(x,t), les termes diagonaux de Cm-1 sont égaux à 1/Om2. Définition du modèle géologique a priori On suppose que le modèle géologique de référence est inconnu et que l'on dispose d'informations géologiques a priori. On considère par hypothèse que la géométrie est bien caractérisée, et que les perméabilités horizontales affectées à chaque couche sont inconnues et que l'on dispose d'informations a priori sur ces perméabilités avec une marge d'incertitude donnée. Les deux perméabilités horizontales ont été choisies comme étant les paramètres à inverser. Les informations disponibles a priori sont modélisées au moyen d'incertitudes gaussiennes (Fig.2). L'écart entre les perméabilités a priori xa et les paramètres x est quantifié dans le second terme de la fonction objective Fx. Les écarts-types des paramètres c(x sont introduits dans l'opérateur de covariance Cx, comme de la même manière que pour la modélisation des erreurs de mesure: les termes diagonaux de Cx 1 sont égaux à 1/Cx2. Les termes extra- diagonaux sont nuls pour cet exemple, ce qui signifie que les informations a priori relatives à une possible corrélation entre les deux paramètres n'ont pas été introduites. Les valeurs numériques choisies a priori sont présentées dans le tableau 1 ci- après: Tableau 1: Modèle géologique a priori Valeur a priori Ecart- type Valeur du cas de Xa x référence K1 160 mD 40 mD 200 mD (couche supérieure) K2 160 mD 40 mD 100 mD (couche inférieure) Comparaison de différentes méthodes d'inversion pour le cas test A) Méthode de calage d'historique automatique On applique une méthode de calage d'historique classique combinée à un algorithme d'optimisation automatique en premier lieu. Les informations géologiques disponibles a priori sont utilisées, dans cette première approche, uniquement pour choisir le point de départ de la procédure. La fonction objective à minimiser inclut uniquement le premier terme Fm, et ne tient pas compte du modèle a priori. Les paramètres initiaux sont fixés à 160 mD (valeurs a priori). On utilise la procédure de calage automatique basée sur une méthode de gradient et l'algorithme d'optimisation classique de Gauss-Newton. On obtient pour ce cas, par la méthode des gradients, la sensibilité de la pression de fond dans le puits producteur aux valeurs de perméabilités. En raison de l'utilisation de mesures synthétiques produites par le simulateur numérique, le modèle adapté ainsi obtenu est très proche du cas de référence. Cinq simulations seulement ont été nécessaires pour obtenir la position du minimum de la fonctionelle. Le processus d'optimisation est résumé, depuis le point de départ jusqu'à un optimum, dans le tableau 2 ci-dessous: Tableau 2: Processus d'optimisation HM Numéro d'itération K1 K2 Fm 0 160 mD 160 mD 0,799 1 172,8 mD 133,6 mD 0,407 2 193,4 mD 104,9 mD 0,042 3 200,8 mD 99,1 mD 0,0008 4 200,5 mD 99,9 mD 0, 0003 Utilisation du formalisme d'inversion de type Bayesien Par la méthode d'inversion bayésienne, on obtient la densité de probabilité a posteriori, cartographiée sur l'ensemble de l'espace des paramètres. Dans ce cas, l'information disponible a priori (Fig.2) est prise en compte dans cette densité de probabilité a posteriori, en la combinant avec la fonction de probabilité (Fig.3) calculée pour chaque groupe de paramètres. L'espace des paramètres a été cartographié avec une large gamme de perméabilités, de 20 mD à 300 mD, avec 35 valeurs pour chaque paramètre (le nombre total de simulations numériques s'élève à 352 = 1225 opérations). Il est avantageux de disposer d'une carte complète de la fonction de probabilité à posteriori dans la mesure o elle permet d'obtenir de nombreuses informations sur le comportement en simulation. Si l'on considère la représentation graphique de la fonction (Fig.4), on peut remarquer les points suivants: - un seul optimum est clairement identifié, qui se situe à 143 mD (perméabilité de la couche supérieure) et à 173 mD (perméabilité de la couche inférieure). La valeur maximale de la densité de probabilité à posteriori est de 1,48.10-7; - une meilleure caractérisation du modèle géologique a pu être obtenue avec une incertitude moindre concernant les paramètres (en comparaison de la densité de probabilité pdf a priori, Fig.2); - une corrélation est établie entre les deux paramètres: la forme de la surface n'est pas symétrique autour de l'optimum, mais déformée sur un axe de corrélation; D'autre part, l'intégration de la surface permet de calculer des paramètres statistiques à posteriori concernant les paramètres. Les tableaux 4 et 5 ci-après montrent les écarts-types a posteriori comparés aux données initiales: Tableau 4: Paramètres initiaux et a posteriori Initial a Postériori K1 160 mD 143 mD K2 160 mD 173 mD Tableau 5: Ecarts-types initiaux et a posteriori Initial Postérieur K1 40 mD 32,8 mD K2 40 mD 19,0 mD Utilisation ddu formalisme de type bayésien combiné avec la méthode des gradients Cette procédure d'inversion modifiée combinant comme on l'a vu le formalisme d'inversion dit de Bayésienne avec une méthode de gradients, a été appliqué au cas expérimental. L'algorithme d'optimisation est mis en oeuvre pour rechercher directement l'optimum de la densité de probabilité pdf à posteriori. On a utilisé quatre points de départ différents pour couvrir la gamme d'incertitudes initiale. Pour chaque processus d'optimisation, l'optimum est atteint avec un nombre d'itérations compris entre 4 et 7 (Fig.5). L'analyse statistique proposée est appliquée aux environs de l'optimum. Les résultats obtenus sont comparés à la méthode classique d'inversion Bayésiennedans le tableau 6 ci-après: Tableau 6: Résultats comparés obtenus suivant l'inversion bayesienne classique et l'inversion modifiée Inversion bayesienne inversion modifiée Nombre d'opérations 1225 20 Optimum (en mD) K1=143 K1=143, 2 K2=173 K2=172,9 Ecarts-types (en mD) K1=32,8 Kl=31,4 K2= 19,0 K2= 18,8 Max. de la densité de probabilité pdf postérieure 1,48.10-7 1,57.10-7 Ce tableau met clairement en évidence la vitesse d'exécution permise par la méthode d'inversion modifiée selon l'invention. Utilisation de la méthode des scénarios Les résultats obtenus pour le cas expérimental sont présentés sur la Fig. 6. Après avoir utilisé la procédure automatique décrite ci-dessus pour rechercher l'optimum de la densité de probabilité a posteriori, et donc pour obtenir le calage avec les données réelles (de t = 0 à t = 300 jours), en prenant en compte le modèle a priori, on a utilisé une simulation numérique directe jusqu'à tf = 600 jours. Cette simulation donne une prévision moyenne de la pression de fond. On choisit ensuite le taux de probabilité a posteriori ppr, soit ppr = 0,1. Les scénarios de production extrêmes sont identifiés, le modèle géologique restant dans ce domaine de probabilité, en utilisant la méthode selon l'invention. Les contraintes d'iso-probabilité sont représentées sur la Fig.7 (en utilisant la cartographie d'inversion de type bayésien) afin de superposer les valeurs correspondant aux trois scénarios (optimum, maxi. et mini.) dans l'espace des paramètres. la localisation des paramètres correspondant aux scénarios mini/maxi est proche de la frontière définie par le taux de probabilité ppr = 0,1. Les résultats numériques sont présentés dans le tableau 7 Tableau 7: Scénarios extrêmes et comparaison avec l'optimum Maxi. Optimum Mini. Pression (psi) 4179 4098 4020 KI (mD) 90,2 143,2 200,8 K2 (mD) 217,4 173,0 136,5 L'incertitude concernant la prévision relative à la pression de fond à tf = 600 jours est donnée par la différence entre la prévision maximale et la prévision minimale. Ceci démontre qu'il est possible de transposer les incertitudes relatives au modèle géologique en incertitudes relatives aux prévisions en matière de production. REVENDICATIONS 1) Méthode pour prédire l'évolution de la production d'une zone souterraine recélant des fluides telle qu'un gisement d'hydrocarbures, caractérisée en ce qu'elle comporte les étapes suivantes: a) on construit un modèle géologique initial défini par un ensemble de paramètres, en intégrant les données disponibles; et l'on applique une méthode d'inversion pour déterminer les valeurs de paramètres qui permettent d'ajuster la réponse du modèle sur un historique de production connue par ladite zone; b) on définit un ou plusieurs scénarios possibles d'évolution de la de production et, pour chacun de ces scénarios, on ajoute des données observées des données fictive, à des instants futurs donnés correspondant à des hypothèses; c) pour chacun desdits scénarios, on ajuste les paramètres dudit modèle pour restituer à la fois des mesures réelles provenant de l'historique de production et des données rajoutées à des temps future et l'on applique une procédure d'inversion pour obtenir, pour chaque scénario, un nouveau jeu de paramètres x caractérisant un modèle géologique modifié et correspondant aux hypothèses du scénario; et d) pour chacun de ces modèles géologiques modifiés, on fait alors des simulations pour établir de nouvelles prévisions de production. 2) Méthode selon la revendication 1, caractérisée en ce que l'on détermine des scénarios correspondant à des prévisions extrêmes, de façon à quantifier les incertitudes de prévisions issues d'une simulation. 3) Méthode selon la revendication 1 ou 2, caractérisée en ce que l'on ajuste le calage des scénarios dans un processus d'optimisation par une combinaison d'un formalisme de type bayésien et une méthode des gradients. 4) Méthode selon la revendication 3, caractérisée en ce que l'on utilise une méthode des gradients pour obtenir une approximation de la densité de probabilité a postenriori.