/home/coin/SVN-release/OS-2.4.0/Couenne/src/cut/sdpcuts/sdpcuts.cpp

Go to the documentation of this file.
00001 /* $Id: sdpcuts.cpp 508 2011-02-15 21:52:44Z pbelotti $
00002  *
00003  * Name:    sdpcuts.cpp
00004  * Author:  Andrea Qualizza
00005  * Purpose: 
00006  *
00007  * This file is licensed under the Eclipse Public License (EPL)
00008  */
00009 
00010 #include <stdio.h>
00011 #include <stdlib.h>
00012 
00013 #include <sdpcuts.hpp>
00014 
00015 #include <misc_util.hpp>
00016 #include <OsiXxxSolverInterface.hpp>
00017 #include <populate.hpp>
00018 #include <CutGen.hpp>
00019 #include <quadratic_cuts_check.hpp>
00020 #include <linquad_cuts.hpp>
00021 #include <orthocut.hpp>
00022 #include <tracer.hpp>
00023 
00024 #ifdef USE_INTERIOR_POINT
00025 #include <OsiCpxSolverInterface.hpp>
00026 #include <cplex.h>
00027 #endif
00028 
00029 
00030 #if 0
00031 int main (int argc, const char **argv) {
00032         if (argc < 2) {
00033                 printf("Missing argument [mps file]\n");
00034                 return 1;
00035         }
00036 
00037         //print current configuration
00038         print_ifdefs(stdout);
00039 
00040         // determine problem name
00041         char name[256];
00042         char *name_pos = strrchr(const_cast <char *> (argv[1]), '/');
00043         if(name_pos != NULL)
00044                 strcpy(name, &(name_pos[1]));
00045         else
00046                 strcpy(name, argv[1]);
00047         char *last_dot_pos = strrchr(name, '.');
00048         if(last_dot_pos !=NULL)
00049                 last_dot_pos[0] = '\0';
00050 
00051         // initialize random generator
00052         struct timeval tv;
00053         gettimeofday (&tv, NULL);
00054         srand48 (tv.tv_usec);
00055 
00056         Timer globaltimer;
00057         globaltimer.start();
00058 
00059 #ifdef F_RES
00060         FILE *f_res = open_f_res();
00061 #endif /* F_RES */
00062 #ifdef SHORT_F_RES
00063         FILE *short_f_res = open_short_f_res();
00064 #endif /* SHORT_F_RES */
00065 
00066         OsiXxxSolverInterface si;
00067 
00068 #if (!defined TRACE_ALL)
00069         si.messageHandler()->setLogLevel(0);
00070 #endif
00071 
00072         si.setDblParam(OsiPrimalTolerance,LP_TOLERANCE);
00073         si.setDblParam(OsiDualTolerance,LP_TOLERANCE);
00074 
00075         int n;                  // number of x_i variables
00076         int t;                  // number of y_i variables
00077         int origCons;           // number of original constraints
00078         double *b;              // obj coefficients for the x_i variables
00079         double *c;              // obj coefficients for the y_i variables
00080         double **Q;             // obj coefficients for the x_i_j variables
00081         double constant;        // obj function additive constant
00082         double **origMat;       // original problem matrix
00083         double *origRhs;        // original constraints' RHS
00084         char *origSense;        // original constraints' sense
00085         double *xlb;            // lower bounds on x_i vars
00086         double *xub;            // upper bounds on x_i vars
00087         double *ylb;            // lower bounds on y_i vars
00088         double *yub;            // upper bounds on y_i vars
00089 
00090         double objValue;
00091 
00092         int status;
00093 
00094         // read problem file
00095         status = populateProblem (argv[1],&n, &t, &origCons, &b, &c, &Q, &constant, &origMat, 
00096                                         &origRhs, &origSense, &xlb, &xub, &ylb, &yub,&si);
00097 
00098         if (status) 
00099                 exit(1);
00100 
00101 
00102         Tracer *tracer = new Tracer(name,n,n*(n+1)/2,t);
00103 
00104 #ifdef LINQUAD_BOUNDS_CUTS
00105         OsiCuts cs_diagbounds;
00106         linQuadCutGenOriginalBounds(xlb , xub , n , cs_diagbounds, tracer);
00107         si.applyCuts(cs_diagbounds);
00108 #endif
00109 
00110 
00111         printf("Problem:                     %-20s\tx_card=%-3d Xext_card=%-6d y_card=%-3d\n",name,n,n*(n+1)/2,t);
00112         printf("Objective function constant: %.8f  (automatically added)\n",constant);
00113 
00114         CutGen  cg (n, t, origCons, constant, 
00115                         b, c, (const double**) Q, 
00116                         (const double**) origMat, origRhs, origSense, xlb, xub, ylb, yub, 
00117                         &si, &globaltimer, tracer);
00118 
00119 #ifdef MAX_CUTS_PER_ITER
00120         cg.set_max_nb_cuts(MAX_CUTS_PER_ITER);
00121 #else
00122         cg.set_max_nb_cuts(100000);
00123 #endif
00124 
00125 #define OBJ_HISTORY_SIZE 10000
00126         double *obj_history = new double[OBJ_HISTORY_SIZE];
00127         for (int i=0;i<OBJ_HISTORY_SIZE;i++)
00128                 obj_history[i] = COIN_DBL_MAX;
00129 
00130 #ifdef DELETE_INACTIVE_CUTS
00131         int *del_row_ind = new int[1000000];
00132 #ifdef FANCY_CUT_ACTIVITY
00133         int card_act = 0;
00134         int *activity = new int[1000000];
00135         int *gen_iter = new int[1000000];
00136 #endif
00137 #endif
00138 
00139         int niter = 0, ncuts = 0, tot_gen_cuts = 0, tot_del_cuts = 0;
00140         bool do_exit = false;
00141 
00142 
00143         printf("Initialization time:         %.2f\n",globaltimer.time());
00144         printf("\n");
00145 
00146         globaltimer.start();
00147 
00148 #ifdef USE_INTERIOR_POINT
00149         si.XxxInitialSolveBaropt();
00150 #else
00151         si.initialSolve();
00152 #endif
00153         solver_status(&si);
00154         cg.updateSol();
00155 
00156         objValue = si.getObjValue() + constant;
00157         obj_history[0] = objValue;
00158 
00159 #ifdef EXIT_IMPROVEMENT_ITER
00160         double *last_lp_val = new double [EXIT_IMPROVEMENT_ITER];
00161 #ifdef EXIT_TAILING_OFF_VALUE
00162         double tailing_off_value = EXIT_TAILING_OFF_VALUE; 
00163         for(int i=0; i<EXIT_IMPROVEMENT_ITER; i++) {
00164                 last_lp_val[i] = objValue + tailing_off_value;
00165         }
00166 #endif
00167 #ifdef EXIT_TAILING_OFF_PERC
00168         for(int i=0; i<EXIT_IMPROVEMENT_ITER; i++) {
00169                 last_lp_val[i] = objValue + (fabs(objValue)*EXIT_TAILING_OFF_PERC/100);
00170         }       
00171 #endif
00172 #endif
00173 
00174 #ifdef DELETE_INACTIVE_CUTS
00175         int del_init_nrows;
00176         del_init_nrows = si.getNumRows(); //never cut these ones (original and rlt (if separation of the latter is not allowed))
00177 #endif
00178         
00179         print_current_sol(niter, globaltimer.time(), tot_gen_cuts, tot_gen_cuts-tot_del_cuts, 
00180                                 objValue, cg.bestObj(),cg.currObj());
00181 #ifdef F_RES
00182         globaltimer.pause();
00183         print_file_current_sol(f_res, name, niter,
00184                                 globaltimer.time(), tot_gen_cuts, tot_gen_cuts-tot_del_cuts, 
00185                                 objValue, cg.bestObj(), cg.currObj());
00186         globaltimer.restore();
00187 #endif
00188 
00189         tracer->setMainBound(objValue);
00190         tracer->setMainIterationTime(globaltimer.time());
00191         tracer->setMainLPTime(globaltimer.time());
00192         double time_at_previous_iteration = globaltimer.time();
00193         globaltimer.pause();
00194 
00195         //determining active rows
00196         const double *dualvars = si.getRowPrice();
00197         int active_rows = 0;
00198         for (int i=0;i<si.getNumRows();i++)
00199                 if (fabs(dualvars[i]) < 2 * LP_TOLERANCE)
00200                         active_rows++;
00201         tracer->setMainActiveCuts(active_rows);
00202         tracer->setMainAddedCuts(0);
00203         tracer->setMainTotalCuts(si.getNumRows());
00204         tracer->setMainDeletedCuts(0);
00205         tracer->setMainTotalEigendecompositions(0);
00206 
00207         tracer->setSDPNumNegativeEV(0);
00208         tracer->setSDPMostNegativeEV(0.0);
00209         tracer->setSDPCutsTime(0.0);
00210         tracer->setSDPCutsTotalCuts(0);
00211 
00212 
00213         tracer->setSparsifyTime(0.0);
00214         tracer->setSparsifyTotalCuts(0);
00215         tracer->setSparsifyDuplicatedCuts(0);
00216         tracer->setSparsifyWiseDecompositions(0);
00217         tracer->addSparsifyNz(0);
00218         tracer->addSparsifySingleColumnSparsity(0);
00219         tracer->addSparsifyColumnPairSparsity(0);
00220         tracer->addSparsifyTop20PercCutsViolation(0.0);
00221 
00222         tracer->setOrthocutTime(0.0);
00223         tracer->setOrthocutTotalCuts(0);
00224 
00225         // setLinquad[...] in linQuadCutGenOriginalBounds()
00226 
00227         tracer->setDisjunctiveCutsTime(0.0);
00228         tracer->setDisjunctiveCutsTotalCuts(0);
00229 
00230         // setHeuristics[...] in cg.updateSol()
00231 
00232         globaltimer.restore();
00233 
00234 #ifdef CHECK_QUADRATIC_CUTS
00235         globaltimer.pause();
00236         QuadraticCuts qc(n,si.getColSolution(),tracer);
00237         globaltimer.restore();
00238 #endif
00239 
00240         do {
00241                 niter++;
00242                 tot_del_cuts = 0;
00243 
00244 
00245 #ifdef EXIT_ON_ITER
00246                 if(niter > EXIT_ON_ITER) {
00247                         do_exit = true;
00248                         printf("Exit: Iteration limit [%d] exceeded\n",EXIT_ON_ITER);
00249                 }
00250 #endif
00251 
00252 #ifdef EXIT_ON_TIME
00253                 if(globaltimer.time() > EXIT_ON_TIME) {
00254                         do_exit = true;
00255                         printf("Exit: Time limit [%g] exceeded\n",EXIT_ON_TIME);
00256                 }
00257 #endif
00258                 
00259                 if (do_exit)
00260                         break;
00261 
00262                 tracer->newIter();
00263 
00264                 cg.setIter(niter);
00265 
00266                 OsiCuts cs;
00267                 cg.generateCuts (si, cs);
00268 
00269                 // the slow way to add cuts:
00270                 // si.applyCuts (cs);
00271 
00272                 // the fast way :)
00273                 int size_cs = cs.sizeRowCuts();
00274                 const OsiRowCut **newRowCuts = new const OsiRowCut * [size_cs];
00275                 for(int i=0; i<size_cs; i++) {
00276                         newRowCuts[i] = &cs.rowCut(i);
00277                 }
00278                 si.applyRowCuts(size_cs, newRowCuts);
00279                 delete[] newRowCuts; 
00280 
00281 
00282 
00283 
00284                 Timer lp_timer;
00285                 lp_timer.start();
00286 #ifdef USE_INTERIOR_POINT
00287                 si.XxxResolveBaropt();
00288 #else
00289                 si.resolve();
00290 #endif
00291                 solver_status(&si);
00292                 lp_timer.pause();
00293 
00294                 objValue = si.getObjValue() + constant;
00295                 cg.updateSol();
00296 
00297 #ifdef WRITE_ITER_LP
00298                 globaltimer.pause();
00299 
00300                 si.writeLp("lastiter");
00301                 globaltimer.restore();
00302 #endif
00303 
00304 #ifdef CHECK
00305                 globaltimer.pause();
00306                 const double *sol  = si.getColSolution ();
00307                 int status = feasibility_check(n,t,origCons,sol,(const double**)origMat,(const double*)origRhs,(const char*)origSense,(const double*)xlb,(const double*)xub,(const double*)ylb,(const double*)yub);
00308                 if (status)
00309                         printf("ERROR: Infeasible Outer-approximation solution [status=%d]\n",status);
00310 
00311                 // check that heuristic value is BELOW the current solution value
00312                 if (cg.currObj() - LP_TOLERANCE > objValue)
00313                         printf("ERROR: heuristic solution value > current solution value\n");
00314                 globaltimer.restore();
00315 #endif
00316 
00317                 globaltimer.pause();
00318 
00319                 tot_gen_cuts += ncuts = cs.sizeRowCuts ();
00320 
00321                 print_current_sol(niter,globaltimer.time(), tot_gen_cuts, ncuts,
00322                                         objValue, cg.bestObj(), cg.currObj());
00323 #ifdef F_RES
00324                 print_file_current_sol(f_res, name, niter,globaltimer.time(), tot_gen_cuts, ncuts,
00325                                         objValue, cg.bestObj(), cg.currObj());
00326 #endif
00327                 tracer->setMainIterationTime(fabs(globaltimer.time() - time_at_previous_iteration));
00328                 time_at_previous_iteration = globaltimer.time();
00329                 tracer->setMainLPTime(lp_timer.time());
00330                 tracer->setMainBound(objValue);
00331                 //determining active rows
00332                 const double *dualvars = si.getRowPrice();
00333                 active_rows = 0;
00334                 for (int i=0;i<si.getNumRows();i++)
00335                         if (fabs(dualvars[i]) < 2 * LP_TOLERANCE)
00336                                 active_rows++;
00337                 tracer->setMainActiveCuts(active_rows);
00338                 tracer->setMainAddedCuts(cs.sizeRowCuts ());
00339                 tracer->setMainTotalCuts(si.getNumRows());
00340                 tracer->setMainDeletedCuts(0); // set later
00341                 globaltimer.restore();
00342 
00343 
00344 #ifdef CHECK_QUADRATIC_CUTS
00345                 globaltimer.pause();
00346                 qc.refresh(si.getColSolution());
00347                 globaltimer.restore();
00348 #endif
00349 
00350                 if (cs.sizeRowCuts () == 0) {
00351                         do_exit = true;
00352                         printf("Exit: No cut generated\n");
00353                 }
00354 
00355 #ifdef EXIT_IMPROVEMENT_ITER
00356                 int ind_iter = niter % EXIT_IMPROVEMENT_ITER;
00357 #ifdef EXIT_TAILING_OFF_VALUE
00358                 if ((last_lp_val[ind_iter] < objValue  + tailing_off_value) 
00359                         && (niter > EXIT_IMPROVEMENT_ITER)) {
00360                         do_exit = true;
00361                         printf("Exit: Solution improvement < %g in the last %d iterations\n"
00362                                 ,EXIT_TAILING_OFF_VALUE,EXIT_IMPROVEMENT_ITER);
00363                 } else
00364                         last_lp_val[ind_iter] = objValue;
00365 #endif
00366 #ifdef EXIT_TAILING_OFF_PERC
00367                 if((last_lp_val[ind_iter] - 
00368                                 (fabs(last_lp_val[ind_iter])*EXIT_TAILING_OFF_PERC/100) < objValue)
00369                         && (niter > EXIT_IMPROVEMENT_ITER)) {
00370                         do_exit = true;
00371                         printf("Exit: Solution improvement < %g%% in the last %d iterations\n"
00372                                 ,EXIT_TAILING_OFF_PERC,EXIT_IMPROVEMENT_ITER);
00373                 } else
00374                         last_lp_val[ind_iter] = objValue;
00375 #endif
00376 #endif
00377                 if (do_exit)
00378                         break;
00379 
00380 #ifdef DELETE_INACTIVE_CUTS
00381                 // consider only cuts from the previous iteration
00382                 int curr_nrows = si.getNumRows() - cs.sizeRowCuts();
00383                 const double *y = si.getRowPrice();
00384 #ifdef FANCY_CUT_ACTIVITY
00385                 int card_old_activity = card_act;
00386                 card_act = 0;
00387 #endif
00388                 int card_del_row_ind = 0;
00389                 for(int del_i = del_init_nrows; del_i<curr_nrows; del_i++) {
00390 
00391 #ifdef FANCY_CUT_ACTIVITY
00392                         int old_act = 0, old_gen = niter;
00393                         int old_pos = del_i - del_init_nrows;
00394                 
00395                         if(old_pos < card_old_activity) {
00396                                 old_act = activity[old_pos];
00397                                 old_gen = gen_iter[old_pos];
00398                         }
00399                                 
00400                         if(fabs(y[del_i]) < DELETE_INACTIVE_CUTS) {
00401                                 if(old_act > FANCY_CUT_ACTIVITY - 2) {
00402                                         del_row_ind[card_del_row_ind] = del_i;
00403                                         card_del_row_ind++;
00404                                 }
00405                                 else {
00406                                         activity[card_act] = old_act + 1;
00407                                         gen_iter[card_act] = old_gen;
00408                                         card_act++;
00409                                 }
00410                         }
00411                         else {
00412                                 activity[card_act] = 0;
00413                                 gen_iter[card_act] = old_gen;
00414                                 card_act++;
00415                         }
00416 #else
00417                         if(fabs(y[del_i]) < DELETE_INACTIVE_CUTS) {
00418                                 del_row_ind[card_del_row_ind] = del_i;
00419                                 card_del_row_ind++;
00420                         }
00421 #endif
00422                 }
00423 
00424 #if (defined CUT_DELETION_OBJ_PERC_IMPROVEMENT_PERC) && (defined CUT_DELETION_OBJ_PERC_IMPROVEMENT_ITER)
00425                 int historyIdx = (niter-CUT_DELETION_OBJ_PERC_IMPROVEMENT_ITER) % OBJ_HISTORY_SIZE;
00426                 if (historyIdx < 0)
00427                         historyIdx = 0;
00428                 double oldObjValue = 
00429                   obj_history[historyIdx];
00430 
00431                 if ((fabs(objValue - oldObjValue)*100/ fabs(oldObjValue))
00432                                 >= CUT_DELETION_OBJ_PERC_IMPROVEMENT_PERC)
00433                 {
00434                         si.deleteRows(card_del_row_ind, del_row_ind);
00435                         tot_del_cuts += card_del_row_ind;
00436                 }
00437 #else
00438                 si.deleteRows(card_del_row_ind, del_row_ind);
00439                 tot_del_cuts += card_del_row_ind;
00440 #endif
00441 #ifdef CPLEX
00442                 si.resolve();
00443 #endif
00444 #endif // DELETE_INACTIVE_CUTS
00445 
00446                 tracer->setMainDeletedCuts(tot_del_cuts);
00447 
00448                 obj_history[niter % OBJ_HISTORY_SIZE] = objValue;
00449 
00450         } while(do_exit == false);
00451 
00452 #ifdef SHORT_F_RES
00453         print_file_short_sol(short_f_res,name,niter,globaltimer.time(),
00454                          tot_gen_cuts,tot_gen_cuts-tot_del_cuts,objValue,cg.bestObj());
00455 #endif
00456 
00457         tracer->detailedReport();
00458         tracer->globalReport();
00459 
00460 
00461 
00462 #ifdef F_RES
00463         fclose(f_res);
00464 #endif
00465 
00466 #ifdef SHORT_F_RES
00467         fclose(short_f_res);
00468 #endif
00469 
00470 #ifdef EXIT_IMPROVEMENT_ITER
00471         delete [] last_lp_val;
00472 #endif
00473 #ifdef DELETE_INACTIVE_CUTS
00474         delete [] del_row_ind;
00475 #ifdef FANCY_CUT_ACTIVITY
00476         delete [] activity;
00477         delete [] gen_iter;
00478 #endif
00479 #endif
00480         delete [] obj_history;
00481 
00482         // freeing structures allocated by populate()
00483         delete [] b;
00484         delete [] c;
00485         for(int i=0;i<n;i++)
00486                 delete [] Q[i];
00487         delete [] Q;
00488         for (int i=0;i<origCons;i++)
00489                 delete [] origMat[i];
00490         delete [] origMat;
00491         delete [] origRhs;
00492         delete [] origSense;
00493         delete [] xlb;
00494         delete [] xub;
00495         delete [] ylb;
00496         delete [] yub;
00497 
00498         delete tracer;
00499 
00500         return 0;
00501 }
00502 #endif
00503 
00504 
00505 
00506 
00507 /***********************************************************************/
00508 void print_ifdefs(FILE *f) {
00509         fprintf(f, "Solver:        ");
00510 #ifdef CPLEX
00511         fprintf(f, " CPLEX");
00512 #endif
00513 #ifdef CLP
00514         fprintf(f, " CLP");
00515 #endif
00516 #ifdef USE_INTERIOR_POINT
00517         fprintf(f, " USE_INTERIOR_POINT");
00518 #endif
00519         
00520 
00521         fprintf(f, "\nTrace:         ");
00522 #ifdef TRACE
00523         fprintf(f, " TRACE");
00524 #endif
00525 #ifdef TRACE_ALL
00526         fprintf(f, " TRACE_ALL");
00527 #endif
00528 #ifdef TRACE_EIGENVALUES
00529         fprintf(f, " TRACE_EIGENVALUES");
00530 #endif
00531 #ifdef TRACE_USED_EIGENVECTORS
00532         fprintf(f, " TRACE_USED_EIGENVECTORS");
00533 #endif
00534 #ifdef TRACE_RLT_MPS
00535         fprintf(f, " TRACE_RLT_MPS");
00536 #endif
00537 #ifdef TRACE_HEUR
00538         fprintf(f, " TRACE_HEUR");
00539 #endif
00540 #ifdef TRACE_CUT_TIME
00541         fprintf(f, " TRACE_CUT_TIME");
00542 #endif
00543 #ifdef TRACE_SPARSIFY_CUT_DEPTH
00544         fprintf(f, " TRACE_SPARSIFY_CUT_DEPTH");
00545 #endif
00546 #ifdef TRACE_DISJUNCTIVE_CUTS
00547         fprintf(f, " TRACE_DISJUNCTIVE_CUTS");
00548 #endif
00549 #ifdef TRACE_ORTHOCUT
00550         fprintf(f, " TRACE_ORTHOCUT");
00551 #endif
00552 #ifdef CHECK_QUADRATIC_CUTS
00553         fprintf(f, " CHECK_QUADRATIC_CUTS");
00554 #endif
00555 #ifdef F_RES
00556         fprintf(f, " F_RES");
00557 #endif
00558 #ifdef SHORT_F_RES
00559         fprintf(f, " SHORT_F_RES");
00560 #endif
00561 
00562 
00563 
00564         fprintf(f, "\nExit:          ");        
00565 #ifdef EXIT_ON_ITER
00566         fprintf(f, " EXIT_ON_ITER=%d",EXIT_ON_ITER);
00567 #endif
00568 #ifdef EXIT_ON_TIME
00569         fprintf(f, " EXIT_ON_TIME=%.2f",EXIT_ON_TIME);
00570 #endif
00571 #if (!defined EXIT_ON_ITER) && (!defined EXIT_ON_TIME)
00572 #error EXIT_ON_TIME or EXIT_ON_ITER must be set
00573 #endif
00574 #ifdef EXIT_IMPROVEMENT_ITER
00575         fprintf(f, " EXIT_IMPROVEMENT_ITER=%d",EXIT_IMPROVEMENT_ITER);
00576 #ifdef EXIT_TAILING_OFF_VALUE
00577         fprintf(f, " EXIT_TAILING_OFF_VALUE=%e",EXIT_TAILING_OFF_VALUE);
00578 #endif
00579 #ifdef EXIT_TAILING_OFF_PERC
00580         fprintf(f, " EXIT_TAILING_OFF_PERC=%g%%",(double)EXIT_TAILING_OFF_PERC);
00581 #endif
00582 #if (!defined EXIT_TAILING_OFF_VALUE) && (!defined EXIT_TAILING_OFF_PERC) 
00583 #error EXIT_TAILING_OFF_VALUE or EXIT_TAILING_OFF_PERC must be set
00584 #endif
00585 #if (defined EXIT_TAILING_OFF_VALUE) && (defined EXIT_TAILING_OFF_PERC) 
00586 #error Only one between EXIT_TAILING_OFF_VALUE or EXIT_TAILING_OFF_PERC can be defined
00587 #endif
00588 #endif
00589 
00590         fprintf(f, "\nAlg. options:  ");
00591 #ifdef RLT_CUTS
00592         fprintf(f, " RLT_CUTS");
00593 #endif
00594 #ifdef RLT_SEPARATION
00595 #ifndef RLT_CUTS
00596 #error RLT_CUTS must be defined in order to have RLT_SEPARATION
00597 #endif
00598         fprintf(f, " RLT_SEPARATION");
00599 #endif
00600 #ifdef INCL_ORIG
00601         fprintf(f, " INCL_ORIG");
00602 #endif
00603 #ifdef SCALE_EIGENV
00604         fprintf(f, " SCALE_EIGENV");
00605 #endif
00606 #ifdef ONLY_NEG_EIGENV
00607         fprintf(f, " ONLY_NEG_EIGENV");
00608 #endif
00609 #ifdef ONLY_MOST_NEG
00610         fprintf(f, " ONLY_MOST_NEG");
00611 #ifdef ONLY_NEG_EIGENV
00612 #error ONLY_NEG_EIGENV and ONLY_MOST_NEG conflicts
00613 #endif
00614 #endif
00615 #ifdef SPARSIFY2
00616         fprintf(f, " SPARSIFY2");
00617 #endif
00618 #ifdef SPARSIFY
00619         fprintf(f, " SPARSIFY");
00620 #endif
00621 #ifdef SPARSIFY_MINOR_SDP_CUTS
00622         fprintf(f," SPARSIFY_MINOR_SDP_CUTS");
00623 #endif
00624 #ifdef WISE_SPARSIFY
00625         fprintf(f, " WISE_SPARSIFY");
00626 #endif
00627 #ifdef SPARSIFY_REMOVE_DUPLICATES
00628         fprintf(f, " SPARSIFY_REMOVE_DUPLICATES");
00629 #endif
00630 #ifdef NORMALIZE_SPARSE_CUTS
00631         fprintf(f, " NORMALIZE_SPARSE_CUTS");
00632 #endif
00633 #ifdef SLIDING
00634         fprintf(f, " SLIDING");
00635 #endif
00636 #ifdef ORTHOCUT
00637         fprintf(f, " ORTHOCUT");
00638 #endif
00639 #ifdef DISJUNCTIVE_CUTS
00640         fprintf(f, " DISJUNCTIVE_CUTS");
00641 #endif
00642 #ifdef ZERO_LB_DIAGONAL
00643         fprintf(f, " ZERO_LB_DIAGONAL");
00644 #endif
00645 #ifdef LINQUAD_CUTS
00646         fprintf(f, " LINQUAD_CUTS");
00647 #endif
00648 #ifdef LINQUAD_BOUNDS_CUTS
00649         fprintf(f, " LINQUAD_BOUNDS_CUTS");
00650 #endif
00651 
00652 
00653         fprintf(f, "\nCut management:");
00654 #ifdef MAX_CUTS_PER_ITER
00655         fprintf(f, " MAX_CUTS_PER_ITER=%d",MAX_CUTS_PER_ITER);
00656 #endif
00657 #ifdef RAND_CUT_ADD
00658         fprintf(f, " RAND_CUT_ADD=%.2f",RAND_CUT_ADD);
00659 #endif
00660 #ifdef DELETE_INACTIVE_CUTS
00661         fprintf(f, " DELETE_INACTIVE_CUTS(y_i<%e)",DELETE_INACTIVE_CUTS);
00662 #endif
00663 #ifdef FANCY_CUT_ACTIVITY
00664 #if (!defined DELETE_INACTIVE_CUTS)
00665 #error FANCY_CUT_ACTIVITY must be defined only with DELETE_INACTIVE_CUTS
00666 #endif
00667         fprintf(f, " FANCY_CUT_ACTIVITY=%d",FANCY_CUT_ACTIVITY);
00668 #endif
00669 
00670 #if (defined CUT_DELETION_OBJ_PERC_IMPROVEMENT_PERC) && (!defined CUT_DELETION_OBJ_PERC_IMPROVEMENT_ITER)
00671 #error CUT_DELETION_OBJ_PERC_IMPROVEMENT_PERC must be defined together with CUT_DELETION_OBJ_PERC_IMPROVEMENT_ITER
00672 #endif
00673 #if (!defined CUT_DELETION_OBJ_PERC_IMPROVEMENT_PERC) && (defined CUT_DELETION_OBJ_PERC_IMPROVEMENT_ITER)
00674 #error CUT_DELETION_OBJ_PERC_IMPROVEMENT_ITER must be defined together with CUT_DELETION_OBJ_PERC_IMPROVEMENT_PERC
00675 #endif
00676 #ifdef CUT_DELETION_OBJ_PERC_IMPROVEMENT_PERC
00677 #if (!defined DELETE_INACTIVE_CUTS)
00678 #error CUT_DELETION_OBJ_PERC_IMPROVEMENT_PERC & CUT_DELETION_OBJ_PERC_IMPROVEMENT_ITER must be defined only with DELETE_INACTIVE_CUTS
00679 #endif
00680         fprintf(f, " CUT_DELETION_OBJ_PERC_IMPROVEMENT=(%.4f,%d)",
00681         CUT_DELETION_OBJ_PERC_IMPROVEMENT_PERC,CUT_DELETION_OBJ_PERC_IMPROVEMENT_ITER);
00682 #endif
00683 #ifdef ADD_ONLY_PERC_DEEPEST_SPARSIFY_CUTS
00684         fprintf(f, " ADD_ONLY_PERC_DEEPEST_SPARSIFY_CUTS=%.2f%%",ADD_ONLY_PERC_DEEPEST_SPARSIFY_CUTS);
00685 #endif
00686 #ifdef ADD_RAND_PERC_SPARSIFY_CUTS
00687         fprintf(f, " ADD_RAND_PERC_SPARSIFY_CUTS=%.2f%%",ADD_RAND_PERC_SPARSIFY_CUTS);
00688 #endif
00689 #if (defined ADD_ONLY_PERC_DEEPEST_SPARSIFY_CUTS) && (defined ADD_RAND_PERC_SPARSIFY_CUTS)
00690 #error either ADD_ONLY_PERC_DEEPEST_SPARSIFY_CUTS or ADD_RAND_PERC_SPARSIFY_CUTS can be defined
00691 #endif
00692 
00693         fprintf(f, "\nHeuristics:    ");
00694 #ifdef HEUR_XXT
00695         fprintf(f, " HEUR_XXT");
00696 #endif
00697 #ifdef HEUR_MIN_MATRIX_NORM
00698         fprintf(f, " HEUR_MIN_MATRIX_NORM");
00699 #endif
00700 #ifdef HEUR_IMPROVE_SOLUTION
00701         fprintf(f, " HEUR_IMPROVE_SOLUTION");
00702 #endif
00703 
00704         fprintf(f, "\nDebug:         ");
00705 #ifdef CHECK
00706         fprintf(f, " CHECK");
00707 #endif
00708 #ifdef WRITE_ITER_LP
00709         fprintf(f, " WRITE_ITER_LP");
00710 #endif
00711 
00712 
00713         fprintf(f, "\n\n");
00714         fflush(f);
00715 } /* print_ifdefs */
00716 
00717 /***********************************************************************/
00718 void solver_status(OsiSolverInterface *solver) {
00719 
00720         if (solver->isAbandoned()) {
00721                 printf("### ERROR: Numerical difficulties in Solver\n");
00722                 exit(1);
00723         }
00724         
00725         if (solver->isProvenPrimalInfeasible()) {
00726                 printf("### WARNING: Problem is infeasible\n");
00727         
00728 #ifdef TRACE_ALL
00729                 solver->writeLp("xxxinfeas.lp");
00730                 printf("Problem written in xxxinfeas.lp\n");
00731 #endif
00732                 exit(1);
00733         }
00734 } /* solver_status */
00735 /***********************************************************************/
00736 FILE *open_f_res() {
00737         FILE *f_res = NULL;
00738         f_res = fopen("f_res.xxx", "r");
00739 
00740         if(f_res == NULL) {
00741                 f_res = fopen("f_res.xxx", "w");
00742                 print_ifdefs(f_res);
00743         }
00744         else {
00745                 fclose(f_res);
00746                 f_res = fopen("f_res.xxx", "a");
00747         }
00748         return f_res;
00749 } /* open_f_res() */
00750 /***********************************************************************/
00751 FILE *open_short_f_res() {
00752         FILE *short_f_res = NULL;
00753         short_f_res = fopen("short_f_res.xxx", "r");
00754 
00755         if(short_f_res == NULL) {
00756                 short_f_res = fopen("short_f_res.xxx", "w");
00757                 print_ifdefs(short_f_res);
00758         }
00759         else {
00760                 fclose(short_f_res);
00761                 short_f_res = fopen("short_f_res.xxx", "a");
00762         }
00763         return short_f_res;
00764 } /* open_short_f_res() */
00765 
00766 /***********************************************************************/
00767 void print_current_sol(int niter, double time, int tot_gen_cuts, int current_cuts, double objValue,double bestHeurObj, double currHeurObj) {
00768 
00769         if (niter == 0) 
00770                 printf("    %4s %8s %5s %5s %12s %12s %12s\n",
00771                         "iter","time","tcuts","gcuts","upper bound","curr heur","best heur");
00772         printf ("::: %4d %8.2f %5d %5d %12.4f ", 
00773                 niter, time, tot_gen_cuts, current_cuts, objValue);
00774         if (currHeurObj > -1e40)
00775                 printf("%12.4f ",currHeurObj);
00776         else
00777                 printf("%12s ","N/A");
00778 
00779         if (bestHeurObj > -1e40)
00780                 printf("%12.4f ",bestHeurObj);
00781         else
00782                 printf("%12s ","N/A");
00783         printf("\n");
00784 } /* print_current_sol() */
00785 /***********************************************************************/
00786 void print_file_current_sol(FILE* file, char* name, int niter, double time, int tot_gen_cuts, int current_gen_cuts, double objValue,double bestHeurObj, double currHeurObj) {
00787         if(niter == 0)
00788                 fprintf(file, "%-16s ", name);
00789         else
00790                 fprintf(file, "                 ");
00791 
00792         fprintf(file, "%12.4f %10d %10.2f %10d %10d ", 
00793                 objValue, niter, time, tot_gen_cuts, current_gen_cuts);
00794 
00795         if (currHeurObj > -1e40)
00796                 fprintf(file,"%12.4f",currHeurObj);
00797         else
00798                 fprintf(file,"%12s ","N/A");
00799 
00800         if (bestHeurObj > -1e40)
00801                 fprintf(file,"%12.4f ",bestHeurObj);
00802         else
00803                 fprintf(file,"%12s ","N/A");
00804         fprintf(file,"\n");
00805         fflush(file);
00806 } /* print_file_current_sol() */
00807 /***********************************************************************/
00808 void print_file_short_sol(FILE* file,char* name, int niter, double time, int tot_gen_cuts, int current_gen_cuts, double objValue,double bestHeurObj) {
00809 
00810         fprintf(file,"%-16s %12.4f ",
00811                 name,
00812                 objValue);
00813         if (bestHeurObj > -1e40)
00814                 fprintf(file,"%12.4f ",bestHeurObj);
00815         else
00816                 fprintf(file,"%12s ","N/A");
00817         fprintf(file,"%10.2f %10d %10d %10d\n",time,niter,tot_gen_cuts,current_gen_cuts);
00818         fflush(file);
00819 } /* print_file_short_sol() */
00820 /***********************************************************************/
00821 int feasibility_check(const int n, const int t, const int cons, const double *sol, const double **origMat , const double *origRhs, const char *origSense, const double *xlb,const double *xub,const double *ylb,const double *yub) {
00822         int N =  n*(n+3)/2;;
00823         bool violation = false;
00824         bool violation_hard = false;
00825 /*
00826 for(int i=0;i<n;i++)
00827         printf("x_%d=%.2f ",i,sol[i]);
00828 for(int i=n;i<N;i++)
00829         printf("X_%d=%.2f ",i,sol[i]);
00830 for(int i=N;i<N+t;i++)
00831         printf("y_%d=%.2f ",i,sol[i]);
00832 printf("\n");
00833 */
00834 
00835         //checking bound feasibility
00836         for(int i=0;i<n;i++) {
00837                 if (( sol[i] < xlb[i] - LP_TOLERANCE )||( sol[i] > xub[i] + LP_TOLERANCE )) {
00838                         printf("feasibility_check:: x_%d=%.8f violates bounds [%.8f,%.8f]\n",i,sol[i],xlb[i],xub[i]);
00839                 }
00840         }
00841         for(int i=0;i<t;i++) {
00842                 if (( sol[N+i] < ylb[i] - LP_TOLERANCE )||( sol[N+i] > yub[i] + LP_TOLERANCE )) {
00843                         printf("feasibility_check:: y_%d=%.8f violates bounds [%.8f,%.8f]\n",i,sol[N+i],ylb[i],yub[i]);
00844                         return FEAS_CHECK_BOUNDS_VIOLATION;
00845                 }
00846         }
00847         
00848         //checking constraint feasibility
00849         for (int k=0;k<cons;k++) {
00850                 double lhs = 0.0; 
00851                 bool y_terms = false;
00852                 for (int i=0;i<N;i++) {
00853                         lhs += origMat[k][i] * sol[i];
00854                 }
00855                 for (int i=N;i<N+t;i++) {
00856                         lhs += origMat[k][i] * sol[i];
00857                         if (origMat[k][i] != 0.0)
00858                                 y_terms = true;
00859                 }
00860 
00861 
00862 /*
00863 printf("\n");
00864 printf("origMat cons %d: ",k);
00865 for(int i=0;i<n;i++)
00866         printf("origMat[%d][%d]x=%.2f ",k,i,origMat[k][i]);
00867 printf("\n");
00868 for(int i=n;i<N;i++)
00869         printf("origMat[%d][%d]X=%.2f ",k,i,origMat[k][i]);
00870 printf("\n");
00871 for(int i=N;i<N+t;i++)
00872         printf("origMat[%d][%d]y=%.2f ",k,i,origMat[k][i]);
00873 printf("\n");
00874 printf("lhs=%.2f rhs=%.2f sense=%c\n",lhs,origRhs[k],origSense[k]);
00875 */
00876                 switch(origSense[k]) {
00877                         case 'L': 
00878                                 if ((lhs - LP_TOLERANCE) > origRhs[k]) {
00879                                         violation = true;
00880                                         if (y_terms == false)
00881                                                 violation_hard = true;
00882 #ifdef TRACE_ALL
00883                                         printf("feasibility_check:: violation in the original 'L' constraint %d     [%d]\n",k,violation_hard);
00884 #endif
00885                                 }
00886                                 break;
00887                         case 'G': 
00888                                 if ((lhs + LP_TOLERANCE) < origRhs[k]) {
00889                                         violation = true;
00890                                         if (y_terms == false)
00891                                                 violation_hard = true;
00892 #ifdef TRACE_ALL
00893                                         printf("feasibility_check:: violation in the original 'G' constraint %d     [%d]\n",k,violation_hard); 
00894 #endif
00895                                 }
00896                                 break;
00897                         case 'E': 
00898                                 if (((lhs + LP_TOLERANCE) < origRhs[k]) || ((lhs - LP_TOLERANCE) > origRhs[k])) {
00899                                         violation = true;
00900                                         if (y_terms == false)
00901                                                 violation_hard = true;
00902 #ifdef TRACE_ALL
00903                                         printf("feasibility_check:: violation in the original 'E' constraint %d     [%d]\n",k,violation_hard);
00904 #endif
00905                                 }
00906                                 break;
00907                         default: 
00908                                 printf("feasibility_check:: unrecognized constraint sense '%c' (constraint index = %d)\n",origSense[k],k);
00909                                 return -999;
00910                 }
00911         }
00912 
00913         if (violation == true)
00914                 if (violation_hard == true)
00915                         return FEAS_CHECK_CONSTRAINT_VIOLATION_NO_RECOVER;
00916                 else
00917                         return FEAS_CHECK_CONSTRAINT_VIOLATION;
00918         else 
00919                 return FEAS_CHECK_NO_VIOLATION;
00920 } //easibility_check()
00921 
00922 /***********************************************************************/
00923 double evaluateSolution(const int n, const int t,const double *heurSol, const double *b, const double *c , const double **Q, const double objConstant) {
00924         int N = n*(n+3)/2;
00925         double heurSolValue = 0.0;
00926         for(int i=0;i<n;i++) {
00927                 heurSolValue += heurSol[i] * b[i];
00928         }
00929         for(int i=0;i<t;i++) {
00930                 heurSolValue += heurSol[i+N] * c[i];
00931         }
00932         for(int i=0;i<n;i++)
00933                 for(int j=i;j<n;j++) {
00934                         heurSolValue += heurSol[indexQ(i,j,n)] * Q[i][j];
00935                         heurSolValue += heurSol[indexQ(i,j,n)] * Q[j][i];
00936                 }
00937 
00938         heurSolValue += objConstant;
00939 
00940         return heurSolValue;
00941 } //evaluateSolution()
00942 /***********************************************************************/

Generated on Thu Sep 22 03:05:57 2011 by  doxygen 1.4.7