MeshKit  1.0
LPSolveClass.cpp
Go to the documentation of this file.
00001 #include "LPSolveClass.hpp"
00002 #include "lp_lib.h"
00003 #include <algorithm>
00004 
00005 namespace MeshKit
00006 {
00007 //---------------------------------------------------------------------------//
00008 // construction function for LPSolveClass class
00009 LPSolveClass::LPSolveClass()
00010 {
00011         
00012 
00013 }
00014 
00015 
00016 
00017 //---------------------------------------------------------------------------//
00018 // deconstruction function for LPSolveClass class
00019 LPSolveClass::~LPSolveClass()
00020 {
00021 
00022 }
00023 
00024 //setup the objective function  : cx
00025 void LPSolveClass::SetupObj(vector<double> left, double const_value)
00026 {
00027         obj_const = const_value;
00028         coefficients.resize(left.size());
00029         for (unsigned int i = 0; i < left.size(); i++)
00030                 coefficients[i] = left[i];
00031 }
00032 
00033 //setup the equal constraints   : Ax = b
00034 void LPSolveClass::SetupEqu(vector<vector<double> > left, vector<double> right){
00035 
00036         assert(left.size()==right.size());
00037         A_equ.resize(left.size());
00038         b_equ.resize(right.size());
00039         for (unsigned int i = 0; i < left.size(); i++){
00040                 A_equ[i].resize(left[i].size());
00041                 b_equ[i] = right[i];
00042                 for (unsigned int j = 0; j < left[i].size(); j++)
00043                         A_equ[i][j] = left[i][j];
00044         }
00045 }
00046 
00047 //setup the inequal constraints: always, Ax <= b
00048 void LPSolveClass::SetupInEqu(vector<vector<double> > left, vector<double> right){
00049 
00050         assert(left.size()==right.size());
00051         A_inequ.resize(left.size());
00052         b_inequ.resize(right.size());
00053         for (unsigned int i = 0; i < left.size(); i++){
00054                 A_inequ[i].resize(left[i].size());
00055                 b_inequ[i] = right[i];
00056                 for (unsigned int j = 0; j < left[i].size(); j++)
00057                         A_inequ[i][j] = left[i][j];
00058         }
00059 }
00060 //setup the const constraint
00061 void LPSolveClass::SetupConst(vector<int> right)
00062 {
00063         b_const.resize(right.size());
00064         for (unsigned int i = 0; i < right.size(); i++)
00065                 b_const[i] =right[i];
00066 }
00067 
00068 //Execute function
00069 int LPSolveClass::Execute()
00070 {
00071         /*
00072         std::cout << "---------------------------------\n";
00073         std::cout << "objective function\n";
00074         for (unsigned int i = 0; i < coefficients.size(); i++)
00075                 std::cout << coefficients[i] << "\t";
00076         std::cout << "\nConstant Value = " << obj_const << std::endl;
00077 
00078         std::cout << "---------------------------------\n";
00079         std::cout << "Equality Constraints\n";  
00080         for (unsigned int i = 0; i < A_equ.size(); i++){
00081                 //std::cout << "Row index = " << i << "\t\t";
00082                 for (unsigned int j = 0; j < A_equ[i].size(); j++)
00083                         std::cout << A_equ[i][j] << "\t";
00084                 std::cout << "\n";
00085         }
00086         std::cout << "b\n";
00087         for (unsigned int i = 0; i < b_equ.size(); i++)
00088                 std::cout << b_equ[i] << "\t";
00089         std::cout << "\n";
00090 
00091 
00092         std::cout << "---------------------------------\n";
00093         std::cout << "InEquality Constraints\n";        
00094         for (unsigned int i = 0; i < A_inequ.size(); i++){
00095                 //std::cout << "Row index = " << i << "\t\t";
00096                 for (unsigned int j = 0; j < A_inequ[i].size(); j++)
00097                         std::cout << A_inequ[i][j] << "\t";
00098                 std::cout << "\n";
00099         }
00100         std::cout << "b\n";
00101         for (unsigned int i = 0; i < b_inequ.size(); i++)
00102                 std::cout << b_inequ[i] << "\t";
00103         std::cout << "\n";
00104         */
00105 
00106         lprec *lp;
00107         int Ncol = coefficients.size(), *colno = NULL, j, ret = 0;
00108         REAL *row = NULL;
00109         
00110         /* We will build the model row by row
00111      So we start with creating a model with 0 rows and n columns */
00112 
00113         lp = make_lp(0, Ncol);
00114         if (lp == NULL)
00115                 ret = 1;/* couldn't construct a new model... */
00116                 
00117         if (ret == 0){
00118                 //let us name our variables
00119                 std::string s = "x";
00120                 for (int i = 0; i < Ncol; i++){
00121                         std::stringstream out;
00122                         out << i;
00123                         s = s + out.str();
00124                         char *cpy = new char[s.size()+1] ;
00125                         strcpy(cpy, s.c_str());                 
00126                         set_col_name(lp, i+1, cpy);
00127                 }
00128 
00129                 /* create space large enough for one row */
00130                 colno = (int *) malloc(Ncol * sizeof(*colno));
00131         row = (REAL *) malloc(Ncol * sizeof(*row));
00132                 if ((colno == NULL) || (row == NULL))
00133                 ret = 2;
00134         }
00135 
00136         set_add_rowmode(lp, TRUE);
00137         //add the equation constraints
00138         if (ret == 0){
00139                 /* makes building the model faster if it is done rows by row */
00140                 if (A_equ.size() > 0){
00141                         for (unsigned int i = 0; i < A_equ.size(); i++){//loop over the rows of equality constraints
00142                                 for (unsigned int j = 0; j < A_equ[i].size(); j++){//loop over the columns of equality constraints
00143                                         colno[j] = j+1;//add the j-th column to lpsolve
00144                                         row[j] = A_equ[i][j];
00145                                 }
00146                                 /* add the row to lpsolve */
00147                                 if(!add_constraintex(lp, A_equ[i].size(), row, colno, EQ, b_equ[i]))
00148                                         ret = 2;
00149                         }
00150                 }
00151         }
00152         
00153         //add the inequality constraints
00154         if (ret == 0){
00155                 /* makes building the model faster if it is done rows by row */
00156                 if (A_inequ.size() > 0){
00157                         for (unsigned int i = 0; i < A_inequ.size(); i++){//loop over the rows of inequality constraints
00158                                 for (unsigned int j = 0; j < A_inequ[i].size(); j++){//loop over the columns of inequality constraints
00159                                         colno[j] = j+1;//add the j-th column to lpsolve
00160                                         row[j] = A_inequ[i][j];
00161                                 }
00162                                 /* add the row to lpsolve */
00163                                 if(!add_constraintex(lp, A_inequ[i].size(), row, colno, LE, b_inequ[i]))
00164                                         ret = 3;
00165                         }
00166                 }
00167         }
00168 
00169         //add the const constraint      
00170         if (ret == 0){
00171                 if (b_const.size()>0){
00172                         for (unsigned int i = 0; i < b_const.size(); i++){
00173                                 if (b_const[i] > 0){
00174                                         for (unsigned int j = 0; j < b_const.size(); j++){
00175                                                 if (i == j){
00176                                                         colno[j] = j+1;//add the j-th column to lpsolve
00177                                                         row[j] = 1.0;                                           
00178                                                 }                               
00179                                                 else{
00180                                                         colno[j] = j+1;//add the j-th column to lpsolve
00181                                                         row[j] = 0.0;
00182                                                 }
00183                                         }
00184                                         if(!add_constraintex(lp, b_const.size(), row, colno, EQ, b_const[i]))
00185                                                 ret = 4;                
00186                                 }
00187                         }
00188                 }
00189         }
00190 
00191         //set the variables to be integer
00192         if (ret == 0){
00193                 for (int i = 0; i < Ncol; i++)
00194                         set_int(lp, i+1, TRUE);
00195         }
00196         
00197         /* rowmode should be turned off again when done building the model */
00198         set_add_rowmode(lp, FALSE);     
00199         //add the objective function
00200         if (ret == 0){
00201                 //set the objective function
00202                 for (unsigned int i = 0; i < coefficients.size(); i++){
00203                         colno[i] = i+1;
00204                         row[i] = coefficients[i];
00205                 }
00206                 //set the objective in lpsolve
00207                 if(!set_obj_fnex(lp, coefficients.size(), row, colno))
00208                 ret = 4;
00209         }
00210 
00211         //set the objective to minimize
00212         if (ret == 0){
00213                 set_minim(lp);
00214 
00215                 /* just out of curioucity, now show the model in lp format on screen */
00216         /* this only works if this is a console application. If not, use write_lp and a filename */
00217         write_LP(lp, stdout);
00218 
00219                 /* I only want to see important messages on screen while solving */
00220         set_verbose(lp, IMPORTANT);
00221 
00222         /* Now let lpsolve calculate a solution */
00223         ret = solve(lp);
00224         if(ret == OPTIMAL)
00225                 ret = 0;
00226         else
00227                 ret = 5;
00228         }
00229 
00230         //get some results
00231         if (ret == 0){
00232                 /* a solution is calculated, now lets get some results */
00233 
00234         /* objective value */
00235         std::cout << "Objective value: " << get_objective(lp) << std::endl;
00236 
00237                 /* variable values */
00238         get_variables(lp, row);
00239 
00240                 /* variable values */
00241                 variables.resize(Ncol);
00242                 for(j = 0; j < Ncol; j++)
00243                         variables[j] = row[j];
00244 
00245                 /* we are done now */
00246         }
00247         else{
00248                 std::cout << "The optimal value can't be solved for linear programming, please check the constraints!!\n";
00249                 exit(1);
00250 
00251         }
00252                 
00253         
00254         std::cout << "print the result\t # of line segments is \n";
00255         for (int i = 0; i < Ncol; i++)
00256                 std::cout << "index = " << i << "\t# = " << variables[i] << std::endl;
00257 
00258         /* free allocated memory */
00259         if(row != NULL)
00260         free(row);
00261         if(colno != NULL)
00262         free(colno);
00263 
00264         /* clean up such that all used memory by lpsolve is freed */
00265         if (lp != NULL)
00266                 delete_lp(lp);
00267 
00268         return ret;
00269 }
00270 
00271 void LPSolveClass::GetVariables(vector<int> &var){
00272         var.resize(variables.size());
00273         for(unsigned int i = 0; i < variables.size(); i++)
00274                 var[i] = variables[i];
00275 
00276 }
00277 
00278 }
00279 
00280 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines