Copyright (c) Kout HASSeNE
Permission is granted to copy, distribute and/or modify this
document under the terms of the GNU Free Documentation License, Version 1.2 or
any later version published by the Free Software Foundation;
with no Invariant Sections, no Front-Cover Texts, and no
Back-Cover Texts. A copy of the license is included in the section entitled
"GNU
Free Documentation License".
Détermination et étude numérique de
conditions aux limites
transparentes. Application à la simulation du
procédé CHFR de
logging de la société SCHL UMBER CER
.
Hassene KOUT Encadrement : DESS Ingenierie Mathematique
option
calcul scientifique et modélisation me-- A Schlumberger
RPC :
canique à l'Université de Savoie. Martin
G.Luling,
année universitaire (2002-2003) Isabelle Dubourg.
A l'INRIA de Rocquencourt : Jacques Henry.
Résumé
Dans le domaine de la prospection pétrolière,
les mesures des propriétés physiques de la roche doivent
être en phase de production c'est le reservoir monitoring.
Malgré les conditions opératoires difficiles (
notamment la présence d'un casing), Schlumberger a
développé un outil efficace pour mesurer la
résistivité en trou fermé le CHFR, un outil d'injection de
courant.
Afin de maItriser et étendre son utilisation, une
simulation fiable est souhaitable. Deux méthodes sont actuellement
utilisées. Le modéle des télégraphistes approche le
problème par une résolution analytique simpliste. Une autre
résolution est implémentée dans le programme Cwnlat.
Elle résout par méthode d'éléments finis les
équations aux dérivées partielles qui gouvernent la
distribution du courant.
Au cours d'un log simulé Cwnlat doit
résoudre plusieurs fois le même problème, ce qui est lourd
en temps de calcul.
L'étude menée dans ce stage a permis d'optimiser
le temps de calcul de Cwnlat en appliquant une méthode de
factorisation.
La méthode de factorisation permet d'effectuer les
calculs sur un sous-domaine fixe qui contient l'ensemble des positions de la
sonde en déterminant des conditions aux limites sur le bord du
sous-domaine qui résument exactement le comportement de la solution sur
le domaine enlevé. Le nombre des inconnues est ainsi réduit.
Une expérimentation sur une configuration
simplifiée a montré une bonne optimisation du temps de calcul.
Le rapport fournit des figures qui illustrent l'étude. Le
rapport a été rédigé en juillet.
Table des matières
1 Introduction 7
2 Le cadre physique 9
2.1 La vie d'un puits de pétrole 9
2.2 Le modèle pour l'interprétation des mesures
10
2.3 Contexte et but du stage 11
2.4 Modéle adopté dans la simulation 13
2.5 Fonctionnement de l'outil 13
3 Le cadre mathématique 15
3.1 Mise en équations 15
3.1.1 Les conditions aux limites 17
3.2 Méthode des éléments finis 22
3.2.1 Formulation variationnelle 22
3.2.2 Forme discrete 23
3.2.3 Discrétisation Q1 24
4 Méthode de Factorisation 29
4.1 Motivation 29
4.2 Calcul formel de la méthode de factorisation 31
4.3 Forme discréte 34 4.4 Prise en compte des
conditions aux limites transparentes dans
Cwnlat-bis 42
4.4.1 Formulation variationnelle 43
4.4.2 Forme discrete 44
4.4.3 Discrétisation Q1 46
4.5 Tests 48
5 Conclusion 59
Bibliographie 60
TABLE DES MATIERES
Chapitre 1
Introduction
Schlumberger Ribout Product Center (SRPC) est un centre
d'études et production de la section Wireline & Testing de la
branche Oilfield Services de la grande multinationale Schlumberger. Le centre
est situé à Clamart au sud de Paris, et emploie environ 500
personnes.
Les activités de SRPC se concentrent autour de la
conception et la simulation des outils de mesures dans les puits de
pétrole. Le centre dispose en outre de deux puits d'essais pour tester
les outils. La chaIne des opérations menées est assez complete on
va de la modélisation par CAO à la réalisation, en passant
par tous les aspects électroniques, mécaniques, matériaux,
etc. Cependant, SRPC a de plus en plus tendance à sous traiter quelques
étape de la production.
SRPC est en liaison étroite avec le centre de
Schlumberger situé à Houston au Texas, de nombreux cadres se
déplacent fréquemment de Paris à Houston et vice versa. Et
comme à Houston, le niveau du prix du pétrole influe
considérablement sur l'activité du centre.
Parmi les unités de SRPC, le projet CHFR travaille sur
l'outil CHFR de mesure de résistivité en puits fermé. Il
est constitué d'une dizaine de personnes qui travaillent entre autre sur
l'éléctronique et la simulation.
A cause d'un problCme de temps d'éxécution de
CWNLAT qui est long, le projet a cherché une collaboration avec le
projet ONDES à l'INRIA Rocquencourt, C'est ainsi que s'est
présenté le stage.
Le projet ONDES à l'INRIA travaile sur le
modélisation et résolution de problCmes de propagation d'ondes
électromagnétiques, acoustiques etc. Il comprend environ six
personnes.
Le stage était donc dans un contrat de collaboration
entre Schlumberger et l'INRIA, il s'est déroulé en partie
à Clamart mais surtout à Rocquencourt, et ce de Avril à
juillet 2003.
Introduction
Chapitre 2
Le cadre physique
2.1 La vie d'un puits de pétrole
Après avoir localisé (par mesures sismiques) une
formation rocheuse susceptible de contenir du pétrole (souvent
en-dessous d'une couche de roches imperméables tel l'argile), on
procède au forage du puits. Le forage creuse jusqu'un peu en-dessous de
la formation fertile, et pendant cette opératoin, une boue
spéciale est déversée dans le puits pour retenir le
pétrole.
Aprés ce forage, le puits est en"trou ouvert" (open
hole). On en profite pour faire une série de mesures des
paramètres de la roche. Ces mesures sont réalisées par des
instruments électriques et/ou nucléaires, elles vont servir
à identifier un modéle physique de la roche et à inspecter
sa contenance effective en hydrocarbures.
Une fois que ces mesures ont montré la
rentabilité du puits, on vient assurer sa tenue mécanique par un
cylindre en acier le casing, qui va jusqu'un peu en-dessous de la formation
fertile. Il est installé par morceaux et fait couramment 3000 m de
longueur. Le casing est fixé par une cimentation autour de lui. Enfin on
procède à des perforations par explosifs au niveau de la
formation pour laisser sortir le pétrole.
A partir de ce moment, le puits peut être mis en
production, les dispositifs nécessaires sont alors mis à sa
surface (well head). Le puits est "fermé", c'est pourquoi on parle de
puits en "trou fermé" (cased hole).
Le suivi régulier de l'état du puits reste
cependant un élément important de l'exploitation. Des mesures
doivent régulièrement être faites pour contrOler la
qualité du pétrole qui est produit et la structure de la roche
(surtout contenance en pétrole et état mécanique). Le but
est de savoir s'il faut arrêter la production dans un puits ou si l'on
peut continuer à produire
et pour combien de temps. Cet ensemble d'opérations de
prospection pendant la production est connu sous le nom de"reservoir
monitoring".
Les phénomènes et grandeurs physiques qui
peuvent déterminer la vie d'un puits sont nombreux. De l'infiltration
d'eau dans la roche jusqu'à l'influence d'autre puits voisins, des
réponses précises et fiables aux questions posées sont
nécessaires pour faire un bon reservoir monitoring.
Les contraintes de coüts et de possibilités
techniques obligent à"faire parler la roche " a travers le casing toutes
les mesures envisageables sont filtrées par cette géante
électrode. Ceci rend les mesures surtout électriques, trés
difficiles.
De plus, il faut arrêter la production pendant la mesure
pour pouvoir opérer à l'intérieur du casing, une boue
spéciale remplit ce dernier pour contrecarrer la pression du
pétrole. Ceci implique que les mesures doivent être faites le plus
rapidement possible. D'autre part, l'usure rapide des instruments de mesures
qui descendent dans le puits est rapide, il est difficile donc de recommencer
la mesure en cas d'erreur.
En résumé, un bon reservoir monitoring
nécessite des outils de mesures fiables et rapides.
Actuellement, pour la mesure de résistivité de
la roche, deux types de mesures sont faites des mesures nucléaires
où l'on envoie des radiations et recueille les échos, des mesures
par injection de courant en profondeur à travers la casing et on le
recueille au well head.
En règle générale, ces trois
procédés fournissent ce que l'on appelle un log le tracé
d'un signal de sortie en fonction de la profondeur dans la zone fertile. En
fait l'outil de mesure fonctionne par station la mesure est
répétée à plusieurs profondeurs du puits.
L'étape suivante est de faire la liaison entre la
resultat de la mesure et les paramètres physiques intéressant de
la roche (tel sa densité, sa résistivité, sa nature, etc.)
qui permettront de tirer des conclusions en matière de réservoire
monitoring. C'est le métier d'interprétation.
2.2 Le modèle pour l'interprétation des mesures
L'interprétation des mesures n'a pas de sens si l'on ne
se donne pas un modèle de la roche, ce qui revient en fait à
faire un choix des zones dans les quelles on suppose les paramètres
physiques constants.
La première hypothèse serait de dire que la
roche est homogène par couches, ainsi un paramètre physique sera
attribué à une profondeur et on a
un modèle unidimensionnel à symétrie de
révolution autour du casing. Cependant, dans une même couche, les
paramètres physiques peuvent varier radialement par rapport au casing.
c'est le cas par exemple des zones limitrophes au casing qui on
été influencées par le contact avec la boue de forage ou
avec le ciment de fixation, on parle couramment de zones envahies. Si l'on veut
rendre compte de ces effets, il faut considérer un modèle
bidimensionnel (dépendance de la profondeur et du rayon), en
première approximation avec touj ours une symétrie de
révolution autour du casing. Le modèle complet sera
tridimensionnel, où l'on ne suppose plus la symétrie de
révolution. ce dernier modéle est indispensable pour rendre
compte par exemple des infiltrations d'eau d'une seule direction.
Le modèle est élaboré en phase de trou
ouvert, où les mesures précises et fiables permettent de
comprendre la roche. Cependant, les mesures ne peuvent être faites
qu'avant l'installation du casing et donc ne peuvent être
répétées pour suivre l'évolution de la roche.
Actuellement, le dernier outil d'injection en trou ouvert ( le HRLA ) permet,
en ayant plusieurs profondeurs d'investigation de constituer un modèle
bidimensionnel de la roche en matière de résistivité.
Une fois constitué, le modèle en trou ouvert
sera une référence pour les mesures ultérieures. Le
reservoir monitoring consiste en fait à déceler les chagements
d'état de la roche par rapport à ce qui était avant
l'exploitation mais en adoptant toujours le même modèle.
2.3 Contexte et but du stage
La simulation du procédé de mesure est un outil
trés apprécié en reservoir monitoring. Elle permet de
prévoir plus ou moins la réponse de l'outil et permet de
décider de son utilité dans une situation précise. Si l'on
veut localiser ou mesurer un phénomène dans la roche ( par
exemple une infiltration d'eau ), il serait profitable de vérifier au
préalable par simulation si l'outil sera sensible à ce
phénomène, et si oui dans quelles mesures.
La simulation d'un procédé de mesure
nécessite en premier lieu une modélisation de la roche, ce qui a
été fait dans les mesures en trou ouvert. Ensuite, il faut
modéliser et mettre en équation intégrales et/ou
différentielles les phénomènes physiques qui entrent en
jeu. Il faut enfin un outil informatique pour résoudre les
équations et calculer les grandeurs qu'aurait sorties l'outil en
réalité. C'est en examinant le résultat des simulations
que
l'on va décider de l'utilisation ou non d'un outil
particulier pour mesurer ou quantifier un certain aspect des
propriétés physiques de la roche.
Ce stage s'inscrit dans les efforts menés à
Schlumberger pour la simulation du CHFR ( Cased Hole Formation Resistivity )
qui est un outil d'injection de courant et qui sert à mesurer la
résistivité dans différentes parties de la roche en trou
fermé.
La mesure de résistivité est l'une des plus
anciennes expertises de Schlumberger dans les service pétroliers. La
résistivité est une information précieuse dans la
prospection pétrolière puisqu'elle indique directement la
saturation de la roche en hydrocarbures. Elle permet donc au départ de
voir la rentabilité du puits et pendant le reservoir monitoring de
suivre à quelle vitesse la roche perd son pétrole, et donc
d'estimer l'espérance de vie du puits.
Le CHFR est une sonde qui descend au fond du puits et puis
remonte en s'arrêtant à des stations pour injecter
latéralement un courant à travers le casing. Le courant
stabilisé à une valeur fixe retourne à la surface du
casing, tandis que des mesures de potentiel en-dessous du point d'injection
permettent d'estimer le courant qui pénètre dans la formation. Le
courant et le potentiel permettent d'estimer le résistivité de la
formation.
Le CHFR est actuellement commercialisé pour le
reservoir monitiring et il a été utilisé plusieurs fois
à travers le monde, cependant on ne dispose pas d'une simulation fiable
de son comportement. Une bonne simulation sera surtout utile pour adapter
l'outil à de nouvelles situations, tel par exemple une avance d'un front
d'eau sur le puits ou une élévation du niveau de la nappe
phréatique en dessous.
Schlumberger dispose d'un code de calcul Cwnlat qui a
été adapté au cas du CHFR, et d'une résolution
analytique du phénomène électromagnétique dans la
roche.
Le but de ce stage est d'essayer d'optimiser le temps de
calcul de ce code de calcul Cwnlat, en appliquant une méthode de
factorisation de domaine qui sera présenter et tester plus loin.
2.4 Modéle adopté dans la simulation
Le modéle adopté pour la simulation de CHFR dans
ce stage sera bidimensionnel, avec symétrie de révolution autour
du casing. Il est résumé dans la figure ci dessous, qui est une
coupe du domaine complet et qui suffit à la simulation puisqu'on a
supposé une symétrie de révolution autour du casing. On
suppose de plus dans le modéle que la formation s'étend à
l'infini, ce qui n'est pas une trés grave approximation.
2.5 Fonctionnement de l'outil
Le courant est injecté à la surface
intérieure du casing, et il retourne à la surface comme le montre
le figure 2. du fait de le forte conductivité de l'acier du casing, la
plus grande partie du courant cirdule dedans et seule une faible partie
traverse la formation avant de regagner le Well head. Les lignes de courant (
tangentes au vecteur densité de courant ) sont donc surtout
concentrées dans le casing, Elles partent toutes du point d'injection et
reviennent au Well head. On note J le courant total dans le casing et j le
courant qui fuit normallement par unité de profondeur à la
surface extérieure de casing. Le contraste de résistivité
entre le casing et le reste du domaine (environ 107) est tellement
grand que le casing est une équipotentielle. De ce faits, les lignes de
courant à sa sortie sont quasi normales et on a j=-dJ/dz. Pour plus de
détaille voir la référence.
Environ cinq métres en dessous de l'injection, trois
potentiométres mesurent le potentiel électrique V aux points A,B
et C. L'outil calcule alors la dérivée seconde du potentiel par
rapport à la profondeur qui, divisée par Rc, n'est autre
que j au point B. En pratique la mesure est un peu plus compliquée vu
qu'il faut estimer Rc par une premiére phase de calibrage.
L'outil fait cette procédure à plusieurs
stations en remontant tout au long de la zone qui nous intéresse de la
formation ( de 300 à 800 m d'epaisseur ). L'hypothése
d'interprétation consiste à dire que la résistivité
de la roche en face du point de mesure B (et non celle du ciment ou de la zone
envahie) est proportionnelle au quotient VB/jB. Le
résultat final des mesures est donc k.VB/jB le long de
la formation, le coefficient k est une constante estimée d'après
la géométrie du casing.
Pour faciliter le traitement électronique, l'outil
injecte un courant
alternatif (AC). Cependant, le procédé est
limité à une trés basse fréquence (quelques Hertz)
de peur que l'épaisseur de peau dans le casing ne soit trop petit et
fasse un écran total entre le courant injecté et le formation.
-*E)=-,i.?
-*
= -i.W.,i.H (3.5)
-* rot(
-* H ?t
Chapitre 3
Le cadre mathématique
3.1 Mise en equations
champ électrique -E,
le champ magnétique
Le phénomène physique de circulation de courant
qui se manifeste dans le fonctionnement de l'outil obéit aux
équations de Maxwell. Avec le **-H,
la permittivité ,i, la conduc-
tivité électrique or, la
permittivité électrique å0 et le
densité volumique de courant *-I,
ces équations sont :
-* rot(
-*E) = -,i. ?-* ?t
(3.1)
H
-* rot(
-* H) =
*-I +E0. ?-* ?t (3.2)
E
div(-* E ) = 0 (3.4)
Le courant injecté étant sinusoIdal à la
fréquence f = 2ð, on se
W
place en régime harmonique où l'opérateur
de dérivation se remplace par la multiplication par i.W (on
adopte cette convention et non -i.W, mais c'est équivalent)
avec i2 = -1. Seules les trois
premières équations sont utiles à la modélisation
et elles donnent :
--* rot(
|
--* H) = or.
|
--*E + C0.?--* E
?t
|
--*
= (or + i.ù.C0).E (3.6)
|
L'équation(3.6) nous donne :
rot( 1
--* . or + i.ù.C0
|
--* rot(
|
--* H)) +
i.ù.ii.(--* H) =
*--0 (3.7)
|
Dans le modèle qu'on a adopté pour la simulation
au chapitre précédent, on suppose une symétrie de
révolution autour de l'axe z (l'axe du casing) en ce qui concerne les
propriétés physiques des différents milieux du domaine
(c'est à dire or et ii). Il est naturel alors de
supposer une telle symétrie dans les champs électriques et
magnétiques.
Dans ce cas (mode transverse magnétique - TM) ces champs
s'expriment par :
*-- H = H(r, z).--*
u9(3.8)
* --E = Er(r, z).--*
ur + Ez.--* uz (3.9)
*--u9 =
|
uz ? *--
|
*--ur
|
(3.10)
|
Q uelques valeurs numériques :
Dans la roche or est comprise entre 0.1
Ù-1m-1 et 100
Ù-1m-1,
ii = 4.ð.10-7
5.I
Dans le casing or
5.106Ù-1m-1,
ii = 100.4.ð.10-7 5.I
Dans le domaine C0 = 36.ð.1095.I,
et enfin f vaut quelques Hertz.
1
On constate déjà que le produit
ù.C0 est complètement négligeable devant
or, ce qui permet de le supprimer. L'équation devient alors
avec la seule inconnue H, la composante orthoradiale du champ
magnétique, dépendant de r et de z :
*--u9 ·
|
rot(1 --* .
or
|
rot(H.--*
--* u9)) + i.ù.ii.H = 0 (3.11)
|
· est le produit scalaire euclidien
La densité volumique du courant est :
|
*--I=
or.
|
--* E = --* rot(H.--*u9).
|
L'équation (3.11) est une équation aux
dérivées partielles en H, les propriétés sur le
courant lui fourniront les conditions aux limites nécessaires
pour avoir un problème bien posé dans un certain
domaine géométrique I. La simulation du CHFR se résume
donc à la résolution de ce problème aux
dérivées partielles avec méthodes numériques
adéquates.
3.1.1 Les conditions aux limites
On résout donc l'équation (3.11) dans un domaine
tridimensionnel I de frontière 8I. Dans le cas du CHFR, I est
le demi espace en dessous de la surface, privé de l'intérieur du
casing (la boue qui le remplit est isolante). On choisit aussi d'enlever toute
la partie dans la continuation du casing, ceci correspond en fait à une
simplification qu'on verra plus loin. Une coupe du domaine par symétrie
de révolution autour de l'axe, est donnée dans la figure 3.1,
comme adopté dans le chapitre précédent, ce domaine
s'étend à +8 en r et à -8 en z. Le bord
8I est donc la réunion des différents segments
Fi,
Fr, F1, F0, F2 et F.
Pour mettre en évidence les propriétés
physiques qui conduisent aux conditions aux limites pour notre équation
(3.11), on utilise le théorème d'Ampère, qui dit que la
circulation du champ magnétique le long d'une courbe fermée vaut
le flux total du courant à travers toute surface limitée par
cette courbe.
On prend un disque horizontal D centré sur l'axe du
casing et de rayon r1 et dont le centre est au-dessus du
casing shoe. Vu que la boue est isolante, le seul courant circulant à
l'intérieur de D est éventuellement celui du cable. Au dessus de
l'électrode d'injection, ce courant vaut -I0 car le
sens conventionnel est celui de l'axe des z. En dessous de l'électrode
d'injection, ce courant est nul. Une première approximation serait de
dire que le courant est nul aussi en dessous du casing shoe bien que, la roche
étant intacte dans cette zone, on puisse avoir des courants. Mais cette
approximation revient à dire que le courant qui y circule est
très faible car de façon prépondérante il remonte
vers la surface pour rejoindre l'électrode de recueil.
?-
Ainsi, commeHH.-? n=
2.ð.r1.H(r1,z), on a comme
conditions
?D
aux limites sur sur F1 :
2.ð.r1.H(r1,z) =
-I0 (3.12)
Le cadre mathematique
z
1'
r
1'
Ù
r
1'
1
r
2
r
1
1'
i
L
1'
1'
0
1'
s
1'
2
1'
Axe de symetrie de revolution
FIG. 3.1 Coupe du domaine de resolution, le casing
hachuré. Les points noirs délimitant les frontières
designees.
sur F0 U F2 on a la condition aux limites suivante :
2.ð.r1.H(r1, z) = 0
(3.13)
En ce qui concerne la frontière F, où z
= z_ ou r = r+ , on a un champ
magnétique nul car très éloigné des sources. Ainsi
on aboutit aux conditions aux limites :
sur F : 2.ð.r.H(r, z surface) = 0
(3.14)
Restent les électrodes d'injection et de retour, qui
sont respectivement un anneau cylindrique petit en épaisseur (1 cm
environ) et l'électrode de retour = anneau à la surface. Sur ces
zones, le plus simple est de sup-poser un potentiel électrique constant,
ou un champ électrique normal à la surface de l'électrode,
c'est équivalent. En notant *--m
la normale extérieure et
*-- r(.) · *--m
la dérivée normale, on a donc sur Fi *--
m= ----* uret
8 8m
8V 8m
|
1
=--Ez =--.
a
|
u9 · --*
|
rot(H.--*
--* u9)
=1.(--*m A *
--u9) · --* u9) = 1 8(r.H)
rot(H.--* .
a a.r8m
(3.15)
|
On va supposer sur les électrodes d'injection et de retour
des conditions aux limites plus générale du type :
1
sur Fi et Fr :
.(--* m A * --u9) ·
--*rot(H.--* u9) -- 'y.H = C
(3.16) a
Ici G est un champ électrique imposé et 'y
est en quelque sorte l'impédance de l'électrode. Cette
condition aux limites est valable sur la partie bidimensionnelle de 8I
correspondante.
Finalement, en utilisant l'expression du rotationnel en
cooordonnées
--* 8H
cylindriques r et z :rot(H.--*uz
+18r.H
u9) = -- 8z .--* 8r .--*
.ur et en posant u=r.H,
r
on obtient de l'équation (3.11) :
8r( 1
8 a.r.8u 8r ) + 8 8z (
1a.r.8u8z ) -- i.ù.u r .u=0 dans I
(3.17)
Et l'equation(3.16) devient :
?u
+ã.ó.r = --G.ó.r
(3.18)
?n
Les conditions aux limites de notre problème sont donc
:
surF1 :u=--
2.ð(3.19)
I0
surF0UF2UF:u=0
(3.20)
?u
sur Fi U Fr: +
ã.ó.r = --G.ó.r (3.21)
?n
Ainsi le problème est bien posé et on pourra
appliquer une méthode de résolution numérique. Pour cette
résolution, on sera obligé de remplacer les frontières
à l'infini par des frontières ayant des coordonnées finies
assez grandes par rapport à la longueur du casing L. La
méthode de résolution utilisée est la méthode des
éléments finis.
Avant d'appliquer cette méthode de résolution,
une approche par perturbations a été appliquée au
problème qui a permis de considérer le casing comme une
frontière (ceci en réduisant son épaisseur) et donc de
créer de nouvelles conditions aux limites sur le casing, cette nouvelle
approche a été appliquée durant le stage de M. El
Moussawi, pour plus de détails sur cette approche veuillez consulter la
référence [1].
Avec les notations de la figure (3.1)
r
S
1
(PV)
|
?
? ????????????????????
?????????????????????
|
a 1 au a 1 au ù.u
Si
S0
S
2
FIG. 3.2 notations dans le domaine I En appliquant cette
méthode le problème devient
sur I, a aa a
( .)+ ( . ) -- i. .u = 0
ró1.rr z ó1.rz r
au
sur Fr, = 0,
az
au
sur S1, -- A.u = --B
ar
au
sur Si, --C.u=0
ar
au
sur S0, --A.u=0
ar
sur S2?F, u=0
A,B et C sont les constantes suivantes :
ó1
A = .
ó0
|
r2 -- r1 ;B=U0
1 ó1.
ó0 r2 -- r1
.C=0. 1
|
Où r1 : le rayon
intérieur du casing
r2 : le rayon éxterieur du casing
a0 : la résistivité du casing
a1 : la résistivité de la
formation
U0 : l'intensité du courant
injecté à travers le casing
3.2 Résolution par la méthode des
éléments fi-
nis
3.2.1 Formulation variationnelle
On a donc le problème (PV) à
résoudre par la méthode des éléments finis.
Soit U = {v E H1(Ù); v
= 0 sur S2 U F} l'espace des fonctions admis-
sibles pour ce problème, (PV) est
équivalent à la formulation variationelle suivante
(
où (,) représente le produit scalaire dans
L2(Ù).
On applique les formules de GREEN on obtient :
Or(1
O a1.r.Ou Or), v)
+(O Oz(1
a1.r.OuOz ),v)-(i. W.u
r.u,v)=0 (3.22)
(i )
a1.r.Ou Oz, Ov ) +
(1 du
Oz ) 0
.u, v . ra1.rdn
(1 .
a1.r
OuOr ,
(1 .
a1.r
dn, v)r =0 (3.23)
du
En substituant on obtient :
(1 .
a1 .r
|
OuOr ,
|
Ov Or ) - ( 1
a1.r.OuOz, OvOz) -
(i.W.ur .u,v) + ( 1
a1.r2.du
dn,v)S1?S0?Si +
|
Alors,
(1 .
ó1.r
|
?u ?r ,
|
?v ?r ) + ( 1 r.u,v)+
ó1.r.?u ?z , ?v
?z )+(i. ù.u
|
( 1
ó1.r2(A.1S1?S0(z)
+ C.1Si(z))u,
v)r=r2 = (B
ó1.r2
.1S1(z),
v)r=r2 (3.24)
Donc le problème précédent admet pour
formulation variationnelle
½ trouver u E U véri~ant
( PV) a(u,v-u)=l(v-u) pour
toutvEU
Où a et l sont respectivement la forme
bilinéaire symétrique et la forme linéaire données
par :
?v ?r ) +(1
ó1.r
.?u?z,?v
?z) + (i. ù.u r.u,v)
a(u,v) = (1
?u ?r ,
.
ó1.r
+( 1
ó1.r2
(A.1S1?S0(z)
+C.1Si(z))u,
v)r=r2
l(v) = ( B 2
.1Si(z), v)S1
.r
ó1
3.2.2 Forme discrete
La méthode consiste à chercher une approximation
uh de u dans un sous espace Uh de dimension finie
mr * mz dont on désignera par ?i,j(i = 1,
..., mr et j = 1, ..., mz) des fonctions de
bases.
On peut donc décomposer uh et vh dans la
base des ?i,j(r, z) :
u(r, z) = Xnr
Xnz
i=1 j=1
v(r,z) = Xnr
Xnz
p=1 q=1
ui,j?i,j(r, z)
vp,q?p,q(r, z)
Le problème discret associé à (PV)
est donc ½ trouver uh E Uh vérifiant
(PV)h a(uh, vh - uh) =
l(vh - uh) pour toutv E Uh
En injectant dans (PVh) les expressions de uh
et de vh, il vient, en
notant wi,j= vi,j-
ui,j
? ??
??
|
trouver u1,1, u1,2, ..,
u2,1, u2,2, .., unr*nz solution
de
Xnr Xnz wi,j
Xnr Xnz a(?i,j,
?p,q)up,q = Xnr Xnz
wi,jl(?i,j) V ù E <nr * <nz i=1
j=1 p=1 q=1 i=1 j=1
|
Les fonctions de bases étant connues, les valeurs
prises par a et l sur celles-ci sont des coefficients connus.
Notons les par Ai,j,p,q et par Fp,q respectivement :
Ai,j,p,q = a(?i,j, ?p,q) Fp,q =
l(?p,q) i,p E (1, ..., r) j, q E (1, ...,
z)
Le problème approché consiste donc à trouver
r * z nombres u solutions des équations
linéaires suivantes :
Xnr Xnz wi,j
Xnr Xnz Ai,j,p,qup,q =
Xnr Xnz wi,jFi,j V ù E <nr
* <nz i=1 j=1 p=1 q=1 i=1
j=1
En optant pour une numérotation des noeuds par ligne :
u 1 = u1,1, u 2 =
u2,1, ..., u nr = u1,2, u
nr+1 = u2,2, ..., u nr*nz = unr,nz
Les équations linéaire précédente
s'écrivent sous forme condensée
trouver le vecteur U E <nr*nz tel que AU =
F
Où la matrice de rigidité (stiffness matrix)
A est de taille ( r * z) * ( r * z) (symétrique
définie positive) de coefficients Ai,j,p,q , et le second
membre F est un vecteur de composante fp,q ·
3.2.3 Discrétisation Q1
On réalise une rectangulation du domaine I comportant
r * z noeuds et e éléments
rectangulaires :
Ù = UneÙe
Notons k le pas de discrétisation en z et h le pas de
discrétisation en r (supposés ici constants).
Les fonctions de bases ?i,j sont choisies affines par
morceaux pour chaque coordonnée séparément
(élément Q1) valant 1 au noeud numéro (i,j) et 0
aux autres :
où
|
?r i (r)=
?z j(z)=
|
? ??
??
? ??
??
|
?i,j(r, z) =
?r i (r)?z
j(z)
(r - ri_1)
hsi r E
[ri_1,ri]
(ri+1 - r)
hsi rE [ri,
ri+1]
(z - zj_1)
ksi z E [zj_1,
zj]
(zj+1 - z)
ksi z E
[zj,zj+1]
Avec h = ri - ri_1 V i = 1, .., nr
et k = zj - zj_1 V j = 1, .., nz
Donc la matrice de rigidité A est construite de
la facon suivante :
Ai,j,p,q = ( 1
w.u
.Vr?i,j, Vr?p,q) + (1 .Vz?i,j,Vz?p,q) +
(i. .?i,j, ?p,q) -
ó1.ró1.rr
( 1
ó1 .r2
(A.1S1?S0 (zj) +
C.1S% (zj))?1,j,
?1,q)r=r2
= ( 1 .Vr?r i , Vr?r
p)(?z j, ?z
q) + (Vz?z
j,Vz?z q)( 1
ó1.r.?r i , ?r
p)
ó1.r
w.u
+(i.r
|
.?i,j, ?p,q) - (1
ó1.r2
(A.1S1?S0(zj) +
C.1S%(zj))?1,j,
?1,q)r=r2
|
Pour w = 0 on a
Ai,j,p,q = ( 1 .Vr(pr i ,
Vr(prp)((pz
j,(pz q) +
(Vz(pz j,Vz(pz q)( 1
ó1.r.(pri,
(prp)
ó1.r
-(1
ó1.r2(A.1S1uS0(zj)
+ C.1S%(zj))(pr
1(pz j,(pr 1(pz
q)r=r2
Le ehoix des fonetions de bases fait que A est
pentadiagonale. Considérons les matrices T, F, Mj de
taille nr * nr définies par :
Mj = [mjt,e],
aveemjt,e= (1
ó1.r2
(A.1S1uS0(zj) +
C.1S%(zj))ä1,e,
seul
mj1,1 est non nulle.
T = [tj,q] = ((pj, (pq), done
tj,q = käj,q
D = [di,p] = ( 1 1
.(pi, (pp), donedi,p =
äi,p
ó1.ró1.ri
F = [fi,p] =
(1.Vr(pi,Vr(pp), done
ó1.r
+
1
ó1
(1
rt+1
1( 1 1
+ ) sii=p+1
ó1 rt rt+1
1(1 1
+ )
ó1 rt rt_1
sii = p-1
1
+ ) sii=p
rt rt_1
2
-
-
1
fi,p = 2
?
? ???????????
????????????
Et B la matriee de taille nz * nz
définies par
B = [bj,q] = (Vz(pj, Vz(pq),
done
bj,q =
|
?
? ???????????????
????????????????
|
k si j = q avec
i=61et p=61
2
k sij=q=1
1
k sij = q+1
1
k sij = q-1
1
Avec ces notations, la matrice A s'écrit par
exemple avec mr = 5 et mz = 3:
A=
|
?
? ? ? ? ? ? ? ? ? ? ? ? ?
|
(F+M1). k 2 + 1
k.D
- 1 k.D
0
- 1
k.D k.(F+M2)+ 2
k.D
-
?
? ? ? ? ? ? ? ? ? ? ? ? ?
1 k.D0
- 1
k.D k.(F+M3)+ 2
k.D
Les 0 représentent des matrices nulles de taille mr *
mr.
Pour le second membre, si on note V j le
vecteur à mr éléments seul le premier coefficient
est non nul et vaut ó1.r2
.1S1 (zj), Alors F
s'écrit :
B
?
? ?
? V 1.k
? 2
F = ? V2.k
V 3.k
Nous utiliserons un gradient conjugué pour
résoudre le systéme linéaire A.U = F.
La discrétisation et le gradient conjugué sont
codés en Scilab sous le nom Cwmlat - bis, ce nom est
inspiré du nom Cwmlat qui est le nom d'un code
éléments finis chez Schlumberger.
Cwnlat résout le même problème en
se basant sur la même formulation variationnelle et la même
discrétisation, seulement le maillage utilisé est à pas
variable contrairement à celui utilisé dans Cwnlat - bis,
pour plus de détails sur le code Cwnlat consulter la
référence [31.
Données numériques utilisées dans Cwnlat
- bis
r1 = 0.078486 m Rayon intérieur du
casing.
r2 = 0.0889 m Rayonextérieurducasing.
L nes de courant c culées pa
rinf 2540 m L'infini en r.
0
ele1 -3000 m Extrémité supérieure de
l'électrode d'injection.
00
ele0 -3500 m Extrémité inférieure de
l'électrode d'injection.
shoe -4500 m Extrémite inférieure du casing.
0
zinf -5300 m L'infini en z.
ó1 = 0.01 Ù
Résistivité de la formation.
E 0.0001 Coefficient de régularisation.
00
FIG. 3.3 - Lignes de courant calculées par Cwnlat-bis dans
une formation homogène à 0.01Ù
Chapitre 4
Méthode de Factorisation
4.1 Motivation
Au cours d'un log simulé le CHFR effectue plusieurs
stations en remontant. A chaque station un courant est injecté à
travers le casing, donc pour chaque station le programme Cwnlat doit
résoudre sur tout le domaine Ù (qui est assez grand) le
même problème (PV). En général le nombre de
stations est assez élevé, c'est pourquoi le temps de calcul du
programme Cwnlat pour un log simulé est très long.
A maillage fixé on se propose d'effectuer les calculs
sur un sous domaine fixe Ù1 qui contient l'ensemble des positions de la
sonde. Le nombre des inconnues est ainsi réduit. Le problème est
alors de déterminer des conditions aux limites sur le bord
t(Ù1 est la partie du domaine qu'on souhaite enlever) qui
résument exactement le comportement de la solution sur Ù1, ou
autrement dit d'éliminer les inconnues sur Ù1.
Les conditions aux limites transparentes seront donc
calculées sur le domaine t, en
utilisant l'opérateur Dirichlet-Neumann :
=Q.u+q. (4.1)
Où Q est un opérateur autoadjoint et
q un vecteur. Ainsi grace à cette opérateur, si on
suppose connaitre une solution u de la famille des solutions de notre
problème sur t de passer à une condition aux
limites de type Neumann sur la frontière transparente,pour plus de
détails sur cette méthode veuillez consulter la
référence [21.
Ces conditions aux limites transparentes une fois
calculées seront inté-
grées dans le programme Cwnlat - bis qui se base
sur le cadre mathématique décrit dans le chapitre 3, et plusieurs
tests seront effectués par la suite.
Dans un premier temps je présente le calcul formel de
la méthode de factorisation, puis le calcul discrét avant de
conclure par différents tests. La figure ci dessous donne les notations
pour le développement de la méthode.
r s
Ù1
z
S
1
zt
t
Si
S0
S
2
Ù = Ù Ù
-
2
1
FIG. 4.1 Notations pour le développement de la
méthode
4.2 Présentation du calcul formel de la méthode
de factorisation
(PV)
|
?
? ????????????????????????
?????????????????????????
|
? 1 ?u ? 1 ?u ù.u
On rappelle le problème à factoriser :
sur Ù, ? ?? ?
( . ) + ( . ) -- i. .u = 0
ró1.rr z ó1.rz r
?u
sur Fr(électrode de recueil), = 0
?z
?u 1
surFs, + .u=0
?z E
?u
sur S1, -- A.u = --B
?r
?u
sur Si, --C.u=0
?r
?u
sur S0, --A.u=0
?r
sur S2?F, u=0
Où A,B et C sont des constantes :
ó1
A = .
ó0
|
r2 -- r1 ;B=U0
1 ó1.
ó0 r2 -- r1
;C=0.
1
|
Avec r1 : le rayon intérieur
du casing
r2 : le rayon extérieur du casing
ó0 : la résistivité du
casing
ó1 : la résistivité de la
formation
?u 1
+
?z
U0 : l'intensité du courant
injecté à travers le casing
.u = 0, le coefficient E est choisi très
petit (10-6).
E
Avant de commencer les calculs, j'aimerais insister sur le
fait que l'opérateur Q et la fonction q
dépendent de z, puisqu'on a choisi de réduire le
domaine parallèlement à la surface. Et nous allons faire jouer
à z un rOle particulier comparable à celui du temps pour
le contrOle de problèmes d'évolut ion.
Soit U = {v E
H1(Fzt);v = 0 sur F}
l'espace des fonctions admis-
sibles pour ce problème.
Remarque Pour les profondeurs z situées
en dessous du pied de casing
on prendrait U = {v E
H1(Fzt); v = 0 sur S2
U F}
Soit donc ? E U,on a
(1 .
ó1.r
|
Ou Oz, ?) = (Q.u, ?) + (q, ?)
(4.2)
|
On dérive cette équation par rapport à la
variable z :
( O Oz ( 1
ó1.r.Ou Oz ), ?) =
(OQ Oz.u, ?) +(Q.Ou Oz , ?)
+(Oq
Oz ,?) (4.3)
De (PV) on a
(iù.u
r
|
u,?) - (O Or(1
ó1.r.Ou Or ), ?) =
(OQ Oz.u,?)+(Q.Ou Oz
,?)+(Oq
Oz ,?) (4.4)
Et de l'équation (4.1) on obtient
|
(iù.u
r
u, ?) - ( O Or ( 1Oz .u, ?) +
(Q.ó1r.Q.u, ?) +
ó1.r.Ou Or ), ?) =
(OQ
(Q.ó1r.q, ?) +(Oq
Oz ,?) (4.5)
On intègre par partie le terme (O
Or( 1
ó1.r.OuOr) par rapport à la
va-
riable r, ce qui nous donne :
(iù.u
r
|
u, ?)+( 1 .
ó1.r
|
OuO?
Or, Or ) +
ó1.r2(A.1S1?S0(z) +
C.1S%(z))u(r2
, z)?(r2) -
1
|
ó1.r2
1S1(z).?(r2) =
(OQ
B Oz .u, ?) + (Q.ó1r.Q.u, ?)
+ (Q.ó1r.q, ?) + (Oq
Oz , ?)
Cette égalité fait intervenir les
différentes conditions aux limites de notre problème. Le calcul
précédent est valide pour tous les u de la famille
vérifiant l'équation aux dérivées partielles de
notre problème et la condition limite sur Fr. Ainsi
u(r, z) est quelconque dans la famille des solutions, on
en déduit de l'identité précédente
:
ù.u (i.r
|
ø, ?) + ( 1 .Vrø, Vr?) -
(Qó1rQø, ?) - (?Q ?z
.ø, ?)+ ó1.r
ó1.r2(A.1S1?S0(z)
+
C.1Si(z))ø(r2)?(r2)
= 0 (4.6)
1
|
Et, (q,?) + (Q.ó1r.q,?) +
ó1.r21S1(z).?(r2)
= 0 (4.7)
B
Avec ? et ø dans U
Il nous reste à déterminer Q et q
pour z = 0 (en surface). sur F5
?u
on a la condition aux limites +
?z
1 ?u
.u = 0, et sur Fr on a = 0,
E ?z
alors d'après la relation Dirichlet-Neumann on a sur
F5 (Q + D.1)u+ q
= 0,
E
comme u est arbitraire alors on a
{ -1E.D
surF5
Q(0)= 0 surFr q(0) = 0
sur F5
1
äi,j,
Ici D est une matrice diagonale de taille
nr'3*nr'3 telle que [di,j]
=ó1.ri
nr'3 est le nombre de
noeuds sur F5.
L'opérateur Q qui dépend de z
et agit sur des fonctions de la variable r vérifie une
équation de Riccati qui a la particularité de ne pas avoir de
terme linéaire. En réunissant les équations (4.6), (4.7)
et les conditions initiales pour Q et q, on obtient le
système
?
? ??????????????
???????????????
ù.u (i. r
ø, ?)+( 1 .Vrø,Vr?) -
(Q.ó1r.Qø, ?) - (?Q ?z
.ø, ?)+ ó1.r
ó1.r2(A.1S1?S0(z)
+
C.1Si(z))ø(r2)?(r2)
= 0
1
(q,?) + (Q.ó1r.q,?) +
B
1 .
ó1 .r
?u + Qu = -q ?z
ó1.r21S1(z).?(r2)
= 0
1
Q(0)= -.D
E
q(0) = 0
Ce système est découplé car on peut
intégrer les deux premières equations en z de 0 à
zinf pour obtenir Q et q, puis u est obtenu
par intégration retrograde de la troisième equation.
4.3 Forme discrete
= Q.u+qpar
1
On approche la relation Dirichlet-Neumann:
?u .
ó1.r?z
un schema differences finies décentré à
gauche :
u j+1 - uj
D.hj+1
|
= Qj+1.uj+1 +
qj+1 (4.8)
|
1
Ici D est une matrice diagonale de taille nr * nr
[di,j] = .äi,j.
ó1.ri
Et hj+1 = zj+1 - zj.
Ona
uj+1 - uj
(D. , ?) =
(Qj+1.uj+1 , ?) +
(qj+1, ?) (4.9)
hj+1
uj - uj_1
(D. ,?)=(Qj.uj,?)+(qj,?) (4.10)
hj
Après une soustraction entre l'equation(4.9) et l'equation
(4.10) on obtient :
hj.u j+1 -
uj(hj+1 + hj) +
hj+1.uj_1
(D., ?) =
(Qj+1.uj+1 - Qj.uj, ?)
+ hj+1.hj
(qj+1 -qj,?) (4.11)
Pour que le terme de gauche de l'equation (4.10) soit bien
symé-
2
trique nous allons multiplier cette equation par :
hj+1+hj
:
(D. hj+1+hj uj+1 -
2uj + 2hj+1
2.hj hj+1+hj
uj_1 2
hj+1hjhj+1 + hj
, ?)=
[(Qj+1uj+1 -
Qjuj,?) +
(qj+1 - qj, ?)] (4.12)
hj
hj+1+hj 2.hj
Pour alléger l'écriture posons hj+1
2 = 2, áj=hj+1
+
et áj =
|
2.hj+1
|
|
. hj+1 +hj
On obtient :
|
áj.u j+1 - 2.uj +
áj.uj-1 1
(D., ?) = [(Qj+1.u j+1 -
Qj.uj, ?)
hj+1.hj hj +1
2
+(qj+1 - qj, ?)] (4.13)
Une semi-discrétisation de l'équation aux
dérivées partielle (3.17) suivant la variable z nous
donne :
áj.u j+1 - 2.uj +
áj.uj-1
D.hj+1.hj
|
ù.u
=(i. r
|
1
Vr. .Vr).uj (4.14)
ó1 .r
|
Où Vr désigne l'opérateur
?r, alors l'équation (4.13)
devient
?
ù.u ((i. r
|
1 .Vr).uj, ?)= 1
Vr. [(Qj+1uj+1
- Qjuj, ?) + (qj+1 - qj, ?)]
ó1.rhj+ 1
2 (4.15)
|
Qu'on peut encore écrire
ù.u ((i .r
|
1
Vr. .Vr).uj, ?) = 1
[((Qj+1 - Qj)uj,
?) ó1.rhj+1
2
|
+(Qj+1(uj+1 -
uj), ?) + (qj+1 - qj, ?)] (4.16)
De l'équation (4.8) on a
uj+1 - uj = hj+1(D -
hj+1Qj+1)-1(Qj+1u j
+qj+1) (4.17)
Alors l'équation (4.16) devient
ù.u ((i. r
|
1 1
Vr. .Vr)uj, ?) =
[hj+1(Qj+1(D -
hj+1Qj+1)-1Qj+1uj,
?)
ó1.r hj+1
2
+((Qj+1 - Qj)uj, ?) +
hj+1(Qj+1(D -
hj+1Qj+1)-1qj+1,
?)
+(qj+1 - qj, ?)](4.18)
|
1
Une intégration par parties sur le terme ((Vr.
.Vr)uj, ?)zt par
ó1.r
rapport à la variable r nous donne
l'équation suivante :
ù.u (i. r
|
duj
uj, ?) + ( 1 .Vruj,Vr?)
+ (1
ó1.r2. dn ,
?)r=r2
ó1.r
|
1
= [((Qj+1 - Qj)uj, ?) +
hj+1(Qj+1(D -
hj+1Qj+1)-1Qj+1uj,
?) hj+ 1
2
+ hj+1(Qj+1(D -
hj+1Qj+1)-1qj+1, ?)
+ (qj+1 - qj, ?)] (4.19)
On obtient donc l'équation :
ù.u (i.r
|
uj,?) + ( 1 B
.Vruj,Vr?) - ó1.r2
1S1(zj).?(r2)
ó1.r
|
|
+ ó1.r2
(A.1S1?S0(zj) +
C.1S%(zj))uj(r2)?(r2)
1
1
= [((Qj+1 - Qj)uj, ?) +
hj+1(Qj+1(D -
hj+1Qj+1)-1Qj+1uj,
?) hj+ 1
2
+ hj+1(Qj+1(D -
hj+1Qj+1)-1qj+1, ?)
+ (qj+1 - qj, ?)] (4.20)
L'équation (4.20) qui fait intervenir les
différentes conditions aux limites de notre problème. Le calcul
précédent est valide pour tous les uj de la
famille vérifiant l'équation aux dérivées
partielles de notre problème et la condition limite sur z = 0.
Ainsi uj(r, zj) est quelconque on en déduit
de l'identité précédente :
ù.u
hj+ 1 2(i.ø, ?) +
hj+ 12( 1
ó1.r.Vrø,Vr?) -
(Qj+1ø, ?) + (Qjø, ?)
r
1
+hj+ 1 2 .
ó.r2(A.1S1?S0(zj) +
C.1S%(zj))ø(r2)?(r2)
- hj+1(Qj+1(D -
hj+1Qj+1)-1Qj+1ø,
?) = 0 (4.21)
Avec ? et ø dans U Et
B
(qj+1
-qj,?)+hj+1(Qj+1(D-hj+1Qj+1)-1qj+1,?)+hj+1
2 . ó1.r2
1S1(zj).?(r2) =
0
(4.22)
Cherchons maintenant une approximation øh de
ø dans un sous espace Uh de dimension finie nr
dont on désignera par îi (i = 1, ...,
nr) les fonctions
de bases, choisies affines par morceaux valant 1 au noeud
numéro i en r et 0 aux autres.
Tout élément Ph de Uh est donc de
la forme
Où les Pri sont des coefficients (projections de
vh sur les fonctions de bases), On a donc
ù.u
hj+ 1 2 (i.îl, îe)
+ hj+12 ( 1
ó1.r.Vrîl,Vrîe)
- (Qj+1îl, îe) + (Qjîl,
îe)
r
1
-hj+1(Qj+1(D -
hj+1Qj+1)-1Qj+1îl,
îe) +hj+ 1 2
.ó.r2(A.1S1?S0(zj)
+C.1Si(zj))î1(r2)81,lî1(r2)81,e
= 0 (4.23)
Et
Xnr
l=1
|
qj+1(rl)(îl,
îe) -
|
Xnr
w=1
|
qj(rw)(îw, îe) +
|
Xnr
l=1
|
hj+1(Qj+1(D -
hj+1Qj+1)-1
îlqj+1(rl) , îe)
|
B
-hj+ 1 2 .ó.r2
1S1
(zj).î1(r2)81,e
= 0
Rappellons que Qj est un opérateur autoadjoint
et que les fonctions de bases îl sont des fonctions de r,
alors Qjîl est aussi une fonction de r
en effet puisqu'on applique à la fonction de bases îl
l'opérateur Qj.
D'aprés le théorème Noyau de Schwartz il
existe un noyau gj
défini sur zt X zt
tel que
Z
Qjîl(r) =
gj(r,x)îl(x)dx
zt
alors
ZZ
(Qjîl, îe) =
gj(r,x)îl(x)îe(r)dx.dr
supp(îl)supp(îe)
Pour réaliser la condensation de masse nous allons
utiliser la méthode de quadrature pour le calcul les intégrales.
Alors,
kt+1 + kt ke+1 +
ke
(Qjît, îe)'2
.2gj(rt,xe).
Posonsgj t,e = ke+1 + ke
2gj(rt,xe)
(Qjît, îe) = kt+1 + kt
2gj t,e
Nous allons donc considérer comme inconnue de
l'équation (4.23) la
matrice Qj = [gj
t,e].
On obtient donc
T. Q j+1 - T . Q
j + h j+1T . Q j+1(
D - hj+1 Qj+1) _1
Qj+1
-hj+ 1 2
(E+F+Mj) = 0 (4.24)
Et
T.qj+1 - T.qj + hj+1
Qj+1(D - hj+1
Qj+1)_1.T.qj+1
- hj+ 1 2 .V j
= 0 (4.25)
Où V j est un vecteur de nr
éléments, seul le premier élément est non nul
et vaut ó1.r21S1(zj),
tous les autres éléments sont nuls, quant à T,
F,
B
Mj, E, ce sont des matrices de taille nr
* nr :
k1
T = [tt,e], avec tt,e=
kt+1 + kt
2 ät,e,si l = e = 1
t1,1 = 2
E = [et,o], avec et,o=
ù.u ät,o
r
t
F = [ft,e], avec
-
-
1( 1 1
+ )sil = e+1
ó1 rt rt+1
1( 1 1
+ )sil = e-1
ó1 rt rt_1
1
+ )sil = e
rt rt_1
2
1
ft,e = 2
?
? ???????????
????????????
En ce qui concerne la matrice Mj seul le
coefficient mj1,1 est non
nul et vaut ó1.r2
(A.1S1uS0(zj) +
C.1Si(zj)), tous les autres coefficients
sont 1
nuls.
Le calcul de l'opérateur Qj et du
résidu qj se fait suivant les indices croissant de zj,
il faut donc exprimer Qj+1 en fonction de Qj et
qj+1 en fonction de qj.
Donc en revenant aux équations (4.24) et (4.25) on a
Qj+1(D -- hj+1 Q j+1)
-1((D -- h j+1 Q j+1) + h
j+1 Qj+1) = Qj +hj+ 2
.T
-1(E+F+Mj) ?
Qj+1(D--hj+1
Qj+1)-1D = Qj+
hj+ 2 .T -1(E + F
+ Mj)
? ( Q-1
j+1 --
hj+1D-1)-1 =
Qj+ hj+ 2 .T -1(E
+ F + Mj)
? Qj+1 =
(hj+1D-1 + ( Qj+
hj+ 2 .T -1(E + F
+ Mj))-1)-1
(4.26)
Et
(I + hj+1 Qj+1(D --
hj+1 Qj+1)-1)qj+1 =
qj + hj+ 2 .T -1.V j
? qj+1 = (I + hj+1
Qj+1(D -- hj+1 Qj+1)-1)-1(qj
+hj+ 2 .T -1.V j)
(4.27)
Nous venons donc de construire un schéma numérique
sur l'opérateur Qj et le résidu qj, il nous
reste à déterminer Q1 et
q1.
Par identification avec (4.8) écrit pour j = 0 en
utilisant le caractère arbitraire de u1 on
a
q1 = 0et (Q1ø, ?) =
--(1 ø, ?)8
donc
²
01
u D0 Q1=01t)
--D1.1
²
Ici D1 est une matrice diagonale de
taille (mr -- m)*(mr -- m) telle que
[d 1 i,k] = = 1
ó1ri
|
äi,k si (mr -- m) < k, i < mr.
|
D0 est une matrice nulle de taille m * m,
où m est le nombre de noeuds de l'emplacement de
l'électrode de recueil.
Et 01 une matrices nulle de taille m *
(mr - m).
Maintenant que Q1 et q1 sont
déterminés et si le courant injecté est continu alors
l'algorithme de calcul de l'opérateur Q3 et le résidu
q3 est
(ALGzt)
|
? ?????????
?????????
|
de j = 1 à mzt faire
calcul de la matrice M3
Q3+1 =
(h3+1D-1 + ( Q3+
h3+ 1 2 T
-1(F +
M3))-1)-1
calcul du vecteur V3
q3+1 = (I + h3+1
Q3+1(D - h3+1 Q3+1)-1)-1(q3
+h3+ 1 2 .T
-1.V 3)
fin de la boucle sur j
L'algorithme (ALGzt) est programmé
en Scilab avec les mêmes données numériques
utilisées dans le programme Cwmlat - bis. L'opérateur
Q et q sont sauvegardé dans des fichiers
différents et seront intégré dans
Cwmlat - bis.
Une façon de vérifier le calcul de
l'opérateur Q et le résidu q est de les
calculer sur tout le domaine Ù suivant les indices croissants en z3
et en calculant la solution grace à la discrétisation de la
relation Dirichlet-Neumann
suivant les indices décroissants en z3 u3
= u3+1 -
h3+1D-1
(Q3+1u3+1 +
q3+1),
avec unzinf = 0 (condition aux limites de type Dirichlet
sur F). Pour résumer
on a l'algorithme :
(ALGzinf)
?
? ???????????????
????????????????
de j = 1 à nzinf faire calcul
de la matrice M3
Q3+1 =
(h3+1D-1 + ( Q3+
h3+ 1 2 T
-1(F +
M3))-1)-1
calcul du vecteur V3
q3+1 = (I +
h3+1
Q3+1(D - h3+1
Q3+1)-1)-1(q3 + h3+ 1
2 .T -1.V 3)
fin de la boucle sur j
de i = nzinf à 1 faire
u -1 = u - h D-1( Q u
+ q )
fin de la boucle sur i
L'algorithme (ALGzinf) est aussi
programmé en Scilab sous le nom cal - Q - q - zinf avec les
mêmes données numériques utilisées dans le programme
Cwnlat - bis. En comparant cette solution avec celle trouvée
avec le programme Cwnlat - bis on a trouvé une erreur moyenne
quadratique de l'ordre de 2.23 10-5 entre les deux
solutions.
La figure ci-dessous représente les lignes de courant
calculées avec Cwnlatbis et cal - Q - q - zinf. Les
deux autre courbes représentes la distribution spatiale des erreurs
absolues et relative entre les deux solutions.
La distribution spatiale des erreurs absolues et relatives
montrent une variation de l'erreur au voisinage des électrodes
d'injection et de recueil, ceci est dû au fort gradient au voisinage des
électrodes, on remarque aussi que la distribution des erreurs relatives
montre une variation de l'erreur relative au voisinage de la frontière
S2, ceci peut étre dû à une mal prise en
compte de la condition aux limites sur cette frontière, dans le
schéma numérique.
Distribution spatiale des erreurs abs Distr bution spatiale
des er eu s ela
-1
-5000
FIG. 4.2 Test de comparaison entre les lignes de courant
calculées par
Cwnlat - bis et les lignes de courant calculées
par cal - Q - q - zinf dans
une formation résistive à 0.01Ù.
4.4 Prise en compte des conditions aux limites
transparentes dans Cwnlat-bis
En réduisant le domaine Ù le problème
devient :
(PV fac)
|
?
? ????????????????????
?????????????????????
|
sur p2= p - p1,
-or(1
o ó1.r.ou
or)-o oz(1
ó1.r.ouoz)+i. ù.u
r.u=0
ou oz
= Q.u+q
1
sur Fzt(frontière transparente),
.
ó1.r
ou
sur S1, - A.u = -B or
ou
sur Si, -C.u=0
or
ou
sur S0, - A.u = 0 or
sur S2UF, u=0
4.4.1 Formulation variationnelle
Soit réduit = {v E
H1(p2); v = 0 sur
S2 U F} l'espace des fonc-
tions admissibles pour ce problème et u, v deux
éléments de réduit , Alors
-( o or( 1
ó1.r.ou or
),v)-( o oz ( 1
ó1.r.ouoz ),v)+(i.
ù.u r.u,v)=0 (4.28)
où (,) représente le produit scalaire dans
L2(p2).
On applique la formule de GREEN, on obtient :
ov or ) +(1
ó1.r.ouoz, ovoz) -
(i.ù.ur.u, v)+( 1
ó1.r.du
(1 .
ó1.r
ou or ,
dn , v)?Ù2=
0
(1 .
ó1.r
dn,v)z=0=0
(4.29)
du
En substituant on obtient :
(1 .
ó1.r
|
ou or ,
|
ov or ) + ( 1
ó1.r.ou oz ,ov
oz ) - (i. ù.ur.u,v)
+ ( 1
ó1.r2.du
dn,v)S1?S0?Si +
|
Alors,
(1 .
ó1.r
?u ?r ,
?v ù./1
+
?r ) +(1
r
ó1.r.?u?z, ?v
?z) + (i u v) +
(Q u v)
( 1
ó1.r2(A.1S1?S0(z)
+ C.1S%(z))u,
v)r=r2 =
--(q, v)z=0
+(B
ó1.r2
.1S1 (z),
v)r=r2 (4.30)
Donc le problème précédent admet pour
formulation variationnelle
(PV)réduit { trouver u E
Uréduit vérifiant
a(u, v -- u) = l(v -- u) pour
tout v E Uréduit
Où a et l sont respectivement la forme
bilinéaire symétrique et la forme linéaire données
par :
?v ù./1
)
?r ) +(1
ó1.r.?u?z, ?v
?z) + (i) +
(Q
r
+( 1
ó1.r2
(A.1S1?S0(z) +
C.1S%(z))u,
v)r=r2
a(u,v) = (1
?u ?r ,
.
ó1.r
l(v) =
--(q,v)z=0 +
(B
ó1.r2.1S1(z),v)S1
4.4.2 Forme discrete
On utilise la même méthode utilisée dans
le chapitre 3, c'est à dire chercher une aproximation uh de
u dans un sous ensemble Uhréduit de
dimension finie mr * mz, avec les mêmes fonctions de bases.
On peut donc décomposer uh et vh dans la
base des ?i,j(r, z) :
u(r, z) = Xnr
Xnz ui,j?i,j(r, z), v(r, z)
= Xnr Xnz vp,q?p,q(r,
z) i=1 j=1 p=1
q=1
Le problème discret associé à (PV)
réduitest donc
(PV)hréduit { trouver uh E
Uhréduit vérifiant
a(uh, vh - uh) = l(vh - uh)
pour toutv E Uhréduit
En injectant dans (PVhréduit)
les expressions de uh et de vh, il vient, en notant
wi,j= vi,j- ui,j, noeuds par ligne :
? ??
??
|
trouver u1,1, u1,2, ..,
u2,1, u2,2, .., unr*nz solution
de
Xnr Xnz wi,j
Xnr Xnz a(?i,j,
?p,q)up,q = Xnr Xnz
wi,jl(?i,j) v ù E <nr * <nz i=1
j=1 p=1 q=1 i=1 j=1
|
Les fonctions de bases étant connues, les valeurs
prises par a et l sur celles-ci sont des coefficients connus.
Notons-les par Ai,j,p,q réduitet par F p,q
réduitrespectivement :
Ai,j,p,q réduit = a(?i,j, ?p,q) F
p,qréduit = l(?p,q) i,p E (1, ..., mr) j, q E
(1, ..., mz)
Le problème approché consiste donc à trouver
mr * mz nombres u solutions des équations
linéaires suivante :
Xnr Xnz wi,j
Xnr Xnz Ai,j,p,q
réduitup,q =
Xnr Xnz wi,jFi,jréduit
v ù E <nr * r:nz i=1
j=1 p=1 q=1 i=1 j=1
Les équations linéaire précédente
s'écrivent sous forme condensée trouver le vecteur
Uréduit E<nr*nztel
que AréduitUréduit =
Fréduit
Où la matrise de rigidité (stiffness matrix)
Aréduit , est de taille
réduit
,j,p,q
(mr *mz) * (mr * mz) (symétrique
définie positive) de coefficients Ai
, et le second membre Fréduit est un
vecteur de composante fp,qréduit.
4.4.3 Discrétisation Q1
Avec le même maillage utilisé dans le chapitre 3, la
matrice de rigidité Aréduit est construite de
la facon suivante :
Ai,j,p,qréduit = ( 1
w.u
.Vr?i,j, Vr?p,q) + ( 1 .Vz?i,j, Vz?p,q) +
(i..?i,j, ?p,q) +
ó1.ró1.rr
(Q.?i,1, ?p,1)z=0 + ( 1
ó1.r2(A.
1S1?S0 (zj)
+C.1Si(zj))?i,j,?p,q)r=r2
= ( 1 .Vr?r i ,
Vr?rp)(?z j,?z
q) + (Vz?z j,
Vz?z q)( 1
ó1.r.?r i , ?r
p) +
ó1.r
(1
ó1.r2
(A.1S1?S0(zj) +
C.1Si(zj))?r
1?z j, ?r 1?z
q)r=r2 +
w.u
(Q.?r i ?z 1, ?r
p?z1)z=0 + (i. .?i,j, ?p,q)
r
Pour w = 0 on a
Ai,j,p,qréduit = ( 1
.Vr?r i , Vr?r p)(?z
j, ?z q) + (Vz?z
j, Vz?z q)( 1
ó1.r.?r i , ?r
p)
ó1.r
+(1
ó1.r2
(A.1S1?S0(zj) +
C.1Si(zj))?z
1?z j, ?r 1?z
q)r=r2 +
(Q.?r i , ?z 1,
?p?z1)z=0
Avec les mêmes notations utilisées dans la chapitre
3 la matrice Aréduit s'écrit par exemple avec mr
= 5 et mz = 3 :
1
(F+M1). k 2 +
k.D+Q*k
- 1 k.D 0
k.(F+M2)+ 2
k.D
- 1 k.D
- 1 k.D
k.(F+M3)+ 2
k.D
0
- 1 k.D
?
? ? ? ? ? ? ? ? ? ? ? ? ?
?
? ? ? ? ? ? ? ? ? ? ? ? ?
Aréduit =
Les 0 représentent des matrices nulles de taille r *
r.
Pour le second membre, si on note V j le
vecteur a r éléments seul le
premier coefficient est non nul et vaut
ó1.r2 .1S1
(zj), Alors Fréduit s'écrit
:
B
? ? ?
? V 1.k 2+
k.q
Fréduit = ? ?
V2.k
V3.k
Le systéme linéaire obtenu est plus rapide à
résoudre car on a diminué le nombre d'inconnues.
4.5 Tests
Pour mettre en évidence le gain en temps de calcul, on
a calculé à plusieurs profondeurs des conditions aux limites
transparentes qu'on a intégrées dans Cwmuat - bis et on
a recupéré à chaque fois le temps que celui ci
1
met à s'exécuter. Ce temps de calcul n'est pas
vraiment objectif car je n'ai pas tenu compte du caractére creux de la
matrice de rigidité dans les deux
160
codes, Cwmuat - bis et Cwmuat - bis - avec - comd -
uim - tram.
La courbe 4.3 représente la variation du temps de calcul
de Cwmuat-bis
4
prenant en compte les conditions aux limites transparentes
calculées à différentes profondeurs en fonction des
profondeurs,
0
FIG. 4.3 Variation du temps de calcul de Cwmuat - bis
avec des conditions aux limites transparentes calculées à
plusieurs profondeurs en fonction de l'emplacement des profondeurs
Attention! Le temps t pour z = 0
représente le temps de calcul de Cwnlat - bis sans conditions
aux limites transparentes.
On remarque que pour une profondeur proche de la surface (un
pas en dessous de la surface) Cwnlat - bis avec conditions aux limites
transparentes calculées à cette profondeur, met plus de temps
à résoudre le problème que Cwnlat - bis sans
conditions aux limites tranparentes. Ceci vient du fait que le nombre
d'inconnues diminues peu, donc la dimension de la matrice de rigidité
Aréduit reste importante de plus
Aréduit est localement pleine car Q est une
matrice pleine . Par contre plus on s'éloigne de la surface plus le
nombre d'inconnue diminues donc la dimension de la matrice de rigidité
Aréduit devient plus petite ce qui rend
l'éxécution du code Cwnlat - bis avec conditions aux
limites tranparentes plus rapide d'ou la décroissance en temps en
fonction de la profondeur sur la courbe.
On s'est intéressé aussi aux erreurs entre les
deux codes Cwnlat - bis avec conditions aux limites transparentes et
Cwnlat - bis sans conditions aux limites transparentes.
Soit U(r, z) la solution calculée par
Cwnlat - bis sur tout le domaine I
et U(r, z)zt la solution
calculée par Cwnlat - bis - avec - cond - aux - limi - trans
sur le domaine réduit I - I1. Pour zt =
{-100m, -500m, -1000m -
1500m, -2000m,
-2500m}, nous avons calculé l'erreur moyenne
quadratique:
U(r,z) - U(r,z)zt
2à chaque emplacement de la frontière transparente
zt.
U(r,z) 2
La figure 4.4 représente la variation de l'erreur moyenne
quadratique en fonc-
tion de zt
FIG. 4.4 - Variation de l'erreur moyenne quadratique en fonction
de zt. Pour savoir comment et où varie l'erreur, nous affichons
pour chaque
zt les solutions U(r, z) et U(r,
z)réduit pour zinf < z < zt,
la distribution
spatiale des erreurs absolues ainsi que la distribution spatiale
des erreurs relatives, sur les figures 4.5, 4.6, 4.7, 4.8, 4.9, 4.10, 4.11.
Ces tests montrent que l'erreur est systématique et
varie au voisinage de l'emplacement de la frontière transparente. Ceci
est dû au fait que dans le calcul de la condition aux limites
transparente, la condition aux limites sur Fr est approchée par
un schéma décentré à gauche, alors que dans le
programme Cwnlat - bis, cette condition aux limites est
approchée par un schéma centré.
Distribution spatiale des erreurs abs es Distribution
spatiale des erreu s relat ves
FIG. 4.5 - La condition aux limites transparente est
calculée à une profondeur de -100m. Les lignes de courant
calculées par Cwnlat--bis sont affichées pour
z <= --100m
Distribution spatiale des erreurs abs es D stribution
spatiale des erreu s rela ves
12
FIG. 4.6 - La condition aux limites transparente est
calculée à une profondeur de -50Dm. Les lignes de courant
calculées par Cwnlat--bis sont affichées pour
z <= --500m
Distribution spatiale des erreurs abs es x 10-3
Distribution spatiale des erreu s relat ves
FIG. 4.7 - La condition aux limites transparente est
calculée à une profondeur de -1000m. Les lignes de courant
calculées par Cwnlat -- bis sont affichées
pour z <= --1000m
Méthode de Factorisation
4000
Distribution spatiale des erreurs abs es
x 10
Distribution spatiale des erreurs abs es x 10-3
Distribution spatiale des erreu s relat ves
FIG. 4.9 - La condition aux limites transparente est
calculée à une profondeur de -2000m. Les lignes de courant
calculées par Cwnlat -- bis sont affichées
pour z <= --2000m
Distribution spatiale des erreurs abs es D stribution
spatiale des erreu s rela ves
x 10
FIG. 4.10 - La condition aux limites transparente est
calculée à une profondeur de -2500m. Les lignes de courant
calculées par Cwnlat -- bis sont
affichées pour z <=
--2500m
-3
10
Distribution spatiale des erreurs abs es
Distribution spatiale des erreu s relat ves
FIG. 4.11 - La condition aux limites transparente est
calculée à une profondeur de -290Dm. Les lignes de courant
calculées par Cwnlat -- bis sont
affichées pour z <=
--2900m
Chapitre 5
Conclusion
Cette étude a montré la faisabilité de la
mise en oeuvre de la méthode de factorisation à la simulation du
procédé CHFR. L'objectif est de diminuer le temps de calcul du
code Cwmuat-bis qui se base sur la même formulation
variationnelle et discrétisation éléments finis que celle
sur laquelle se base le code Cwmuat, en réduisant le domaine de
la simulation. On a ici étudié le cas d'un domaine limité
supérieurement par une surface parallèle au sol. La
méthode de factorisation permet de calculer l'opérateur
DirichletNeumann sur cette surface qui fournit la condition aux limites
transparentes. On a calculé cet opérateur théoriquement et
on a présenté un schéma numérique pour obtenir une
approximation numérique.
Cependant il est nécessaire de regarder de prêt
certaint points :
· Tenir compte du caractère creux de la matrice de
rigidité dans
les programmes Cwmuat - bis et Cwmuat - bis -
avec - comd - uim - trams,
pour une meilleur estimation du temps de calcul.
· Améliorer le schéma numérique
construit dans cette étude pour diminuer l'erreur sur le
frontière transparente.
Conclusion
60
Bibliographie
[1 Ali EL Moussawi. Rapport de stage : Simulation du CHFR
Amélioration du programme CWNLAT et proposition d'une nouvelle approche
par perturbation. Septembre 2000.
[21 Jacques Henry and Angel M.Ramos. Article : Factorization
of Second Order Elliptic Boundary Value Problems by Dynamic Programming.
December 2001.
[31 John Lovell. Book : Modelling Frequency Effects on
Laterologs.July 1990.
BIBLIOGRAPHIE
Remerciement s
Je tiens à remercier toutes les personnes qui ont
contribué au bon déroulement de mon stage, Martin G.Luling, et
Isabelle Dubourg de la société Schlumberger RPC, ainsi que toute
l'équipe du projet Ondes de l'INRIA de Rocquencourt Jacques Henry,
Houssem Haddar, Michel Kern, Gilles Scarella, Abdelaaziz Ezziani...
GNU Free Documentation License Ve sion 12November 2002
Copyright (C) 2000,2001,2002 Free Software Foundation, Inc.
59 Temple Place, Suite 330, Boston, MA 02111-1307 USA Everyone is
permitted to copy and distribute verbatim copies of this license document, but
changing it is not allowed.
0. PREAMBLE
The purpose of this License is to make a manual, textbook, or
other functional and useful document "free" in the sense of freedom: to assure
everyone the effective freedom to copy and redistribute it, with or without
modifying it, either commercially or noncommercially. Secondarily, this License
preserves for the author and publisher a way to get credit for their work,
while not being considered responsible for modifications made by others.
This License is a kind of "copyleft", which means that
derivative
works of the document must themselves be free in the same sense.
It complements the GNU General Public License, which is a copyleft license
designed for free software.
We have designed this License in order to use it for manuals for
free software, because free software needs free documentation: a free program
should come with manuals providing the same freedoms that the software does.
But this License is not limited to software manuals; it can be used for any
textual work, regardless of subject matter or whether it is published as a
printed book. We recommend this License principally for works whose purpose is
instruction or reference.
1. APPLICABILITY AND DEFINITIONS
This License applies to any manual or other work, in any medium,
that contains a notice placed by the copyright holder saying it can be
distributed under the terms of this License. Such a notice grants a world-wide,
royalty-free license, unlimited in duration, to use that work under the
conditions stated herein. The "Document", below, refers to any such manual or
work. Any member of the public is a licensee, and is addressed as "you". You
accept the license if you
copy, modify or distribute the work in a way requiring permission
under copyright law.
A "Modified Version" of the Document means any work containing
the Document or a portion of it, either copied verbatim, or with
modifications and/or translated into another language.
A "Secondary Section" is a named appendix or a front-matter
section of the Document that deals exclusively with the relationship of the
publishers or authors of the Document to the Document's overall subject (or to
related matters) and contains nothing that could fall directly within that
overall subject. (Thus, if the Document is in part a textbook of mathematics, a
Secondary Section may not explain any mathematics.) The relationship could be a
matter of historical connection with the subject or with related matters, or of
legal, commercial, philosophical, ethical or political position regarding
them.
The "Invariant Sections" are certain Secondary Sections whose
titles are designated, as being those of Invariant Sections, in the notice that
says that the Document is released under this License. If a section does not
fit the above definition of Secondary then it is not allowed to be designated
as Invariant. The Document may contain zero Invariant Sections. If the Document
does not identify any Invariant Sections then there are none.
The "Cover Texts" are certain short passages of text that are
listed,
as Front-Cover Texts or Back-Cover Texts, in the notice that
says that the Document is released under this License. A Front-Cover Text may
be at most 5 words, and a Back-Cover Text may be at most 25 words.
A "Transparent" copy of the Document means a machine-readable
copy, represented in a format whose specification is available to the general
public, that is suitable for revising the document
straightforwardly with generic text editors or (for images
composed of pixels) generic paint programs or (for drawings) some widely
available drawing editor, and that is suitable for input to text formatters or
for automatic translation to a variety of formats suitable for input to text
formatters. A copy made in an otherwise Transparent file
format whose markup, or absence of markup, has been arranged to
thwart or discourage subsequent modification by readers is not Transparent. An
image format is not Transparent if used for any substantial amount of text. A
copy that is not "Transparent" is called "Opaque".
Examples of suitable formats for Transparent copies include plain
ASCII without markup, Texinfo input format, LaTeX input format, SGML or XML
using a publicly available DTD, and standard-conforming simple HTML, PostScript
or PDF designed for human modification. Examples of transparent image formats
include PNG, XCF and JPG. Opaque formats include proprietary formats that can
be read and edited only by proprietary word processors, SGML or XML for which
the DTD and/or processing tools are not generally available, and the
machine-generated HTML, PostScript or PDF produced by some word
processors for output purposes only.
The "Title Page" means, for a printed book, the title page
itself,
plus such following pages as are needed to hold, legibly, the
material this License requires to appear in the title page. For works in
formats which do not have any title page as such, "Title Page"
means the text near the most prominent appearance of the work's title,
preceding the beginning of the body of the text.
A section "Entitled XYZ" means a named subunit of the Document
whose title either is precisely XYZ or contains XYZ in parentheses following
text that translates XYZ in another language. (Here XYZ stands for a specific
section name mentioned below, such as "Acknowledgements", "Dedications",
"Endorsements", or "History".) To "Preserve the Title"
of such a section when you modify the Document means that it
remains a section "Entitled XYZ" according to this definition.
The Document may include Warranty Disclaimers next to the notice
which states that this License applies to the Document. These Warranty
Disclaimers are considered to be included by reference in this License, but
only as regards disclaiming warranties: any other implication that these
Warranty Disclaimers may have is void and has no effect on the meaning of this
License.
2. VERBATIM COPYING
You may copy and distribute the Document in any medium, either
commercially or noncommercially, provided that this License, the copyright
notices, and the license notice saying this License applies
to the Document are reproduced in all copies, and that you add no
other conditions whatsoever to those of this License. You may not use technical
measures to obstruct or control the reading or further
copying of the copies you make or distribute. However, you may
accept compensation in exchange for copies. If you distribute a large enough
number of copies you must also follow the conditions in section 3.
You may also lend copies, under the same conditions stated above,
and you may publicly display copies.
3. COPYING IN QUANTITY
If you publish printed copies (or copies in media that commonly
have printed covers) of the Document, numbering more than 100, and the
Document's license notice requires Cover Texts, you must enclose the copies in
covers that carry, clearly and legibly, all these Cover Texts: Front-Cover
Texts on the front cover, and Back-Cover Texts on the back cover. Both covers
must also clearly and legibly identify you as the publisher of these copies.
The front cover must present the full title with all words of the title equally
prominent and
visible. You may add other material on the covers in addition.
Copying with changes limited to the covers, as long as they preserve the title
of the Document and satisfy these conditions, can be treated as verbatim
copying in other respects.
If the required texts for either cover are too voluminous to fit
legibly, you should put the first ones listed (as many as fit reasonably) on
the actual cover, and continue the rest onto adjacent pages.
If you publish or distribute Opaque copies of the Document
numbering more than 100, you must either include a machine-readable Transparent
copy along with each Opaque copy, or state in or with each Opaque copy a
computer-network location from which the general network-using public has
access to download using public-standard network protocols a complete
Transparent copy of the Document, free of added material. If you use the latter
option, you must take reasonably prudent steps, when you begin distribution of
Opaque copies in quantity, to ensure that this Transparent copy will remain
thus accessible at the stated location until at least one year after the last
time you distribute an Opaque copy (directly or through your agents or
retailers) of that edition to the public.
It is requested, but not required, that you contact the authors
of the Document well before redistributing any large number of copies, to
give
them a chance to provide you with an updated version of the
Document.
4. MODIFICATIONS
You may copy and distribute a Modified Version of the Document
under the conditions of sections 2 and 3 above, provided that you release the
Modified Version under precisely this License, with the Modified Version
filling the role of the Document, thus licensing distribution and modification
of the Modified Version to whoever possesses a copy of it. In addition, you
must do these things in the Modified Version:
A. Use in the Title Page (and on the covers, if any) a title
distinct from that of the Document, and from those of previous versions (which
should, if there were any, be listed in the History section of the Document).
You may use the same title as a previous version if the original publisher of
that version gives permission.
B. List on the Title Page, as authors, one or more persons or
entities responsible for authorship of the modifications in the Modified
Version, together with at least five of the principal authors of the Document
(all of its principal authors, if it has fewer than five), unless they release
you from this requirement.
C. State on the Title page the name of the publisher of the
Modified Version, as the publisher.
D. Preserve all the copyright notices of the Document.
E. Add an appropriate copyright notice for your modifications
adjacent to the other copyright notices.
F. Include, immediately after the copyright notices, a license
notice giving the public permission to use the Modified Version under the terms
of this License, in the form shown in the Addendum below.
G. Preserve in that license notice the full lists of Invariant
Sections and required Cover Texts given in the Document's license notice.
H. Include an unaltered copy of this License.
I. Preserve the section Entitled "History", Preserve its Title,
and add to it an item stating at least the title, year, new authors, and
publisher of the Modified Version as given on the Title Page. If there is no
section Entitled "History" in the Document, create one stating the title, year,
authors, and publisher of the Document as given on its Title Page, then add an
item describing the Modified Version as stated in the previous sentence.
J. Preserve the network location, if any, given in the Document
for public access to a Transparent copy of the Document, and likewise the
network locations given in the Document for previous versions
it was based on. These may be placed in the "History" section.
You may omit a network location for a work that was published at least four
years before the Document itself, or if the original publisher of the version
it refers to gives permission.
K. For any section Entitled "Acknowledgements" or "Dedications",
Preserve the Title of the section, and preserve in the section all the
substance and tone of each of the contributor acknowledgements and/or
dedications given therein.
L. Preserve all the Invariant Sections of the Document,
unaltered in their text and in their titles. Section numbers or the equivalent
are not considered part of the section titles.
M. Delete any section Entitled "Endorsements". Such a section
may not be included in the Modified Version.
N. Do not retitle any existing section to be Entitled
"Endorsements" or to conflict in title with any Invariant Section.
O. Preserve any Warranty Disclaimers.
If the Modified Version includes new front-matter sections or
appendices that qualify as Secondary Sections and contain no material copied
from the Document, you may at your option designate some or all of these
sections as invariant. To do this, add their titles to the list of Invariant
Sections in the Modified Version's license notice. These titles must be
distinct from any other section titles.
You may add a section Entitled "Endorsements", provided it
contains nothing but endorsements of your Modified Version by various
parties--for example, statements of peer review or that the text has been
approved by an organization as the authoritative definition of a standard.
You may add a passage of up to five words as a Front-Cover Text,
and a passage of up to 25 words as a Back-Cover Text, to the end of the list of
Cover Texts in the Modified Version. Only one passage of Front-Cover Text and
one of Back-Cover Text may be added by (or through arrangements made by) any
one entity. If the Document already includes a cover text for the same cover,
previously added by you or by arrangement made by the same entity you are
acting on behalf of, you may not add another; but you may replace the old one,
on explicit permission from the previous publisher that added the old one.
The author(s) and publisher(s) of the Document do not by this
License give permission to use their names for publicity for or to assert or
imply endorsement of any Modified Version.
5. COMBINING DOCUMENTS
You may combine the Document with other documents released under
this License, under the terms defined in section 4 above for modified versions,
provided that you include in the combination all of the Invariant Sections of
all of the original documents, unmodified, and list them all as Invariant
Sections of your combined work in its
license notice, and that you preserve all their Warranty
Disclaimers.
The combined work need only contain one copy of this License,
and multiple identical Invariant Sections may be replaced with a single copy.
If there are multiple Invariant Sections with the same name but different
contents, make the title of each such section unique by adding at the end of
it, in parentheses, the name of the original author or publisher of that
section if known, or else a unique number. Make the same adjustment to the
section titles in the list of Invariant Sections in the license notice of the
combined work.
In the combination, you must combine any sections Entitled
"History" in the various original documents, forming one section Entitled
"History"; likewise combine any sections Entitled "Acknowledgements", and any
sections Entitled "Dedications". You must delete all sections Entitled
"Endorsements".
6. COLLECTIONS OF DOCUMENTS
You may make a collection consisting of the Document and other
documents released under this License, and replace the individual copies of
this License in the various documents with a single copy that is included in
the collection, provided that you follow the rules of this License for verbatim
copying of each of the documents in all other respects.
You may extract a single document from such a collection, and
distribute it individually under this License, provided you insert a copy of
this License into the extracted document, and follow this License in all other
respects regarding verbatim copying of that document.
A compilation of the Document or its derivatives with other
separate and independent documents or works, in or on a volume of a storage or
distribution medium, is called an "aggregate" if the copyright resulting from
the compilation is not used to limit the legal rights of the compilation's
users beyond what the individual works permit. When the Document is included in
an aggregate, this License does not apply to the other works in the aggregate
which are not themselves derivative works of the Document.
If the Cover Text requirement of section 3 is applicable to these
copies of the Document, then if the Document is less than one half of the
entire aggregate, the Document's Cover Texts may be placed on covers that
bracket the Document within the aggregate, or the electronic equivalent of
covers if the Document is in electronic form. Otherwise they must appear on
printed covers that bracket the whole aggregate.
8. TRANSLATION
Translation is considered a kind of modification, so you may
distribute translations of the Document under the terms of section 4. Replacing
Invariant Sections with translations requires special permission from their
copyright holders, but you may include translations of some or all Invariant
Sections in addition to the original versions of these Invariant Sections. You
may include a translation of this License, and all the license notices in
the
Document, and any Warranty Disclaimers, provided that you also
include the original English version of this License and the original versions
of those notices and disclaimers. In case of a disagreement between the
translation and the original version of this License or a notice or disclaimer,
the original version will prevail.
If a section in the Document is Entitled "Acknowledgements",
"Dedications", or "History", the requirement (section 4) to Preserve its Title
(section 1) will typically require changing the actual title.
9. TERMINATION
You may not copy, modify, sublicense, or distribute the Document
except as expressly provided for under this License. Any other attempt to
copy, modify, sublicense or distribute the Document is void, and
will automatically terminate your rights under this License. However, parties
who have received copies, or rights, from you under this License will not have
their licenses terminated so long as such parties remain in full compliance.
10. FUTURE REVISIONS OF THIS LICENSE
The Free Software Foundation may publish new, revised versions
of the GNU Free Documentation License from time to time. Such new
versions will be similar in spirit to the present version, but may differ in
detail to address new problems or concerns. See
http://www.gnu.org/copyleft/.
Each version of the License is given a distinguishing version
number. If the Document specifies that a particular numbered version of this
License "or any later version" applies to it, you have the option of following
the terms and conditions either of that specified version or of any later
version that has been published (not as a draft) by the Free Software
Foundation. If the Document does not specify a version number of this License,
you may choose any version ever published (not as a draft) by the Free Software
Foundation.
|