[Top][All Lists]
[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
[Help-glpk] oops... mistake !
|
From: |
Nicoco BERGER |
|
Subject: |
[Help-glpk] oops... mistake ! |
|
Date: |
Mon, 10 Oct 2005 12:57:48 +0200 (CEST) |
I have sent this message to the wrong adress... I have two problems with two
different
solvers ;) ! Here my problem still is that glpk does not solve anything, I
attach my code
to this mail. See the function named ..._glpk
___________________________________________________________________________
Appel audio GRATUIT partout dans le monde avec le nouveau Yahoo! Messenger
Téléchargez cette version sur http://fr.messenger.yahoo.comextern "C" {
#include "lp_solve/lp_lib.h"
#include "glpk.h"
}
#include <iostream>
#include <vector>
#include <fstream>
#include <iostream>
#include <string>
using namespace std ;
long evaluation(vector<int>* coef_obj, vector<int>* sol_coef_a_un) {
if((coef_obj ==NULL) || (sol_coef_a_un==NULL)) {
cout << "pointeur nul en entrée de evaluation\n" ;
exit(-1) ;
}
vector<int>::iterator it ;
long val = 0 ;
for(it = sol_coef_a_un->begin() ; it != sol_coef_a_un->end() ; it ++) {
if((*it) < coef_obj->size()) {
//cout << "variable " << *it <<" coef = "<< (*coef_obj)[*it] <<"\n";
val += (*coef_obj)[*it] ;
}
else {
cout << "Attention : utilisation d'une variable hors limites dans
evaluation" ;
exit(-1) ;
}
}
return val ;
}
void affichage(std::vector<int>& vector, std::string nom){
int i = vector.size();
int j = 1;
cout << '\n';
cout << nom << '\n';
cout << "taille : ";
cout << i << '\n';
cout << '\n';
cout << '\n';
while (i != 0)
{
cout << vector[vector.size() - i ] ;
if (j == 10)
{
cout << '\n';
j = 0;
}
else {cout << ' ';}
j = j + 1;
i = i - 1;
}
cout << '\n';
}//fin de l'affichage
bool is_readable( const std::string & file )
{
std::ifstream fichier( file.c_str() );
return (fichier != 0);
}
bool existe(std::string nom)
{
using std::cout;
if ( is_readable( nom ) )
{
//cout << "Fichier existant et lisible.\n";
return 1;
}
else
{
cout << "Fichier inexistant ou non lisible.\n";
return 0;
}
}
bool charger_DAT(ifstream* fichier, vector<int>* vector_obj1, vector<int>*
vector_obj2,
vector<int>* vector_cont, int* N, int* W )
{
bool ret = false ;
int p;
int k;
int i;
std::string ligne; // variable contenant chaque ligne lue
// on cherche à lire le n
getline( *fichier, ligne ) ;
//cout << ligne << "\n" ;
getline( *fichier, ligne ) ;
//cout << ligne << "\n" ;
if(ligne[0] != '#') {
ret = false ;
}
else {
string n_bogue ;
if((n_bogue = ligne.substr(ligne.find("N") + 1)).size() > 1) {
(*N) = atoi(n_bogue.c_str()) ;
}
else {
getline( *fichier, ligne ) ;
(*N) = atoi(ligne.substr(0).c_str()) ;
}
//cout << ligne << "\n" ;
//initialise les vecteurs objectifs
vector_obj1->resize(*N);
vector_obj2->resize(*N);
vector_cont->resize(*N);
//cout << "a priori il y a " << (*N) << " variables\n" ;
// on cherche le premier objectif
while(std::getline( *fichier, ligne ) && (ligne.substr(0,3) != "# O")) ;
//cout << ligne << "\n" ;
int n_obj1 = 0 ;
// on lit le premier objectif
while(getline( *fichier, ligne ) && (ligne[0] != '#')) {
//cout << (int)ligne[0] << "\n" ;
if((ligne.size() != 0) && ((int)ligne[0] != 13)) {
(*vector_obj1)[n_obj1] = atoi(ligne.c_str()) ;
n_obj1 ++ ;
}
}
//cout << "il y a en réalité " << n_obj1 << "variables\n" ;
int n_obj2 =0 ;
// on lit le deuxième objectif
while(getline( *fichier, ligne ) && (ligne[0] != '#')) {
//cout << ligne << "\n" ;
if((ligne.size() != 0)&& ((int)ligne[0] != 13)) {
(*vector_obj2)[n_obj2] = atoi(ligne.c_str()) ;
n_obj2 ++ ;
}
}
//cout << "il y a en réalité " << n_obj2 << "variables\n" ;
int n_w = 0 ;
// on lit le deuxième objectif
while(getline( *fichier, ligne ) && (n_w < (*N)) ) {
if((ligne.size() != 0)&& ((int)ligne[0] != 13)) {
(*vector_cont)[n_w] = atoi(ligne.c_str()) ;
n_w ++ ;
}
}
//cout << "il y a en réalité " << n_w << "variables\n" ;
if((n_obj1 != (*N))||(n_obj1 != (*N))) {
ret = false ;
}
else {
//cout << "il y a en réalité " << n_obj1 << "variables\n" ;
getline( *fichier, ligne ) ;
//cout << "derniere ligne : " << ligne << "\n";
(*W) = atoi(ligne.c_str()) ;
//cout << " w = " << *W << "\n" ;
ret = true ;
}
}
return ret ;
}
bool charger_KNP(ifstream* fichier, vector<int>* vector_obj1, vector<int>*
vector_obj2,
vector<int>* vector_cont, int* N, int* W ) {
bool ret = false ;
int p;
int k;
int i;
std::string ligne; // variable contenant chaque ligne lue
while ( std::getline( *fichier, ligne ) && (ligne[0] == 'c')) ; // les
commentaires de début de fichier
if(ligne[0] != 'n') {
ret = false ;
}
else {
(*N) = atoi(ligne.substr(2).c_str()) ;
//initialise les vecteurs objectifs
vector_obj1->resize(*N);
vector_obj2->resize(*N);
vector_cont->resize(*N);
//cout << "a priori il y a " << (*N) << " variables\n" ;
std::getline( *fichier, ligne ) ;
int ind1 = 0, ind2 = 0, ind3 = 0, ind4 = 0, n = 0 ;
while(getline( *fichier, ligne ) && (ligne[0] == 'i')) {
ind1 = ligne.find_first_not_of(" ", 1) ;
ind2 = ligne.find(" ", ind1) ;
ind3 = ligne.find_first_not_of(" ", ind2) ;
ind4 = ligne.find(" ", ind3) ;
//cout << ligne.substr(ind1, ind2 - ind1 + 1) << " "<< ligne.substr(ind3,
ind4 - ind3 + 1) ;
//cout << ligne.substr(ind4) << "\n" ;
(*vector_cont)[n] = atoi(ligne.substr(ind1, ind2 - ind1 + 1).c_str()) ;
(*vector_obj1)[n] = atoi(ligne.substr(ind3, ind4 - ind3 + 1).c_str()) ;
(*vector_obj2)[n] = atoi(ligne.substr(ind4).c_str()) ;
n++ ;
}
if(n != (*N)) {
ret = false ;
}
else {
//cout << "il y a en réalité " << n << "variables\n" ;
//getline( *fichier, ligne ) ;
if(getline( *fichier, ligne ) && (ligne[0] == 'W')) {
//cout << "derniere ligne : " << ligne << "\n";
(*W) = atoi(ligne.substr(2).c_str()) ;
ret = true ;
}
else
ret = false ;
//cout << "W = " << (*W) << "\n" ;
}
}
return ret ;
}
bool charger_fichier(string nom, vector<int>* vector_obj1, vector<int>*
vector_obj2,
vector<int>* vector_cont, int* N, int* W )
{
bool ret = false ;
if(!existe(nom)) {
ret = false ;
}
else {
//le fichier existe
int p;
int k;
int i;
// le constructeur de ifstream permet d'ouvrir un fichier en lecture
ifstream fichier(nom.c_str());
if( fichier ) {// ce test échoue si le fichier n'est pas ouvert
string ligne; // variable contenant chaque ligne lue
if(getline( fichier, ligne )) {
// on teste le type de syntaxe du fichier
if(ligne[0] == '#') {
ret = charger_DAT(&fichier, vector_obj1, vector_obj2, vector_cont, N,
W) ;
}
else {
if(ligne[0] == 'c') {
ret = charger_KNP(&fichier, vector_obj1, vector_obj2, vector_cont,
N, W) ;
}
else {
ret = false ;
}
}
}
else {
ret = false ;
}
}
else {
ret = false ;
}
}
return ret ;
}
int demo_lpsolve()
{
//cout << "debut demo\n" ;
lprec *lp;
int Ncol, *colno = NULL, j, ret = 0;
REAL *row = NULL;
/* We will build the model row by row
So we start with creating a model with 0 rows and 2 columns */
Ncol = 2; /* there are two variables in the model */
lp = make_lp(0, Ncol);
if(lp == NULL)
ret = 1; /* couldn't construct a new model... */
if(ret == 0) {
/* let us name our variables. Not required, but can be useful for debugging
*/
set_col_name(lp, 1, "x");
set_col_name(lp, 2, "y");
/* create space large enough for one row */
colno = (int *) malloc(Ncol * sizeof(*colno));
row = (REAL *) malloc(Ncol * sizeof(*row));
if((colno == NULL) || (row == NULL))
ret = 2;
}
if(ret == 0) {
set_add_rowmode(lp, TRUE); /* makes building the model faster if it is
done rows by row */
/* construct first row (120 x + 210 y <= 15000) */
j = 0;
colno[j] = 1; /* first column */
row[j++] = 120;
colno[j] = 2; /* second column */
row[j++] = 210;
/* add the row to lpsolve */
if(!add_constraintex(lp, j, row, colno, LE, 15000))
ret = 3;
}
if(ret == 0) {
/* construct second row (110 x + 30 y <= 4000) */
j = 0;
colno[j] = 1; /* first column */
row[j++] = 110;
colno[j] = 2; /* second column */
row[j++] = 30;
/* add the row to lpsolve */
if(!add_constraintex(lp, j, row, colno, LE, 4000))
ret = 3;
}
if(ret == 0) {
/* construct third row (x + y <= 75) */
j = 0;
colno[j] = 1; /* first column */
row[j++] = 1;
colno[j] = 2; /* second column */
row[j++] = 1;
/* add the row to lpsolve */
if(!add_constraintex(lp, j, row, colno, LE, 75))
ret = 3;
}
if(ret == 0) {
set_add_rowmode(lp, FALSE); /* rowmode should be turned off again when done
building the model */
/* set the objective function (143 x + 60 y) */
j = 0;
colno[j] = 1; /* first column */
row[j++] = 143;
colno[j] = 2; /* second column */
row[j++] = 60;
/* set the objective in lpsolve */
if(!set_obj_fnex(lp, j, row, colno))
ret = 4;
}
if(ret == 0) {
/* set the object direction to maximize */
set_maxim(lp);
/* just out of curioucity, now show the model in lp format on screen */
/* this only works if this is a console application. If not, use write_lp
and a filename */
write_LP(lp, stdout);
/* write_lp(lp, "model.lp"); */
/* I only want to see important messages on screen while solving */
set_verbose(lp, IMPORTANT);
/* Now let lpsolve calculate a solution */
ret = solve(lp);
if(ret == OPTIMAL)
ret = 0;
else
ret = 5;
}
if(ret == 0) {
/* a solution is calculated, now lets get some results */
/* objective value */
printf("Objective value: %f\n", get_objective(lp));
/* variable values */
get_variables(lp, row);
for(j = 0; j < Ncol; j++)
printf("%s: %f\n", get_col_name(lp, j + 1), row[j]);
/* we are done now */
}
/* free allocated memory */
if(row != NULL)
free(row);
if(colno != NULL)
free(colno);
if(lp != NULL) {
/* clean up such that all used memory by lpsolve is freed */
delete_lp(lp);
}
return(ret);
}
int demo_glpk() {
LPX *lp ;
int ia[1+1000], ja[1+1000];
double ar[1+1000], Z, x1, x2, x3;
lp = lpx_create_prob();
lpx_set_prob_name(lp, "sample");
lpx_set_obj_dir(lp, LPX_MAX);
lpx_add_rows(lp, 3);
lpx_set_row_name(lp, 1, "p");
lpx_set_row_bnds(lp, 1, LPX_UP, 0.0, 100.0);
lpx_set_row_name(lp, 2, "q");
lpx_set_row_bnds(lp, 2, LPX_UP, 0.0, 600.0);
lpx_set_row_name(lp, 3, "r");
lpx_set_row_bnds(lp, 3, LPX_UP, 0.0, 300.0);
lpx_add_cols(lp, 3);
lpx_set_col_name(lp, 1, "x1");
lpx_set_col_bnds(lp, 1, LPX_LO, 0.0, 0.0);
lpx_set_obj_coef(lp, 1, 10.0);
lpx_set_col_name(lp, 2, "x2");
lpx_set_col_bnds(lp, 2, LPX_LO, 0.0, 0.0);
lpx_set_obj_coef(lp, 2, 6.0);
lpx_set_col_name(lp, 3, "x3");
lpx_set_col_bnds(lp, 3, LPX_LO, 0.0, 0.0);
lpx_set_obj_coef(lp, 3, 4.0);
ia[1] = 1, ja[1] = 1, ar[1] = 1.0; /* a[1,1] = 1 */
ia[2] = 1, ja[2] = 2, ar[2] = 1.0; /* a[1,2] = 1 */
ia[3] = 1, ja[3] = 3, ar[3] = 1.0; /* a[1,3] = 1 */
ia[4] = 2, ja[4] = 1, ar[4] = 10.0; /* a[2,1] = 10 */
ia[5] = 3, ja[5] = 1, ar[5] = 2.0; /* a[3,1] = 2 */
ia[6] = 2, ja[6] = 2, ar[6] = 4.0; /* a[2,2] = 4 */
ia[7] = 3, ja[7] = 2, ar[7] = 2.0; /* a[3,2] = 2 */
ia[8] = 2, ja[8] = 3, ar[8] = 5.0; /* a[2,3] = 5 */
ia[9] = 3, ja[9] = 3, ar[9] = 6.0; /* a[3,3] = 6 */
lpx_load_matrix(lp, 9, ia, ja, ar);
lpx_set_int_parm(lp, LPX_K_MSGLEV, 0) ;// pas de messages
lpx_simplex(lp);
Z = lpx_get_obj_val(lp);
x1 = lpx_get_col_prim(lp, 1);
x2 = lpx_get_col_prim(lp, 2);
x3 = lpx_get_col_prim(lp, 3);
printf("\nZ = %g; x1 = %g; x2 = %g; x3 = %g\n", Z, x1, x2, x3);
lpx_delete_prob(lp);
return 0 ;
}
int ajoute_contrainte(lprec* lp, int Ncol, vector<int>* coef, int type_cont,
double borne) {
if((lp==NULL) || (coef==NULL)) {
cout << "pointeur nul en entrée de ajoute_contrainte\n" ;
exit(-1) ;
}
int *colno = NULL, j, ret = 0;
REAL *row = NULL;
// on cree de l'espace pour une ligne
colno = (int *) malloc(Ncol * sizeof(*colno));
row = (REAL *) malloc(Ncol * sizeof(*row));
if((colno == NULL) || (row == NULL))
ret = 2;
if(ret == 0) {
set_add_rowmode(lp, TRUE); /* makes building the model faster if it is
done rows by row */
vector<int>::iterator it ;
it = coef->begin() ;
for(j = 0 ; j < Ncol ; j ++) {
if(it == coef->end()) {
cout << "Problème dans les données : moins de " << Ncol << " entiers
dans w\n" ;
exit(-1) ;
}
colno[j] = j+1; /* column par colonne*/
row[j] = (*it) ;
it ++ ;
}
// on ajoute la ligne au solveur
if(!add_constraintex(lp, j, row, colno, type_cont, borne))
ret = 3;
set_add_rowmode(lp, FALSE); /* rowmode should be turned off again when done
building the model */
}
/* free allocated memory */
if(row != NULL)
free(row);
if(colno != NULL)
free(colno);
return ret ;
}
lprec* init_modele_lpsolve(int n, vector<int>* w, int cap) {
if((w==NULL)) {
cout << "pointeur nul en entrée de init_model_lpsolve\n" ;
exit(-1) ;
}
lprec *lp; // variables utiles pour lpsolve
int Ncol, *colno = NULL, j, ret = 0;
REAL *row = NULL;
vector<int>::iterator it ;
// Construction du modèle sans objectif
Ncol = n; // il y a n variables
lp = make_lp(0, Ncol);
set_verbose(lp, IMPORTANT); // afficher uniquement les messages importants
set_maxim(lp); // le but est de trouver un maximum
if(lp == NULL)
ret = 1; // quand on n'a pas pas pu construire de modèle
// On construit le modèle ligne par ligne
// Au départ il y a n colonnes mais 0 ligne
if(ret == 0) {
// on nomme les variables, ce n'est pas obligatoire mais ca peut etre utile
pour le débogage
char* var_name;
for(int i = 0 ; i < Ncol ; i++) {
var_name = (char*) malloc(5) ;
sprintf(var_name, "x%d", i) ;
set_col_name(lp, i+1, var_name);
// nos variables sont binaires
set_binary(lp, i+1, TRUE) ;
}
}
// on construit la contrainte des poids des objets
ret = ajoute_contrainte(lp, n, w, LE, cap) ;
if(ret == 0)
return lp ;
else {
cout << "echec de l'initialisation du solveur\n" ;
exit(-1) ;
}
}
int set_objectif(lprec* lp, vector<int>* coef, int Ncol) {
if((lp==NULL)||(coef==NULL)) {
cout << "pointeur nul en entrée de set_objectif\n" ;
exit(-1) ;
}
int *colno = NULL, j, ret = 0;
REAL *row = NULL;
// on cree de l'espace pour une ligne
colno = (int *) malloc(Ncol * sizeof(*colno));
row = (REAL *) malloc(Ncol * sizeof(*row));
if((colno == NULL) || (row == NULL))
ret = 2;
if(ret == 0) {
vector<int>::iterator it = coef->begin() ;
for(j = 0 ; j < Ncol ; j ++) {
if(it == coef->end()) {
cout << "Problème dans les données : moins de " << Ncol << " entiers
dans p_prim\n" ;
exit(-1) ;
}
colno[j] = j+1; /* first column */
row[j] = (*it);
it ++ ;
}
// on mémorise
if(!set_obj_fnex(lp, j, row, colno))
ret = 4;
}
/* free allocated memory */
if(row != NULL)
free(row);
if(colno != NULL)
free(colno);
return ret ;
}
int get_solution(lprec* lp, int n, vector<int>* sol, double* max){
if((lp==NULL)||(sol==NULL)||(max==NULL)) {
cout << "pointeur nul en entrée de get_solution\n" ;
exit(-1) ;
}
int Ncol = n, *colno = NULL, j, ret = 0;
REAL *row = NULL;
if(ret == 0) {
if(get_solutioncount(lp) > 1) {
cout << "il y a plus d'une solution !" ;
exit(-1) ;
}
sol->clear() ;
(*max) = get_objective(lp) ;
// on cree de l'espace pour une ligne
colno = (int *) malloc(Ncol * sizeof(*colno));
row = (REAL *) malloc(Ncol * sizeof(*row));
if((colno == NULL) || (row == NULL))
ret = 2;
get_variables(lp, row);
for(j = 0; j < Ncol; j++) {
if(row[j] != 0) {
// mémoriser la solution : les indices des variables à un
sol->push_back(j) ;
}
}
}
/* free allocated memory */
if(row != NULL)
free(row);
if(colno != NULL)
free(colno);
return ret ;
}
int resolution_lexicographique_lpsolve(int n, int cap, vector<int>* p_prim,
vector<int>* p_sec,
vector<int>* w, vector<vector<int>*>* E)
{
if((p_prim==NULL)||(p_sec==NULL)||(w==NULL)||(E==NULL)) {
cout << "pointeur nul en entrée de resolution_lexicographique_lpsolve\n" ;
exit(-1) ;
}
lprec *lp; // variables utiles pour lpsolve
int Ncol = n, *colno = NULL, j, ret = 0;
REAL *row = NULL;
vector<int>::iterator it ;
vector<int>* sol_coef_a_un = new vector<int>() ;
lp = init_modele_lpsolve(n, w, cap) ;
// a) obtenir la solution optimale selon le premier objectif à maximiser
if(ret == 0) {
// on définit l'objectif primaire
ret = set_objectif(lp, p_prim, Ncol) ;
}
//cout << "après le premier objectif ret vaut "<< ret <<"\n" ;
if(ret == 0) {
//cout << "avant z_prim\n";
//write_LP(lp, stdout);
// on calcule la solution optimale (pour le premier objectif seul)
ret = solve(lp);
if(ret == OPTIMAL)
ret = 0;
else
ret = 5;
}
double max_prim ;
// on récupère la solution et le maximum atteint par l'objectif primaire
ret = get_solution(lp, n, sol_coef_a_un, &max_prim) ;
// libère la mémoire
delete_lp(lp) ;
// b) chercher à améliorer l'autre objectif SANS perdre l'optimalité sur
le premier
// i) ajouter une contrainte imposant de ne pas descendre sur le premier
objectif
lp = init_modele_lpsolve(n, w, cap) ; // on réinitialise un modèle (bug
sinon)
ret = ajoute_contrainte(lp, n, p_prim, GE, max_prim) ;
// ii) définir l'objectif secondaire
ret = set_objectif(lp, p_sec, Ncol) ;
if(ret == 0) {
//cout << "avant z_sec\n";
//write_LP(lp, stdout);
// on calcule la solution optimale (pour le deuxième objectif sous le
contrainte de ne pas faire chuter le premier)
ret = solve(lp);
if(ret == OPTIMAL)
ret = 0;
else
ret = 5;
}
// on sauvegarde la solution
ret = get_solution(lp, n, sol_coef_a_un, &max_prim) ;
if(ret == 0) {
E->push_back(sol_coef_a_un) ;
}
delete_lp(lp) ;
return ret ;
}
int resolution_recursive_lpsolve(vector<int>* xA, int z1_A, int z2_A,
vector<int>* xB, int z1_B, int z2_B,
vector<vector<int>*>* E,
vector<int>* p1, vector<int>* p2, vector<int>*
w, int cap, int n) {
//cout << "\nresol_recurs avec \n" ;
//affichage(*xA, "xA") ;
//cout << "z1(xA) = " << z1_A << " z2(xA) = " << z2_A << "\n" ;
//affichage(*xB, "xB") ;
//cout << "z1(xB) = " << z1_B << " z2(xB) = " << z2_B << "\n" ;
if((xA==NULL)||(xB==NULL)||(E==NULL)||(p1==NULL)||(p2==NULL)||(w==NULL)) {
cout << "pointeur nul en entrée de \n" ;
exit(-1) ;
}
int ret ;
lprec* lp = init_modele_lpsolve(n, w, cap) ;
// on ajoute les contraintes de ne pas aller en dessous de z2(xA) ni en
dessous de z1(xB) => éviter de boucler
ret = ajoute_contrainte(lp, n, p1, GE, z1_B + 1) ; // attention : si pas +1
on boucle
ret = ajoute_contrainte(lp, n, p2, GE, z2_A + 1) ;
// on pose : l1 = z2(xB) - z2(xA) et l2 = z1(xA) - z1(xB) et on calcule la
solution optimale
// de max{l1 * z1(x) + l2 * z2(x)} => direction de recherche othogonale au
vecteur xAxB
int l1 = z2_B - z2_A ;
int l2 = z1_A - z1_B ;
// on crée le vecteur des coefs du nouveau polynome objectif
vector<int>::iterator it1 = p1->begin(), it2 = p2->begin() ;
vector<int>* p_ortho = new vector<int>() ;
for(int j = 0 ; j < n ; j++) {
if((it1 == p1->end()) || (it2 == p2->end())){
cout << "Problème dans les données : moins de " << n << " entiers dans
p1 ou p2\n" ;
exit(-1) ;
}
p_ortho->push_back(l1 * (*it1) + l2 * (*it2)) ;
it1 ++ ;
it2 ++ ;
}
// on met en place l'objectif
ret = set_objectif(lp, p_ortho, n) ;
// on résoud
if(ret == 0) {
ret = solve(lp);
if(ret == OPTIMAL)
ret = 0;
else
ret = 5;
}
if((ret == 0) && (get_solutioncount(lp) != 0)){
double max ;
vector<int>* xC = new vector<int>() ;
// on sauvegarde la solution
ret = get_solution(lp, n, xC, &max) ;
int z1_C = evaluation(p1, xC),
z2_C = evaluation(p2, xC);
//cout << "nouvelle solution :\n" ;
//affichage(*xC, "xC") ;
//cout << "z1(xC) = " << z1_C << " z2(xC) = " << z2_C << "\n" ;
// on mémorise la nouvelle solution
E->push_back(xC) ;
// on efface avant le rappel pour ne pas saturer la mémoire
delete_lp(lp) ;
// on relance le processus sur chacun des deux segments qu'on obtient grace
a la nouvelle solution
resolution_recursive_lpsolve(xA, z1_A, z2_A,
xC, z1_C, z2_C,
E, p1, p2, w, cap, n) ;
resolution_recursive_lpsolve(xC, z1_C, z2_C,
xB, z1_B, z2_B,
E, p1, p2, w, cap, n) ;
}
return ret ;
}
int algo_lpsolve(int n, int cap, vector<int>* p1, vector<int>* p2, vector<int>*
w){
if((p1==NULL)||(p2==NULL)||(w==NULL)) {
cout << "pointeur nul en entrée de \n" ;
exit(-1) ;
}
vector<vector<int>*> E ; // pour stocker les solutions
// étape 1 : résolution lexicographique (z1,z2) puis (z2,z1) => deux
solutions x1 et x2 dans E
resolution_lexicographique_lpsolve(n, cap, p1, p2, w, &E) ; // calcul de x1
resolution_lexicographique_lpsolve(n, cap, p2, p1, w, &E) ; // calcul de x2
//affichage(*(E.front()), "x1") ;
//cout << "z1(x1) = " << evaluation(p1, E.front()) << " z2(x1) = " <<
evaluation(p2, E.front()) << "\n" ;
//affichage(*(E.back()), "x2") ;
//cout << "z1(x2) = " << evaluation(p1, E.back()) << " z2(x2) = " <<
evaluation(p2, E.back()) << "\n" ;
// étape 2 : calcul de toutes les solutions efficaces entre x(1) et x(2) =>
N solutions dans E
resolution_recursive_lpsolve(E.front(), evaluation(p1, E.front()),
evaluation(p2, E.front()),
E.back(), evaluation(p1, E.back()),
evaluation(p2, E.back()),
&E, p1, p2, w, cap, n) ;
cout << "# il y a " << E.size() << " solutions\n" ;
vector<vector<int>*>::iterator it ;
for(it = E.begin() ; it != E.end() ; it ++) {
cout << evaluation(p1, (*it)) << " " << evaluation(p2, (*it)) << "\n" ;
}
}
int ajoute_contrainte_glpk(LPX* lp, int Ncol, vector<int>* coef, int type_cont,
double borne_inf,double borne_sup) {
if((lp==NULL) || (coef==NULL)) {
cout << "pointeur nul en entrée de ajoute_contrainte\n" ;
exit(-1) ;
}
int *colno = NULL, j, i, ret = 0;
double *row = NULL;
// on cree de l'espace pour une ligne
colno = (int *) malloc((Ncol+1) * sizeof(*colno));
row = (double *) malloc((Ncol + 1) *sizeof(*row));
if((colno == NULL) || (row == NULL))
ret = 2;
if(ret == 0) {
//set_add_rowmode(lp, TRUE); /* makes building the model faster if it is
done rows by row */
vector<int>::iterator it ;
it = coef->begin() ;
for(j = 1; j <= Ncol; j ++) {
if(it == coef->end()) {
cout << "Problème dans les données : moins de " << Ncol << "
entiers dans w\n" ;
exit(-1) ;
}
colno[j] = j; /* column par colonne*/
row[j] = (*it) ;
it ++ ;
}
// on ajoute la ligne au solveur
//if(!add_constraintex(lp, j, row, colno, type_cont, borne))
// ret = 3;
i = lpx_add_rows(lp,1) ;
lpx_set_mat_row(lp,i,Ncol,colno,row);
lpx_set_row_bnds(lp,i,type_cont,borne_inf,borne_sup) ;
//set_add_rowmode(lp, FALSE); /* rowmode should be turned off again when
done building the model */
}
/* free allocated memory */
if(row != NULL)
free(row);
if(colno != NULL)
free(colno);
return ret ;
}
LPX* init_modele_glpk(int n, vector<int>* w, int cap) {
if((w==NULL)) {
cout << "pointeur nul en entrée de init_model_lpsolve\n" ;
exit(-1) ;
}
LPX *lp;
int Ncol, *colno = NULL, j, ret = 0;
double *row = NULL;
vector<int>::iterator it ;
// Construction du modèle sans objectif
Ncol = n; // il y a n variables
//lp = make_lp(0, Ncol);
lp = lpx_create_prob() ;
lpx_add_cols(lp,n) ;
//set_verbose(lp, IMPORTANT); // afficher uniquement les messages importants
//set_maxim(lp); // le but est de trouver un maximum
lpx_set_obj_dir(lp,LPX_MAX) ;
lpx_set_class(lp,LPX_MIP) ;
lpx_set_int_parm(lp, LPX_K_MSGLEV, 0) ;// pas de massages
lpx_set_int_parm(lp, LPX_K_PRESOL, 1) ;// prétraitement
if(lp == NULL)
ret = 1; // quand on n'a pas pas pu construire de modèle
// On construit le modèle ligne par ligne
// Au départ il y a n colonnes mais 0 ligne
if(ret == 0) {
// on nomme les variables, ce n'est pas obligatoire mais ca peut etre utile
pour le débogage
char* var_name;
for(int i = 0 ; i < Ncol ; i++) {
var_name = (char*) malloc(5) ;
sprintf(var_name, "x%d", i) ;
lpx_set_col_name(lp, i+1, var_name);
// nos variables sont binaires
lpx_set_col_kind(lp,i+1,LPX_IV) ;
lpx_set_col_bnds(lp, i+1, LPX_DB, 0.0, 1.0);
//set_binary(lp, i+1, TRUE) ;
}
//cout << "il y a " << lpx_get_num_bin(lp) << " variables binaires\n" ;
}
// on construit la contrainte des poids des objets
ret = ajoute_contrainte_glpk(lp, n, w, LPX_UP, 0.0, cap) ;
// lpx_adv_basis(lp) ;
//printf("contrainte ajoutee\n") ;
if(ret == 0)
return lp ;
else {
cout << "echec de l'initialisation du solveur\n" ;
exit(-1) ;
}
}
int set_objectif_glpk(LPX* lp, vector<int>* coef, int Ncol) {
if((lp==NULL)||(coef==NULL)) {
cout << "pointeur nul en entrée de set_objectif\n" ;
exit(-1) ;
}
int *colno = NULL, j, ret = 0,k = 0;
double *row = NULL;
// on cree de l'espace pour une ligne
colno = (int *) malloc(Ncol * sizeof(*colno));
row = (double *) malloc(Ncol * sizeof(*row));
if((colno == NULL) || (row == NULL))
ret = 2;
if(ret == 0) {
vector<int>::iterator it = coef->begin() ;
for(j = 0 ; j < Ncol ; j ++) {
if(it == coef->end()) {
cout << "Problème dans les données : moins de " << Ncol << "
entiers dans p_prim\n" ;
exit(-1) ;
}
colno[j] = j+1; /* first column */
row[j] = (*it);
it ++ ;
}
// on mémorise
//if(!set_obj_fnex(lp, j, row, colno))
// ret = 4;
for (k=0;k<Ncol;k++)
{
lpx_set_obj_coef(lp,k+1,row[k]);
}
}
/* free allocated memory */
if(row != NULL)
free(row);
if(colno != NULL)
free(colno);
return ret ;
}
int get_solution_glpk(LPX* lp, int n, vector<int>* sol, double* max){
if((lp==NULL)||(sol==NULL)||(max==NULL)) {
cout << "pointeur nul en entrée de get_solution\n" ;
exit(-1) ;
}
int Ncol = n, *colno = NULL, j, k,ret = 0;
double *row = NULL;
/*
if(get_solutioncount(lp) > 1) {
cout << "il y a plus d'une solution !" ;
exit(-1) ;
}
*/
cout << "dans get_solution_glpk\n" ;
sol->clear() ;
(*max) = lpx_mip_obj_val(lp);
// on cree de l'espace pour une ligne
colno = (int *) malloc(Ncol * sizeof(*colno));
row = (double *) malloc(Ncol * sizeof(*row));
if((colno == NULL) || (row == NULL))
ret = 2;
for(k = 0; k < Ncol; k++) {
row[k] = lpx_mip_col_val(lp,k+1) ;
}
for(j = 0; j < Ncol; j++) {
if(row[j] != 0) {
// mémoriser la solution : les indices des variables à un
sol->push_back(j) ;
}
}
/* free allocated memory */
if(row != NULL)
free(row);
if(colno != NULL)
free(colno);
return ret ;
}
int resolution_lexicographique_glpk(int n, int cap, vector<int>* p_prim,
vector<int>* p_sec,
vector<int>* w, vector<vector<int>*>* E)
{
if((p_prim==NULL)||(p_sec==NULL)||(w==NULL)||(E==NULL)) {
cout << "pointeur nul en entrée de resolution_lexicographique_lpsolve\n"
;
exit(-1) ;
}
LPX *lp; // variables utiles pour lpsolve
int Ncol = n, *colno = NULL, j, ret = 0;
double *row = NULL;
vector<int>::iterator it ;
vector<int>* sol_coef_a_un = new vector<int>() ;
//printf("initiation du modele\n");
lp = init_modele_glpk(n, w, cap) ;
// a) obtenir la solution optimale selon le premier objectif à maximiser
// on définit l'objectif primaire
ret = set_objectif_glpk(lp, p_prim, Ncol) ;
lpx_print_prob(lp, "prob.txt") ;
cout << "avant la résolution ret vaut " << ret << "\n" ;
if(ret == 0) {
cout << "avant z_prim\n";
//write_LP(lp, stdout);
// on calcule la solution optimale (pour le premier objectif seul)
//cout << "Calcul sol optimale\n" ;
if(lpx_get_obj_dir(lp) != LPX_MAX) {
cout << "houston ? we have a problem\n" ;
}
else {
cout << "on est sur un max\n" ;
}
ret = lpx_simplex(lp) ;
cout << "le solveur a renvoye " << ret << "\n" ;
//ret = lpx_integer(lp);
if(ret == LPX_E_OK)
ret = 0;
else
ret = 5;
}
double max_prim = -1 ;
cout << "après la résolution ret vaut " << ret << "\n" ;
if(ret == 0) {
// on récupère la solution et le maximum atteint par l'objectif
primaire
ret = get_solution_glpk(lp, n, sol_coef_a_un, &max_prim) ;
}
cout << "max_prim = " <<max_prim << "\n" ;
affichage(*sol_coef_a_un, "sol prim") ;
// libère la mémoire
lpx_delete_prob(lp) ;
// b) chercher à améliorer l'autre objectif SANS perdre l'optimalité
sur le premier
// i) ajouter une contrainte imposant de ne pas descendre sur le premier
objectif
lp = init_modele_glpk(n, w, cap) ; // on réinitialise un modèle (bug
sinon)
ret = ajoute_contrainte_glpk(lp, n, p_prim, LPX_UP,0.0, max_prim) ;
// ii) définir l'objectif secondaire
ret = set_objectif_glpk(lp, p_sec, Ncol) ;
if(ret == 0) {
cout << "avant z_sec\n";
//write_LP(lp, stdout);
// on calcule la solution optimale (pour le deuxième objectif sous le
contrainte de ne pas faire chuter le premier)
ret = lpx_simplex(lp) ;
//ret = lpx_integer(lp);
if(ret == LPX_E_OK)
ret = 0;
else
ret = 5;
}
// on sauvegarde la solution
ret = get_solution_glpk(lp, n, sol_coef_a_un, &max_prim) ;
if(ret == 0) {
E->push_back(sol_coef_a_un) ;
}
lpx_delete_prob(lp) ;
return ret ;
}
int resolution_recursive_glpk(vector<int>* xA, int z1_A, int z2_A,
vector<int>* xB, int z1_B, int z2_B,
vector<vector<int>*>* E,
vector<int>* p1, vector<int>* p2, vector<int>*
w, int cap, int n) {
//cout << "\nresol_recurs avec \n" ;
//affichage(*xA, "xA") ;
//cout << "z1(xA) = " << z1_A << " z2(xA) = " << z2_A << "\n" ;
//affichage(*xB, "xB") ;
//cout << "z1(xB) = " << z1_B << " z2(xB) = " << z2_B << "\n" ;
if((xA==NULL)||(xB==NULL)||(E==NULL)||(p1==NULL)||(p2==NULL)||(w==NULL)) {
cout << "pointeur nul en entrée de \n" ;
exit(-1) ;
}
int ret ;
LPX* lp = init_modele_glpk(n, w, cap) ;
// on ajoute les contraintes de ne pas aller en dessous de z2(xA) ni en
dessous de z1(xB) => éviter de boucler
ret = ajoute_contrainte_glpk(lp, n, p1, LPX_UP, 0.0,z1_B + 1) ; // attention
: si pas +1 on boucle
ret = ajoute_contrainte_glpk(lp, n, p2, LPX_UP,0.0, z2_A + 1) ;
// on pose : l1 = z2(xB) - z2(xA) et l2 = z1(xA) - z1(xB) et on calcule la
solution optimale
// de max{l1 * z1(x) + l2 * z2(x)} => direction de recherche othogonale au
vecteur xAxB
int l1 = z2_B - z2_A ;
int l2 = z1_A - z1_B ;
// on crée le vecteur des coefs du nouveau polynome objectif
vector<int>::iterator it1 = p1->begin(), it2 = p2->begin() ;
vector<int>* p_ortho = new vector<int>() ;
for(int j = 0 ; j < n ; j++) {
if((it1 == p1->end()) || (it2 == p2->end())){
cout << "Problème dans les données : moins de " << n << " entiers
dans p1 ou p2\n" ;
exit(-1) ;
}
p_ortho->push_back(l1 * (*it1) + l2 * (*it2)) ;
it1 ++ ;
it2 ++ ;
}
// on met en place l'objectif
ret = set_objectif_glpk(lp, p_ortho, n) ;
// on résoud
if(ret == 0) {
ret = lpx_integer(lp);
if(ret == LPX_E_OK)
ret = 0;
else
ret = 5;
}
if((ret == 0) /*&& (get_solutioncount(lp) != 0)*/){
double max ;
vector<int>* xC = new vector<int>() ;
// on sauvegarde la solution
ret = get_solution_glpk(lp, n, xC, &max) ;
int z1_C = evaluation(p1, xC),
z2_C = evaluation(p2, xC);
//cout << "nouvelle solution :\n" ;
//affichage(*xC, "xC") ;
//cout << "z1(xC) = " << z1_C << " z2(xC) = " << z2_C << "\n" ;
// on mémorise la nouvelle solution
E->push_back(xC) ;
// on efface avant le rappel pour ne pas saturer la mémoire
lpx_delete_prob(lp) ;
// on relance le processus sur chacun des deux segments qu'on obtient grace
a la nouvelle solution
resolution_recursive_glpk(xA, z1_A, z2_A,
xC, z1_C, z2_C,
E, p1, p2, w, cap, n) ;
resolution_recursive_glpk(xC, z1_C, z2_C,
xB, z1_B, z2_B,
E, p1, p2, w, cap, n) ;
}
return ret ;
}
int algo_glpk(int n, int cap, vector<int>* p1, vector<int>* p2, vector<int>* w){
if((p1==NULL)||(p2==NULL)||(w==NULL)) {
cout << "pointeur nul en entrée de \n" ;
exit(-1) ;
}
vector<vector<int>*> E ; // pour stocker les solutions
// étape 1 : résolution lexicographique (z1,z2) puis (z2,z1) => deux
solutions x1 et x2 dans E
resolution_lexicographique_glpk(n, cap, p1, p2, w, &E) ; // calcul de x1
//resolution_lexicographique_glpk(n, cap, p2, p1, w, &E) ; // calcul de x2
//affichage(*(E.front()), "x1") ;
//cout << "z1(x1) = " << evaluation(p1, E.front()) << " z2(x1) = " <<
evaluation(p2, E.front()) << "\n" ;
//affichage(*(E.back()), "x2") ;
//cout << "z1(x2) = " << evaluation(p1, E.back()) << " z2(x2) = " <<
evaluation(p2, E.back()) << "\n" ;
// étape 2 : calcul de toutes les solutions efficaces entre x(1) et x(2)
=> N solutions dans E
/*
resolution_recursive_glpk(E.front(), evaluation(p1, E.front()),
evaluation(p2,E.front()),
E.back(), evaluation(p1, E.back()), evaluation(p2,
E.back()),
&E, p1, p2, w, cap, n) ;*/
cout << "il y a " << E.size() << " solutions\n" ;
vector<vector<int>*>::iterator it ;
for(it = E.begin() ; it != E.end() ; it ++) {
cout << evaluation(p1, (*it)) << " " << evaluation(p2, (*it)) << "\n" ;
}
}
int main( int argc, char ** argv ) {
// créer trois vecteur de donnees d'entiers vide
vector<int> vector_obj1;
vector<int> vector_obj2;
vector<int> vector_cont;
int n ;
int w ;
string fic ;
// cout << "argc = " << argc << "\n" ;
if(argc <= 1) {
cout << "# Pas de fichier en paramètre, utilisant par défaut
2KP50-11.DAT\n" ;
fic = "2KP50-11.DAT" ;
}
else {
fic = string(argv[1]) ;
}
if(charger_fichier(fic, &vector_obj1, &vector_obj2, &vector_cont, &n, &w)) {
/*
cout << "n = " << n << '\n';
cout << "w = " << w << '\n';
affichage(vector_obj1, "Objectif 1");
affichage(vector_obj2, "Objectif 2");
affichage(vector_cont, "contrainte 1");
*/
cout << "# Avec glpk :\n" ;
//demo_glpk() ;
algo_glpk(n, w,&vector_obj1, &vector_obj2, &vector_cont) ;
//cout << "# " << fic << " avec lp_solve :\n" ;
//algo_lpsolve(n, w, &vector_obj1, &vector_obj2, &vector_cont) ;
//demo_lpsolve() ;
}
else {
cout << fic << " n'a pas pu etre ouvert\n" ;
}
}
| [Prev in Thread] |
Current Thread |
[Next in Thread] |
- [Help-glpk] oops... mistake !,
Nicoco BERGER <=