Programmation lineaire
Description
Voici un programme que j'ai concu il y a quelques mois.
Il resout le probleme de la programmation lineaire par
l'algorithme du simplexe. Tout est tres bien explique
dans numerical recipes au chapitre 10, vous expliquez
serait seulement de la paraphrase...
L'atout majeur du programme est l'apercu graphique du
probleme, ce qui permet d'avoir une bonne visualisation
du deroulement de l'algorithme, et en meme permet de
justifer la validite de l'implementation (pas ou moins
de bug)
La visualisation se fait en toute dimension, alors ames
sensibles s'abstenir pour ceux qui pensent qu'il n'est
pas possible d'afficher un espace N dimensions sur un
ecran 2 dimensions. Les axes sont marques en bas a droite,
la rotation des axes s'effectue avec les touches 'P' et 'M',
et le choix des deux axes par les nombres 1,2,3...
et SHIFT+1,2,3... au dessus du clavier (pas le pave numerique)
La touche espace pour faire une etape de l'algorithme, les
aretes visitees se colorient au fur et a mesure en violet
(ou orange), et le polytope est en vers.
En haut a gauche il y a le tableau de l'algorithme similaire
a celui decrit dans le numerical recipes (LA reference...)
On possede n hyperplans et il y a p variables (on vit
donc dans un espace de dimension p)
Le calcul du polytope est simple et bourrin : on calcule tous
les sommets (on calcule les intersections de p hyperplans
choisis parmi n, ce qui en fait C(p,n), ce qui est grand,
l'intersection se calcule en inversant une matrice...
pour ce faire l'algorithme du pivot de Gauss est le
bien-venu)
Ensuite il suffit de calculer les aretes a partir de
l'ensemble des sommets : deux sommets qui differs qu'un
hyperplan sont les sommets d'une aretes du polytope,
c'est logique ! facile.
L'exemple du polytope de gauche est un hypercube en 4D.
L'exemple du polytope de droite est genere aleatoirement.
En pratique l'algorithme du simplexe est polynomial,
mais en theorie il est exponentiel (voir contre-exemple
de Klee-Minty)
Petit rappel :
La programmation lineaire consiste a maximiser une fonction
linaire (gradient constant) des variables x1,x2...xp sous
contraintes lineaires de la forme "dans un demi-espace"
(i.e. a1.x1+a2.x2+...ap.xp<=A).
Ce probleme ce rencontre dans les problemes de flots
maximals, ou dans des problemes d'optimisations de productions
en economie... les applications ne manquent pas.
L'ensemble des solutions de ce probleme est un polytope
(parti convexe dont la frontiere est des bouts d'hyperplans)
L'idee du simplexe est d'abord de dire que la solution
optimale est forcement sur un sommet, ensuite de dire que
l'on part initialement d'un sommet et qu'a chaque iteration
de l'algorithme on change de sommet en se deplacant suivant
(ou "le long de") d'une arete, et en se deplacant ainsi,
on a augmente la fonction objective (celle que l'on cherche
a maximiser). Le choix de l'arete est bien explique dans
numerial recipes.
... arretons-la la paraphrase.
N'hesite pas a poser des questions !
NB: le code qui suit, et l'algorithme sur une etape,
il gere seulement un seul pivot. Il est imcomprehensible
sans la connaissance de numerical recipes...basic vars
et slack vars...
Il resout le probleme de la programmation lineaire par
l'algorithme du simplexe. Tout est tres bien explique
dans numerical recipes au chapitre 10, vous expliquez
serait seulement de la paraphrase...
L'atout majeur du programme est l'apercu graphique du
probleme, ce qui permet d'avoir une bonne visualisation
du deroulement de l'algorithme, et en meme permet de
justifer la validite de l'implementation (pas ou moins
de bug)
La visualisation se fait en toute dimension, alors ames
sensibles s'abstenir pour ceux qui pensent qu'il n'est
pas possible d'afficher un espace N dimensions sur un
ecran 2 dimensions. Les axes sont marques en bas a droite,
la rotation des axes s'effectue avec les touches 'P' et 'M',
et le choix des deux axes par les nombres 1,2,3...
et SHIFT+1,2,3... au dessus du clavier (pas le pave numerique)
La touche espace pour faire une etape de l'algorithme, les
aretes visitees se colorient au fur et a mesure en violet
(ou orange), et le polytope est en vers.
En haut a gauche il y a le tableau de l'algorithme similaire
a celui decrit dans le numerical recipes (LA reference...)
On possede n hyperplans et il y a p variables (on vit
donc dans un espace de dimension p)
Le calcul du polytope est simple et bourrin : on calcule tous
les sommets (on calcule les intersections de p hyperplans
choisis parmi n, ce qui en fait C(p,n), ce qui est grand,
l'intersection se calcule en inversant une matrice...
pour ce faire l'algorithme du pivot de Gauss est le
bien-venu)
Ensuite il suffit de calculer les aretes a partir de
l'ensemble des sommets : deux sommets qui differs qu'un
hyperplan sont les sommets d'une aretes du polytope,
c'est logique ! facile.
L'exemple du polytope de gauche est un hypercube en 4D.
L'exemple du polytope de droite est genere aleatoirement.
En pratique l'algorithme du simplexe est polynomial,
mais en theorie il est exponentiel (voir contre-exemple
de Klee-Minty)
Petit rappel :
La programmation lineaire consiste a maximiser une fonction
linaire (gradient constant) des variables x1,x2...xp sous
contraintes lineaires de la forme "dans un demi-espace"
(i.e. a1.x1+a2.x2+...ap.xp<=A).
Ce probleme ce rencontre dans les problemes de flots
maximals, ou dans des problemes d'optimisations de productions
en economie... les applications ne manquent pas.
L'ensemble des solutions de ce probleme est un polytope
(parti convexe dont la frontiere est des bouts d'hyperplans)
L'idee du simplexe est d'abord de dire que la solution
optimale est forcement sur un sommet, ensuite de dire que
l'on part initialement d'un sommet et qu'a chaque iteration
de l'algorithme on change de sommet en se deplacant suivant
(ou "le long de") d'une arete, et en se deplacant ainsi,
on a augmente la fonction objective (celle que l'on cherche
a maximiser). Le choix de l'arete est bien explique dans
numerial recipes.
... arretons-la la paraphrase.
N'hesite pas a poser des questions !
NB: le code qui suit, et l'algorithme sur une etape,
il gere seulement un seul pivot. Il est imcomprehensible
sans la connaissance de numerical recipes...basic vars
et slack vars...
Source / Exemple :
//-------------------------------------------------
BOOL oneStep(P_DATA_WINDOW data)
{
int i,j;
double deltaZj,deltaRVj;
double deltaZMax,deltaRVjMax;
int jMax,iLimjMax;
deltaZMax = 0.;
for(j=1;j<=data->nBasicVar;j++)
{
int iLim;
double deltaMax;
iLim = 0;
for(i=1;i<=data->nSlackVar;i++)
{
double C,beta;
C = GetItemMatrix(data->mat,i,0);
beta = GetItemMatrix(data->mat,i,j);
if(beta < 0)
{
double delta;
delta = - C / beta;
if(0 == iLim || delta < deltaMax)
{
iLim = i;
deltaMax = delta;
}
}
}
ErrorIf((0 == iLim,"le polytope n'est pas borné !!"));
deltaRVj = deltaMax;
deltaZj = GetItemMatrix(data->mat,0,j) * deltaRVj;
if(deltaZj > deltaZMax)
{
deltaZMax = deltaZj;
deltaRVjMax = deltaRVj;
jMax = j;
iLimjMax = iLim;
}
}
if(deltaZMax <= 0.)
{
//message(("\n\nFIN DE L'ALGORITHME DU SIMPLEXE\n\n"));
return FALSE;
}
// on echange les deux variables
{
int i0,j0;
double betai0j0;
i0 = iLimjMax;
j0 = jMax;
betai0j0 = GetItemMatrix(data->mat,i0,j0);
for(i=0;i<=data->nSlackVar;i++)
{
if(i != i0)
{
SetItemMatrix(data->matTmp,i,0,
GetItemMatrix(data->mat,i,0)
-
GetItemMatrix(data->mat,i,j0)*GetItemMatrix(data->mat,i0,0)/betai0j0
);
for(j=1;j<=data->nBasicVar;j++)
{
if(j != j0)
{
SetItemMatrix(data->matTmp,i,j,
GetItemMatrix(data->mat,i,j)
-
GetItemMatrix(data->mat,i0,j)*GetItemMatrix(data->mat,i,j0)/betai0j0
);
}
else
{
SetItemMatrix(data->matTmp,i,j,
GetItemMatrix(data->mat,i,j0)/betai0j0
);
}
}
}
else // i == i0
{
SetItemMatrix(data->matTmp,i,0,
-
GetItemMatrix(data->mat,i0,0)/betai0j0
);
for(j=1;j<=data->nBasicVar;j++)
{
if(j != j0)
{
SetItemMatrix(data->matTmp,i,j,
-
GetItemMatrix(data->mat,i0,j)/betai0j0
);
}
else
{
SetItemMatrix(data->matTmp,i,j,
1./betai0j0
);
}
}
}
}
}
SwapPointer(&data->mat,&data->matTmp);
SwapPointer(&data->RVar[jMax-1],&data->LVar[iLimjMax-1]);
// on attribue les valeurs
for(i=0;i<data->nSlackVar;i++)
{
data->LVar[i]->val = GetItemMatrix(data->mat,1+i,0);
}
for(j=0;j<data->nBasicVar;j++)
{
data->RVar[j]->val = 0.;
}
return TRUE;
} // oneStep()
Vous n'êtes pas encore membre ?
inscrivez-vous, c'est gratuit et ça prend moins d'une minute !
Les membres obtiennent plus de réponses que les utilisateurs anonymes.
Le fait d'être membre vous permet d'avoir un suivi détaillé de vos demandes et codes sources.
Le fait d'être membre vous permet d'avoir des options supplémentaires.