A.

Jouaiti

I: Méthodes numériques :
En chimie physique on est souvent confronté à la résolution d’équations sous la forme
f(x) =0 ; la solution ou les solutions de cette équation sont appelées racine ou zéro de la
fonction f(x).

Exemple :

Prenons le cas du polynôme d’ordre 2, f(x)=ax²+bx+c=0 la résolution de cette équation

abouti à deux racines :

−�± �2 −4��
X1,2= 2�

La résolution de cette équation est donc simple, pour les équations d’ordre trois et
quatre on arrive encore à trouver une formule permettant de trouver les racines. Mais,
pour un ordre supérieur à 4, on a pas de formule, donc on doit faire recours aux
méthodes numériques.

I.a: Méthode de la bissection.

La méthode de la bissection repose sur une idée simple, à savoir que la fonction f(x)
change de signe lorsqu’elle passe par un zéro.

Fig 1 : Méthode de la bissection : f(x)=e-x-x

 A. Jouaiti

Algorithme de la méthode de la bissection :

1-Etant donné un intervalle [x1,x2] pour lequel f(x) possède un changement de signe.
2-Etant donné , le critère d’arrêt, et N le nombre maximal d’itérations.
� +�
3-Poser �� = 1 2 2
�2 −�1
4-Si 2��
<�
Convergence atteinte
Ecrire la racine de xm
Ecrire f(xm)
Arrêt
5-Ecrire x1,x2,xm, f(x1), f(x2) et f(xm)
6-Si f(x1)f(xm)<0 alors x2=xm
7-Si f(x2)f(xm)<0 alors x1=xm
8-si le nombre maximal d’itération est N est atteint :
Convergence non atteinte
Arrêt.
9- Retour à l’étape 3

Travail à faire :

Faire un programme avec le langage python, en utilisant l’algorithme précédent, pour une
résolution numérique par la méthode de la bissection.

Définir une fonction en python

def fonction( argument1, argument2,…):

Appel d’une fonction en python

fonction(argument1, argument2,…)

Exemple (addition de deux valeurs)

def addition(x1, x2):

return x1+x2

print(addition(2,4))
 A. Jouaiti

Cet exemple à copier dans l’IDE Thonny

Résultat :

A faire avant de passer à la suite.

On peut faire de même pour la soustraction (x1-x2), la multiplication (x1*x2), la division x1/x2
et la puissance (x1**x2)
 A. Jouaiti

Passons à l’algorithme de la méthode de la bissection.

Notre fonction est : f(x)=e-x-x
On va définir les conditions d’arrêt :
le nombre d’itération N=100 (par exemple)
epsilon ()=0.001
intervalle de travail [a,b]=[0.5,1]

Code python de la bissection :

#ceci est un commentaire

#Méthode de la bissection de la fonction f(x)=exp(-x)-x

import math #on import la bibliothèque math pour pouvoir utilser exp, sin,
cos....

def f(x):
return [Link](-x)-x

#déclaration des différents paramètres

a=0.5
b=1
N=100
epsilon=0.001

#boucle de calcul

for i in range(N): #signifie que i varire de 0 à (N-1) on peut aussi préciser le

début et la fin
xm=0.5*(a+b) # xm c'est le centre de l'intervalle
if f(a)*f(xm)<0:
b=xm
if f(b)*f(xm)<0:
a=xm
if abs(b-a)/(2*xm)<epsilon:
print('la convergence est atteinte pour un nombre d"itérations égal à: ', i)
print('la fonction f(x) est nulle pour x= ', xm)
break
 A. Jouaiti

Résultat:

L’intervalle [a,b] se rétrécit à chaque itération

i=0 [a,b]=[0.5,1]
a = 0.5 b = 1 Xm = 0.75
f(a)= 0.1065 f(b)= -0.6321 f(xm) -0.2776
f(a)*f(xm)<0 donc b=xm

i=1 [a,b]=[0.5,75]
a = 0.5 b = 0.75 Xm = 0.625
f(a)= 0.1065 f(b)= -0.2776 f(xm) -0.0897
f(a)*f(xm)<0 donc b=xm

i=2 [a,b]=[0.5,0.625]
a = 0.5 b = 0.625 Xm = 0.5625
f(a)= 0.1065 f(b)= -0.0897 f(xm) 0.0073
f(b)*f(xm)<0 donc a=xm

i=3 [a,b]=[0.5625,0.625]
a = 0.5625 b = 0.625 Xm = 0.5938
f(a)= 0.0073 f(b)= -0.0897 f(xm) -0.0415
f(a)*f(xm)<0 donc b=xm

i=4 [a,b]=[0.5625,0.5938]
a = 0.5625 b = 0.5938 Xm = 0.5781
f(a)= 0.0073 f(b)= -0.0415 f(xm) -0.0172
f(a)*f(xm)<0 donc b=xm

i=5 [a,b]=[0.5625,0.5781]
a = 0.5625 b = 0.5781 Xm = 0.5703
f(a)= 0.0073 f(b)= -0.0172 f(xm) -0.005
f(a)*f(xm)<0 donc b=xm

i=6 [a,b]=[0.5625,0.5703]
a = 0.5625 b = 0.5703 Xm = 0.5664
f(a)= 0.0073 f(b)= -0.005 f(xm) 0.0012
f(b)*f(xm)<0 donc a=xm

i=7 [a,b]=[0.5664,0.5703]
a = 0.5664 b = 0.5703 Xm = 0.5684
f(a)= 0.0012 f(b)= -0.005 f(xm) -0.0019
f(a)*f(xm)<0 donc b=xm

i=8 [a,b]=[0.5664,0.5684]
a = 0.5664 b = 0.5684 Xm = 0.5674
 A. Jouaiti

f(a)= 0.0012 f(b)= -0.0019 f(xm) -0.0004

f(a)*f(xm)<0 donc b=xm
la convergence est atteinte pour un nombre d"itérations égal à: 8
la fonction f(x) est nulle pour x= 0.5673828125

I.b : Méthodes des Points fixes.

Un point fixe d’une fonction g(x) est une valeur de x qui reste invariante pour cette
fonction, c'est-à-dire toute solution de x=g(x) est un point fixe de la fonction g(x).

Pour déterminer les points fixes il suffit d’effectuer les itérations de la manière suivante
à partir d’une valeur initiale X0:
X0 donné
Xn+1=g(Xn)

Algorithme des points fixes

1. Etant donné , un critère d’arrêt
2. Etant donné N, nombre maximal d’itérations
3. Etant donné X0, valeur initiale estimée du point fixe
4. Effectuer Xn+1=g(Xn)
��+1−��
5. Si ��+1 < �
Convergence atteinte

Ecrire la solution Xn+1

Arrêt
6. Si le nombre maximal d’itérations N est atteint :
Convergence non atteinte
Arrêt.
7. Retour à l’étape 4.

Revenons à la fonction précédente f(x)=e-x-x.

f(x)=e-x-x=0  x = e-x=g(x)
 A. Jouaiti

x0=0.5
N=100,
=0.001
On donne à x la valeur x0 et on calcule g(x0)=e-x0 on obtient une nouvelle valeur de x0 et
ainsi de suite jusqu’à ce que le système converge.

Code python de la méthode du point fixe:

#ceci est un commentaire
#Méthode des points fixes de la fonction f(x)=exp(-x)-x

import math #on import la bibliothèque math pour pouvoir utilser exp, sin, cos....

def g(x):
return [Link](-x)

#déclaration des différents paramètres

x0=0.5
N=100
epsilon=0.001

#boucle de calcul

for i in range(N):
x=g(x0) # f(x)=exp(-x)-x=0 ---> x=g(x)=exp(-x)

if i>N-2:
print('la convergence n"est pas atteinte ')
if abs(x-x0)/x<epsilon:
print('la convergence est atteinte pour un nombre d"itérations égal à: ', i)
print('la fonction f(x) est nulle pour x= ', x)
break
x0=x

Résulat:
la convergence est atteinte pour un nombre d"itérations égal à: 10
la fonction f(x) est nulle pour x= 0.5672771959707785
 A. Jouaiti

Détail du calcul

x0=0.5

i= 0
x= g(x0)=g( 0.5 )= 0.60653
i= 1
x= g(x0)=g( 0.60653 )= 0.54524
i= 2
x= g(x0)=g( 0.54524 )= 0.5797
i= 3
x= g(x0)=g( 0.5797 )= 0.56006
i= 4
x= g(x0)=g( 0.56006 )= 0.57117
i= 5
x= g(x0)=g( 0.57117 )= 0.56486
i= 6
x= g(x0)=g( 0.56486 )= 0.56844
i= 7
x= g(x0)=g( 0.56844 )= 0.56641
i= 8
x= g(x0)=g( 0.56641 )= 0.56756
i= 9
x= g(x0)=g( 0.56756 )= 0.56691
i = 10
x= g(x0)=g( 0.56691 )= 0.56728
la convergence est atteinte pour un nombre d"itérations égal à: 10
la fonction f(x) est nulle pour x= 0.5672771959707785
 A. Jouaiti

Si on change de fonction et on prend

f(x)=x2+x-1

f(x)=x2+x-1=0 x=g(x)=1-x2

Valeurs de x
0.6636
0.5596350400000001
0.6868086220041982
0.5282939167406944
0.7209055375347763
0.4802952059516953
0.7693165151398186
0.40815209953312526
0.8334118636467018
0.30542466553293124
0.9067157736840971
0.17786650575244933
0.9683635061314139
0.062272119992875186
0.996122183071593
 A. Jouaiti

0.007740596392683874
0.9999400831674856
0.00011983007500193654
0.9999999856407531
2.871849358321299e-08
0.9999999999999992
1.5543122344752192e-15
1.0
0.0
Erreur division par zéro (abs(x-x0)/x<epsilon

Dans cet exemple la méthode du point fixe le système diverge.

I.c. Méthode de Newton

La méthode de Newton est largement utilisée pour la résolution d’équations non
linéaires, elle est basée sur le développement de Taylor.
Soit une fonction à résoudre f(x)=0
A partir d’une valeur X0 initiale de la solution, on cherche une correction X telle que :
0=f(X0+X)
En faisant un développement de Taylor autour de X=X0 on trouve.

��2
0 = � �0 + �' �0 �� + �'' �0 +…
2
Il suffit de négliger les termes d’ordres supérieurs à 1 en X pour obtenir
0 f(X0)+ f ‘(X0) X
On peut alors obtenir X=- f(X0)/ f’(X0)
La correction X est la correction à ajouter à X0 pour annuler la fonction f(X).
On pose X1= X0+ X
On recommence le processus en cherchant à corriger X1 d’une nouvelle quantité X, on
obtient l’algorithme suivant.
Algorithme de Newton
1. Etant donné , un critère d’arrêt
2. Etant donné N, nombre maximal d’itérations
3. Etant donné X0, valeur initiale de la solution
�(� )
4. Effectuer ��+1 = �� − �'(�� )
�
��+1−��
5. Si ��+1 < �
Convergence atteinte
Ecrire la solution Xn+1
Arrêt
6. Si le nombre maximal d’itérations N est atteint :
Convergence non atteinte
Arrêt.
7. Retour à l’étape 4.
 A. Jouaiti

Réaliser le programme permettant de calculer le zéro de la fonction

f(x)=x2+x-1 (en utilisant la méthode du point fixe on a trouvé
précédemment que le système diverge)

Code python de la méthode dNewton:

#ceci est un commentaire
#Méthode de Newton f(x)=x**2+x-1

def f(x):
return x**2+x-1

def df(x): # la dérivée de f(x)

return 2*x+1

#déclaration des différents paramètres

x0=0.5
N=100
epsilon=0.001

#boucle de calcul

for i in range(N): #signifie que i varier de 0 à (N-1) on peut aussi préciser le début et la
fin
x=x0-f(x0)/df(x0)
if i>N-2:
print('la convergence n"est pas atteinte ')
if abs((x-x0)/x)<epsilon:
print('la convergence est atteinte pour un nombre d"itérations égal à: ', i)
print('la fonction f(x) est nulle pour x= ', x)
break
x0=x

Résultat

la convergence est atteinte pour un nombre d’itérations égal à: 2

la fonction f(x) est nulle pour x= 0.618033988957902
 A. Jouaiti

Détail du calcul:
i= 0
x0 = 0.5 x= x0-f(x0)/df(x0)= 0.625
i= 1
x0 = 0.625 x= x0-f(x0)/df(x0)= 0.61806
i= 2
x0 = 0.61806 x= x0-f(x0)/df(x0)= 0.61803
la convergence est atteinte pour un nombre d’itérations égal à: 2
la fonction f(x) est nulle pour x= 0.618033988957902
 A. Jouaiti

II- Systèmes linéaires

De façon générale, la résolution d’un système d’équations linéaires consiste à trouver un

vecteur � = �1 , �2 , … solution de l’équation suivante.

�11 �1 + �12 �2 + … + �1� �� = �1

�21 �1 + �22 �2 + … + �2� �� = �2
.
.
.
��1 �1 + ��2 �2 + … + ��� �� = ��

On écrit le système précédent sous la forme :

�� = �

�11 ⋯ �1�
Où A est la matrice: ⋮ ⋱ ⋮
��1 ⋯ ���
�1
Et � = .
��
Les coefficients de la matrice A et du vecteur � sont connus. Il reste à déterminer le vecteur �.
La solution de l’équation
�� = �

Elimination de Gauss.
Pour transformer un système quelconque en système diagonal, il suffit de faire les opérations
suivantes qu’on va appliquer à un exemple.

2�1 + �2 + 2�3 = 10 L1 (ligne 1)

6�1 + 4�2 + 0�3 = 26 L2
8�1 + 5�2 + �3 = 35 L3
Et sous forme matricielle :
2 1 2 �1 10
6 4 0 �2 = 26
8 5 1 �3 35

2 1 2 10
6 4 0 26
8 5 1 35
#Système
A=[[2,1,2],[6,4,0],[8,5,1]]
B=[10,26,35]
 A. Jouaiti

�1 �1
On divise la première ligne L1 et b1 par a11=2 (a11 est le pivot). �11
�� �11
N=3 # matrice carrée de degré 3
Pivot= A[0][0]
B[0]= B[0]/Pivot
for j in range(3):
A[0][j]= A[0][j]/Pivot

On aura
1 1/2 1 5
6 4 0 26
8 5 1 35

On fait la transformation suivante de la matrice :

�2 − �21 �1 �� �2 − �21 �1
�3 − �31 �1 �� �3 − �31 �1
#ligne 2
Coef= A[1][0]
B[1]=B[1]-Coef*B[0]
for k in range(N):
A[1][k]= A[1][k]-Coef*A[0][k]

#ligne 3
Coef= A[2][0]
B[2]=B[2]-Coef*B[0]
for k in range(N):
A[2][k]= A[2][k]-Coef*A[0][k]

On obtient :
1 1/2 1 5
0 1 −6 −4
0 1 −7 −5
Le niveau pivot est a22=1 on divise la ligne L2 par le pivot de a22 (ici pas de changement
� �
puisque a22=1) � 2 �� � 2
22 22
Pivot= A[1][1]
B[1]= B[1]/Pivot
for j in range(3):
A[1][j]= A[1][j]/Pivot

On obtient :
1 1/2 1 5
0 1 −6 −4
0 1 −7 −5

On fait la transformation suivante de la matrice :

�1 − �12 �2 �� �1 − �12 �2
�3 − �32 �2 �� �3 − �32 �2
 A. Jouaiti

#ligne 1
Coef= A[0][1]
B[0]=B[0]-Coef*B[1]
for k in range(N):
A[0][k]= A[0][k]-Coef*A[1][k]
#ligne 3
Coef= A[2][1]
B[2]=B[2]-Coef*B[1]
for k in range(N):
A[2][k]= A[2][k]-Coef*A[1][k]

On obtient :
1 0 4 7
0 1 −6 −4
0 0 −1 −1
�3 �3
Le niveau pivot est a33=-1 on divise la ligne L3 par le pivot de a33. �33
�� �33
Pivot= A[2][2]
B[2]= B[2]/Pivot
for j in range(3):
A[2][j]= A[2][j]/Pivot

On obtient
1 0 4 7
0 1 −6 −4
0 0 1 1
On fait la transformation suivante de la matrice :
�1 − �13 �3 �� �1 − �13 �3
�2 − �23 �3 �� �2 − �23 �3
#ligne 1
Coef= A[0][2]
B[0]=B[0]-Coef*B[2]
for k in range(N):
A[0][k]= A[0][k]-Coef*A[2][k]
#ligne 2
Coef= A[1][2]
B[1]=B[1]-Coef*B[2]
for k in range(N):
A[1][k]= A[1][k]-Coef*A[2][k]

print(B)
On obtient :
1 0 0 3
0 1 0 2
0 0 1 1

3
Donc la solution est � = 2
1
 A. Jouaiti

Ecriture d’un vecteur ou d’une matrice en code python

Vecteur :
3
�= 2
1
x=[3,2,1]
Matrrice:
2 1 2
�= 6 4 0
8 5 1
Exemple:
A=[[2,1,2],[6,4,0],[8,5,1]]
print('afiche la marice ', A )
print('afiche la ligne 1 ', A[0] )
print('afiche la ligne 2 ', A[1] )
print('afiche la ligne 3 ', A[2] )
print('affiche le 3ème élément de la ligne2 :', A[1][2])

# Code Python
#Système
A=[[2,1,2],[6,4,0],[8,5,1]]
B=[10,26,35]
N=3 # matrice carrée de degré 3
Pivot= A[0][0]
B[0]= B[0]/Pivot
for j in range(3):
A[0][j]= A[0][j]/Pivot
#ligne 2
Coef= A[1][0]
B[1]=B[1]-Coef*B[0]
for k in range(N):
A[1][k]= A[1][k]-Coef*A[0][k]

#ligne 3
Coef= A[2][0]
B[2]=B[2]-Coef*B[0]
 A. Jouaiti

for k in range(N):
A[2][k]= A[2][k]-Coef*A[0][k]
Pivot= A[1][1]
B[1]= B[1]/Pivot
for j in range(3):
A[1][j]= A[1][j]/Pivot

Pivot= A[1][1]
B[1]= B[1]/Pivot
for j in range(3):
A[1][j]= A[1][j]/Pivot

#ligne 1
Coef= A[0][1]
B[0]=B[0]-Coef*B[1]
for k in range(N):
A[0][k]= A[0][k]-Coef*A[1][k]
#ligne 3
Coef= A[2][1]
B[2]=B[2]-Coef*B[1]
for k in range(N):
A[2][k]= A[2][k]-Coef*A[1][k]
Pivot= A[2][2]
B[2]= B[2]/Pivot
for j in range(3):
A[2][j]= A[2][j]/Pivot
#ligne 1
Coef= A[0][2]
B[0]=B[0]-Coef*B[2]
for k in range(N):
A[0][k]= A[0][k]-Coef*A[2][k]
#ligne 2
Coef= A[1][2]
B[1]=B[1]-Coef*B[2]
for k in range(N):
A[1][k]= A[1][k]-Coef*A[2][k]

print(B)

# Code Python pour la méthode de gauss

# Pour une valeur quelconque de N
print('+++++++++++++++')
print('Système AX=B ')
#Données
A=[[2,1,2],[6,4,0],[8,5,1]]
B=[10,26,35]
N=3
print('A ',A )
print('B ',B)
 A. Jouaiti

#calcul
for i in range(N):
Pivot=A[i][i]
B[i]=B[i]/Pivot
for j in range (N):
A[i][j]=A[i][j]/Pivot #On divise par le pivot

for k in range(N):
if k!=i : #if k not equal to i
B[k]=B[k]-A[k][i]*B[i]
Coef=A[k][i]
for l in range (N):
A[k][l]=A[k][l]-Coef*A[i][l]
print('-----------Solution---------- ')
print('B ',B )
 A. Jouaiti

III . METHODE DU SIMPLEXE

Position du problème :

Trouver les x1,x2,…,xn qui optimise la fonction objective(ou économique)

�(�1 , �2 , ………, ��) = �1 �1 + �2 �2 ……… + �� ��

et vérifient les contraintes :

�11 �1 + �12 �2 + … + �1� �� < �1

�21 �1 + �22 �2 + … + �2� �� < �2

��1 �1 + ��2 �2 + … + ��� �� < ��

Avec c1,c2,…,cn, amn, b1, b2,…,bm connus.

En notation condensée :

Max Z=CTX sous A.X≤b et x≥0

Avec C (n*1), A(m*n) et b(m*1) connus et X(n*1)

La forme standard est un problème d’optimisation sous la forme

Max Z=CTX sous A.X=b et x≥0

Avec C(nx1), A(mxn) et b(mx1) connus et X(nx1) inconnus

 A. Jouaiti

On peut toujours remplacer les inégalités :

�k1 �1 + �k2 �2 + … + ��� �� < �k 1≤�≤�

Par des égalités en introduisant des variables d’écart.

�k1 �1 + �k2 �2 + … + ��� �� + �� = �k

�k1 �1 + �k2 �2 + … + ��� �� + �� = �k

�k1 �1 + �k2 �2 + … + ��� �� + �� = �k

En forme condensée

��1 �1 + ��2 �2 + … + ��� �� + �� = �� 1≤�≤�

Le problème devient :

Maximiser Z’ : fonction objective ou économique

�' �1 , �2 , ………, �� = �1 �1 + �2 �2 ……… + �� �� + 0�1 + 0�2 + …. + 0��

Sous les contraintes suivantes :

��1 �1 + ��2 �2 + … + ��� �� + �� = �� 1≤�≤�

Ou sous la forme matricielle:

� � �
Max Z’=[CT 0T] sous [A 1] =b et ≥0
� � �

avec C, A et b connus.
 A. Jouaiti

Définitions :

Déf 1 : On appelle solution tout vecteur X tel que A.X=b

Déf 2 : Une solution X qui vérifie A.X≤b X≥0 est une solution admissible.

Déf 3 : A étant une matrice (m*n) de rang m. Une base est une matrice B (m*m) extraite de
A dont le déterminant est non nul.

On remarque que le système à résoudre possède m égalité et n+m inconnues.

Donc la valeur de n variables peut être fixée arbitrairement, par exemple à zéro. On dit
qu’elles sont mises hors base.

Déf 4 : On appelle variables hors base les n variables de n+m fixées à zéro. Les m
variables restantes sont appelées variables de base.

Déf 5 : On appelle solution de base une solution où en ayant choisi n variables hors
base, on obtient une solution unique en résolvant les m contraintes d’égalité obtenues
en ajoutant les variables d’écart.

Notons xB les variables de base correspondant aux m variables positives ou nulles et xH

les variables hors base, correspondant aux n variables nulles.

Déf 6 : Une solution de base telle que xB>0 est dite admissible.

Principe de la méthode simplexe.

L’algorithme procède de la manière suivante :

1- On recherche un sommet de départ

2- On teste si ce sommet est l’optimum ou si la fonction objective n’est pas bornée
inférieurement ; dans ce cas le problème n’a pas de solution finie.
3- Si le sommet que l’on vient d’examiner n’est pas optimal, on se déplace sur un
sommet voisin pour lequel la fonction économique diminue et on repasse à
l’étape précédente.
Le nombre de sommets étant fini et tout minimum local étant absolu, le sommet optimal
est atteint lorsqu’aucun des sommets voisin ne permet plus de diminution du critère.
 A. Jouaiti

Première méthode de résolution (graphique)

Méthode graphique:

 A partir des contraintes déterminer l'ensemble des solutions admissibles (ensemble

des points qui vérifient toutes les contraintes simultanément)
 A partir de la fonction économique tracer la courbe d'indifférence (ensemble des
points qui donnent le même profit(ou le même coût)).
 Identifier la direction de déplacement de la courbe d'indifférence pour l'optimiser.
 Déterminer le dernier point en contact avec la courbe d'indifférence, c'est le point
optimal.

 Un point de plan cartésien (x1, x2) est dit admissible s'il satisfait toutes les
contraintes.
 L'ensemble des points admissibles forme la région admissible.
 Un point (x1, x2) de l'ensemble des solutions admissibles est aussi appelé
une solution possible.

Exemple 1:

Soit le programme linéaire suivant:

���� = 2. � + 2. �
��� ��� 2. � + � ≤ 10
2. � + 3. � ≤ 18
� ≤ 4 � + 3� ≥ 3
�≥0 �≥0
 A. Jouaiti

La zone hachurée représente l'ensemble des points admissibles comme solution. Afin
de trouver la solution optimale on trace la droite de la fonction économique pour un
profit minimal et on la déplace jusqu'à atteindre une valeur optimale de Z (un profit
maximal qui vérifient toutes les contraintes).

Le point M est le dernier point de contact entre la droite de la fonction économique et

la zone hachurée ce qui fait de ses coordonnées les valeurs optimal pour avoir un profit
maximal.

Il suffit de résoudre le système � = 10 − 2. � pour trouver les coordonnées du

point M. 2. � + 3. � = 18

Deuxième méthode de résolution

Exemple 2

Maximiser la fonction � = 5�1 + 4�2 + 3�3

2�1 + 3�2 + �3 ≤ 5
Sous les contraintes 4�1 + �2 + 2�3 ≤ 11
3�1 + 4�2 + 2�3 ≤ 8
�1, �2, �3 ≥ 0
Introduisant des variables d’écart.

� = 5�1 + 4�2 + 3�3

2�1 + 3�2 + �3 + �1 = 5
4�1 + �2 + 2�3 + �2 = 11
3�1 + 4�2 + 2�3 + �3 = 8
La contrainte devient �1, �2, �3, �1, �2, �3 ≥ 0

Partons d’une solution évidente x1=x2=x3=0. On trouve e1=5, e2=11, e3=8 et z=0

(S0) x1=x2=x3=0, e1=5, e2=11, e3=8

Les variables x1, x2 et x3 sont les variables hors base (qu’on fixe à zéro).
Les variables e1, e2 et e3 sont les variables de base.
 A. Jouaiti

Le but de la méthode est de faire croître Z, pour cela on va faire croître qu’une variable,
en choisissant celle qui a le plus grand coefficient dans Z=5x1+4x2+3x3
C'est-à-dire x1, on garde x2=x3=0.

De combien on va faire croître x1 (en gardant x2=x3=0) tout en maintenant e1, e2, e3
≥0

La condition e1=5- (2x1+3x2+x3) implique x1≤5/2

La condition e2= 11-(4x1+x2+2x3) implique x1≤11/4
La condition e3=8-(3x1+4x2+2x3) implique x1≤8/3

Pour que ces trois contraintes soient respectées, on doit choisir la condition sur x1 la
plus petite, c’est x1≤5/2. Faisons croitre x1 jusqu’à 5/2.
On obtient alors la solution :

(S1) x1=5/2, x2=0, x3=0, e1=0, e2=1, e3=1/2.

Ce qui donne Z=25/2 (ce qui mieux que Z=0)

x1 entre à la base et e1 sort de la base.

On exprime Z en fonction des nouvelles variables hors base x2, x3 et e1

Z=25/2 -7/2x2+1/2x3-5/2 e1

Le nouveau système est alors : (exprimer les variables de base en fonction des variables
hors base)

x1=5/2-3/2x2-1/2x3-1/2 e1
e2=1+5x2+0x3+2 e1
e3=1/2+1/2x2-1/2x3+3/2 e1

On remplace x1 par son expression dans la fonction Z.

Z =25/2 -7/2x2+1/2x3-5/2 e1

On recommence le processus en prenant x2=e1=0 et en faisant croitre x3 (signe positif).

x1= 5/2-1/2x3
e2=1+0 x3
e3=1/2-1/2x3
La contrainte �1, �2, �3 ≥ 0 implique: �3 ≤ 5 �� �3 ≤ 1

On prend la plus petite valeur de x3 pour que toutes les contraintes soient satisfaites
(�1, �2, �3 ≥ 0)
 A. Jouaiti

Faisons croitre x3 jusqu’à 1. On obtient alors la solution :

(S1) x1=2 x2=0 x3=1 e1=0 e2=1 e3=0

Ce qui donne Z=13 (ce qui est mieux que Z=25/2)

e3 sort de la base et x3 entre à la base

Exprimer les nouvelles variables de base en fonction des nouvelles variables hors bases

x1= 5/2-1/2(1-2 e3)=5/2-1/2+e3=2+e3

e2=1
x3=1-2e3

On remplace x1 par son expression dans la fonction Z.

Z=13-7/2 x2-1/2 e3-5/2 e1

Le coefficient de x2 est négatif donc on maintient la valeur de x2 à 0.

La solution est alors :

x1=2, x2=0 et x3=1
Z=13

Troisième méthode de résolution (tableau)

Exemple 3 :

Soit la fonction économique suivante :

Maximiser Z � = 3. � + 2. �
�+�≤7
2. � + � ≤ 9
On ajoute des variables d'écart positif pour transformer le système de deux inéquations à
2inconnues à un système de deux équations à 4 inconnues x, y, e1, e2.
On aura donc : � = 3. � + 2. �
� + � + �1 = 7
2. � + � + �2 = 9
La contrainte est : � ≥ 0; � ≥ 0; �1 ≥ 0; �2 ≥ 0

Nombre d'équation est inférieur au nombre d'inconnues : une infinité de solution. Alors
on doit fixer deux variables pour et trouves les autres variables.
 On démarre avec la solution de base S0=( 0, 0, 7, 9) x et y sont des variables hors
base, e1 et e2 sont des variables de base.
Z= 0
On aboutit à une analyse par tableaux ou matrices
 A. Jouaiti

Soit le tableau suivant:

� �1 , �2 , �1 , �2 = 3�1 + 2�2 + 0�1 + 0�2

C1 C2 C3 C4
Coeff 3 2 0 0
dans Z
Base x y e1 e2 bi (contrainte) �
R( � )
���
Coef.Z Var.
base
0 e1 �11 = �12 = �13 = �14 = � + � + �1 + 0�2 =7
1 1 1 0
0 e2 �21 = �21 = �31 = �41 = 2� + � + 0�1 + �2
2 1 0 1 =9
Zj 0 0 0 0 0
Cj-Zj 3 2 0 0

L’encadré vert correspond aux zj : c’est-à-dire à la somme �� ∗ ���

Choix de la variable entrante (dans la base)

Maximum des Cj-Zj pour des problèmes de max.

(‘Minimum des Cj-Zj pour des problèmes de min.)

La variable x donne une valeur Cj-Zj la plus grande (3) donc la variable entrant dans
la base est x.
k=1
 A. Jouaiti

Choix de la variable sortante

�
Dans un problème de min ou de max, la variable sortante sera le minimum des� �
��

C1 C2 C3 C4
Coeff 3 2 0 0
dans Z
Base x y e1 e2 bi (contrainte) R
Coef.Z Var. base
0 e1 1 1 1 0 7 7/1
0 e2 2 1 0 1 9 9/2
Zj 0 0 0 0 0
Cj-Zj 3 2 0 0

9/2 est le plus petit R, donc la variable sortante est e2.

� �1 , �2 , �1 , �2 = 3�1 + 2�2 + 0�1 + 0�2

C1 C2 C3 C4
Coeff 3 2 0 0
dans Z
Base x y e1 e2 bi (contrainte) R
Coef.Z Var. base
0 e1 1 1 1 0 7
3 x 2 1 0 1 9
Zj 0
Cj-Zj

2 est le pivot : c’est l’intersection de la colonne de la variable entrante et la ligne de la variable sortante

On divise la ligne du pivot par la valeur du pivot (ici égale à 2)

Pour calculer les coefficients du tableau A(aij) on utilise la formule suivante :

é�é���� �� �� ����� �� ����� � é�é���� �� �� ������� �� �����

��� −
�����
 A. Jouaiti

Puis pour la colonne du pivot on met des zéro sauf pour la valeur du pivot qui égale à
1.

� �1 , �2 , �1 , �2 = 3�1 + 2�2 + 0�1 + 0�2

C1 C2 C3 C4
Coeff 0 2 0 3
dans Z
Base e2 y e1 X bi (contrainte) R
Coef.Z Var. base
0 e1 0 0.5 1 -0.5 5/2
3 x 1 0.5 0 0.5 9/2
Zj 3 3/2 0 3/2 0
Cj-Zj -3 1/2 0 -3/2

La variable y donne une valeur Cj-Zj la plus grande (1/2) donc la variable entrant
dans la base est y. k=2

Recherche de la variable sortante

�
La variable sortante sera le minimum des� �
��

� �1 , �2 , �1 , �2 = 3�1 + 2�2 + 0�1 + 0�2

C1 C2 C3 C4
Coeff 3 2 0 0
dans Z
Base x y e1 e2 bi (contrainte) R
Coef.Z Var. base
0 e1 0 0.5 1 -0.5 5/2 5
3 x 1 0.5 0 0.5 9/2 9
Zj 0
Cj-Zj

R=5 est le rapport le plus petit donc la variable sortante est e1.
 A. Jouaiti

� �1 , �2 , �1 , �2 = 3�1 + 2�2 + 0�1 + 0�2

C1 C2 C3 C4
Coeff 3 2 0 0
dans Z
Base x y e1 e2 bi (contrainte) R
Coef.Z Var. base
2 y 0 1 2 -1 5
3 x 1 0.5 0 0.5 9/2
Zj 0
Cj-Zj

Transformation du tableau en procédant de la même manière que précédemment

� �1 , �2 , �1 , �2 = 3�1 + 2�2 + 0�1 + 0�2

C1 C2 C3 C4
Coeff 3 2 0 0
dans Z
Base x y e1 e2 bi (contrainte) R
Coef.Z Var. base
2 Y 0 1 2 -1 5
3 X 1 0 -1 1 2
Zj 3 2 1 1 0
Cj-Zj 0 0 -1 -1

x=2 et y=5 et Z=16

Le critère d’arrêt
Nous arrêtons lorsque nous obtenons le critère d'optimalité. L'algorithme du simplexe
s'arrête lorsque:
Cj-Zj ≤0 pour un problème de max
Cj-Zj ≥0 pour un problème de min
 A. Jouaiti

##############Code python de la méthode des simplex#################

import numpy as np

#Problème de maximisation.

#génère un vecteur de dimension var avec des coefficients égaux à

#zéros

def gen_vect(var):

tab = [Link](var)

return tab

#permet de déterminer les indices des variables de base

x=[]

def base_index(cons): #indices du système de base

for i in range(var,var+cons):

[Link](i)

###########Système à étudier ####################################

#fonction objective

#Z =3x1+2x2

Z_coef=[3,2]

#Contraintes

#x + y =7

#2x +y =9

Cont_Coef=[[1,1],[2,1]]

b=[7,9]

###################################################################

cons=len(b) #dimension du vecteur b

b=np.array_split(b, cons) # b=[[7],[9]]

#construction des matrices et vecteurs

var=len(Z_coef) #dimension du vecteur des variables hors base

Id=[Link](cons) #permet de créer une matrice Identité de

#dimension cons x cons

 A. Jouaiti

base_index(cons) #permet de déterminer les indices des variable de la base

Zj=gen_vect(var+cons+1) #[0,0,0,0,0]

Coef_Z=gen_vect(cons) # [0,0]

R=gen_vect(cons) # [0,0]

Cj=[Link](Z_coef,gen_vect(cons), axis=0) # permet de fusionner

#les deux vecteurs [3,2] & [0,0]= [3,2,0,0]

R=b #permet d initialser R

A=[Link](Cont_Coef, Id, axis=1) #permet de fusioner la matrice

# Cont_Coef et la matrice Identité

A=[Link](A, b, axis=1) #permet de fusioner la matrice A et la

# matrice b

#A = Cont_Coef + Id + b ce n est pas une sommation c est une fusion

print('A ',A)

index_c=[Link](Cj) #permet de déterminer l'indice dde la valeur

#maximale d'un vecteur

#construction des matrices et vecteurs

def CoefZ():

for i in range(cons):

Coef_Z[i]=Cj[x[i]]

def Z_j():

for i in range(var+cons+1):

Zj[i]=0

for j in range(cons):

Zj[i]=Zj[i]+Coef_Z[j]*A[j][i]

Cj_Zj=gen_vect(var+cons)

def CjZj():

for i in range(var+cons):

Cj_Zj[i]=Cj[i]-Zj[i]
 A. Jouaiti

#Calcul de R

def R_cal():

for i in range(cons):

R[i]=A[i][var+cons]/A[i][index_c]

#division/pivot

def Div_pivot():

pivot=A[index_l][index_c]

for j in range(var+cons+1):

A[index_l][j]=A[index_l][j]/pivot

# permet d introduire une nouvelle variable de base

def permutation():

x[index_l]=index_c

# permet de calculer les nouveaux coefficient de la matrice A

def calcul_A():

for i in range(cons):

for j in range(var+cons+1):

if (i!=index_l)and (j!=index_c):

A[i][j]=A[i][j]-A[i][index_c]*A[index_l][j]

for i in range(cons):

if i!=index_l:

A[i][index_c]=0

############### CALCUL ##################

N=0 # condition d arrêt sur le nombre d itération

Zmax=10000

while (Zmax>=0):

N=N+1
 A. Jouaiti

Z_j()

CjZj()

Zmax=[Link](Cj_Zj)

if Zmax<=0:

print('le résultat obtenu:')

for i in range(var):

print('x',x[i],' = ',A[i][var+cons])

print('z = ',Zj[var+cons])

break

index_c=[Link](Cj_Zj)

R_cal()

index_l=[Link](R)

Div_pivot()

calcul_A()

permutation()

CoefZ()

if N==var+cons:

print('le système diverge:')

break

Problème
Afin d’assurer une bonne santé de cobayes, il faut les nourrir en leur donnant un
minimum de 24 grammes de lipides, 36 grammes de glucides et 4 grammes de protéines
par jour. Il ne faut pas leur donner plus de 5 unités de nourriture. Nous disposions de
deux sources d’alimentation. Les croquettes royal cobaye qui contiennent 8 grammes de
lipides, 12 grammes de glucides et 2 grammes de protéines par unité et qui coutent 3
euros les 10 unités, et les boulettes « vite se casse » qui contiennent 12 grammes de
lipides, 12 grammes de glucides et 1 gramme de protéines par unité et qui coutent 2
euros les 10 unités.