Brevet FR2757947A1 — publié le 19980703

L'invention concerne une niéthode destinée à déterminer la perméabilité de fissure équivalente d'un réseau fisurè dans un milieu multicouche souterrain t'issOk. permettant de réaliser une modélisation plus réaliste d'une structure géologique souterraine fissurée. La méthode peut être mise en oeuvle par exemple par des ingénieurs de gisement en ue d'obtenir des prédictions d'écoulement d'huile fiables. Les réservoirs fissurés constituent un type extrême de réservoirs hétérogènes 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 de fissure représentant moins de I ic de l'huile en place et hautement conducteur. Le milieu de fissure lui-même peut être complexe, avec différents ensembles de fissures caractérisés par leur densité, longueur, orientåtion, inclinaison et ouverture respectives. Les images en 3D de réservoirs fissurés ne peuvent pas être utilisées directement sous forme de données d'entrée de simulation de réservoir. La représentation du réseau de fissures au sein de simulateurs d'écoulement de réservoir a été longtemps considérée comme irréaliste car la configuration du réseau est en partie inconnue et à cause des limitations numériques liées à la juxtaposition de nombreuses cellules présentant des dimensions et des propriétés extrêmement contrastées. C'est pourquoi la modélisation simplifiée mais réaliste de tels milieux présente un grand intérêt pour les ingénieurs de gisement. L"'approche double porosité" telle qu'elle est enseignée, par exemple, par Warren J.E. et al dans "The Behavior of Naturally Fractured Reservoirs", SPE Journal (Septembre 1963), 245-255, est connue dans l'art pour interpréter le comportement d'un écoulement monophasique observé en testant un réservoir fissuré. Selon ce modèle de base. tout volume élémentaire du réservoir fissuré est modélisé sous la forme d'un ensemble de blocs parallélépipédiques identiques limités par un système orthogonal de fissures uniformes continues orientées dans la direction de l'un des trois principaux sens de l'écoulement. L'écoulement des fluides à l'échelle du réservoir effectue à travers le milieu de fissure seulement et des échanges de fluide surviennent localement entre les fissures et les blocs matriciels. De nombreux simulateurs de réservoirs fissurés ont été développés en utilisant un tel modèle. avec des amélioratlons spécifiques relatives à la modélisation des échanges de fluide matrice-fissure régis par des forces capillaires. gravitationnelles. visqueuses et des mécanismes compositionnels. ainsi que des échanges matrice-matrice (simulateurs double perméabilité double porosité). De nombreux exemples de techniques de l'art antérieur sont cités dans les références indiquées ci-après. Thomas. L.K. et al . "Fractured Reservoir Simulation". SPE Journal (Fé'rier 1953) 42-54: Quandalle, P. et aI : "Typical Features of a New Multipurpose Reservoir Simulatol-". SPE 16007. présenté au 9ème Symposium sur la simulation de réservoir de San Antonio. Texas, 1-4 février 1987. Coats. K.H. "Implicit Compositional Simulation of Single-Porosity and Dual Porosity Reservoirs", SPE 18427 présenté au Symposium SPE sur la simulation de réservoir de Houston, Texas, 6-8 février 1989. L'un des problèmes que rencontrent les ingénieurs de gisement consiste à doter ce modèle de base de paramètres en vue d'obtenir des prédictions d'écoulement fiables. En particulier, les propriétés pétrophysiques de base des fissures et des matrices, ainsi que la taille des blocs matriciels, doivent être connues pour chaque cellule du simulateur d'écoulement. Alors que la perméabilité matricielle peut être évaluée à partir de carottes. la perméabilité du réseau de fissures contenu dans la cellule, c'est-à-dire la perméabilité de fissure équivalente, ne peut pas être évaluée simplement et nécessite la prise en compte de la géométrie et des propriétés du réseau de fissures. Une méthode directe connue permet de déterminer l'écoulement en régime stable dans un réseau de fissures. Elle comprend l'utilisation de mailles fines et régulières conventionnelles discrétisant à la fois les fissures et les blocs matriciels du volume de roche fissurée parallélépipédique considéré. Pour différentes raisons. cette méthode connue ne fournit pas de résultats fiables, sauf si le volume de roche fissurée est discrétisé en utilisant une maille comportant un nombre extrêmement élevé de cellules. ce qui nécessitc des moyens de calcul très importants. D'autres modèles spécifiques servant à calculer les perméabilités équivalentes de réseaux fissurés en 2D ou 3D sont également connus. par exemple par - Odling. N.E. ''Perme;lbility of Natural and Simulated Fracture Patterns". Structural and Tectonic Modelling and its Application to Petroleum Geology. NPF Special Publication 1,365-380. Elsevier. Norwegian Petroleum Society (N'PF) 199': - Long. J.C.S. et al : "A Modes for Steady Fluid Flou in Random Three-Dimensional Networks of Disc-Shaped Fractures", Water Resources Research (août 1985 > . vol.9 I. n 8. 1105-1115: - Cacas. M.C. et al : ",NIodeling Fracture Flow with a Stochastic Discrete Fracture Network : Calibration and Validation. 1. The Flow Model". Water Resources Research (mars 1990) vol.26, n 3; - Billaux, D. : Hydrogéologie des milieux fracturés. Géométrie, connectivité et comportement hydraulique". Thèse de doctorat préentée à l'école Nationale Supérieure des Mines de Paris, document du BRGM n 186, Editions du BRGM, 1990; - Robinson, P.C. "Connectivity, Flow and Transport in networks Models of Fractured Media", Thèse de doctorat, St Catherine's College, Oxford University, Réf. TP1072, mai 1984. La présente invention concerne une méthode de détermination de la perméabilité de fissure équivalente d'un réseau de fissures dans un milieu multicouche souterrain. La méthode est caractérisée en ce qu'elle comporte les étapes suivantes - discrétisation du réseau de fissures en éléments de fissure (tels que des rectangles par exemple) et définition de noeuds représentant les éléments de fissure interconnectés dans chaque couche du milieu, et - détermination des écoulements de fluides à travers le réseau discrétisé en imposant des conditions de pression limites et des transmissivités de fluide à chaque couple de noeuds voisins. Plus précisément, la méthode est caractérisée en ce que: - I'on divise le milieu cn une série de couches parallèles s'étendant chacune dais un plan de référence perpendiculaire à l'axe de référence et définies chacune par une coordonnée le long dudit axe, - on divise chaque fissure en une série de rectangles limités le long dudit axe de rétëlence par deux couches adjacentes et on catalogue les rectangles en leur associant des caractères géométriques et physiques tels que coordonnées et tailles des rectangles et conductivités hydrauliques des fissures, - on positionne des noeuds dans chaque couche pour toutes les fissures interconnectées, et - pour tous les couples de noeuds voisins, on calcule des facteurs de transmissivité et résolution des équations d'écoulement en vue de déterminer les perméabilités équivalentes du milieu dans trois directions orthogonales. Dans un mode de réalisation préféré, la perméabilité équivalente du milieu comprend ia détermination directe du tenseur d'anisotropie de perméabilité équivalente et l'étalonnage de valeurs absolues de perméabilité à partir de résultats d'essais de puits. La méthode ainsi résumée permet de relier de manière systématique des modèles de caractérisation de réservoirs fissurés et des simulateurs double porosité en vue de réaliser une modélisation plus réaliste d'une structure géologique souterraine fissurée. La méthode peut être mise en oeuvre par exemple en production pétrolière par des ingénieurs de gisement en vue d'obtenir des prédictions d'écoulement fiables. D'autres caractéristiques et avantages de la méthode selon l'invention apparaîtront plus clairement à la lecture de la description ci-après de modes de réalisation donnés à titre d'exemples non limitatifs, en se référant aux schémas figurant en annexe parmi lesquels - la figure I montre par exemple une image en 3D d'un réseau de fissures stochastiquement obtenue à partir d'observations et de mesures prises sur un affleurement de gres. - les figures 7. 3 montrent une fissure divisée en une série de rectangles R. - la figure 4 montre la structure des données d'entrée spécifiant les caractères de la fissure. - la figure 5 montre le mode préféré de discrétisation d'un plan de fissure. - les figures 6. 7. 8 et 9 illustrent schématiquement le calcul des facteurs de transmissivité pour différentes positions de noeuds les uns par rapport aux autres ou avec une limite. Les perméabilités équivalentes d'un réseau de fissures 3D sont déterminées ci-après en utilisant une technique numérique basée sur la méthode connue de "réseau de résistance" décrit par exemple dans le document d'Odling N.E. cité ci-avant. Dans la présente méthode, on admet que la matrice est imperméable afin qu'elle soit compatible avec l'approche double porosité. Dans les simulateurs de réservoir, les écoulements matricefissure et matrice-matrice sont en fait calculés séparément par rapport aux écoulements au sein des fissures. On admet que le réseau de fissures en 3D considéré représente. dans un volume égal à une cellule de réservoir, la distribution réelle des fissures donnée par intégration des caractères de fissure du champ dans un modèle de caractérisation. Le principal objectif de calculs d'écoulements monophasiques dans le réseau de fissures en 3D est d'évaluer l'anisotropie de perméabilité équivalente (Kv/Kh et Ky/Kx) de la cellule de fissure considérée, qui constitue un paramètre important régissant le comportement des écoulements polyphasiques dans les réservoirs. Les valeurs de perméabilité équivalente obtenues par ces calculs pourraient en pratique être comparées aux résultats d'essais de puits en vue d'étalonner les caractères de fissures tels que les conductivités hydrauliques de fissure (ou les ouvertures hydrauliques équivalentes). qui ne peuvent être définies que de manière médiocre a priori. De plus. les résultats dc perméabilité équivalente peuvent être utilisés en vue dc déterminer un tenseur de perméabilité dont les principales directions permettent tine orientation optimale de la maille du modèle de réservoir. Cependant. des conditions limites spécifiqucs sont necessaires pour obtenir ces informations Des contlittons limites sans écoulement latéral imposées sur les quatre faces latérales du volume parallèlépipédique étudié ne donnent pas accès aux termes non diagonaux du tenseur de perméabilité équivalente. alors que des potentiels (ou des pressions; à 'aviation linéale SUI les faces latérales permettent d'imposer la direction du gradient de potentiel au selon du milieu anisotrope et de déduire diirctement les termes de perméabilité non diagonaux des débits latéraux. Les techniques d'intégration de données de fissuration naturelle dans des modèles de réservoirs fissurés sont bien connues. Les données de fissuration sont essentiellement dc nature géométrique et comportant des mesures de densité. de longueur. d'azimut et d'inclinaison des plans de fissure soit observés sur des affleurements. dans des galeries de mines, sur des carottes. soit déduits à partir de diagraphies. Différents ensembles de fissures peuvent être distingués et caractérisés par différentes distributions statistiques de leurs caractères de fissure. Une fois les schémas de fissuration caractérisés. on peut créer des réseaux numériques de ces ensembles de fissures en utilisant un procédé stochastique respectant les distributions statistiques des paramètres de fissures. De tels procédés sont décrits par exemple dans les brevets FR-A-2.725.814; 2.725.794 ou 2.733.du demandeur. La méthode selon l'invention s'applique à des images de structures géologiques fissurées présentant différentes tailles ou volumes et/ou situées en différents endroits. créées par un générateur de modèle de fissure. La figure 1 présente une telle image. DONNEES D'ENTREE Avant de développer les procédures recommandées en vue de déterminer les paramètres hydrauliques équivalents d'images de fissures en 3D, une étape importante consiste à définir une structure commune de données d'entrée pqur ces images afin qu'elles puissent être traitées indépendamment de l'outil de traitement utilisé pour leur création. Comme le montrent les figures 2. 3, on admet que les fissures F sont sensiblement verticales (c'est-à-dire perpendiculaires aux limites de la couche). Cependant. une nlcllle structure de données peut être appliquée aux fissures s'écartant légèrement de la verticale. L'image 3D est discrCtisée verticalement en respectant la stratification géologique réelle lorsque cette information est disponible. Dans le cas cintraire. on applique une discrétisation arbitraire à l'image. Chaque couche horizontale L est caractérisée par sa coordonnée verticale zL dans le système de référence de coordonnées > ON. OY. OZ). Une série de rectangles R doit être définie pour chaque couche L. Chaque rectangle consiste en un élément de plan de fissure compris entre les limites d'une couche donnée. Ainsi. chaque fissure naturelle consiste en une série de rectangles superposés R et on lui attribue une origine (origine de la fissure). Chaque rectangle est défini par: - les trois coordonnées (xO, yO, zO) de l'origine du rectangle O. Pour une fissure naturelle donnée, tous les points d'origine des rectangles constitutifs sont situés sLir la même ligne verticale (ou plus fort pendage) tracée depuis l'origine de la fissure: - les coordonnées du vecteur d'unité horizontal i (xH, yH) et du vecteur d'unité vertical j (xV, yV) définissant l'orientation du rectangle dans le système de référence de coordonnées, avec xVertical et yVertical égaux à zéro en cas de fissures verticales mais considérés comme des données d'entrée afin de pouvoir traiter des fissures non verticales; - les deux longueurs horizontales algébriques 1- et 1+ séparant l'origine du rectangle et les deux limites latérales (verticales) de ce rectangle; - la hauteur h du rectangle, c'est-à-dire la longueur du rectangle dans le sens j. qui est l'épaisseur de la couche si la discrétisation en sens j correspond à la géologie; - la conductivité hydraulique c dérivée de l'application de la loi de Darcy relative à AP l'écoulement dans les fissures (pour un gradient de pression , le débit au sein de la ch A P fissure de hauteur h est -. , tt étant la viscosité du fluide). La conductivité c est donnée par la relation c=k.a. où k = a2/l 2 (en utilisant la représentation idéalisée clefs fissures de Poiseuille) est la perméabilité intrinsèque de la fissure et a son oLIvcltulc hydraulique équivalente. La conductivité hydraulique c est une valeur dc rêlèrence donnée pour une direction de la contrainte horizontale maximale parallèle au sens de la fisstire: - les deux rectangles sLipérieLir et inférieur voisins l UR. LR; - la série de fissures FS à laquelle appartient le rectangle considéré: - l'angle d'orientation ao de la direction de la contrainte horizontale maximale prise à partir de l'axe (OX) dans le système de référence de coordonnées: - pour chaque série de fissures, une table de corrélation corrélant I ) l'angle compris entre la direction de la contrainte horizontale maximale et la direction de la fissure (azimut) avec 2) la conductivité hydraulique c ou l'ouverture hydraulique équivalente a précédemment définies. Les termes "horizontal" et "vertical" utilisés dans ce contexte se rapportent aux directions respectivement parallèles et perpendiculaires aux limites de couches considérées ici comme étant horizontales. Les limites de couche discrétisent les plans de fissure dans le sens "vertical" Il convient de noter que les données d'entrée mentionnées ci-dessus I ) conviennent pour tous les outils logiciels existants utilisés pour caractériser et générer des fissures et 9) qu'elles pourraient ête utilisées pour discrétiser un réseau de fissures légèrement non verticales, c'est-à-dire qui ne sont pas perpendiculaires aux limites de couche. MODES OPERATOIRES Les modes opératoires et les tests de validation de la méthode permettant de calculer l'anisotropie de perméabilité du réseau de fissures dans son ensemble sont décrits ci-après à partir de l'image en 3D ainsi codifiée. La procédure numérique destinée à calculer les perméabilités équivalentes d'un réseau de fissures 3D est également décrite. Le problème consiste à trouver la distribution des débits dans le réseau pour les conditions frontières suivantes aux limites du volume parallélépipédique étudié, c'est-à- dire les pressions fixes imposées sur deux faces opposées et les pressions variant de manière linéaire sur les quatre faces latérales (entre les valeurs inlposees sur les deux autres faces) Les principales étapes sont résumées ci-dessous 1) Discrétisation du réseau En utilisant les définitions données pour la structure des données d'entrée. le résealLI cie tissure est discrétlsé sous la forme d'une série de "noeuds" N. chaque noeud étant placé au centre des segments d'intersection IS de deux rectangles R (c'est-à-dire de deux plans de fissure au sein d'une couche donnée). Comme le montre la figure 5. des noeuds supplémentaires AN sont placés au-dessus et au-dessous des noeuds précédents N afin de représenter d'autres rectangles discrétisant les fissures et de minimiser les longueurs d'écoulement au sein d'une fissure donnée. Sur la figure 5, BL est une limite latérale de deux cellules de fissure voisines. Une fois le réseau discrétisé, une procédure de tri est appliquée à ce réseau de fissures en vue d'éliminer les noeuds isolés ou des groupes de noeuds sans lien avec l'une des limites latérales FL du volume 3D étudié, car les fissures ainsi "triées" ne contribuent pas au transport du fluide et peuvent entraver les procédures de résolution mises en oeuvre en vue de trouver les pressions aux noeuds de fissures lors d'un écoulement en régime stable à travers le réseau. 2) Calcul des transmissivités Un facteur de transmissivité T est calculé pour chaque paire de noeuds associés en utilisant la relation ch kah T = - = -- i i où c est la conductivité hydraulique de la fissure, k la perméabilité intrinsèque de la fissure, a l'ouverture de la fissure, h la hauteur de la fissure et I la distance entre deux noeuds d'une fissure. Différentes situations doivent être envisagées selon la position respective des deux noeuds. Pour des noeuds contenus au sein de la même couche (fig.6 > . le facteurs' dc transmissivité horizontale T est obtenu directement sous la forme de la distance (11+12 > séparant les deux noeuds dans le sens de l'écoulement (fi.7). Pour des noeuds SitLiés dans deux couches difi;e'rentes (Fie ())* le facteur de transrllissivité horizonLale est la somme arithmétique des facteurs de transmissivité (T'+T") relatifs aux deux éléments du plan de fissure des cellules de fissure superposées. II comporte une longueur d'écoulement égale à la demi-somme des deux épaisseurs de couche hl et hl. Pour les noeuds supplémentaires tels que dètinîs précédemment. reliés pal un seul plan de fissure. un seul facteur de transmissivité est calculé pour cet élément du plan de fissure. Le facteur de transmissivité T entre un noeud et une limite du volume 7D étudié est exprimé comme dans le cas de deux noeuds. selon les deux cas suivants. Pour une limite verticale latérale, le facteur de transmissivité T peut être exprimé directement pour un seul élément du plan de fissure (fig.8), et sous la forme de la somme de deux transmissivités si deux plans de fissure relient le noeud et la limite. Pour une limite inférieure ou supérieure horizontale. le facteur de transmissivité verticale peut être exprimé en considérant une longueur d'écoulement égale à la moitié de l'épaisseur de la couche (fig.9). 3) Equations d'écoulement En régime stable, un écoulement monophasique incompressible à travers le réseau de fissure est déterminé en résolvant un série de n équations, une pour chaque noeud, comme il est connu dans l'art. Chaque équation exprime le fait que le débit total est nul à chaque noeud de fissure. Pour calculer un tenseur de perméabilité, on considère qu'une pression constante est imposée à chacune des limites aval et amont. Une pression variant de manière linéaire en fonction de la position entre les limites aval et amont est imposée. La matrice de perméabilité équivalente (Kij) déterminée précédemment est diagonalisée pour calculer les directions principales de l'écoulement avec les perméabilités équivalentes respectives dans ces directions. En pratique, le problème se limite souvent à la découverte des principales directions horizontales de l'écoulement U et V dans la mesure où la direction perpendiculaire aux limites de couche t ènéralement verticales) est toujours prise comme l'axe des z. Dans un tel cas, seuls les termes extra-diagonaux Kxy et Kyx doivent être calculés: ils peuvent être obtenus avec les conditions limites mixtes suivantes - les écoulements horizontaux sont calculés avec des faces Ints-rieui-es et supérieures mperméahles et des pressions à variation linéaire sur les faces verticale parallèles au .sens de l'écoulement. - I'écoulement vertical est calculé avec toutes les faces latérales imperméables. Ainsi. un tenseur de perméabilité simplifié est obtenu. à partir duquel on peut facilement déduire les principales directions horizontales de l'écoulement U et V Validation La méthode a été validée avec succès par rapport aux calculs d'écoulement monophasique de référence mentionnés ci-avant, réalisés avec un simulateur de gisement conventionnel. Les calculs de référence ont été obtenus sur de fines mailles régulières discrétisant les fissures, ainsi que les blocs matriciels du volume de roche fissurée parallélépipédique considéré. Pour un sens d'écoulement donné, une pression d'injection et une pression de production fixes ont été imposées sur les faces d'entrée et de sortie et le débit résultant a été calculé pour des conditions exemptes d'écoulement latéral. On a réalisé trois étapes pour valider le calcul de - la perméabilité verticale équivalente d'un volume de roche traversé par une seule fissure, cette dernière étant représentée par plusieurs noeuds correspondant aux intersections avec de petites fissures isolées; - les perméabilités horizontales équivalentes (dans une géométrie d'écoulement en 2D > et les principaux sens d'écoulement: - les perméabilités équivalentes et l'anisotropie de perméabilité dans Lin seul éseaLi présentant une géométrie d'écoulement en 3D. Les résultats obtenus pour la troisième étape (pour une géométrie d'écoulement en 3D > sont donnés dans le tableau ci-après. Pour les directions horizontales d'écoulement. une solution analytique de référence peut également être calculée puisque la géométrie de l'écoulement est un écoulement 2D dans ces directions (la géométrie d'écoulement en 3D concerne la direction z). Perméabilités Simulation METHODE Solution équivalentes (md) MAILLE FINE PRESENTEE ANALYTIQUE Kx 0,119 0,120 0,120 Ky 0,2240,227 0,226 Kz 0,255 0,267 Anisotropie Kz/(KxKy)0-5 1,56 1,62 Il apparaît clairement que les résultats obtenus avec la méthode présentée sont très proches des valeurs correspondantes obtenues avec la solution analytique et la simulation à maille fine pour les directions X et Y. De plus, la différence entre les valeurs de perméabilité équivalente verticale concernant l'écoulement 3D reste acceptable. Ainsi, le rapport d'anisotropie, égal à 1,6, est prédit de manière satisfaisante par la méthode avec un nombre très limité de cellules. La méthode selon l'invention, qui fournit une représentation facilement transposable d'un réseau de fissures naturel, se prête bien au calcul d'écoulements au sein de fissures. Elle peut également être utile pour améliorer l'image d'origine du réseau de fissures. Une telle image est en fait obtenue à partir d'un générateur de fissure stochastique utilisant comme données d'entrée les résultats de l'intégration de données de fissuration archivées dans un modèle de caractértsation de fissure comme il est décrit dans les brevets FR-A2.725.814: 2.725.794 ou ' 733.073 du demandeur déja mentionnés. De telles Images. une fois discrétisées par la méthode selon l'invention. peuvent être facilement modifiées afin d'être adaptées aux règles géologiques. Par exemple. l'interruption systématiqLie d'une fissure donnée contre une autre série de fissures peut être prise en compte dans l'image d'origine en annulant les éléments du plan de fissure d'une série donnée 'étendant au-delà des fissures entrecoupées de l'autre série. REVENDICATiONS 1) Méthode permettant de déterminer la perméabilité de fissure êqui'alente d'lll réseau de fissures dans un milieu multicouche souterrain fissuré. à partir diane représentation prédéterminée dudit réseau. comportant les étapes de - discrétisation de chaque fissure (F) du réseau de fissures en éléments de fissure et de définition de noeuds (N) représentant des éléments de fissure interconnectés dans chaque couche du milieu, et - détermination d'écoulements de fluides à travers le réseau de fissures en imposant des conditions limites de pression et de transmissivités de fluides pour chaque couple de noeuds voisins. 2) Méthode selon la revendication I, caractérisée en ce que: - I'on divise le milieu en une série de couches parallèles s'étendant chacune dans un plan de référence perpendiculaire à l'axe de référence et définies chacune par une coordonnée le long dudit axe, - on divise chaque fissure en une série de rectangles limités le long dudit axe de référence par deux couches adjacentes et on catalogue les rectangles en leur associant des caractères géométriques et physiques tels que coordonnées et tailles des rectangles et conductivités hydrauliques des fissures, - on positionne des noeuds dans chaque couche pour toutes les fissures interconnectées, et - pour tous les couples de noeuds voisins, on calcule des facteurs de transmissivité et résolution des équations d'écoulement en vue de déterminer les perméabilités équivalentes du milieu dans trois directions orthogonales. 3) Méthode selon l'une des revendications précédentes, caractérisée en ce que la perméabilité équivalente du milieu comprend la détermination directe d'un tenseur de perméabilité équivalent et l'étalonnage de valeurs absolues de perméabilité à partir de résultats d'essais de puits. A) Méthode selon l'une des revendications précédentes. caractérisée cn ce que lesdits caractères géométriques et phvsiques sont les coordonnées et les dimensions des rectangles et les conductivités hydrauliques des fissures.