Brevet FR3021141A1 — publié le 20151120

PROCEDE DE DETERMINATION D'UNE CARTE DE HAUTEUR D'HYDROCARBURE LIQUIDE DANS UN RESERVOIR La présente invention concerne le domaine de la détermination des réserves en 5 hydrocarbures liquides dans un réservoir et, notamment, le domaine de la détermination de cartes de hauteur d'hydrocarbure liquide dans un réservoir. Afin de déterminer facilement les réserves en hydrocarbures (gaz, pétrole ou huiles) des réservoirs en exploitation ou à exploiter, les géologues manipulent et 10 étudient le plus souvent des cartes en deux dimensions représentant le « HuPhiSo » (dites cartes « HuPhiSo »). En effet, ces cartes sont synthétiques et renseignent rapidement le géologue quant à la valeur économique du réservoir (3D). Pour un volume rocheux en trois dimensions, la carte « HuPhiSo » (en deux dimensions) comprend une pluralité de valeurs, chaque valeur étant le produit de 15 plusieurs paramètres : HuPhiSo = h - NTG - cf) - So le paramètre h étant la hauteur de la colonne (selon l'axe 2> si la carte est selon le plan(5>c,9)), le paramètre NTG (ou « Net to Gross » en anglais) étant la fraction du volume rocheux favorable à la présence d'hydrocarbures, le paramètre cf) étant la 20 porosité de la roche et le paramètre So étant la saturation en huile des porosités. Ainsi, si V est le volume rocheux considéré, V - NTG - Phi - So est le volume d'huile contenue dans le volume rocheux. Une valeur « HuPhiSo » représente la « hauteur de colonne d'huile » dans une colonne considérée. Pour autant, cette valeur est reliée indirectement au volume de l'huile et peut représenter la « densité 25 surfacique d'huile ». Pour un modèle maillé en cellules parallélépipédiques, il est possible de calculer une valeur « HuPhiSo » pour la carte finale en sommant la valeur « HuPhiSo » de chaque cellule d'une colonne de ce modèle. Bien entendu, la détermination des valeurs « HuPhiSo » peut être simple si les paramètres h, NTG, Phi et So sont connues de manière certaine pour tout point/cellule du modèle. Néanmoins, en pratique, ces paramètres ne sont parfaitement connus qu'en quelques points isolés : les points des puits de forage. Bien entendu, il est possible de simuler ces paramètres pour les autres cellules grâce au croisement de diverses analyses de données (ex. sismiques) et de l'utilisation d'algorithmes de simulation stochastique basés sur des méthodes géostatistiques afin de « compléter » le modèle. A cela s'ajoute la distribution naturelle des propriétés pétro-physiques (impédance acoustique, porosité effectives, porosité des roches dans le NTG, etc.). Chaque grandeur physique peut donc être considérée comme une variable aléatoire, éventuellement non stationnaire (ex. la porosité diminue le plus souvent avec la profondeur). Il existe un grand nombre d'incertitudes qui pèsent sur les différents paramètres dont dépend la carte de HuPhiSo, de sorte que cette dernière peut varier considérablement d'une carte possible à une autre. Pour les prises de décision, les ingénieurs réservoir apprécient particulièrement les cartes d'espérance du HuPhiSo et les cartes d'écart-type (ou de variance) autour de l'espérance du HuPhiSo. Habituellement, ces cartes sont obtenues en simulant un grand nombre de cartes HuPhiSo. Pour chaque cellule de la carte, est alors calculée une moyenne et une variance en fonction de cellules correspondantes dans les cartes HuPhiSo simulées. Néanmoins, chaque simulation peut être coûteuse en ressources et en temps de calcul et il n'est pas rare que la détermination des cartes d'espérance et d'écart-type prenne plusieurs centaines d'heures. Il existe donc un besoin pour la détermination de carte de HuPhiSo moyen (espérance) et/ou de l'intervalle d'incertitudes correspondant à cette carte (variance ou écart-type) de manière plus rapide, sans réaliser un grand nombre de simulations de cartes HuPhiSo. La présente invention vient améliorer la situation. A cet effet, la présente invention propose une méthode analytique permettant, après certaines simplifications, de déterminer directement l'espérance de la hauteur d'hydrocarbure liquide dans un modèle géologique. La présente invention vise alors un procédé de détermination d'une carte d'espérance de hauteur d'hydrocarbure liquide dans un modèle géologique, le modèle géologique comportant des mailles pouvant être associée à un faciès parmi une pluralité de faciès et à un élément architectural parmi une pluralité d'éléments architecturaux. Le procédé comporte, pour au moins une colonne du modèle géologique : /a/ réception d'un cube de répartitions pour chaque élément architectural AE parmi la pluralité d'éléments architecturaux, le cube de répartitions décrivant pour chaque maille c du modèle la probabilité pAE(c) que ladite maille appartienne audit élément architectural, /b/ réception d'un cube de proportions pour chaque faciès parmi la pluralité de faciès et pour chaque élément architectural parmi la pluralité d'éléments architecturaux, le cube de proportions associé à un faciès A et à un élément architectural AE décrivant pour chaque maille c du modèle la proportion dudit faciès dans ladite maille si ladite maille appartient audit élément architectural AE, ladite proportion étant une variable aléatoire de distribution triangulaire définie par une valeur minimale pA,AE,min(c), une valeur maximale pA,AE,max(c), et un mode PA,AE,mode(C), /c/ détermination, pour chaque maille c de ladite colonne et pour chaque faciès A de la pluralité de faciès, de la somme sur la pluralité d'éléments 1 architecturaux de n c AE (l - n - A,AE,min(c) PA,AE,max(c) PA,AE,mode(c)). 3 /d/ détermination d'une valeur d'espérance de hauteur d'hydrocarbure liquide pour ladite colonne de manière analytique en fonction de la somme déterminée à l'étape /c/ pour chaque maille c de ladite colonne et pour chaque faciès A. Normalement, l'espérance de la valeur HuPhiSo est très difficile à calculer directement Esp(HuPhiSo) = Esp(h - NT G - (p - So) car l'ensemble des variables NTG, Phi Sa dépendent notamment des faciès et ainsi de la localisation des faciès dans le modèle. Le mode d'une distribution triangulaire peut également être appelé base. En faisant l'hypothèse de l'indépendance des paramètres NTG, Phi et Sa, il est possible de simplifier grandement le calcul de cette espérance. Expérimentalement cette hypothèse a été vérifiée et cette simplification, qui n'est pas la simple application de principes mathématiques, permet de mutualiser le calcul de la la somme (sur la pluralité d'éléments architecturaux) de pAE(c)-Î (PA,AE,min(c) + PA,AE,max(C) PA,AE,mode(C)) et éventuellement de calculer séparément les espérances de Esp(cpA(c)), Esp(NTGA(c)), Esp(SoA(c)). Cette simplification permet une meilleure efficacité lors des calculs à réaliser de l'espérance et ainsi améliorer sensiblement les temps de calculs nécessaires. En outre, chaque maille c du modèle pouvant être associée à une proportion de porosité pour chacun des faciès de ladite pluralité de faciès, ladite proportion de porosité pouvant être une variable aléatoire de distribution triangulaire définie par une valeur minimale cpA,min(c), une valeur maximale w A,max(C) et un mode (PA,mode(C) 5 le procédé peut comprendre en autre: /c-1/ détermination, pour chaque maille c de ladite colonne et pour chaque faciès A parmi la pluralité de faciès, de la valeur î ((PA,min(c) + (PA,max(c) + (PA,mode(c)), La détermination de l'étape /d/ peut alors être fonction de la valeur déterminée à l'étape /c-1/ pour chaque maille c de ladite colonne et pour chaque faciès A parmi la pluralité de faciès. Ainsi, la valeur Esp(cpA(c)) peut être aisément calculée, indépendamment des autres paramètres, si la variable cpA(c) comporte une incertitude. Dans un mode de réalisation possible, chaque maille c du modèle pouvant être associée à une proportion de volume rocheux favorable à la présence d'hydrocarbures pour chacun des faciès de ladite pluralité de faciès, ladite proportion de volume rocheux favorable à la présence d'hydrocarbures pouvant être une variable aléatoire de distribution triangulaire définie par une valeur minimale NTGA,min(c), une valeur maximale NTGA,max(c), et un mode NTGAmode(c), le procédé peut comprendre en outre : /c-2/ détermination, pour chaque maille c de ladite colonne et pour chaque faciès A parmi la pluralité de faciès, de la valeur î (NTGA,min(c) + NTGA,m'(c) + NTGA,mode(c)), La détermination de l'étape /d/ peut alors être fonction de la valeur déterminée à l'étape /c-2/ pour chaque maille c de ladite colonne et pour chaque faciès A parmi la pluralité de faciès. Ainsi, la valeur Esp(NTGA(c)) peut être aisément calculée, indépendamment des autres paramètres, si la variable NTGA(c) comporte une incertitude. Avantageusement, chaque maille c du modèle pouvant être associée à une saturation en hydrocarbure liquide, ladite saturation en hydrocarbure liquide pouvant être caractérisée par une cote d'interface hydrocarbure liquide-eau z, ayant une probabilité p(zw) et par une distribution triangulaire de la saturation en hydrocarbures pour chaque cote h au-dessus de cette cote d'interface, ladite distribution de saturation en hydrocarbure pouvant être définie par une valeur minimale CSHc,A,min(h), une valeur maximale CSHcAm'(h), et un mode C SHC,A,mode( h), le procédé peut comprendre en outre : /c-3/ détermination, pour chaque maille c de ladite colonne ayant une cote z et pour chaque faciès A parmi la pluralité de faciès, de la valeur f zw,min csHc,max (z - zw) + CSFic,mode(Z - Zw) CSFic,min(Z - Zw) Zw,max p(Zw) - dzw 3 La détermination de l'étape /d/ peut être fonction de la valeur déterminée à l'étape /c-3/ pour chaque maille c de ladite colonne et pour chaque faciès A parmi la pluralité de faciès. Par notation, il est choisi que la notation zw,min correspond à la cote la plus haute (avec un axe 2> vers le haut) en comparaison à la cote zw,m'. Ainsi, nous avons Zw,min > Zw,max en prenant cette notation. En effet, la cote de l'interface zw,min (respectivementzw,max) correspond à la cote pour laquelle la quantité d'hydrocarbure est la moins importante (respectivement la plus importante). Toutes les notations de cote zx,min et zxmax respectent le même principe (avec X pouvant être toute chaine de caractères). Dans un mode de réalisation particulier de l'invention, ladite saturation en hydrocarbure liquide pouvant être caractérisée par une cote d'interface hydrocarbure liquide-gaz zg ayant une probabilité p(z9), de distribution comprise entre une valeur minimale zemin, et une valeur maximale zemax, le procédé peut comprendre en outre : /c-4/ détermination, pour chaque maille c de ladite colonne ayant une cote z, de la valeur OSZ zemin Z Zemin izz p(za)dza si zemax emax 1 si z zemax La détermination de l'étape /d/ peut être fonction de la valeur déterminée à l'étape /c-3/ pour chaque maille c de cote z de ladite colonne et pour chaque faciès A parmi la pluralité de faciès multiplié par la valeur déterminée à l'étape /c-4/ pour chaque maille c de cote z de ladite colonne. En outre, la présente invention peut viser également un procédé de détermination d'une carte de variance de hauteur d'hydrocarbure liquide dans un modèle 25 géologique. le procédé peut comporter, pour au moins une colonne du modèle géologique : /e/ détermination d'une valeur d'espérance de hauteur d'hydrocarbure liquide dans le modèle géologique pour ladite colonne, comme exposée plus haut ; /f/ détermination d'une valeur d'espérance de hauteur d'hydrocarbure liquide pour ladite colonne de manière analytique en fonction de la somme déterminée à l'étape /e/. Dans un mode de réalisation, les valeurs des cubes de proportions pour deux cellules différentes du modèle et pour des éléments architecturaux différents peuvent être considérées indépendantes entre elles. En outre, le procédé peut comporter en outre : /g/ réception d'un corrélogramme p(A, Az) selon une direction de ladite colonne pour ledit faciès A dans un élément architectural donné ; /h/ détermination de l'espérance du produit de la présence du faciès A dans ledit élément architectural dans une maille cl par la présence du faciès A dans ledit élément architectural dans une maille c2, la distance entre cl et c2 selon ladite direction du corrélogramme étant Az, la probabilité de présence du faciès A dudit élément architectural pour la maille cl étant n i A,AE,c15 la probabilité de présence du faciès A dudit élément architectural pour la maille c2 étant n i A,AE,c2 , en fonction de : (1 - p(A, Az)) - (Esp(n \i- A,AE,c1). ESP(PA,AE,c2) + _\II 1 ar (pA,AE xi) ' V ar (p A,AE,c2) 1 \ 1 + p (A, Az) - (-2E sp ( n \i- A,AE,c1 j + - ESP (PA,AE,c2) 2 La détermination de l'étape /f/ peut alors être fonction de la détermination de l'étape /h/. Un dispositif destiné à la détermination d'une carte d'espérance de hauteur d'hydrocarbure liquide dans un modèle géologique peut être avantageux, en lui-même. Ainsi, la présente invention vise également un dispositif destiné la détermination d'une carte d'espérance de hauteur d'hydrocarbure liquide dans un modèle géologique, le modèle géologique comportant des mailles pouvant être associée à un faciès parmi une pluralité de faciès et à un élément architectural parmi une pluralité d'éléments architecturaux. Le dispositif comporte, pour au moins une colonne du modèle géologique : /a/ une interface pour la réception d'un cube de répartitions d'élément architectural pour chaque élément architectural AE parmi la pluralité d'éléments architecturaux, le cube de répartitions décrivant pour chaque maille c du modèle la probabilité pAE(c) que ladite maille appartienne audit élément architectural, /b/ une interface pour la réception d'un cube de proportions pour chaque faciès parmi la pluralité de faciès et pour chaque élément architectural parmi la pluralité d'éléments architecturaux, le cube de proportions associé à un faciès A décrivant pour chaque maille c du modèle la proportion dudit faciès dans ladite maille si ladite maille appartient audit élément architectural AE, ladite proportion étant une variable aléatoire de distribution triangulaire définie par une valeur minimale pA,AE,min(c), une valeur maximale pA,AE,max(c), et un mode PA,AE,mode(C), /c/ un circuit pour la détermination, pour chaque maille c de ladite colonne et pour chaque faciès A de la pluralité de faciès, de la somme sur la pluralité d'éléments architecturaux de P AE(C) - 7(P A,AE,min(c) PA,AE,max(c) PA,AE ,mode (0). /d/ un circuit pour la détermination d'une valeur d'espérance de hauteur d'hydrocarbure liquide pour ladite colonne de manière analytique en fonction de la somme déterminée par le circuit /c/ pour chaque maille c de ladite colonne et pour chaque faciès A. Un programme informatique, mettant en oeuvre tout ou partie du procédé décrit ci-avant, installé sur un équipement préexistant, est en lui-même avantageux. Ainsi, la présente invention vise également un programme informatique comportant des instructions pour la mise en oeuvre du procédé précédemment décrit, lorsque ce programme est exécuté par un processeur. Ce programme peut utiliser n'importe quel langage de programmation (par exemple, un langage objet ou autre), et être sous la forme d'un code source interprétable, d'un code partiellement compilé ou d'un code totalement compilé. La figure 3 décrite en détails ci-après, peut former l'organigramme de l'algorithme général d'un tel programme informatique. D'autres caractéristiques et avantages de l'invention apparaîtront encore à la lecture de la description qui va suivre. Celle-ci est purement illustrative et doit être lue en regard des dessins annexés sur lesquels : - la figure la illustre un exemple possible de modèle géologique ; - la figure lb décrit un exemple possible de la loi de probabilité peut être modélisée par une distribution triangulaire ; - la figure 2a décrit un exemple de coupe transversale d'un réservoir comportant de l'eau, du gaz et de l'huile ainsi que le taux de saturation en eau de la roche de ce réservoir. - la figure 2b décrit un exemple de courbe de saturation en eau de la roche d'un réservoir modélisé comme étant continue par morceau ; - la figure 3 est un exemple d'ordinogramme d'un mode de réalisation de l'invention ; - la figure 4 est un exemple de dispositif permettant la mise en oeuvre d'un mode de réalisation de l'invention ; 302 1 1 4 1 10 - Les figures 5a et 5b sont respectivement des exemples de cartes d'espérance et de variance de la hauteur d'hydrocarbure liquide dans le modèle déterminé à l'aide des outils de multi-réalisations ; - Les figures 5c et 5d sont respectivement des exemples de cartes 5 d'espérance et de variance de la hauteur d'hydrocarbure liquide dans le modèle déterminé à l'aide des procédés selon un mode de réalisation possible de l'invention. La figure la illustre un exemple possible d'une coupe transversale d'un modèle 10 géologique 100. Ce modèle géologique 100 est constitué d'un grand nombre de pixels (ex. 104, 105) décrivant chacun un faciès présent dans le sous-sol géologique. Les données pétro-physiques NTG, Phi et So sont souvent déterminées pour chaque faciès. Pour chaque pixel c (de dimensions hpxhpxhp), il est possible de 15 calculer une valeur de HuPhiSo « pixel » définie comme suit : hp. NTGc,Ai. SOc,Ai Ai où lAi(x,y,z) est la fonction indicatrice du faciès Ai (égale à 1 si le faciès pour le pixel de coordonnées (x,y,z) est Ah égale à 0 sinon). Bien entendu, chaque pixel n'est associé qu'à un unique faciès et : lAi(x,y,z) =1 Ai 20 Par ailleurs, selon l'histoire géologique du sous-sol, le remplissage sédimentaire et le contexte géodynamique qui ont donné naissance aux roches du réservoir, la sédimentation peut s'être faite de manière hétérogène mais structurée. Les environnements sédimentaires ainsi créés sont plus ou moins riches en certains faciès. 25 Par exemple, dans un chenal, les faciès les plus présents seront les sables grossiers. Plus la roche est éloignée du chenal et plus elle est argileuse. Ces structures sont généralement appelées les éléments architecturaux (ou AE en anglais pour « Architectural Element »). Pour un pixel donné du modèle, on considère qu'il n'existe qu'un seul élément architectural. Le modèle géologique 100 comprend trois éléments architecturaux 101, 102, et 103. Les interfaces entre deux éléments architecturaux sont souvent repérées sur des cartes par les géophysiciens. Par exemple, l'interface 102,in peut représenter une interface possible entre l'élément architectural 102 et l'élément architectural 103. Il est ainsi possible d'introduire, comme pour les faciès, une indicatrice d'éléments architecturaux 1AEi(x,y,z) indiquant la présence ou non de l'élément architectural AEi pour le pixel de coordonnées (x,y,z). Les proportions de faciès en un pixel donné c sont, le plus souvent, fonction de l'élément architectural associé à ce pixel c. Ainsi l'indicatrice de faciès Ai, pour un élément architectural donné AEi, s'exprime sous la forme lAilmi(x,y,z) (également notée 1Ai,AE1 (x, y, z) ou lAi,AE j(c)). Ainsi, la valeur de HuPhiSo « pixel » est définie comme suit : 1AE i(x, y, z) - lAivE i(x, y, z) - h,. - cp"i - NTG"i - SO"i Ai AE avec hc la hauteur du pixel c considéré. Il est également possible de déterminer les valeurs HuPhiSo « colonne » en sommant les valeurs HuPhiSo « pixel » pour toutes les colonnes du modèle (par exemple, colonne106). Cette valeur HuPhiSo « colonne » est donnée par la formule suivante : 1AE j(C) lAilAE j(C) - h,. - (Pc,Ai - NTG c,Ai - Soc,Ai c Ai \AE avec c les pixels de la colonne considérée et avec hc la hauteur du pixel courant. Dans le cas général, l'indicatrice lAivEi(x,y,z) n'est pas connue de manière précise, pixel par pixel : seules les cubes de proportions des différents faciès Ai est connu. Afin de palier à ce problème, les solutions antérieures proposaient : - de déterminer les limites des éléments architecturaux de manière stochastique ; - de déterminer une carte de faciès en simulant de manière stochastique une valeur possible de faciès pour chaque pixel grâce aux cubes de proportions donnée pour chaque facies de chaque élément architectural, - de déterminer la valeur HuPhiSo « pixel » de chaque pixel ; - de sommer les valeurs HuPhiSo « pixel » pour chaque colonne du modèle (i.e. carte HuPhiSo) ; - de réitérer les quatre premières étapes un grand nombre de fois afin, d'obtenir une moyenne et une variance convergente des cartes HuPhiSo. La présente invention ne choisit pas cette voie stochastique. En effet, il est possible de déterminer pour chaque colonne du modèle, la moyenne des valeurs HuPhiSo colonne et la variance des valeurs HuPhiSo colonne de manière analytique, après certaines hypothèses et simplifications géologiques présentées ci-après. DETERMINATION DE L'ESPERANCE L'espérance des valeurs HuPhiSo colonne peut s'exprimer sous la forme : Esp(HuPhiS0) = h,. - (Esp[lAEi(c)1- E sp [1A,IAE j(c)1) c=1 A, \AE - (Esp(cpc,Ai) - Esp(NTGc,Ai) - Esp(S0c,A)) en supposant que les pixels c de la colonne sont indexés de 1 à n (n étant un entier positif supérieur à 2). L'expression de cette moyenne fait certaines hypothèses géologiques afin de simplifier l'expression de cette moyenne : ex. l'indépendance des paramètres cf) NTG,,Ai et SO,,Ai (ainsi, l'espérance du produit de ces variables est également au produit de leur espérance). Ainsi, cette expression n'est pas la simple résolution mathématique du problème posé. Si l'indicatrice de faciès 1A(c) est une variable aléatoire suivant une loi de Bernouilli dont le paramètre est la proportion locale du faciès, alors l'espérance de cette variable sachant pAi(c) peut s'exprimer sous la forme suivante : Esp(lAi(c) I pAi(c)) = Esp(lAi(c) = 1 I pAi(c)) = PA(c) En conclusion : Esp(lAi(c) ) = Esp(pAi(c)) Les proportions pAi(c) ne sont déterminées avec certitude qu'au niveau des puits. Il existe ainsi une incertitude sur la proportion de faciès dans le reste du réservoir. A cause de cette incertitude, la proportion de faciès peut être, elle aussi, une variable aléatoire, dont la loi de probabilité peut être modélisée par une distribution triangulaire 150 décrite dans la figure lb (pour un faciès Ai donné et pour un pixel c) ayant comme valeur minimale pmin,Ai, comme valeur de mode Pmode,Ai et comme valeur maximale n r-max,Ai- Ainsi : Pmin,AilAEi(c) Pmode,AilAEi(c) Pmax,AilAEi(c) E sp (lAilAE j(c)) = 3 Bien entendu, le réservoir n'est pas nécessairement stationnaire, i.e. ses proportions pAi(c) (et n min,Ai(c),Pmode,Ai(C),pmax,Ai(C)) peuvent varier en fonction de la cellule c considérée : par exemple, en se rapprochant d'une limite ou d'un point particulier, un faciès peut tendre à disparaître ou tendre à prédominer. Le même raisonnement est valable pour calculer l'espérance de l'indicatrice d'éléments architecturaux. Néanmoins, la proportion d'éléments architecturaux ne connait pas d'incertitude même si plusieurs positions d'interfaces entre éléments architecturaux sont possibles (il est souvent demandé à l'opérateur de pointer une interface minimum 102min, et une interface maximum 102max (éventuellement une interface probable 102made, si cette dernière interface n'est pas déterminé, il est considéré que cette « interface probable » se situe à équidistance de l'interface minimum et maximum). Ainsi : Esp (1AEi(c)) - pAEi(c) En conclusion, et par analogie, l'expression EAEjEsp[lAEi(c)1- Esp[lAivEi(c)1 Pmin,AilAE j(C)+Pmode,AilAE j(c)+Pmax,AmE ,(C) peut s'exprimer sous la forme n (cl AE; AE; , ' qui 3 est en fait la proportion de faciès Ai pour le pixel c (qui prend en compte la possible présence de plusieurs AE et donc plusieurs cubes de proportions pour le même faciès). Une roche est constituée de grains et de vides entre les grains. La porosité de la roche est le rapport entre le volume des vides et le volume de la roche : Vvide = Vtotal Cette porosité est le plus souvent répartie selon une distribution naturelle comme une loi normale. Sa moyenne est, elle aussi, incertaine et suit une loi de probabilité triangulaire. Phi est donc géologiquement la somme de deux variables aléatoires : - (Pmoy (Pstoch OÙ cpmoy est une variable aléatoire qui peut suivre une loi triangulaire (semblable à la figure lb), de valeur minimale pmin, de valeur de mode (n T mode et de valeur maximale cpmax, et où co stoch est une variable aléatoire qui suit une loi normale centrée (de moyenne nulle). Ainsi : Esp(cp) = Esp('Pm0y) + Esp(cpstoch) = Esp('Pm0y) = (Pmin (Pmode (Pmax 3 La variable NTG représente la fraction d'un volume rocheux favorable à la présence d'hydrocarbures. Si V est un volume rocheux, alors NTG*V est le volume favorable à la présence d'hydrocarbures, c'est le volume utile de la roche. La valeur NTG peut être répartie selon une distribution naturelle qui suit en général (mais non nécessairement) une loi normale de moyenne NTGmoy. La moyenne NTGmoy peut être, elle-même, incertaine et peut suivre une loi de probabilité triangulaire (semblable à la figure lb), de valeur minimale NTGmin, de valeur de mode NTG mode et de valeur maximale NTGmax. La distribution NTG est donc géologiquement la somme de deux variables aléatoires : NTG = NTGmoy + NTG,t'h OÙ NTGmoy est une variable aléatoire qui suit la loi triangulaire décrite précédemment et où NTGstoch est une variable aléatoire qui suit une loi normale centrée (de moyenne nulle). NTGmin+ NTGmode + NTGm', Esp(NTG) = Esp(NTGmoy) = 3 La variable SO représente la saturation en huile de la roche. Les pores des roches peuvent être remplis de gaz, d'eau ou d'huile (pétrole). La saturation en huile So, la saturation en eau Sw et la saturation en gaz Sg se complètent donc : So + Sw + Sg = 1 De plus, on note SHc la saturation en hydrocarbures (i.e. gaz et huile) : SHC - So + S9 On a alors : SHc - 1 Sw Dans les réservoirs (roches poreuses et perméables), les fluides sont soumis à la gravité et répartis selon leur densité. Ainsi, sous une roche perméable 201, il est possible de trouver (voir la figure 2a) en haut le gaz (zone 202), au milieu l'huile/pétrole (zone 203) et enfin l'eau (zone 204). En-dessous du contact huile/eau 205 (cote Z205), les pores ne sont remplis que par l'eau et, au-dessus du contact 205 (zone 203), ils peuvent être remplis d'huile mais de l'eau résiduelle peut également être adsorbée aux grains de la roche. Ainsi, en dessous de la cote Z205, la saturation en eau Sw est de 1 (courbe 211). 5 Au-dessus de la cote Z205, la saturation en eau Sw diminue progressivement (courbe 212) jusqu'à une cote (zR) à laquelle la saturation en eau Sw reste sensiblement constante (i.e. la valeur résiduelle de saturation en eau) et égale à VR (courbe 213). Les courbes de saturation (211, 212, 213) peuvent être, elles-mêmes, incertaines (pour une cote Z205 fixée et connue de manière certaine). Il est ainsi possible de 10 déterminer trois courbes de saturation en eau : - une courbe présentant un minimum de saturation en eau (214), - une courbe présentant un maximum de saturation en eau (215), - une courbe de mode présentant le cas le plus probable (212 et 213). Ainsi, pour chaque hauteur au-dessus du contact 205, il est possible de définir une 15 loi de probabilité triangulaire pour la saturation en eau à l'aide de trois points correspondants des courbes ainsi définies. Par ailleurs, la cote zw d'interface entre l'huile et l'eau (Z205 dans l'exemple de la figure 2a) peut ne pas être connue de manière certaine. Cette cote zw peut également être représentée par une variable aléatoire qui suit une loi de probabilité 20 triangulaire de valeur maximale zw,min, de valeur minimale zw,max et de valeur de mode Zw,mode- On note SHc(z) = CSHc(z - zw) pour mettre en évidence le fait que la courbe de saturation en hydrocarbures peut être translatée verticalement du fait de l'incertitude sur zw : 25 - la courbe présentant un minimum de saturation en hydrocarbures est notée CSHc,min(z - zw), - la courbe présentant un maximum de saturation en hydrocarbures est notée CSHc,m'(z - zw), - la courbe de mode présentant le cas le plus probable est notée CSHc,mode(z - zw). De même, la cote zg représentant la cote de l'interface gaz/huile (z206 dans l'exemple de la figure 2a pour l'interface 206) peut ne pas être connue de manière certaine. Cette cote zg peut également être représentée par une variable aléatoire qui suit une loi de probabilité triangulaire de valeur maximale zg,mm, de valeur minimale zg,max et de valeur de mode Zg,mode. Ainsi, l'espérance de la saturation en hydrocarbures SHc peut s'exprimer sous la forme : Esp(SHc(z)) = Esp(1- Sw(z)) = Esp(CSHc(z - zw)) = Esp(Esp(CSHc(z - zw)I zw),zw) Les variables CSFic et Zw pouvant être considérées, géologiquement, comme deux variables aléatoires indépendantes, il est possible d'écrire : Zw,minr CSHC,max Esp(SHc(z)) = CSFIc(Z Zw) - p(CSFIc(z - zw) zw) - p(zw) - dcsFic Zw,max CSHC,min - dzw Zw,min Esp(SHc(z)) =p(zw) - Esp(CSHc(z - zw) sachant zw) - dzw zw,max zw'rnin csnc,min(z zw) + CSHc,mode(Z - Zw) + CS Hc,max (Z - zw) Esp(SHc(z)) = 3 P(zw) - dzw zwmia csHc,xxx(z-zw), Ainsi, il existe trois termes à calculer de type j zwmax 3 p(zw) - dzw (avec XXX parmi min, max et mode). Si la distribution de probabilité de zw est supposée triangulaire, il est possible d'établir deux cas : Si zw E [zw,max, zw,model Zw,max P(zw) = Plp zw + Pic 2 - zw (zw,min zw,max) (zw,made zw,max) -2zw,max (zw,min zw,max). (zw,mode zw,max) Si Zw E [Zw,mode,Zw,minl P(zw) = P2p zw + P2c -2 - zw (zw,min zw,max). (zw,min zw,mode) 2 - zw,min (zw,min zw,max). (zw,min zw,mode) Ainsi : I z p(zw) - dzw Zw,max w'rnin csHc,xxx(z - zw) 3 = j-Zw,mode rc zw,max I 3 -C,XXX(z zw) -FI 1 zw + Plc) zw'rnin csHc,xxx (z zw) f 3 U12p zw + P2c) ' dzw Zw,mo de Puis en réalisant un changement de variable h = z - zw, il est possible d'écrire : f zw'rnin csHc,xxx(z - zw) - dzw 3 zw,max p(zw) Iz3 z - d-zwmax csHc,xxx(h) ,(h) + pic )1h Z-Zw,mo de fZ-Zw,m°de [CS IIC,XXx (h) f Z-Zw,min Vi2p - h) + P2c)1' dh 3 dzw En pratique, la courbe CSHc,xxx (et donc la courbe SWyyy avec YYY=min si XXX=max, YYY=mode si XXX=mode et YYY=max si XXX=min, car Swmin = 1- CSHc,max) est souvent linéaire par morceaux, par exemple, entre deux cotes successives des listes suivantes : - hl 1, hl 2, hl ni -1, hl ni (ordonnées par ordre croissant dans l'intervalle [z - Zw,max,Z - Zw,model , les cotes extrêmes correspondant aux bornes de l'intervalle, ces cotes extrêmes pouvant ne pas correspondre à un point de rupture de la linéarité de la courbe, voir ci-dessous) et - h21, h22, h2n2-1, h2n2 (ordonnées par ordre croissant dans l'intervalle [z - Zw,mode, Z les cotes extrêmes correspondant aux bornes de l'intervalle, ainsi h21= hi ni, ces cotes extrêmes pouvant ne pas correspondre à un point de rupture de la linéarité de la courbe, voir ci-dessous). Par exemple, la figure 2b représente un exemple de courbe de saturation en eau 1 - CS HC,XXX zw) (250), la cote zw étant représentée par la ligne h=0. Cette courbe 250 est linéaire par morceau entre les points successifs 251, 252, 253, 254, 255, 256 et 257. Si le point de cote (z - zw,max) ne correspond pas à un des points précédents, il est possible de définir un nouveau point 271 sur la partie linéaire de la courbe de saturation en eau entre les points 252 et 253. Il en va de même pour la cote (z - zw,mode) (correspondant alors au nouveau point 272) et pour la cote (z - zw,min) (correspondant alors au nouveau point 273). Ainsi : 302 1 1 4 1 20 Z-Zw,max c r ----'11C,XXX(Z zw) p(Zw) - dzw n1-1 h - [CSHC,XXx(h) ( - h) + pic)1 - dhl i=1 [ 3 f Z-Zw,min + f sc Hc,xxx (h ) (P2p(z - h) P2c)1 dhl n2-1 h2. h2 i+ i=1 Dans l'exemple de la figure 2b, chacune des deux sommes précédentes comprend deux termes : 5 _ fhl, rcsHc,xxx(h) , u,p(z - h) + Pic)] - dhl (intervalle 281) comprend un 3 terme correspondant à l'intégral entre les cotes du point 271 et du point 253 et un terme à l'intégral entre les cotes du point 253 et du point 272 ; 2,»C.1t 3 rh2i rcsHc,xxx(h) (P2p - h) + Pu)]. dhl (intervalle 282) comprend un terme correspondant à l'intégral entre les cotes du point 272 et du point 254 et un terme à l'intégral entre les cotes du point 254 et du point 273. Il est considéré que l'équation décrivant la courbe CSHC,XXX entre deux points hl; et 10 hl 41 est Al i,xxx h + Bli,xxx et que l'équation décrivant la courbe CSHC,XXX entre deux points h2; et h24-1 est A2i,xxx h + B2i,xxxx. Autrement dit, la courbe CSHC,XXX est linéaire par morceau entre deux points consécutifs suivants : (Cl txxx, hl 1), (C12,xxx, h12), (C1 n1-1,XXX, hl n1-1), (C1n1,XXX, hl ni)=(C2i,xxx, h21), (C22,xxx, h22), (C2n2-1,XXX, h2n2-1), (C2n2,XXX, h2n2). 15 Alors : cl i+i,xxx -cl i,xxx h1i+1-h1i C2i+i,xxx-C2i,xxx h2i+1-h2i = et C1 Cl i+i,xxx-Cli,xxx t = Bli,XXX et hi 1-hl Al i,xxx = et C2i+1,xxx- u C2i,xxx 7',-)L t = B2i,XXX et A2xxx C2,,xxx h2+1-h2i Ainsi : f zw,max zw'min CSHC,XXX(Z zw) p(zw) - dzw 3 121-1 hl- = [f E+1 [Ali,XXX ' h + (Pip(z - h) + Pic)] hli i=1 n2-1 h2i+i [A2i xxx ' h + B2 ] dh i,xxx 3Pu) i=1 h2 dh] ou encore : Izw'rnin csHc,xxx(z - zw) p(zw) - dzw n1-1 = -P1p 3 22 zw,max 3 11 3 3 hli+12 ^ (Ali,XXX - (Pic + Pip - z) P ip- Bli,XXX ) - 2 + Bli,xxx ' (Pic ^ Pip z) - 111t+1) h1i3 , hli2 ( Pip - Ali,XXX - 3 + Ylli,XXX - (Pic + Pip - z) - Pip ' Bli,XXX) 2 [( h2j+1 Bli,XXX - (Pic Pip z) - hli)1 n2-1 3 + Li V _ P2p A2t,XXX 3 3 i=i h2i+12 ^ (A2i,XXX - (P2c P2p - Z) P2p - B2i,xxx) - + B2i,xxx - (P2c 2 ^ P2p z) - h2t+1) ( h2i3 P2p - A2i,XXX - -3 + (A2i,xxx ' (P2c + P2p - z) - P2p ' B2i,XXX) 2 ^ B2i,XXX (P2c P2p z) - h2i)1 Ainsi, chacun des termes de EsKSHc(z)) peut être aisément calculé (i.e. rzw,min csHC,XXX(z-zw) p(Zw) - dzw avec XXX parmi min, max et mode). ^ Jzw,max 3 Par ailleurs, il est possible de définir l'indicatrice « huile + eau » (ou encore indicatrice « liquide-gaz », le liquide étant alors ici l'huile et l'eau) low par : {1 si z Il est ainsi possible de déterminer (étape 305) l'espérance de la valeur HuPhiSo pour la cellule sélectionnée : (pAE;(c) Pmin,AilAE j(c) Pmode,AilAE j(c) Pmax,AilAE j(c) \AE; 3 - (Esp (q)c,Ai) - Esp(NTG"i) - Esp(S0c,A3) Cette formule fait un certain nombre de simplifications fondées sur des considérations géologiques évoquées ci-avant. S'il est considéré que les variables cp"i, NTG"i, ou SO"i sont certaines (i.e. valeurs fixées sans incertitude) pour la cellules c et pour le faciès Ai, alors l'espérance de ces valeurs est égale à ces valeurs certaines et la valeur déterminée à l'étape 305 peut être associée directement à la cellule sélectionnée dans le modèle. Si jamais la variable Cpc,Ai est considérée comme incertaine, il est fourni trois valeurs représentatives d'une distribution triangulaire de probabilité pour cette cellule sélectionnée et pour chaque faciès Ai : une valeur minimale Cpc,Ai,min, une valeur maximale w et une valeur de mode w Dès lors, l'espérance de la T c,A i,max T c,Ai,mode - porosité pour la cellule sélectionnée et pour le faciès Ai est (étape 306a) : ESp(Cpc,A) = (Pc,Ai,min (Pc,Ai,mode (Pc,Ai,max 3 De même, si la variable NTG"i est considérée comme incertaine, il est fourni trois valeurs représentatives d'une distribution triangulaire de probabilité pour cette cellule sélectionnée et pour chaque faciès Ai : une valeur minimale NTG'i,min, une valeur maximale NTG'i,max, et une valeur de mode NTG"i,mode - Dès lors, l'espérance de la valeur NTG pour la cellule sélectionnée et pour le faciès A est (étape 306b) : Esp(NTG"i) = NTG"i,min+ NT G c,Ai,mode + NTG"i,m' 3 Enfin, si la variable SOc,Ai est considérée comme incertaine, il est fourni trois courbes de saturation en eau fonction de la profondeur z de la maille c, une courbe minimale, maximale et de mode, linéaires par morceaux et respectivement définies par une pluralité de points et par trois valeurs représentatives d'une distribution triangulaire de probabilité pour la valeur de la cote d'une interface eau-huile : une valeur minimale zw,m', une valeur maximale zw,min, et une valeur de mode z -w,mode - Dès lors, l'espérance de la valeur SO pour la cellule sélectionnée est (étape 306c) calculée comme indiquée ci-avant. Une distribution de probabilité p(zg) de la cote zg huile-gaz peut également être fournie. La valeur déterminée à l'étape 305 et à l'étape 306a et/ou 306b et/ou 306c peut être utilisée pour calculer une valeur d'espérance de la variable HuPhiSo. Cette espérance peut alors être associée à cette cellule. Si toute les cellules du modèle ont une valeur d'espérance HuPhiSo associée (test 304, sortie KO), il est possible de sommer (étape 307) les espérances d'HuPhiSo des cellules appartenant à une même colonne afin d'obtenir une espérance d'HuPhiSo colonne pour ladite colonne considérée. Cette carte 308 est alors fournie à un opérateur pour visualisation et/ou traitement ultérieur. DETERMINATION DE LA VARIANCE Par ailleurs, la variance des valeurs HuPhiSo colonne peut s'exprimer sous la forme : Var (HuPhiSocolonne) = Esp(HuPhiSocolonne 2 ) Esp(HuPhiSocolonne) 2 Le calcul du terme Esp(HuPhiSocoionne) et donc de Esp(HuPhiSocoionne)2 est détaillé ci-dessus et celle-ci peut être calculée pour le calcul de la variance. La formule ci-dessous fait l'hypothèse que, géologiquement, les variables cpc,Ai, NTGc,Ai et SOc,Ai sont des variables indépendantes : Esp(HuPhiSo2) hclhc2. EspOLA(c1)1B(c2)) ESP(q) cl,A(1) c2,B) A B=A cl c2=c1 - Esp(NTGc1,ANTGc2,B)- EsP(Soci,ASoc2,B) hclhc2. EspOLA(c1)1B(c2)) ESP(q) cl,A(1) c2,B) A BA cl c2=c1 - Esp(NTGc1,ANTGc2,B)- EsP(Soci,ASoc2,B) hclhc2. EspOLA(c1)1B(c2)) ESP(çPc1,A(Pc2,B) A B=A cl c2c1 - Esp(NTGc1,ANTGc2,B)- EsP(Soci,ASoc2,B) hclhc2 EspOLA(c1)1B(c2)) ESP(çPc1,A(Pc2,B) A BA cl c2c1 - Esp(NTGc1,ANTGc2,B)- EsP(Soci,ASoc2,B) avec A, B un faciès (éventuellement identique) parmi les différents faciès possible, et cl, c2 une cellule (éventuellement identique) de la colonne considérée pour le calcul du HuPhuiS0 colonne. La première somme quadruple correspond aux termes où les faciès A et B sont identiques et les cellules cl et c2 sont identiques, la deuxième somme correspond aux termes où les faciès A et B sont différents mais les cellules cl et c2 identiques, la troisième somme correspond aux cas où les faciès A et B sont les mêmes mais les cellules cl et c2 différentes et enfin la dernière somme correspond aux cas où les faciès A et B sont différents et les cellules cl et c2 aussi. Il est fait l'hypothèse géologique qu'il ne peut pas y avoir deux faciès différents dans la même cellule, soit : Esp(1A(c1) - 1B(c2)) = 0 si cl = c2 et A # B Ainsi : Esp(HuPhiSo2) =11c.2 ' ESP(1A(C)2) - ESP(Pc,A2) Esp(NTGc,A2) - EsP(Soc,A2) c A -F hciha ' ESP(1A(C1)1A(C2)) - ESK(Pc1,A(Pc2,A) A cl c2c1 - Esp(NTGc1,ANTGc2,A). EsP(Soci,ASoc2,A) -F liciha - Esp(1A(c1)1B(c2)) - E SP(CP cl,A(Pc2,B) A 13,11 cl c2c1 - ESKNT G ciANT G c2,B) . ESKSOci,AS0c2,B) Cette somme fait apparaitre un certain nombre de termes identiques qu'il est possible de regrouper afin de diminuer le temps de calcul nécessaire. Si la colonne considérée comporte n cellules et s'il existe N faciès possible (ceux- ci étant alors indexés), il est possible d'écrire : Esp(HuPhiSo2) [ = lic2 Esp(cpc,A2) - Esp(NT G c,A2) - Esp(Soc,A2) c A -Esp(PAE(c)) - Esp (PA,AE (C)1 AE N n-1 n +2 [hc1hc2Esp(1A(c1)1A(c2)) ESP(Pc1,A(P c2,A) A=1 c1=1 c2=c1+1 - E sp (NT G ciANT G c2,A) - EsP(Soci,AS0c2,A) -Esp(1-AEi(C1)1AE2 (C2)) ' Esp (1A,AE1(c1)1A,AE2(c2))1 AE1 AE2 n-1 n N N +2 [hc1hc2Esp(1A(c1)1B(c2))'EsP(Pc1,A(Pc2,B) c1=1 c2=c1+1 A=1 13,11 - Esp(NTGc1,ANTGc2,B) - EsP(Soci,AS0c2,B) - Esp(1AE1(c1)1AE2(C2)) ' Esp (1A,AE1(C1)1B,AE2(C2))1 AE1 AE2 Le terme EAE ESKPAE(C)) . Esp (pA,AE(c)) a été calculé précédemment et, suite aux simplifications géologiques apportées, vaut : P AE (c.) . (P A,AE.min(C) + P A,AE.mode(C) + P A,AE.max(C) 3 AE Le terme Esp(1AE1(c1)1AE2(c2)) peut s'exprimer comme suit s'il est fait l'hypothèse géologique selon laquelle les incertitudes concernant les éléments architecturaux ne sont que horizontales (il n'existe pas d'incertitudes verticales) : Esp(1AE1(c1) . 14E2(c2)) PAE1(c1) - P AE2(c2) si au moins l'un des AE est certain = pAE (ci) si les AE sont tous incertains et AE1 = AE2 1 0 si les AE sont tous incertains et AE1 # AE2 C'est-à-dire : Esp(1AE1(c1) . 14E2(c2) PAE2(C2) Si (PAE1(C1) = 1) PAE1(C1) Si (PAE2(C1) = 1) 0 Si (13 AE1(C1) = 0) ou si (PAE2 (a) = 0) cl) si AE1 = AE2 sinon fPAE ( L 0 si AE1 # AE2 Dès lors, nous avons : ESP(1AE1(C1-) . 1AE2(C2)) . Esp (1A,AE1(c1) ' 1B,AE2(C2)) AE1 AE2 =Esp(1-AEi(c1) ' 1,4E1(C2)) . Esp (1A,AE1(c1) ' 1B,AE1(C2)) AE1 -F ESP(1AE1(C1) . 1AE2(C2)) Esp (1A,AE1(c1) ' 1B,AE2(C2)) AE1 AE2AE1 Il est supposé que, géologiquement, les indicatrices de faciès de deux cellules situées dans des éléments architecturaux différents sont indépendantes entre elles. Ainsi, pour AE1 #AE2 : Esp (1A,AE1(c1) ' 1B,AE2(C2)) = Esp (1A,AE1(C1))' Esp (1B,AE2(C2)) - Esp (p"E1(c1)). Esp (pB,AE2(c2)) Par ailleurs, après certaines considérations géologiques, pour AE1=AE2, il est possible d'écrire, si les faciès sont identiques (i.e. A=B) : Esp (1A,AE(c1) - 1A,AE (C2)) = (1 - p(A, Az)) . ( (PA,AE,ci) Esp ' EsP(PA,AE,c2) .\1+ var(pA,AE,c1)-Var(PA,AE,c2) 1 \ 1 f + p (A, Az) - (-2Esp ( n \i- A,AE,c1) + -2 Esp VJA,AE,c2) avec p(A, Az) la fonction « corrélogramme » du faciès A (notée également pA(Az) ou pA(c1, c2) et définie ci-après) et Az la différence de cotes verticales selon l'axe vertical z (i.e. l'axe de la colonne) entre les cellules cl et c2. De plus, pour des faciès différents (i.e. A#I3), il est possible d'écrire : Esp (1A,AE(c1) - 1B,AE(c2)) - (Esp(pA,AE). Esp(pB,AE) + _\IV ar(pA,AE (ci)) . V ar(pB,AE(c2))) - (1 - p(A, Az)) On appelle variogramme d'une variable aléatoire X, la fonction y: 1 y(X , h) = 7. Esp ((X(u + h) - X(u))2) Ce variogramme y représente la structuration spatiale de la covariance entre la variable X en un point de l'espace et la même variable X en un autre point. Le plus souvent, les variogrammes y d'une variable aléatoire réelle respectent les conditions suivantes : - Le variogramme y est une fonction croissante de la distance, ce qui signifie que plus deux points sont éloignés et moins les données en ces points sont corrélées. - Au-delà d'une distance Ry (également appelée « range » en anglais), le variogramme y a atteint sa valeur asymptotique (ou presque, soit 95% de cette valeur), ce qui signifie que la structuration spatiale des données n'existe que jusqu'à cette distance caractéristique Ry; au-delà de cette distance Ry, les données ne sont plus corrélées. - Pour une distance nulle, le variogramme y est nul, ce qui veut dire que la corrélation entre les données est parfaite pour deux points très proches l'un de l'autre. On appelle le corrélogramme p la covariance normalisée : Cov(X(u + h), X (u)) y(X, h) p(X, h) = = 1 C ov (X (u), X (u)) Cov(X(u), X (u)) Le plus souvent, pour une variable aléatoire réelle X donnée (comme le faciès), le 10 corrélogramme peut être modélisé sous la forme (Ry étant une valeur prédéterminée) : - d'un corrélogramme exponentiel : 3h) p(h) = exp (--Rv - d'un corrélogramme sphérique : h G Ry 1 3 h 1 i h \3 si Y Y 0 sinon - d'un corrélogramme gaussien : 3h)2 p(h) = texP (- (Fy ) si h G Ry 0 sinon 15 Comme mentionné ci-dessus, géologiquement, les cellules de même faciès peuvent se répartir par blocs, qui sont souvent le résultat de l'histoire de l'environnement sédimentaire auquel ils appartiennent. Par exemple, pour un environnement sédimentaire de type chenal en milieu marin (sables, turbidites), les coulées de sable peuvent avoir provoqué la formation de lobes de sable. Ces lobes seront modélisés par le géologue comme un groupe de cellules de faciès sableux. Cette tendance qu'ont les cellules de même faciès à se regrouper par groupes, de même forme et de même taille moyenne, définit une corrélation spatiale des faciès. Une cellule proche d'une autre cellule a ainsi plus de chance de faire partie du même faciès que cette dernière. Il est possible de modéliser cette corrélation spatiale avec les variogrammes ou un corrélogramme : Cov(1A(c1),1A(c2)) = pA(cl, c2) . var(1A(c1)) - Var(1A(c2)) Les valeurs Ry ou « ranges » (voir ci-dessus) des faciès (dans les trois dimensions de l'espace) sont les dimensions des formes géologiques observées sur le terrain, dans un élément architectural donné. Un variogramme ou un corrélogramme des faciès peut donc être défini par élément architectural. Si les indicatrices de faciès sont corrélées spatialement via les variogrammes, les proportions de faciès peuvent être, quant à elles, parfaitement corrélées entre elles. Pour le calcul du HuPhiSo, les variogrammes peuvent être des données d'entrée. On utilise notamment les variogramme pour corréler spatialement les variables NTG et Phi. Géologiquement, les corrélations spatiales peuvent s'appliquer qu'à la partie stochastique (noté Xstoch) des variables aléatoires (X), les incertitudes sur la moyenne étant alors parfaitement corrélées entre elles. De plus, les corrélations entre les porosités ou les variables NTG pour deux cellules différentes se font dans un même faciès. Ainsi, il n'y a corrélation entre deux porosités ou variables NTG, dans deux cellules différentes, que pour les mêmes faciès. Pour deux faciès différents, les parties stochastiques des variables NTG et Phi peuvent être indépendantes. Il est alors possible de calculer les covariances entre les variables stochastiques pour deux cellules différentes et le même faciès : C OV (NT G stoch(C1), NT G stoch(C2)) = PNTG (cl, c2) - Var(NTGstoch (cl)) - Var(NTGstoch(c2)) Cov@Pstoch (cl), (Pstoch (a)) = (cl, c2) - Var(cpstoch (cl)) - Var(cpst0ch(c2)) ou PNTG est le corrélogramme vertical de la variable NTG et p9 est le corrélogramme vertical de la variable Phi. Concernant le terme Esp(cpc,A2), il est possible d'écrire, après certaines considérations géologiques, que Esp(cpc,A2) = Esp(cpc,A)2 + Var (0 C c,A,moy) V ar(( P c,A,stoch)- Il en va de même pour la variable NTG et ainsi Esp(NTGc,A2) = Esp(NTGc,A)2 + Var(NTGc,A,m0y) + Var(NTGc,A,stoch)- Concernant l'espérance du produit de la porosité pour un même faciès et deux cellules différentes, il est possible d'écrire, après certaines considérations géologiques (notamment en considérant qu'il existe une corrélation parfaite entre les incertitudes sur la moyenne la partie moyenne du corrélogramme de la variable Phi) : Esp((pci,A - (pc2,A) = EsP((Pci,A) ' EsK(Pc2,A) + var ((pci,A,inc) - var((pc2,A,inc) P(Pstoch Var4c1,A,stoch Var4c2,A,stoch ou Apstoch est la partie stochastique du corrélogramme de la variable Phi (voir plus haut) et (D cl,A,inc est l'incertitude sur la moyenne (i.e. si la distribution de (D ci,A est triangulaire de borne minimal a, de borne maximale c et de mode b, T c 1,A, inc 2 7-2 2 7- a - v -ab-bc-ca 18 Il en va de même pour la variable NTG et ainsi : ). Esp(NT G ci,A NT G cz,A) - V ar (NT G c2,A,inc) + ar (NT G cl,A,inc) = E Sp(NT G ci,A) - E Sp(NT G c.2,A) PNTGstoch ar (NT G ci,A,stoch) ' V ar(NT G c2,A,stoch) ou PNTGstoch est la partie stochastique du corrélogramme de la variable NTG (voir plus haut). De même, si A#I3, alors il est possible d'écrire : E sK(P ci,A - (f) cz,B) = EsP(q)ci,A) E sP(q) cz,B) + V ar4 - V ar(q) c2,B,inc) Esp (NTGcl,A ' NTGc2,B) = Esp (NTGctA ) - Esp (NTGc2,B) + -\/Var(NTGcl,A,inc) - V ar (NT G c2,B,inc) Ainsi, en réutilisant les notations ci-dessus : 1 E sp (S (z)2) = z(E sP Hc ,min(z)2) + E sP (S. HC ,mode (Z)2) E sP (S. Hc ,max(z)2) + Esp (SHC,min(Z) SHC,mode(Z)) Esp (SHC,min(Z) SHC,max(Z)) + Esp (SHC,mode(Z) SHC,max(Z))) A titre d'illustration, la valeur Esp HC,XX1(z) SHC,XX2(Z)) peut être calculée (avec XX1 et XX2 parmi mode, min et max et avec XX1 et XX2 différents) comme suit : E spS hli+14 (- - - Hc,xxi(z) ' SHC,XX2(z)) n1-1 {= A 1i,XX1 . A li,XX2 i=1 4 + ' A1i,xx2 ' (- Plp z- Pic) Plp Bli,xx2 3 + Ali,xx2 B1i,xx1)J hli+i ^ {(Ali,XX1 Bli,XX2 Ali,XX2 Bli,XX1) (Plp z Pic)} + 2 I' Bli,xxi ' B1i,xx2 ' hli-Ei h1i4 - -Plp . Ali,XX1. Ali,xx2-I' [ - + {Ali,xxi ' A1i,xx2 ' (- Plp . z- + Pic) - Plp . (Ali,xxi Bli,xx2 hli3 + Ali,xx2 B1i,xx1).1 3 \.) h1i2 ^ {(Ali,XX1 Bli,XX2 Ali,XX2 Bli,XX1) (Plp z Pic 11 2 Pic 2 ' ' B1i,xx2 ' hli n2-1 h2i+14 {[-P2p - A2i,xxi - A2i,xx2 4 + {A2i,xx1 - A2i,xx2 ' (- P2p z- P2c) P2p h2i+1 + A2i,xx2 B2i,xxi)} 3 3 (A2,1 B2 i,XX2 2 h2i+i + {(A 2 i,xxi B2 i,xx2 A2i,XX2 B2i,XX1) - (P2p - z P2c)} - 2 P2c I- B2i,xxi - B2i,xx2 - h2i-Ei [h2i4 - -P2p - A2i,xxi - A2i,xx2 4 + {A2i,xxi ' A2i,xx2 ' (P2p . z + P2c) P2p . (A2,1 B2 i,XX2 h2i3 + A2i,XX2 B2i,XX1)} 3 h2i2 + {(A2i,xx1 B2i,xx2 + A2i,xx2 B2i,xx1) - (P2p - z + P2c)} - 2 P2c - B2i,xxi - B2i,xx2 - h2i1} Par ailleurs, avec z1 différent de z2 Esp (SHc,A(z1) - SH-c,B(Z2)) 1 = (E (SHC,A,min(z1) SHC,B,min(z2)) E (SHC,A,mode(z1) SHC,B,mode(z2)) + E (SFic ,A,max(z1) SHC,B,max(z2))) 1 ( + 12 (SHC,A,min(z1) SHC,B,mode(z2)) E (SHC,A,mode(z1). SHC,B,min(z2)) + E (SHC,A,min(z1) SHC,B,max(z2)) E (SHC,A,max(z1). SHC,B,min(z2)) + E (SHC,A,mode(z1) SHC,B,max(z2)) E (SHC,A,max(z1) SHC,B,mode(z2))) Chacun des termes de cette formule peut se calculer comme suit (avec XX1 et XX2 parmi mode, min et max) : E (SH",xxi(z1) ' 5HC,B,XX2(z2)) i=1 n1 rw,i+i (Plp zw pic) - (A111 (Z1 - Zw) B11) - (A2i,xx2 - (Z2 - zw) B2i,xx2) - dzw n2 z zw,i+i (P2p Zw P2c) - (1111,XX1 (Z1 - zw) B11) i=2 Zwi - (i121,xx2 - (Z2 - zw) B21,XX2) - dz, De plus : Eseoc,A) =- EsP(Slic,A(c)) EsP(low(c)) Eseo,,A2) = EsP(SHc,A(02)' EsP(low(c)) et Esp(Soc" - S0c2,B) = Esp(SH"(c1) - SHc,B(C2)) - ESp(10w(C1) - low(C2)) sachant que : Esp(low(zi) - low(z2)) = Esp(low(max(zi,z2))) En outre, il peut exister certaines incertitudes quant à la position des structures 10 dans le réservoir, pouvant ainsi modifier les cartes « HuPhiSo » (de variance ou d'espérance). Ces incertitudes peuvent être liées - aux erreurs et à l'imprécision des mesures, par exemple concernant la position et l'orientation des surfaces. Cette incertitude concerne 15 principalement la géométrie des structures ; - à la variabilité naturelle des objets géologiques, plus ou moins importante selon les objets géologiques (ex. la variabilité naturelle d'une surface de faille est beaucoup plus faible que celle d'une surface de minéralisation d'un gisement d'uranium) ; - au manque de connaissances, sur l'existence d'une structure par exemple. Il est possible de prendre en compte ces incertitudes en les intégrant à l'incertitude sur la position de l'interface huile-eau mentionné plus haut. En reprenant les notations précédentes, il est possible de considérer que l'incertitude de l'interface huile-eau suit une distribution triangulaire (incertitude mentionnée précédemment) mais auxquelles s'ajoute une incertitude de distribution normale de déviation o-(z) (liée à l'incertitude de la structure). Sous hypothèses géologiques, le terme : fzw'min CSHC,XXX(Z - zw) p(zw) - dzw zw,max 3 (mentionné ci-dessus) peut être remplacé par, un changement de variable h=z-z, ayant été effectuée : z+ 30-(z) z-3c)-(z) j- Z-Zw,mode Z-Zw,max (z' - h) + pic) csHc'xxx(h) dh 3 (Plp dz' z+3 cr(z) f Z-Zw,mode csFicxxx (h) (P2p (zi h) + P2c) 3 dh Z-Zw,min z-30-(z) dz' Ainsi, en reprenant les notations précédentes, il est possible de remplacer : fzw,max zw'min CSHC,XXX(Z - zw) p(zw) - dzw 3 par : z+30-(z) n1-1 Plp Ali,XXX h1i+13 f 131 3l (pic + pip - Z1) - 10pB1 z-3o-(z) i=1 hli+12 Bli,XXX (Pic Plp Z1) - hli+1 2 hli3 . Ali,xxx ' -3 (Plp + (Ali,xxx ' (Pic + PlpZ1) - Pip hli2 - - Bli,XXX) 2 Bli,XXX (Pic Plp Z1) - 1111 n2-1 h2i+13 2 h2i+1 1 3 A2i,XXX P2 3 p i=1 - B2i,XXX). 2 1- D L'i,XXX (P2C - (A2i,XXX (P2C P2p Zi) P2p P2p zi) - h2t+1) h2i2 h2.3 \ ( P 2p i,XXX (A2i,xxx ' (P2c + P2p . zi) . A2 + . B2i,XXX) ; P2p 2 7 a( z, \ )2 2a(z) B2i,XXX - h2i)11 2n-u(z) (P2c P2p Z1) Par ailleurs, il est possible de calculer la variance de la valeur HuPhiSo prenant en compte les incertitudes sur la structure (ou Var(HuPhiSostr inc), à partir de la valeur de la variance de cette valeur HuPhiSo sans incertitude (ou Var(HuPhiSo) calculée précédemment) : \ [Esp(HuPhiSostr inc )12 ( Var(HUPhiSOstr inc) = Var (HuPhiSo) + Esp(Hstr inc) v aru iu str inc) avec Hstr inc la hauteur de la colonne au-dessus de l'interface eau-huile. La figure 4 représente un exemple de dispositif de détermination d'une carte d'espérance de hauteur d'hydrocarbure liquide dans un mode de réalisation de l'invention. Dans ce mode de réalisation, le dispositif comporte un ordinateur 400, comprenant une mémoire 405 pour stocker des instructions permettant la mise en oeuvre du procédé, les données de mesures reçues, et des données temporaires pour réaliser les différentes étapes du procédé tel que décrit précédemment. L'ordinateur comporte en outre un circuit 404. Ce circuit peut être, par exemple : - un processeur apte à interpréter des instructions sous la forme de programme informatique, ou - une carte électronique dont les étapes du procédé de l'invention sont décrites dans le silicium, ou encore - une puce électronique programmable comme une puce FPGA (pour « Field- Programmable Gate Array » en anglais). Cet ordinateur comporte une interface d'entrée 403 pour la réception de données notamment les cubes de proportions nécessaires aux calculs de l'invention, et une interface de sortie 406 pour la fourniture des cartes d'espérance de la hauteur d'hydrocarbure liquide. Enfin, l'ordinateur peut comporter, pour permettre une interaction aisée avec un utilisateur, un écran 401 et un clavier 402. Bien entendu, le clavier est facultatif, notamment dans le cadre d'un ordinateur ayant la forme d'une tablette tactile, par exemple. Par ailleurs, le schéma fonctionnel présenté sur la figure 3 est un exemple typique d'un programme dont certaines instructions peuvent être réalisées auprès de l'équipement décrit. A ce titre, la figure 3 peut correspondre à l'organigramme de l'algorithme général d'un programme informatique au sens de l'invention. La figure 5a est un exemple de carte de l'espérance de la variable HuPhiSo pour un modèle maillé calculé à l'aide des méthodes de l'état de la technique connues (i.e. calcul à l'aide d'une pluralité de réalisations ou de « multi-réalisations », l'espérance de ces réalisations étant alors calculée). La figure 5b est un exemple de carte de la variance de la variable HuPhiSo pour un modèle maillé calculé à l'aide des méthodes de l'état de la technique connues (multi-réalisations). Le modèle maillé, considéré dans ces exemples, comporte 80000 mailles. L'incertitude concernant notamment sur les cubes de proportions, la porosité, les saturations en eau est prise en compte. Le calcul d'une centaine de réalisations prend environ 20h. Néanmoins, il est habituel de déterminer plutôt 2000 réalisations lors d'un tel calcul en « multi-réalisation » afin d'assurer une convergence correcte des résultats. Ainsi, le temps de calcul peut largement excéder 20h. Pour ce même modèle et les mêmes prises en compte d'incertitudes, la figure 5c est un exemple de carte de l'espérance de la variable HuPhiSo déterminée à l'aide d'une méthode mettant en oeuvre un des modes de réalisation de l'invention. La figure 5d est un exemple de carte de la variance de la variable HuPhiSo déterminée à l'aide d'une méthode mettant en oeuvre un des modes de réalisation de l'invention. Le calcul des cartes 5c et 5d à l'aide d'un des modes de réalisation de l'invention prend environ 2h pour un résultat tout à fait comparable aux cartes 5a et 5b (respectivement 2% de différences en moyenne et 0,5% de différences en moyenne). Bien entendu, la présente invention ne se limite pas aux formes de réalisation décrites ci-avant à titre d'exemples ; elle s'étend à d'autres variantes. D'autres réalisations sont possibles. REVENDICATIONS1. Procédé de détermination d'une carte d'espérance de hauteur d'hydrocarbure liquide dans un modèle géologique (100), le modèle géologique comportant des mailles (104, 105) pouvant être associée à un faciès parmi une pluralité de faciès et à un élément architectural parmi une pluralité d'éléments architecturaux (101, 102, 103), le procédé comporte, pour au moins une colonne du modèle géologique : /a/ réception d'un cube de répartitions pour chaque élément architectural AE (3021, ..., 302m) parmi la pluralité d'éléments architecturaux, le cube de répartitions décrivant pour chaque maille c du modèle la probabilité pAE(c) que ladite maille appartienne audit élément architectural, /b/ réception d'un cube de proportions pour chaque faciès parmi la pluralité de faciès et pour chaque élément architectural parmi la pluralité d'éléments architecturaux (3011,1, 3012,1, ... 301,,m), le cube de proportions associé à un faciès A et à un élément architectural AE décrivant pour chaque maille c du modèle la proportion dudit faciès dans ladite maille si ladite maille appartient audit élément architectural AE, ladite proportion étant une variable aléatoire de distribution triangulaire (150) définie par une valeur minimale pA,AE,mm(c), une valeur maximale pA,AE,max(c), et un mode PA,AE, mode (C) , /c/ détermination (305), pour chaque maille c de ladite colonne et pour chaque faciès A de la pluralité de faciès, de la somme sur la pluralité d'éléments architecturaux de n c AE ()- -3 (i n A,AE,min (c.) PA,AE,max(c) PA,AE,mode(c)). /d/ détermination d'une valeur d'espérance de hauteur d'hydrocarbure liquide pour ladite colonne (307) de manière analytique en fonction de la somme déterminée à l'étape /c/ pour chaque maille c de ladite colonne et pour chaque faciès A. 2. Procédé selon la revendication 1, dans lequel, chaque maille c du modèle étantassociée à une proportion de porosité pour chacun des faciès de ladite pluralité de faciès, ladite proportion de porosité étant une variable aléatoire de distribution triangulaire définie par une valeur minimale cf) Amin(c), une valeur maximale P A ,max (C) et un mode cpAmode(c), le procédé comprend en outre : /c-1/ détermination (306a), pour chaque maille c de ladite colonne et pour chaque faciès A parmi la pluralité de faciès, de la valeur î ((PA,min(c) + P A,max(C) P A,mo de (C)), dans lequel, la détermination de l'étape /d/ est fonction de la valeur déterminée à l'étape /c-1/ pour chaque maille c de ladite colonne et pour chaque faciès A parmi la pluralité de faciès. 3. Procédé selon l'une des revendications précédentes, dans lequel, chaque maille c du modèle étant associée à une proportion de volume rocheux favorable à la présence d'hydrocarbures pour chacun des faciès de ladite pluralité de faciès, ladite proportion de volume rocheux favorable à la présence d'hydrocarbures étant une variable aléatoire de distribution triangulaire définie par une valeur minimale NTGAmin(c), une valeur maximale NTGAmax(c), et un mode NTGAmode(c), le procédé comprend en outre : /c-2/ détermination (306b), pour chaque maille c de ladite colonne et pour chaque faciès A parmi la pluralité de faciès, de la valeur î(NTGA,min(c) + NTGA,max (c) + NTGA,mode (c)), dans lequel, la détermination de l'étape /d/ est fonction de la valeur déterminée à l'étape /c-2/ pour chaque maille c de ladite colonne et pour chaque faciès A parmi la pluralité de faciès. 4. Procédé selon l'une des revendications précédentes, dans lequel, chaque maille c du modèle étant associée à une saturation en hydrocarbure liquide, ladite saturation en hydrocarbure liquide étant caractérisée par une cote d'interfacehydrocarbure liquide-eau z, ayant une probabilité p(zw) de distribution comprise entre une valeur maximale zw,min, et une valeur minimale zw,m', et par une distribution triangulaire de la saturation en hydrocarbures pour chaque cote h au-dessus de cette cote d'interface, ladite distribution de saturation en hydrocarbure étant définie par une valeur minimale CSFic,A,min(h), une valeur maximale CSHc,A,m'(h), et un mode CSHcA,mode(h), le procédé comprend en outre : /c-3/ détermination (306c), pour chaque maille c de ladite colonne ayant une cote z et pour chaque faciès A parmi la pluralité de faciès, de la valeur zw,min csHc,miii(z - zw) + CSFic,mode(Z - Zw) + CSFic,max(Z - Zw) f p(Zw) - dzw dans lequel, la détermination de l'étape /d/ est fonction de la valeur déterminée à l'étape /c-3/ pour chaque maille c de ladite colonne et pour chaque faciès A parmi la pluralité de faciès. 5. Procédé selon la revendication 4, dans lequel, ladite saturation en hydrocarbure liquide étant caractérisée par une cote d'interface hydrocarbure liquide-gaz zg ayant une probabilité p(z9), de distribution comprise entre une valeur maximale Zemin, et une valeur minimale zemax, le procédé comprend en outre : /c-4/ détermination, pour chaque maille c de ladite colonne ayant une cote z, de la valeur 0 si z Z Zemin z p(za)dza si zemax z emax 1 si z zemax dans lequel, la détermination de l'étape /d/ est fonction de la valeur déterminée à l'étape /c-3/ pour chaque maille c de cote z de ladite colonne et pour chaque faciès A parmi la pluralité de faciès multiplié par la valeur déterminée à l'étape /c-4/ pour chaque maille c de cote z de ladite colonne. 6. Procédé de détermination d'une carte de variance de hauteur d'hydrocarbure zw,max 3liquide dans un modèle géologique, le procédé comporte, pour au moins une colonne du modèle géologique : /e/ détermination d'une valeur d'espérance de hauteur d'hydrocarbure liquide dans le modèle géologique pour ladite colonne, selon la revendication 1 ; /f/ détermination d'une valeur d'espérance de hauteur d'hydrocarbure liquide pour ladite colonne de manière analytique en fonction de la somme déterminée à l'étape /e/. 7. Procédé selon la revendication 6, dans lequel, les valeurs des cubes de proportions pour deux cellules différentes du modèle et pour des éléments architecturaux différents sont considérées indépendantes entre elles. 8. Procédé selon l'une des revendications 6 ou 7, dans lequel, le procédé comporte en outre : /g/ réception d'un corrélogramme p(A, Az) selon une direction de ladite colonne pour ledit faciès A dans un élément architectural donné ; /h/ détermination de l'espérance du produit de la présence du faciès A dans ledit élément architectural dans une maille cl par la présence du faciès A dans ledit élément architectural dans une maille c2, la distance entre cl et c2 selon ladite direction du corrélogramme étant Az, la probabilité de présence du faciès A dudit élément architectural pour la maille cl étant n i A,AE,c15 la probabilité de présence du faciès A dudit élément architectural pour la maille c2 étant n i A,AE,c2 , en fonction de : (1 - p(A, Az)) - (Esp(n \i- A,AE,c1). ESP(PA,AE,c2) + _\II 1 ar(pA,AE xi) ' V ar(p A,AE,c2) 1 \ 1 + p (A, Az) - (-2E sp ( n \i- A,AE,c1 j + - ESP (PA,AE,c2)) 2 dans lequel, la détermination de l'étape /f/ est fonction de la détermination del'étape /h/. 9. Dispositif de détermination d'une carte d'espérance de hauteur d'hydrocarbure liquide dans un modèle géologique, le modèle géologique comportant des mailles pouvant être associée à un faciès parmi une pluralité de faciès et à un élément architectural parmi une pluralité d'éléments architecturaux, le dispositif comporte, pour au moins une colonne du modèle géologique : /a/ une interface pour la réception d'un cube de répartitions d'élément architectural pour chaque élément architectural AE parmi la pluralité d'éléments architecturaux, le cube de répartitions décrivant pour chaque maille c du modèle la probabilité pAE(c) que ladite maille appartienne audit élément architectural, /b/ une interface pour la réception d'un cube de proportions pour chaque faciès parmi la pluralité de faciès et pour chaque élément architectural parmi la pluralité d'éléments architecturaux, le cube de proportions associé à un faciès A décrivant pour chaque maille c du modèle la proportion dudit faciès dans ladite maille si ladite maille appartient audit élément architectural AE, ladite proportion étant une variable aléatoire de distribution triangulaire définie par une valeur minimale pA,As,min(c), une valeur maximale pA,AE,max(c), et un mode PA,AE,mode(C), /c/ un circuit pour la détermination, pour chaque maille c de ladite colonne et pour chaque faciès A de la pluralité de faciès, de la somme sur la pluralité d'éléments architecturaux de PAE ) - -3(P A,AE,min(c) PA,AE,max(c) PA,AE,mode(c)). /d/ un circuit pour la détermination d'une valeur d'espérance de hauteurd'hydrocarbure liquide pour ladite colonne de manière analytique en fonction de la somme déterminée par le circuit /c/ pour chaque maille c de ladite colonne et pour chaque faciès A. 10. Produit programme informatique comportant des instructions pour la mise en oeuvre du procédé selon l'une des revendications 1 à 8, lorsque ce programme est exécuté par un processeur.