REPUBLIQUE ALGERIENNE DEMOCRATIQUE ET POPULAIRE

MINISTERE DE L’ENSEIGNEMENT SUPERIEUR ET DE LA RECHERCHE SCIENTIFIQUE

UNIVERSITE DES SCIENCES ET DE LA KASDI MERBAH-OUARGLA

FACULTE DES SCIENCES APPLIQUEES
Département de Génie Civil et Hydraulique

Polycopié
Modélisation hydraulique

Présenté par :
Dr : Kateb Samir

Ce cours est destiné à l’étudiante 1ère année doctorat hydraulique

Domaine Science et Technologie
Année Universitaire 2020
 Modélisation hydraulique

1.1. Méthodes de résolution des systèmes d'équations algébriques.

L’objet de l’analyse numérique est de concevoir et d’étudier des méthodes de résolution de
certains problèmes mathématiques, en général issus de la modélisation de problèmes “réels",
et dont on cherche à calculer la solution à l’aide d’un ordinateur

1.2. Méthodes itératives.

1.2.1. méthodes des substitutions successives.
Une équation de type F(x) peut toujours s'écrire sous la forme équivalent f (x) où f est une
nouvelle fonction de x.
x (0) la valeur initiale.
x (1) = f (x (0) )
x (2) = f (x (1) )
: :
: :
x (i) = f (x (i-1) )
: :
: :
x (n) = f (x (n-1) )

Le problème est ce que in converge vers la solution ?

 y y
Y2(x)=x Y2(x)=x y
Y1=f(x)
Y2(x)=x

Y1=f(x)
Y1=f(x)
x x x

On diverge On converge
Sans racine
Pour que la méthode des substitutions successives converge il faut que │f'' (x) │ soit
inférieur à 1. Cette condition n'est pas nécessaire mais elle est suffisante. Plus f'(x) est
petite plus la convergence est rapide.
1.2.2. méthode du promoteur de convergence de WEGSTEIN.
La méthode du promoteur de convergence de Wegstein consiste à modifier la méthode des
substitutions successives de façon a accéléré ou forcer systématiquement sa convergence
dans la méthode précédente.
x (n+1) =x (n) + ∆x
∆x = f(x (n)) – x (n)

Wegstein suggère de modifier l'équation comme suit :

x (n+1) =x (n) +α ∆x
1- x (n) = f(x (n-1))
f x n   x n 
2-   n 
x  x n 1
1
3-  
1 
 
4- x n1  x n   f x n   x n  
n1 n 
5- s'arrête si x x   , ε est la précision.
 y

Y2(x)=x

Y1=f(x)
x
x(2) x(1) x(0)

1. Application 1.
On veut résoudre l'équation
X2 – 100 X + 1 = 0 (1)
A l'aide de la méthode des substitutions successives l'équation (1) peut s'écrire :

X
1
100
 X 2 1

On obtient la formule de récurrence

X k 1 
1
100
X k2 1
On propose X0 =1 (valeur initial)
On calcul

X1 
1
100
1 1 0.02
2

X2 
1
100
0.02 1 0.04
2

X3 
1
100
0.04 1 0.01001
2

X4 
1
100
0.01001 1 0.01001
2

La solution exacte X = 0.01001

2. Application 2.
On concéder un écoulement turbulent (IR ≥ 2300)
On calcul le coefficient de perte de charge linéaire λ par

1  D 2.51 
  2 log   
  3.71 Re  
Où :
ε : rugosité relative de la conduite.
D : diamètre.
IRe : nombre de Rynolds.

Pour obtenir une estime de la solution, on pourra utiliser l'une des relations empiriques :
1- Blasius :   0.316 R 0.25 2300 <Re< 105 avec ε = 0
2- Leeds :   0.00714  0.6104 R 0.35 Re< 4.7 105
3- Hermar :   0.0032  0.221 R 0.237 106 <Re< 108
On pose :
1
 X , on aura :

  D 2.51 
X   2 log   X 
 3.71 Re 
Équation implicite de Coolebrok
On appliquant la méthode des substitutions successives
1
X 0  
X0

0  0.0032  0.221 106 

0.237
 0.011
Donc :
1
X 0    9.53
0
X    f  X   4.181939692
1
0

X    f  X   4.183378536
2
1

X    f  X   4.1833785149
3
2

X    f  X   4.1833785149
4
3

Donc la solution est : X  4.183378149

D'où
1
 0.057
X2

1.3. méthode directe.

1.3.1. systèmes à matrice A quelconque méthode de GAUSS.
Théorème de GAUSS : étant donnée une matrice carré A, il existe des matrices
inversible S tell que S A = A' où A' est une matrice triangulaire supérieurs.
Principe : la méthode de GAUSS consiste a transformé le système Ax =b a matrice
quelconque A ont un système équivalent A'x = b', où A' est une matrice triangulaire
supérieurs.

La résolution de ce dernier système étant immédiate avec l'algorithme suivant :

1  n 
xi   bi   aij x j  i = n, n-1, ……1
 
aii  j i 1 
a11 a12 a13 
A   0 a22 a23 
 0 0 a33 

a. Algorithme de GAUSS.
On doit remarque que toute les transformations suppose que les termes akk appelle pivot
sont nommé, il faudra donc ajouter une phase de vérification de la non nullité du pivot akk
et l'algorithme générale devient.

b. Première phase :
Si akkk 1  0 et alkk 1  0 avec l  k

On permute les lignes l et k de la matrice A b

k

c. Deuxième phase :
La matrice A, b
k 1
étant éventuellement modifier par la première phase on opère les

calcule des éléments non nul de A b

k

 k 1
k  k 1 a 
aij  aij   ik  akjk 1
 akk 
 k 1
k   k 1 a 
bi  bi   ik  bkk 1
 akk 
i = k+1, n
j = k+1, n+1
k = 1, n-1
 3. Application:
Lors de traitement d'un problème d'écoulement par la méthode des différences finies on a
obtenu le système suivant :
1 3 3 1   0 
3 2 0     2 
  2  
3 2 6 3  11

Où Φ potentiel hydraulique.

Solution :
n nombre d'inconnues = 3
Donc : k = 1, n – 1 k = 1 ,2

L'étape k =1:
i=k+1,n i = 2, 3 a22 a23
j=k+1,n j = 2, 3 a32 a33
0 
i 2 1 0  a  0  3
 a22  a22   21  a12  2  3   7
 j 2  a11  1
0 
i  2 1 0  a  0  3
 a23  a23   21  a13  0  3   9
 j 3  a11  1
0 
 i 3 1 0  a  0  3
 a32  a32   31  a12  2  3   7
 j 2  a11  1
0 
i 3 1 0  a  0  3
 a33  a33   31  a13  6  3   3
 j 3  a11  1

Pour les valeurs de (b) (i = 2, 3)

0 
1 0  a  2
b2  b2   21  b10   2   0  2
 a11  1
0 
1 0  a  3
b3  b3   31  b10   11  0  11
 a11  1

1 3 3  1   0 
   
A  0  7  9 2   20

0  7  3 3  11
l'étape k =2:
i=3 a33
j=3
1
i 3 a   7
   9  6
2  1 1
 a33  a33   32  a23   3
 j 3  a22   7
1
2  1 a   7
b3  b3   32  b21  11   2  9
 a22   7

1 3 3  1  0
   
A  0  7  9 2   2

0 0 6  3  9

Résolution :
9 3
3 6  9    
6 2

 3  1   3    31
 a b2  a23 3    7  2   9  2    14
1 1
2   b2   a2 j  j  

a22  j 3  22  

 3  1   31  3  30
1 
1
 b1   a1 j  j  
  a
1
b1  a12  2  a13 3   
 0  3    3   
a11  j  2  11 1   14  2  14

1.4. Système à matrice A symétrique définie positive.

1.4.1. Méthode de CHOLOVSKI.
Soit A une matrice symétrique définie positive (les termes de la diagonale sont positive)
A = At et aii > 0  i
On décompose la matrice A comme suit :
A = L D Lt
Avec :
L : matrice triangulaire inférieur unitaire tell que Lii = 1  i
D : matrice diagonale Dii ≠ 0 Dij = 0
Lt : matrice transposé de L donc triangulaire supérieur unitaire.
Les éléments non nul de D et L sont donnée par l'algorithme suivant :
Lii = 1 i=1,n
r 1

d rr  arr   Lrk d kk Lrk 
k 1 
 r 1, n
1  r 1  
L jr  arj  Lrk d kk L jk  j  r 1, n
d rr  k 1  

soit L y =b
Dz=y
Lt x = z
Ces systèmes étant respectivement une matrice triangulaire inférieure, diagonale et
matrice supérieur sont des solutions immédiates :
1- pour la matrice triangulaire inférieure on a.
1  i 1

xi   i  aij x j  i 1, n
b 
aii  j 1 
2- la matrice diagonale :
b
xi  i i 1, n
aii
3- pour la matrice triangulaire supérieur on a :
1  n 
xi   i  aij x j  i  n ,1
b 
aii  j 1  1 
Application :
4 1  x1   6 
A    
1 9  x2  19
Calcule les éléments des matrices L, D, Lt :
Lii = 1, i = 1, n (n nombre d'inconnue) i = 1, 2 donc L11=L22 = 1
r = 1, 2
r=1:
d11 = a11 = 4
j = r +1 , n j=2
1 1 1
L21  a12  1
d11 4 4
r=2:
 21
d 22  a22   L2 k d kk L2 k
k 1

1 1 35
d 22  a22  L21 d11 L21  9   4  
4 4 4
La matrice triangulaire inférieure L est :
 1 0
L 
1 4 1

La matrice diagonale est :

4 0 
D  
0 35 4
La matrice triangulaire supérieur Lt est :
1 1 4
Lt   
0 1 
Déterminer x1 et x2 :
 1 0  y1   6   y1   6 
1 4 1  y   19   y   70 4
   2    2  
4 0   z1   6   z1  6 4
0 35 4  z   70 4   z    2 
   2    2  
1 1 4  x1  6 4  x1  1
 0 1   x    2    x   2 
   2    2  

Contrainte totale contrainte effective et pression interstitielle.

q
hw
h h h

σvT=γ z σvT = γw hw+ γ h σvT = γh +q

 γ Poids volumique du sol (grains + eau).
1. Contrainte verticale totale
Dans un volume de sol la contrainte verticale totale à une profondeur donnée est due au poids
de tous ce qui est ou dessus grains, eau, fondation.

zw
zw
zw

u u u

zw1
zw2
zw3

u1 u2 u3

2. Pression interstitielle u
L'eau dans les pores d'un sol saturé exerce une pression connue par pression interstitielle elle
est conventionnellement représenté par la hauteur d'eau.
Si l'eau dans le sol est stationnaire le niveau d'eau est horizontal cependant si la surface d'eau
n'est pas niveau il y a une infiltration (percolation) et l'eau s'écoule dans les pores des grains
du sol saturé vers le niveau le plus bas.
 Zone sèche
U =0

Partiellement saturé
U =?
Niveau de Saturé (tension) U négative
la nappe
Sol saturé hw U positive

hw

Cette figure schématise la variation de la pression interstitielle dans une région ou voisinage
de la surface libre dans des conditions hydrostatique.
À la surface de la nappe la pression interstitielle est nulle u=0 au dessous de la nappe u est
proportionnelle à la hauteur zw u = γw zw par convention la pression est positive.
Immédiatement au dessus du niveau de la nappe le sol demeure saturé par la phénomène de
capillarité l'eau dans cette zone est on tension ainsi la pression est négative u = - γw zw , la
hauteur de la zone saturé par capillarité dépend des dimension des grains et plus
particulièrement des dimension des pores entre la zone saturé qui contient les grains eau et
gaz (vapeur et l'air) c'est un milieu tri phasique , la pression dans ce domaine partiellement
saturé est encore mal connus.
Contrainte effective σ' il est évident que les instabilités peuvent être posé par le changement
de contrainte totale généré par le chargement d'une fondation ou excavation ou remblaiement
mais peut être n'est pas évident que l'instabilité des terrains peut résulte d'un changement de
pression généré sou changement significative d'une contrainte totale (talus stable peut être
instable après une pluie)
Le concept de contrainte effective à été postulé par Terzagé 1936 elle est défini comme la
déférence entre contrainte totale et la pression interstitielle б ' =б – u .
La repense d'un sol à des modifications des contraintes dépend presque exclusivement du
niveau des contraintes effective à l'intérieur de ce sol.
Le principe de contrainte effective σ ' est probablement le plus important on géotechnique.
 Après soulèvement du
γsat
niveau de la nappe

γw

Un changement important de б' sons changement significative de бv le soulèvement de niveau

de la nappe
U augmente   dim inuer   r dim inuer  M réssistent dim inuer  Fs dim inuer

 r  c    tg  (Lois de coulomb).

Exemple 1.
Soit le récipient de la figure si dessous dont le poids volumique du sol saturé γsat = 20 KN/ m3
Calculer les contrainte totale, effective et la pression interstitielle au niveau A lorsque :
- la nappe phréatique est au niveau A.
- la nappe phréatique est élevée au niveau B.

Niveau B
Zw=2m

H=5m
sol sol

Niveau A

1- saturé 2- déjaugé

Solution :
- 1er cas au niveau A.
Contrainte totale :
σt = γsat H = 20*5 =100 σt = 100 KN/m2

Pression interstitielle :
U=0
Contrainte effective :
σ' = σt –U = 100 – 0 = 100 σ' = 100 Kpa
 - 2ème cas au niveau B.
Contrainte totale :
σt = γsat H + γw Z = 20*5+10*2 =120 σt = 120 KN/m2
Pression interstitielle :
U = γw ( H +Z ) = 10*(5+2) = 70 U= 70 Kpa
Contrainte effective :
σ' = σt –U = 100 – 70 = 50 σ' = 50 Kpa

Dans le cas du rabattement du niveau de la nappe on est dans le cas de tassement et pour le
cas contraire dans le cas de soulèvement du niveau de la nappe on a un risque de glissement.

Exemple 2.
Donner l'expression de σ' au point M en fonction de γ' σ' (M) =f(γ' ) avec γ' = γsat - γw

γ1 Z1

γ2=γsat Z2

.M
Solution
σt = σ' +U
σ' = σt –U
σ' =( γsat Z2 + γ1 Z1) - γw Z2
σ' = γsat Z2 + γ1 Z1 - γw Z2
σ' = (γsat - γw ) Z2 + γ1 Z1
σ' = γ' Z2 + γ1 Z1

Dans le cas générale la formule du contrainte effective en fonction du γ' d'écrit comme suit :
 '   'j Z j   i Z i
γ'j Zj les couches qui se trouvent au dessous de la nappe d'eau déjaugé.
γ1 Z1 les couches qui se trouvent au dessus de la surface de la nappe d'eau.

Exemple 3.
Soit la coupe verticale d'un terrain stratifie (figure si dessous).
Calculer les contraintes totales, effectives et la pression interstitielle.

h1= 2m γ =13.5 KN/ m2

Sable
h2= 2m γsat 1 =18.5 KN/ m2

h3= 4m γsat 2 =20 KN/ m2 Argile

.
A

Solution.
- 1ère méthode :
Contrainte totale :
σt = γ h1 + γsat 1 h2 + γsat 2 h3 = 13.5*2 + 18.5*2 + 20*4 = 144 σt = 144 KN/m2

Pression interstitielle :
U = γw ( h2 +h3 ) = 10*(2+4) = 60 U= 60 Kpa

Contrainte effective :
σ' = σt –U = 144 – 60 = 84 σ' = 84 Kpa

- 2ème méthode :
 '   'j Z j   i Z i

σ' = γ h1 + (γsat 1 - γw ) h2 + (γsat 2 - γw ) h3

σ' = 13.5*2 + (18.5-10) *2 +(20-10) *4
σ' = 84 Kpa
 Chapitre VII
Méthodes des différences finis

2.1. Introduction.
La méthode des différences finies est la plus ancienne technique numérique utilisée pour
résoudre des systèmes d'équations différentielle, bien sur avec des conditions initiales ou
des conditions limite. L'idée principale de cette méthode est d'approximé les dérivés. Toute
dérivé présente dans un système d'équation est remplacer directement par une expression
algébrique écrite en terme de variation intervenant dans le système d'équation on dit "discret
l'espace"

Φ
2.2. Présentation des dérivées par la MDF.
2.2.1. problème à une dimension.
Soit Φ =f(x) une fonction quelconque.
  f x  x   f ( x)
 lim 
x x x x
x  0
     f x  x   f ( x)   f x   f ( x  x) 
    
2
 x   x x  ∆x ∆x
 lim 
x 2 x x
x  0

Les différences f(x+∆x)-f(x) et f(x)-f(x-∆x) sont les différences entre les valeurs de la
fonction f à un intervalle finie ∆x. les différences sont appeler "différence finie".
 f(x+∆x)-f(x) est appeler différence avant au point (i) d'abscisse x.
 f(x)-f(x-∆x) est appeler différence arrière au point (i) d'abscisse x.
 x   x 
 f  x    f  x   est appeler différence centrale au point (i) d'abscisse x.
 2   2 

La différence centrale est généralement la plus précise.

1. formulation par différence arrière.

   i  i 1
  
 x  i h

     
   
     x  i  x  i 1
2
h
 2   h
 x  i h
h
h

  2  i  i 1   i 1  i  2 
 2  
 x i h2 i+2 i+1 i i-1 i-2

  2  i  2i 1  i  2
 2  
 x  i h2

  2    2 
 2   2 
  3   x  i  x  i 1 i  2i 1  i  2    i 1  2i  2  i  3 
 3   
 x  i h h3

  3  i  3i 1  3i  2  i  3
 3  
 x  i h3

  4  i  3i 1  3i  2  i  3  i 1  3i  2  3i  3  i  4 

 4  
 x  i h4
  4  i  4i 1  6i  2  4i  3  i  4
 4  
 x  i h4
2. formulation par différence avant.

   i  1  i
  
 x  i  1 h

     
   
  
2
 x  i 1  x  i
 2  
 x  i 1 h

  2          2i i1
 2   i 1 i 2 i i 1  i 1
 x  i 1 h h2
  2    2 
 2    2 
  3   x  i 1  x  i i 1  2i  i 1  i  2i 1  i  2 
 3   
 x  i 1 h h3

  3    3i  3i 1 i  2

 3   i 1
 x  i 1 h3

  4    3i  3i 1 i  2  i  3i1  3i  2 i  3 

 4   i 1
 x  i 1 h4

  4    4i  6i 1  4i  2  i  3

 4   i 1
 x  i 1 h4

2.2.2. formulation par différence centrale.

La première centrale de f(x) au point xi est défini par :

  
f ( x) f x  2  f x  2

h h 
x h
Il est plus pratique de travailler avec des fonctions à des intervalles complets que des demi
intervalle (h/2). Ainsi il est préférable d'utiliser la moyenne de la première différence centrale.
 f ( x)  f ( x  h)  f ( x  h) i 1  i 1
   
 x  i 2h 2h

 f ( x) 
  
 f ( x) 

i1  i   i  i1  i1  i   i  i1 
  2 f ( x)   x  i 1  x  i 1 h h h
    
 x  i
2
h h h
  2 f ( x)    2i  i 1
   i 1
 x  i
2
h2

  2 f ( x)    2 f ( x) 
    
  3 f ( x)   x  i 1  x  i 1 i 2  2i 1  i  i  2i 1  i 2 
2 2

   
 x  i
3
2h 3 2h 3

  3 f ( x)  i 2  2i 1  2i 1  i 2
  
 x 
3
2h 3

  4 f ( x)  i 1  3i  3i 1  i 2  i  3i 1  3i 2  i 3 

  
 x 
4
h4

  4 f ( x)  i 1  4i  6i 1  4i 2  i 3

  
  x 4
 h4

2.2.3. calcul de l'erreur.

L'expression de Taylor pour f (xi+h) est :
h f ( xi ) h 2  2 f ( xi ) h 3  3 f ( xi )
f ( xi  h)  f ( xi )     .....
1! x 2! x 2 3! x 3
h i h 2 i h 3 i
f ( xi  h)  i 1  i    .....
1! x 2! x 2 3! x 3

h j  j i
i 1  
j 0 j ! x
j

i 2  f ( xi  2h)  

2h j  ji
j 0 j ! x j

i 1  f ( xi  h)  

h j  ji
j 0 j ! x j

i 1  f ( xi  h)  

 h j  j i
j 0 j! x j

i 21  f ( xi  2 h)  

 2h j  ji
j 0 j! x j
On utilisant ces expressions la formulation de la première dérivée par la différence finie
centrale est :
i 1  i 1

1
i1  i1 
2h 2h
 1  i h 2  2i h 3  3i h 4  4i i h 2  2i h 3  3i h 4  4i 

  i  h     .....  i  h     ..... 
2h  x 2 x 2 6 x 3 24 x 4 x 2 x 2 6 x 3 24 x 4 
i  h 2  3i h 4  5i 
     ........
x  6 x 3
12 x 5


La formulation de la première dérivée par la différence finie avant est :

i 1  i

1
i1  i 
h h
i 1  i 1  h 2  2i h 3  3i 
  i  h i    .......  i 
h h x 2 x 2
6 x 3

i  h  3i 
    .......... 
x  2 x 3


i i 1  i 1
  i
x 2h
C'est l'intervalle (h) est choisie l'erreur dans la formulation par différence centrale est plus
petite que celle des deux autres formulation (différence avant et différence arrière) ainsi :
  2  h 2  4i h 4  6i
 2   2 i 1  2i  i 1    2
1
2   
 x i h 12 x 4 360 x 6

Application 1.
Calculer le moment de l'exemple suivant :
M0=0 M4=0
4 y q

x 4
EI
2 y M 2 y 0 1 2 3 4
  M   EI 2
x 2 EI x L
Y0=0 Y4=0
 M 2
 y q4
  EI 4   EI   q0 h= L/4
x 2
x EI

2M
  q0
x 2
1- Calcule les moments.
 2i i 1  2i  i 1

x 2 h2
  2 M i M i 1  2M i  M i 1

x 2 h2
Noeud 0 : M0 = 0 (appuis simple).
Noeud 1 : M2 – 2 M1 + M0 = - q0 h2
Noeud 2 : M3 – 2 M2 + M1 = - q0 h2
Noeud 3 : M4 – 2 M3 + M2 = - q0 h2
Noeud 4 : M4 = 0 (appuis simple)

- 2 M1 + M3 = - q0 h2
M1– 2 M2 + M3 = - q0 h2
M2 – 2 M3 = - q0 h2

La résolution de ce système d'équation nous donne :

 3 q0 L2
M 1 
 32
 q0 L2
 2
M 
 8
M 3  M 1



2- calcule les déplacements.

4 y q 2 y M
 
x 4
EI x 2
EI
Noeud 0 : y0 = 0 (déplacement selon y n'est pas libre).
Noeud 1 : y2 – 2 y1 + y0 = - M1/EI h2
Noeud 2 : – 2 y3 – 2 y2 + y1 = - M2/EI h2
Noeud 3 : y4 – 2 y3 + y2 = - M3/EI h2
Noeud 4 : y4 = 0 (déplacement selon y n'est pas libre)

- 2 y1 + y2 = - M1/EI
y1– 2 y2 + y3 = - M1/EI
y2 – 2 y3 = - M1/EI

Application 2.
Calculer le moment de l'exemple suivant :
 M0=0 M4=0

0 1 2 3 4
L
Y0=0 Y4=0
h= L/4

2M
  q0 (q0 la charge au niveau du noeud)
x 2
Noeud 0 : M0 = 0 (appuis simple).
Noeud 1 : M2 – 2 M1 + M0 = - (q0/4) h2 = - q0 (L2/64)
Noeud 2 : – 2 M2 – 2 M2 + M1 = - (q0/2) h2 = - q0 (L2/32)
Noeud 3 : M4 – 2 M3 + M2 = - q0 (3/4) h2 = - q0 (3/64) L2
Noeud 4 : M4 = 0 (appuis simple)

III- Application de la MDF à l'écoulement bidimensionnel.

Considérons le maillage indiquer sur la figure ci-dessous où i, j les coordonnée ont un point
nodale, utilisons le théorème de Taylor.

hi, j+1

hi-1, j hi,j hi+1, j

h
hi, j-1

1- pour un écoulement bidimensionnel permanent anisotrope.

2 h  2h
kx 2  k y 0
x y

kx
x 2
hi 1, j  2hi , j  hi 1, j 
ky
 y2
hi, j1  2hi, j  hi, j1  0
Soit α = kx / ky , ∆x = ∆y
 hi1, j  2hi , j  hi1, j  hi , j 1  2hi , j  hi , j 1  0

 hi1, j  2  hi , j   hi1, j  hi , j 1  2hi , j  hi , j 1  0

 2 ( 1) hi , j   hi1, j   hi1, j  hi , j 1  hi , j 1

hi , j 
1
2 ( 1)
 hi1, j   hi1, j  hi, j1  hi, j1 

2- pour un milieu isotrope α = 1

hi , j   hi 1, j  hi 1, j  hi , j 1  hi , j 1 
1
4

a- Pour les nœuds appartenant à une limite imperméable horizontal.

Condition de Newman i,j+1
i-1,j i,j i+1,j
h hi , j 1  hi , j 1
 0 
M EF
0
x 2y
hi , j 1  hi , j 1

hi , j 
1
 hi1, j  hi1, j  2 hi, j1 
4 .
Nœud fictif i,j-1

i,j+1
b- Pour les nœuds appartenant à une limite imperméable vertical..
hi , j   2 hi 1, j  hi , j 1  hi , j 1 
1
4
Nœud fictif i-1,j . i,j i+1,j
Condition de Newman
h hi 1, j  hi 1, j
 0  
M EF
0 i,j-1
x 2x
hi 1, j  hi 1, j
i,j+1 x
c- Pour les nœuds appartenant à une coint.
Condition de Newman
y
h hi , j 1  hi , j 1
 0  
M EF
0
x 2y .
Nœud fictif i-1,j
i,j i+1,j

hi , j 1  hi , j 1
Nœud fictif x
h
 0  
M EF
hi 1, j  hi 1, j
0
. i,j-1
x 2x
hi 1, j  hi 1, j

hi , j 
1
 hi1, j  hi, j1 
2
 La pression.

 U i , j   w hi , j  yi , j 
U i, j
hi , j  yi , j 
w
Le gradient hydraulique.
hi 1, j  hi 1, j hi , j 1  hi , j 1
ii , j ( x )  ii , j ( y ) 
2 x 2 y

2.2.4. méthode de Relaxation en différence finie.

Cette méthode permet de déterminer les variables nodales sans construction de la matrice
globale du système d'équation.
Nous présentons aussi l'application de cette méthode pour un écoulement bidimensionnel.
L'écoulement dans un milieu perméable par l'équation de Laplace :
2
φ=0
L'utilisation d'une approximation en différence finie donne l'équation suivante :
1   2  3   4
 0
4
1 2 3 4  40  0
Dont laquelle :
1 , 2 ,3 , 4 représente les valeurs de la fonction φ au point voisine. Comme il est indiqué
sur la figure.

φ4
F0(4)

φ3 φ0 φ1
F0(3) F0(0) F0(1)

φ2
F0(2)

Soit dans la zone concernée, les valeurs de φ dans le nœud déterminé par l'application de
quelque règle simple (valeur grossièrement approché) on peut voir sur l'équation précédente
que si les valeurs de φ0 à φ4 ne sont pas correctes, l'équation peut être écrite sous la forme
suivante :
1 2 3 4  40  F0
 F0 F0 F0 F0 F0
1 , 1 , 1 , 1 ,  4
1  2 3  4  0
Pour un nœud à l'intérieur du domaine, une variation de + 1 de la valeur, φ 0 modifiera la
valeur de la fonction résiduelle 1 , 2 ,3 , 4 de +1 et modifiera la valeur de F0 de φ0 de -4,
ceci peut être utilisé au moyen de la molécule de relaxation qui représente l'effet d'une
variation de +1 de la valeur de φ en chaque nœud sur la valeur de la fonction résiduelle ont
ce point et au 4 point adjacent.
L'utilisation de cette méthode pour la détermination de la distribution de la charge
hydraulique d'écoulement et de la pression interstitielle est comme suit :
1- diviser le domaine d'écoulement en maillage carré (profité de la symétrie si elle
existe).
2- Identifier les conditions aux limites (ligne de courant, ligne équipotentielle).
3- Appliquer la charge connus correspondant et proposé des valeurs raisonnables au
nœud intermédiaire.
4- Appliquer les équations convenables à chaque nœud.
a- à l'intérieur du domaine.
b- Limite imperméable.
c- A un coint.
d- A un nœud où la charge est connue.
On détermine les résidus F0 de chaque nœud.
5- procéder ensuite à la détermination de la valeur F0 jusqu'à obtenir des valeur presque
nul à l'intérieur de tous les domaines.
A la fin la valeur de φ dans chaque nœud représente la solution de ce problème.

Application 1 :

φ5 =10 φ1 φ3 =0
1m
0
2m

10

φ6 =10 φ2 φ 4 =0

Nœud 1 : le nœud 1 appartient à une limite imperméable horizontale.

i1, j  i1, j  2i , j 1  4i , j  0
φi-1,j φi,j φi+1,j
4i , j i 1, j  i 1, j  2i , j 1
φ1
41 5  3  2 2
41 10  0  2 2
10  22  41  0
φi,j-1
Nœud 2 : le nœud 2 appartient à une limite imperméable horizontale
φi,j+1
4i , j i 1, j  i 1, j  2i , j 1
4 2 6   4  21
42 10  0  21
φ2
10  21  4 2  0
φi-1,j φi,j φi+1,j

Donc on a :
 41  2 2  10

21  4  2  10

1- résoudre ce système d'équations par la méthode du déterminant.

4 2
16  4 12
2 4

 10 2 60
 40  20  60  1  5
 10  4 12

 4  10 60
 40  20  60   2  5
2  10 12
2- méthode de relaxation.
On donne des valeurs pour φ1 φ2, (φ1 = 10 φ2 = 10) et on calcule F0 -4

10 + 2(10) – 4 (10) = -10 = F0(1)

Nœud 1
10 + 2(10) – 4 (10) = -10 = F0(2)
2

valeur donnée
2
F0(1)

Nœud 2
Valeur donnée
F0(2) -4

10 (1) 10-2 (1) 8 8 8

-10 -10+ (-2*(- -2 -2+(-3*2) -8
→ 4)) → → →
10 10 10-3 7
10
-10 -10+(-2*2) -14 (2) -14+(-3*(- (2) -2
 4))

(1) 8-2 (1) 6 6 6 (1) 6-1

-8+(-2*(-4)) 0 0+(-1*2) -2 -2+(-1*(-4))
→ → → → →
7 7 7-1 6 6
-2+(-2*2) -6 (2) -6+(-1*(-4)) (2) -2 -2*(-1*2)

(1) 5 5 5
2 2+ (-1*2) 0
→ → →
6 6-1 5
-4 (2) -4+(-1*(-4)) 0

Application 2.

φ9= 10 φ8= 10 φ7= 10

10

0
φ1 φ6=5
φ2

φ3 φ4 φ5=5

Nœuds 1 : appartient à une limite imperméable verticale.

2i 1, j  i , j 1 i , j 1  4i , j  0

4i , j  2i 1, j  i , j 1  i , j 1

41  2 2  3 9

41  2 2  3 10

10  2 2  3  41  0
Nœud 2 : au milieux du domaine.
4i , j  i 1, j  i 1, j i , j 1  i , j 1

i1, j  i1, j i , j 1  i , j 1  4i , j  0

1  6 4  8  42  0
1  5 4  10  42  0
Nœud 3 : appartient au coint.
2i , j i 1, j i , j 1

i1, j i , j 1  2i , j  0

 4 1  23  0
Nœud 4 : limite horizontal
4i , j i 1, j  i 1, j  2i , j 1

i1, j  i1, j  2i , j 1  4i , j  0

3  5  2 2  4 4  0
3  5  22  44  0

 41   2   3  10 (1)

  4    15
 1 2 4 (2)

1  2 3   4  0 (3)
2 2   3  4 4   5 (4)

(3) → 1  23  4 (3)

On remplace (3)’ dans (1), (2) et (4) on auras :
 4 (2 3  4 )   2   3  10

 2 3  4  4 2   4  15
2    4   5
 2 3 4

 2 2  7 3  4 4  10

  4 2  2 3  15
2     4    5
 2 3 4

2 7 4
4 2 0
2 1  4   16  (16)  0  (16)  0  (112)  64
2 7 4
4 2 0
 10  7 4
 15 2 0
480 15
5 1  4  80  (60)  0  (40)  0  (420)  480  2    7.5
64 2
 10  7 4
5 2 0

2  10 4
 4  15 0
480 15
2 5  4 120  80  0  (120)  0  (160)  480  3    7.5
64 2
2  10 4
4 5 0

2  7  10
4 2  15
440 55
2 1  5   20  40  210  (40)  (30)  (140)  440  4  
64 8
2  7  10
4 2  15
15 55 65
1  23  4  2  
2 8 8
1 méthode de relaxation.
Calcule F0 : on donne φ1 = 8, φ2 = 8, φ3 =8, φ4 =8 -4 2 1 -4
noeud 1: F0(1)   4 (8)  2(8)  8  10  2

noeud 2: F2( 2)  8  4 (8)  8  15  1 1 1

noeud 3: F 2
( 3)
 8  2 (8)  8  0
noeud 4: F2( 4)  2 (8)  8  4 (8)  5   3
1 2


-2 1 1 -4

(4) 8 8 (4) 8 8 (2) 8 8-1

2 -1+(-1*2) 2 -3 2+(-1*1) -3+(-1*(-
→ → 4)) →
8 8-1 8 7 8
7
0+(- -3+(-1*(- -1 1 -1 1+(-1*1)
1*1) 4))
(2) 8 7 (3) 8 7 (3) 8 7
1 1 1+(-1/2*1) 1 1/2 1
→ → →
8 7 8-1/2 7 7.5 7
-1 0 -1+(-1/2*(-2)) 0+(-1/2*1) 0 -1/2

(2) 8 7+1/2 (2) 8 7.5 (1) 8+1/4 7.5

1/2+(1/2*1) -1+(1/2*(- 1 -1 1+(1/4*(- -1+(-
4)) → → 4)) 1/4*2)
7.5 7.5 7
7 7.5 7
0 -1/2+(1/2*1) 0 0 0+(1/4*1) 0

(1) 8.25 7.5 (2) 8.25 7.5-1/8 (2) 8.25 7.375

0 -1/2 0+(-1/8*1) -1/2+(-1/8*(-4)) -1/8 0
→ → →
7.5 7 7.5 7 7.5 7
1/4 0 1/4 0+(-1/8*1) 1/4 -1/8

(3) 8.25 7.375 8.25 7375

-1/8+(1/8*1) 0 0 0
→ →
7.5+1/8 7 7.625 7
1/4+(1/8*(-2)) -1/8+(1/8*1) 0 0