FR2765708A1 - Methode pour determiner des parametres hydrauliques representatifs a grande echelle d'un milieu fissure - Google Patents

Abstract

- Méthode pour déterminer, à partir d'une image 3D, les valeurs de paramètres hydrauliques tels que la perméabilité équivalente à grande échelle, la perméabilité des blocs et le coefficient d'échange alpha entre matrice et fissures d'un milieu poreux fissuré tel qu'une formation géologique.- La méthode comporte essentiellement la discrétisation du milieu par un maillage et la résolution rapide et approchée dans ce maillage d'équations modélisant la diffusion des fluides dans le milieu, la détermination de la variation au cours du temps d'une fonction de transfert à grande échelle (f (t), f (s) ) caractérisant les flux de fluide de la matrice vers les fissures, par simulation du mouvement de particules effectuant des marches aléatoires en temps continu sur ledit maillage et le traitement approprié d'une fonction d'état (epsilon(tau)) indicatrice de leur présence soit dans la matrice soit dans une fissure.- Application à la modélisation à grande échelle de réservoirs pétroliers fissurés, permettant l'interprétation d'essais de puits.

Description

La présente invention concerne une méthode pour déterminer, à partir d'une

image d'un milieu fissuré tel par exemple qu'une formation géologique traversée par un réseau de fissures, des paramètres hydrauliques qui le caractérisent au mieux à grande échelle, à savoir: la perméabilité équivalente à grande échelle, la perméabilité des blocs et le coefficient d'échange a entre matrice et fissures. La méthode selon l'invention trouve des applications notamment dans le domaine de la production pétrolière o l'on cherche à modéliser le comportement d'un gisement en vue d'obtenir des prédictions fiables de production d'hydrocarbures. La connaissance de ces paramètres est importante car elle permet de procéder à des simulations d'écoulements complexes sur le milieu et à des inversions de données obtenues dans des essais de puits qui, en retour, permettent d'améliorer l'image géologique/pétrophysique du milieu, et à terme, par simulation, l'optimisation de l'exploitation ultérieure du

gisement. La méthode selon l'invention peut convenir également en hydrogéologie.

ART ANTERIEUR

Il est fréquent de décrire les gisements ou réservoirs fissurés comme comportant deux milieux contrastés, un milieu matriciel contenant la plus grande part de l'huile en place et présentant une faible perméabilité, et un milieu plus ou moins complexe, en général hautement perméable de fissures caractérisées par leurs densité, longueur,

Orientation, inclinaison et ouverture respectives.

Dans un contexte monophasique, le problème de la simulation des essais de puits à petite échelle s'écrit sous la forme d'une équation aux dérivées partielles: )p(r;t) (r)) 99,1 = 17.(k (r)VP(r, t" (1) -t o p(r, t) désigne la valeur de la pression au point r à l'instant t, p. est la viscosité du fluide et (p, la porosité, c,, la compressibilité totale (roche+fluide) et k(r), la perméabilité. Cette équation est à résoudre avec des conditions aux limites et des conditions initiales convenables au N puits perforés et aux frontières du gisement de façon à assurer l'existence et l'unicité de sa solution, en utilisant comme données de base, une carte de perméabilité k(r) établie au préalable. Les mesures obtenues sur le terrain sont les N enregistrements des pressions p(ri, t), i = 1, N, et des débits de fluide en ces même N puits. Tout le problème pratique consiste à retrouver des informations sur la carte de perméabilité k(r) en fonction de ces N mesures. C'est l'ensemble des techniques

employées dans ce but que l'on appelle l'interprétation des essais de puits.

Pour établir une telle carte dans des cas de réservoirs fissurés, on utilise par exemple des logiciels de modélisation qui, à partir des connaissances disponibles sur le milieu: photographies. données sismiques interprétées, études géologiques, mesures dans des puits forés au travers, analyses de carottes. etc., génèrent des représentations 1 5 spatiales du réseau de fissures, inclus dans une matrice elle-même perméable selon des règles bien précises décrites par exemple dans les brevets FR 2 725 794, FR 2 725 814, FR 2 733 073 du demandeur. Comme les zones du réservoir entre les puits sont imparfaitement connues, ces logiciels de modélisation utilisent eux-mêmes des paramètres mal connus. On souhaiterait être capable d'affiner cette connaissance à partir des essais de puits. Dans ce contexte fissuré, c'est l'ensemble des techniques employées

dans ce but que l'on appelle l'interprétation des essais de puits en milieu fissuré.

Après utilisation de ce type de logiciel de modélisation, le réservoir pétrolier Q2 est partitionné en au moins deux sous-domaines Q2f et Qm représentant respectivement

les fissures et la matrice. La perméabilité k(r) vaut kf(r) si la position r est dans le sous-

domaine des fissures Qf, et km(r), sinon, c 'est à dire si r est dans la matrice 2m. En général, les perméabilités des fissures kf(r) sont très supérieures à celles, km(r), de la

matrice: kf(r) >> km(r).

Les dimensions typiques des fissures (épaisseur de quelques microns ou plus, longueurs très variables) et leur grande densité (plusieurs fissures par mètre cube de roche) interdisant en pratique l'utilisation d'un modèle finement maillé pour résoudre l'équation aux dérivées partielles (1), il est d'usage dans la profession pétrolière de considérer une modélisation à grande échelle telle que décrite par exemple par - J. E. Warren & P. J. Root, The Behavior of Natural Fracture Reservoirs, The Society

of Petroleum Engineers Journal 3, n 3, 245--255, 1963.

Selon ce modèle, qui est une idéalisation de la réalité, on représente un volume élémentaire sous la forme d'un ensemble de blocs parallelépipédiques identiques limités par un réseau de fissures orthogonales uniformes continues orientées dans la direction de l'un des trois principaux axes. L'écoulement des fluides au travers du réservoir s'effectue au travers du milieu de fissures seulement et des échanges de fluides interviennent localement entre les fissures et les blocs matriciels. On résoud les équations suivantes P, (r, t) Kf " ó0/pc, I V.- Pf(r,t) +qm(r,t) t,. -. P. (r,t) -q(2) Les pressions Pm et Pf représentent des moyennes à grande échelle des pressions dans la matrice et dans la fissure. Les coefficients Km et Kf décrivent les perméabilités équivalentes à grande échelle de la matrice et du réseau de fissures. Les quantités Of et om désignent les proportions volumiques de fissures et de matrice. On a évidemment la relation Of + Om = 1. Le terme qm(r, t) représente un terme source d'échange entre les deux pressions qui, dans le cas "pseudo permanent", s'écrit sous la forme suivante c,,(r, t)= (P., (r, t) - P, (r, t)) (3) Ce terme représente l'échange de matière entre les deux milieux, ce qui justifie l'appellation "terme d'échange", et l'appellation "coefficient d'échange" donnée à cc. On

remarque que la quantité (xKm/ ópct est homogène à l'inverse d'un temps.

Il est fondamental de comprendre qu'à l'échelle considérée, les fissures

individuelles ont disparu de la description. En particulier, les paramètres Km, Kf et a

sont en général indépendants du point r considéré, et les détails des propriétés du réseau

de fissures n'interviennent que via la valeur des dits paramètres.

L'idée physique de base de la présente modélisation est de considérer que le fluide stocké dans la matrice se déplace vers les puits via le réseau de fissures. Dans le cas d'un milieu constitué de blocs matriciels cubiques de côté L séparés par des fissures, et est de l'ordre de 1I/L2, d'o l'expression "taille de blocs équivalente" parfois rencontrée pour désigner ce coefficient. D'autres représentations plus complexes du terme qm(r, t) existent, sous forme de convolution temporelle et portent le nom de "modèles transitoires". Cette formulation connue notamment par Warren et al déjà cités, quand elle est

étendue aux écoulements polyphasiques, constitue ce que l'on appelle la description

"double milieu". Elle sert de base à toutes les simulations d'écoulements en réservoir fissurés d'usage dans le milieu pétrolier. Dans le cas particulier des essais de puits monophasiques, la connaissance des solutions analytiques des équations 2, permet de caler les paramètres à partir des variations de pression observées au(x) puits. En particulier, si le puits est de rayon rw, il est d'usage d'introduire le coefficient k = cKnKm(L/rw)2. En pratique, on arrive à des calages corrects pour des essais de puits, mais

dès lors que l'on désire caler la description à l'échelle du champ, on peut trouver d'autres

valeurs pour les coefficients.

Ce modèle, dont le grand mérite est sa simplicité, possède une grande lacune: les paramètres Km, Kf et oc y sont a priori phénoménologiques et ajustables; leur relation avec la structure détaillée du milieu et le processus d'écoulement est très mal connue, et

on dispose de peu d'outils de calcul efficaces.

Pour le calcul de ces paramètres, il faut une remise à l'échelle de la description

détaillée du milieu. Par: - Quintard, M. et al, One- and Two-Equation Models for Transient Diffusion Processes in Two-Phase Svstems, Advances in Heat Transfer, 23, 369--464, 1993, par exemple, on connaît une méthode pour obtenir la valeur de ces paramètres qui comporte la résolution de trois "problèmes de fermeture" définis sur une portion représentative du milieu. Le calcul cependant nécessite la résolution de problèmes aux limites stationnaires exigeant des représentations finement maillées du réservoir. La notion de volume représentatif est délicate, et sa détermination pratique est inconnue. De plus, le coût informatique de ces méthodes devient colossal dès que le réseau de fissures est complexe, car elles n'évitent pas l'obligation de résoudre des systèmes linéaires de grande taille (plusieurs millions d'inconnues) sitôt que l'on veut représenter

convenablement quelques milliers de fissures.

Il existe une technique dite des "marches aléatoires", qui est utilisée pour déterminer les variations de la perméabilité du milieu homogène et hétérogène. On la trouve décrite notamment par: - McCarthy, J. F., Continuous--time Random Walks on Random Media}, J. Phys. A

Math. Gen. 26, 2495--2503, 1993.

On la rappelle ci-après pour la bonne compréhension. Le principe en est le

suivant. On considère un certain nombre de particules de fluide indépendantes placées à.

lI'origine sur certaines mailles d'un maillage régulier discrétisant le milieu étudié, et l'on étudie leur déplacement sur les mailles voisines avec une équiprobabilité pondérée pour tenir compte des hétérogénéités du milieu. Le retour à une vision globale du milieu est réalisé en calculant la moyenne des carrés des déplacements de toutes ces particules. On suppose que le problème de la résolution de l'équation (1) a été préalablement discrétisé par une méthode de différence finies ou d'éléments finis avec un maillage régulier en espace, de N mailles. Les équations discrétisées s'écrivent sous la forme suivante Vaincp/, (tlf =:s Tj,(Pj,(t)- pi(t)) dt (4)

I 1 1

avec - - ±

Tj k1 kj

o mj désigne les mailles voisines de la maille i.

Les N quantités Pi(t) représentent donc la pression évaluée au centre des N noeuds. Les coefficients Tij sont les transmissivités entre les mailles i et j. et sont données par les moyennes harmoniques perméabilités ki et kj des deux mailles considérées. Le caractère hétérogène ou fissuré du milieu est donc pris en compte via les valeurs de ces transmissivités. Le coefficient Vi = (AX) D est le volume des mailles,

supposé ne pas dépendre de i (hypothèse de maillage régulier).

On considère maintenant la marche aléatoire construite de la façon suivante. A l'instant initial t = 0 (itération "zero", k = 0), on tire au hasard une maille d'indice i 0 et on suppose qu'une particule y est lâchée. A l'itération k de l'algorithme, on suppose que

la particule est parvenue sur la maille repérée par l'indice ik, et que l'on est à l'instant tk.

Pour trouver sa position ik + 1 à l'itération k + 1, on choisit au hasard une maille voisine de ik. indexée par j avec la probabilité suivante T, ,, - Zj,< (5) La somme sur l'indice j portant sur les mêmes mailles voisines que dans (4). Les probabilités b sont normalisées: leur somme valant 1 par construction; la particule est obligée de sauter d'un site à l'autre, sans possibilité pour elle de rester sur place. Dans le cas particulier o la particule est sur un site au bord du domaine, on emploie les conditions aux limites périodiques: tout se passe comme si la particule pouvait "sortir" du milieu en y retournant de l'autre côté (sa position absolue étant tout de même retenue). Ce choix étant fait, la position ik + 1 à l'itération k + 1 est connue, et le compteur de temps est alors incrémenté de la façon suivante tkit = tk _ c Ln (ran(k + 1)) (6) O ran (k + 1) est un nombre aléatoire de loi uniforme sur l'intervalle (0, 1). Le logarithme en étant toujours négatif, le compteur de temps croît bien avec l'indice k. AX représente le pas spatial du maillage. La quantité en facteur est bien homogène à un temps, et représente un temps de diffusion typique sur la distance AX. L'accroissement du compteur de temps correspond à la distribution des temps de séjour d'une particule

sur la maille considérée.

Le processus peut être maintenant itéré, et il est arrêté lorsque le temps tk dépasse

une valeur de consigne Tmax fixée par l'utilisateur.

Un parcours ou marche aléatoire est donc constitué d'un tableau de (déplacements) de taille K(Tmax) (taille dépendante de la réalisation) comportant les valeurs de K(Tmax) positions ik et de K(Tmax) instants tk associés. On considère maintenant un grand nombre de telles marches (parcours) aléatoires indépendantes, donc Q listes de positions et d'instants. La probabilité pour une particule d'être sur la maille i entre t et t+dt définit une densité de probabilité Pi(t) qui serait

calculable en effectuant un très grand nombre de marches.

On peut démontrer que l'ensemble des densités Pi(t) obéit aux équations (5) (5?).

Dans le cas considéré, o les particules partent à t=0 d'un site aléatoire du réseau, la condition initiale est Pi(t=0) =I/N. Si l'on choisissait de partir toujours du même site io, à t = 0 on aurait Pi(t=0) = Si io, o 8 est le symbole de Kroenecker. On obtiendrait ainsi la fonction de Green discrète du problème, i.e. la probabilité d'être sur la maille "i" à

l'instant t sachant que l'on était sur la maille io, à t = 0.

On désigne ainsi par "marches aléatoires à temps continu" (CRTW en anglais) ce type de déplacement. C'est la gestion apparemment complexe du pas de temps (équation 6) qui permet aux particules de se déplacer à chaque pas sans jamais faire de pause sur une maille et in fine d'obtenir de bonnes performances numériques (car chaque pose correspondrait à des itérations numériques, et donc à du temps-machine inutilement employé). La méthode selon l'invention évite les inconvénients signalés ci-dessus des méthodes antérieures, car elle permet d'éviter la résolution complète de systèmes linéaires de très grande taille, au prix d'une définition statistique des paramètres à grande échelle. A partir d'une image d'une formation géologique traversée par au moins un réseau de fissures, la méthode permet de déterminer des paramètres hydrauliques de la formation géologique tels notamment que le coefficient d'échange (o) entre matrice et le

réseau de fissures, dans le but d'interpréter des essais de puits.

Malgré les incertitudes sur la structure exacte du milieu fissuré et le caractère

incomplet de la description à grande échelle des écoulements (équations (2)). qui rendent

impossible la détermination de la solution complète sur toutes les mailles du modèle fin.

la méthode permet pourtant de disposer d'une valeur correcte pour les perméabilités équivalentes (Km, Kf) à grande échelle de la matrice et des fissures ainsi que le

coefficient d'échange (x). De plus, elle fournit des informations sur le régime transitoire.

La méthode comporte la discrétisation du milieu par un maillage et la résolution probabiliste (de type Monte-Carlo) dans ce maillage d'équations modélisant la diffusion de la pression dans le milieu. Elle est caractérisée principalement en ce qu'elle comporte la détermination de la variation au cours du temps d'une fonction de transfert à grande 1 5 échelle (f(t)) ou (f(s), s étant le paramètre de Laplace) caractérisant les flux de fluide de la matrice vers le réseau de fissures, par application aux particules de fluide de déplacements aléatoires en temps continu, en leur associant une fonction d'état (e(r))

indicative de leur présence soit dans la matrice soit dans le dit réseau de fissures.

La méthode peut comporter par exemple la détermination d'une fonction (C(T)) caractérisant la décroissance au cours du temps de la pression moyenne de fluide dans les fissures liée à la probabilité de présence des particules. par une auto-corrélation temporelle de ladite fonction d'état, de laquelle on déduit ladite fonction de transfert

(f(s)) ainsi que le coefficient d'échange (x) par une formule de quadrature numérique.

La connaissance de ces paramètres obtenus par la méthode selon l'invention, permet de procéder alors à des simulations d'écoulements complexes sur ce milieu et à des inversions de données d'essais de puits qui donneront en retour une meilleure image géologique/pétrophysique. La technique utilisée permet d'éviter toute résolution de système linéaire, et donc de procéder à des calculs sur des maillages à ultra-haute résolution ingérables à l'aide de méthodes classiques en l'état actuel des capacités de calcul. On obtient de surcroît le moyen de préciser les régimes d'écoulement dans le milieu: régime pseudo-permanent ou transitoire en obtenant dans ces deux cas la forme explicite des fonctions de transfert. La méthode permet en effet d'obtenir - au delà du coefficient d'échange ot - les variations dans le temps de la fonction de transfert f(t) (ou f(s) en représentation de Laplace) caractérisant les flux de fluide de la matrice vers les

fissures, et de traiter ainsi les cas transitoires.

D'autres caractéristiques et avantages de la méthode selon l'invention,

apparaîtront à la lecture de la description ci-après d'un exemple non limitatif de

réalisation, en se référant aux dessins annexés o: - la Fig. l montre un exemple de variation en fonction du temps d'une fonction ú(t) rendant compte de la position d'une particule de fluide soit dans la matrice, soit dans une fissure; - la Fig.2 montre les variations comparées d'une fonction normalisée C(t), les résultats

étant obtenus avec la méthode selon l'invention et avec une méthode analytique.

- la Fig.3 schématise un milieu fissuré de type "feuilleté".

Description de la méthode selon l'invention

La méthode comporte une discrétisation du milieu hétérogène fissuré par un maillage cartésien suffisamment fin. On considère qu'une maille "i" quelconque est une maille "fissure" si son centre est dans une fissure, et une maille "matrice", sinon. A toute maille, on associe donc un scalaire F(i) tel que E(i) = I si la maille "i" est sur une fissure,

et ú(i) = O, sinon.

Si le milieu est suffisamment discrétisé, la proportion de mailles fissure vaut volume occupé par les fissures .t =, soil t: volume total sot N e(i) / N = (7) *=1 La technique rappelée de "marches aléatoires", est appliquée à une particule de fluide. A tout instant tk, on connaît la position ik de la particule et on peut donc associer à tk,la valeur E(tk) = s(ik). Ce nombre indique si, à l'instant tk, la particule est ou non dans une fissure. Par construction, cette fonction vaut zéro ou 1, et elle est

artificiellement prolongée sur l'intervalle entre tk et tk + 1 par sa valeur à l'instant tk.

Comme le montre l'exemple de la Fig. l, on obtient ainsi une fonction définie sur la demi droite réelle que l'on peut facilement définir de façon univoque et coder en retenant l'état

initial de la particule, et les dates de ses changements de milieu.

On considère maintenant la fonction suivante C(X) = < ú(t) E(t + t) >, o les termes entre crochets correspondent à une moyenne sur toutes les marches possibles. On peut remarquer que l'argument t n'intervient pas dans la fonction C. Ceci tient à ce que la position de la particule à t + t, une fois connue à l'instant t, ne dépend pas de son

parcours aux instants antérieurs à t.

La fonction C(t) n'est autre que la fonction d'autocorrélation de présence de la particule dans les fissures: elle décrit la densité de probabilité d'être à la fois dans une

fissure à l'instant t et à l'instant (t + t) suivant.

Elle possède les propriétés simples suivantes C('r = 0) = < a(t) a(t) > = < ú(t) > (car ú(t) 2 = E(t)) Cette moyenne n'est autre que la probabilité pour la particule d'être dans une fissure à l'instant t, soit la proportion of, donc on a: C(T = 0) = of (8) C(' = -o) = Lim (ô -->) (< s(t) ú(t + ') > = < e(t) > <s(t + Yc) > = Of 2 (car les positions de la particule sont indépendantes si T est assez grand), soit C(c =o) = Of 2 (9) Il est plus simple de considérer la fonction C(Tr) normalisée variant de 1 à O, que l'on note: C('r) = (C(T) - f2)/of ( -_ f) =(C(z) - Of 2)f m (1 0) En utilisant maintenant les probabilités conditionnelles, C(T) s'écrit sous la forme: I 5 C(Tr) = <(t s(t) (t +) > ú(t) = I < s(t) >, o l'indice ú(t) = 1 indique que l'on sait que e(t) = I à t. Comrnme < ú(t) > = of correspond à la probabilité pour que la particule soit dans une fissure à l'instant t, la fonction C(t) / Of n'est autre que la probabilité de trouver la particule dans une fissure à l'instant t sachant qu'elle y était à l'instant 0 (cette nouvelle interprétation permettant d'ailleurs de retrouver les deux limites précédentes

très simplement).

Cette dernière fonction possède une interprétation en termes de solution de l'équation aux dérivées partielles (1). Il s'agit de trouver la pression P(r, t) avec pour condition initiale P(r, t = 0) = 1 si r est dans une fissure, ou P(r, t = 0) =0, sinon. Ces deux premières conditions traduisent la condition: la particule est dans les

fissures à l'instant initial 0.

Alors on a C(t) = I P(r,-)drl) fissures Physiquement, ce problème revient à mettre l'ensemble de la zone fissurée au potentiel unité à t = 0, la matrice étant au potentiel nul, et à laisser diffuser la pression ensuite suivant l'équation (1). La fonction C caractérise la décroissance de la pression moyenne dans les fissures, ce qui explique pourquoi la connaissance de la fonction C donne des informations sur la fonction de transfert (matrice vers fissures), et correspond dans l'image particulaire à la fraction du nombre total de particules dans les fissures à

l'instant t.

Il est possible comme on va le voir, de calculer la fonction C(T) à partir d'une seule réalisation du processus aléatoire. Pour cela, on utilise l'ergodicité des trajectoires aléatoires qui permet d'écrire l'égalité suivante: f E(t)e(t + r)dt C(t) = Lim(t, --oo) o (12) max -

Evaluation sur une seule trajectoire de la fonction ú(t).

L'évaluation de l'intégrale est rendue aisée si l'on se souvient que E(t) vaut zéro ou

1, et correspond à l'évaluation du produit de la fonction ú(t) et de sa translatée de t.

Comme la fonction s(t) est définie par le tableau des dates o la particule a changé de milieu, si ti ti+l correspond à un intervalle o e(t) vaut I (un intervalle sur deux), il suffit de calculer l'intégrale T e(t +r)dt ce qui revient à évaluer les instants de changement d'état de la

fonction ú(t + r) dans l'intervalle considéré.

L'ensemble du calcul pour un nombre arbitraire de valeurs de t peut être effectué

rapidement à l'aide d'un programme judicieux d'un type connu.

Un exemple de fonction normalisée C(t) est représentée sur la Fig 2. Le milieu considéré est un milieu en feuillets constitué de fissures planes espacées les unes des

autres d'une distance L (Fig. 3).

Il est possible à présent, de faire le lien avec le système d'équations (2) de Warren 1 5 & Root transitoire, avec comme condition initiale: Pf(r, t = 0) =1 Pm(r, t = 0) = 0, sinon, o Pm et Pf correspondent à des pressions "à grande échelle" qui sont les

moyennes des pressions dans la matrice et dans les fissures.

En intégrant ce système par rapport à la variable r sur tout le milieu, on obtient un système différentiel décrivant les variations des pressions moyennes dans les fissures et la matrice qui s'écrit: 0 dPf (t) _ dr q,(rt) 0cdP (13) Pa" cf ()_ |drq.,(r.t) En termes particulaires, ces grandeurs Pf(t) et Pm(t) correspondent au nombre de

particules dans chacune des deux régions.

Dans le modèle transitoire, le terme source qm(r, t) est relié à Pf(r, t) par une convolution temporelle qm(r, t) = Gm(t) * Pf(r, t), l'hypothèse sous-jacente étant le le

flux d'échange en un point ne dépend que de la valeur de Pf(r, t) en ce même point.

Après intégration, on a donc 0Q1., ' = G m, (t) *Pr ( t) dt c dP t) =-G (t)* Pf (t) dt (14) Ce système peut maintenant être intégré par une transformation de Laplace g(s) =e 'tg(t)dt (15) On trouve: p. (s) = (-' Yf(s) (16) avec s f(s) = s f - Gm(s) Pour le modèle Warren & Root quasi-permanent, l'expression connue de la fonction f(s) est: s + cD,,, (c)ms + oD,,, (17) Avec Dm = Km/ opct Dans ce cas particulier, la fonction Pf(t) est une exponentielle décroissante avec un temps de relaxation égal à. Ce temps correspond au temps typique d'équilibrage des pressions dans le milieu, et est de l'ordre du temps de diffusion typique

d'une particule sur un bloc matriciel d'après les estimations de o(.

Cette fonction qui varie de I à Of quand s parcourt l'intervalle (0, o) (ou de Of à 1 quand t varie de 0 à l'infini) représente une "porosité effective". On peut démontrer de façon générale que la solution dans l'espace de Laplace des essais de puits en milieu fissuré, se déduit des solutions "en simple milieu" par le remplacement de l'argument s

par s f(s).

1 5 On montre que pour un modèle plus général de terme d'échange, non réductible à la simple différence de pression (3), la seule connaissance de la fonction f(s) suffit. Cette fonction a été calculée analytiquement dans quelques géométries simples (feuillets, sphères), comme décrit par: A. de Swaan O., Analytic solutions for determrnining naturallv fractured reservoir

properties bv well testing, trans. AIME, 1976.

En utilisant les équations (10), (11) avec le résultat (16), on a: C(s) =,(/ 1 -- s(l,,(f--])(8 Cela signifie que la connaissance de C(t) (donc de C(s)) permet de déterminer f(s). Détermination du coefficient cf. Pour déterminer ar, on va utiliser un développement limité de C(s) au voisinage de s = 0 en utilisant la forme (17) de f(s). On a: ( = o) -f 0ón, (19) Lg0),, Cette propriété, ici vérifiée dans le cas pseudo-permanent est aussi exacte dans les cas transitoires en utilisant les formules données par A. de Swaan. Elle correspond à une approximation "aux temps longs"; correspondant à la limite s petit. Or aux temps longs, le modèle quasi-permanent est asymptotiquement exact, d'o la généralité de ce résultat. Comme, par définition de la transformation de Laplace on a: C(s) = C(t)dt, on on obtient finalement: .() re(. 1 (D", (D l"" D,,, C(t)dt (20) Comme la méthode des marches aléatoires permet de calculer la fonction C(t),

une quadrature numérique permet de remonter au coefficient cc par cette relation.

Validation et tests On considère un milieu "feuilleté" (Fig. 3), o les fissures F sont supposées être des plans parallèles régulièrement espacés d'une distance L dans une matrice M. Dans ce cas, on montre que a = Pour mettre en oeuvre la méthode des marches aléatoires, on a considéré une cellule élémentaire de taille L = 1, supposée maillée 10 par 10, les mailles de la première rangée de mailles étant supposées être des mailles de fissure. On a donc Of = 0,1 et on a

pris un contraste de perméabilité de 100 entre les fissures et la matrice.

Sur un même graphe ont été représentées les variations de la fonction C(t) obtenue par un logiciel de mise en ceuvre de la méthode (Courbe RW) et parune inversion de Laplace de la relation (18) (Courbe LA) utilisant la formule analytique définie dans l'article de Swann déjà cité, valable pour les feuillets, en prenant a c'est-à-dire égal à 12 ici. On observe un excellent accord entre les deux courbes. Les coefficients ca ont été calculés à l'aide de la formule (20) par une quadrature numérique, et les valeurs en sont indiquées sur le graphe. L'erreur que l'on a pu noter sur la valeur de ces coefficients obtenue par la méthode, est du même ordre que celle provenant de la solution exacte; elle est donc a priori imputable à la procédure de quadrature numérique

elle-même.

On ne sortirait pas du cadre de l'invention en appliquant la méthode qui vient d'être décrite à un milieu traversé par un nombre N (N >1) de réseaux de fissures

différents auxquels on associerait respectivement N fonctions d'état sN(T).

Claims (3)

REVENDICATIONS

1) Méthode pour modéliser, à partir d'une image 3D d'un milieu fissuré tel qu'une formation géologique constituée d'une matrice poreuse traversée par au moins un réseau de fissures, des paramètres hydrauliques caractérisant au mieux à grande échelle le milieu fissuré, dans le but d'interpréter des essais de puits, comportant la discrétisation du milieu par un maillage et la résolution dans ce maillage d'équations modélisant la diffusion de la pression dans le milieu, caractérisée en ce qu'elle comporte la détermination de la variation au cours du temps d'une fonction de transfert à grande échelle (f(t), f(s)) caractérisant les flux de fluide de la matrice vers le dit réseau de fissures, par application aux particules de fluide de déplacements aléatoires en temps continu, en leur associant une fonction d'état (s(T)) indicatrice de leur présence soit dans

la matrice soit dans le dit réseau de fissures.

2) Méthode selon la revendication 1, caractérisée en ce que l'on détermine une fonction (C(t)) caractérisant la décroissance au cours du temps de la pression moyenne de fluide dans les fissures du dit réseau, liée à la probabilité de présence des particules, par le calcul de la fonction d'auto-corrélation temporelle de ladite fonction d'état (E(t)),

et on en déduit ladite fonction de transfert (f(s)) ainsi que le coefficient d'échange (c).

3) Méthode selon l'une des revendications précédentes, caractérisée en ce que, le

milieu comportant N (N >l) réseaux de fissures, elle comporte la détermination de la variation au cours du temps de fonctions de transfert à grande échelle (f(t), f(s)) caractérisant les flux de fluide de la matrice vers les dits réseaux de fissures, par application aux particules de fluide de déplacements aléatoires en temps continu, en leur associant N fonctions d'état (SN(T)) indicative de leur présence soit dans la matrice soit

dans les dits N réseaux de fissures.