Hierarchical Index

00001 #line 2 "cp_cuts.lw"
00002 #include<scil/scil.h>
00003 
00004 #define EPS 0.001
00005 #define MIN_MAX_VIOLATION 0.00001
00006 #define MIN_ADD_VIOLATION 0.0000001
00007 #define MAX_SIZE 100
00008 
00009 #include<LEDA/p_queue.h>
00010 #include<LEDA/array.h>
00011 
00012 class infer_status;
00013 class cp_sym_constraint : public sym_constraint {
00014 
00015  public:
00016   virtual
00017     void infer(infer_status& s) {
00018   };
00019 
00020   virtual
00021     void init(infer_status& s) {
00022   };
00023 
00024   virtual
00025     void reinit() {
00026   };
00027 
00028   virtual
00029     void notify(var v, infer_status& s) {
00030   };
00031 
00032   virtual
00033     ~cp_sym_constraint() {
00034   };
00035 };
00036 
00037 class infer_status {
00038 
00039   p_queue<int, cp_sym_constraint*> Act;
00040   list<cp_sym_constraint*>* NV;
00041   h_array<var,int>& VarIndex;
00042 
00043   private:
00044 
00045   // etablish an index for all variables
00046   var_map<int> VM;
00047 
00048   // an int for all variables 0: fixed to 0; 1 fixed to 1; 2 free
00049   int* varstat;
00050 
00051   // for every fixed variable, list the fixed variables to derive the fixation
00052   // empty list means fixed by lp value.
00053   list<var>* usedvars;
00054 
00055   // number of fixed variables; -1 means unsolvable system.
00056   int number_of_fixed;
00057 
00058   // vars used to derive infeasbibility
00059   list<var> infeasvars;
00060 
00061   public:
00062 
00063   infer_status(h_array<var,int>& VarIndex_) : VarIndex(VarIndex_) {
00064   };
00065 
00066   // A new infer status for the variables in s. The boolean indecates if
00067   // the variables with integral value should be fixed.
00068   void init(subproblem& s, bool fix) {
00069     // TODO only use binary variables
00070     NV=new list<cp_sym_constraint*>[s.nVar()];
00071     int i;
00072     var v;
00073     forall_variables(v, s) {
00074       VM[VarIndex[v]]=v;
00075       NV[VarIndex[v]].clear();
00076       i++;
00077     }
00078     varstat=new int[i];
00079     for(int j=0; j<i; j++) varstat[j]=2;
00080 
00081     usedvars=new list<var>[i];
00082     for(int j=0; j<i; j++) usedvars[j].clear();
00083     number_of_fixed=0;
00084 
00085     if(fix) {
00086       for(int j=0; j<i; j++) {
00087         if(s.value(VM[j])<EPS) { 
00088           varstat[j]=0;
00089           number_of_fixed++;
00090           if(usedvars[j].size()!=0) cout<<"errrer\n";
00091         };
00092         if(s.value(VM[j])>1-EPS) {
00093           varstat[j]=1;
00094           number_of_fixed++;
00095         };
00096       };
00097     };
00098   };
00099 
00100   void reinit(subproblem& s, bool fix) {
00101     int i=s.nVar();
00102 
00103     for(int j=0; j<i; j++) varstat[j]=2;
00104 
00105     for(int j=0; j<i; j++) usedvars[j].clear();
00106     infeasvars.clear();
00107     number_of_fixed=0;
00108 
00109     if(fix) {
00110       for(int j=0; j<i; j++) {
00111         if(s.value(VM[j])<EPS) { 
00112           varstat[j]=0;
00113           number_of_fixed++;
00114           if(usedvars[j].size()!=0) cout<<"errrer\n";
00115         };
00116         if(s.value(VM[j])>1-EPS) {
00117           varstat[j]=1;
00118           number_of_fixed++;
00119         };
00120       };
00121     };
00122   };
00123 
00124   ~infer_status() {
00125   };
00126 
00127   void remove() {
00128     //delete[] usedvars;
00129     //delete[] varstat;
00130     //delete[] NV;
00131     Act.clear();
00132     infeasvars.clear();
00133     //VM.init(-1);
00134     number_of_fixed=0;
00135   };
00136 
00137   void make_active(int i,cp_sym_constraint* c) {
00138     Act.insert(i, c);
00139   };
00140 
00141   bool has_active() {
00142     if(unsolvable()) return false;
00143     return Act.size()>0;
00144   };
00145 
00146   cp_sym_constraint* get_active() {
00147     cp_sym_constraint* c=Act.inf(Act.find_min());
00148     Act.del_min();
00149     return c;
00150   };
00151 
00152   void notify(cp_sym_constraint* c, var v) {
00153     NV[VarIndex[v]].append(c);
00154   };
00155 
00156   // The following function can be used to access the status of a variable
00157   bool fixed(var v) {
00158     if(unsolvable()) {
00159       //cout<<"Warning: modifying unsolvable system not possible\n";
00160       return true;
00161     };
00162     return varstat[VarIndex[v]]!=2;
00163   };
00164 
00165   bool zero_fixed(var v) {
00166     if(unsolvable()) {
00167       //cout<<"Warning: modifying unsolvable system not possible\n";
00168       return true;
00169     };
00170     return varstat[VarIndex[v]]==0;
00171   };
00172 
00173   bool one_fixed(var v) {
00174     if(unsolvable()) {
00175       //cout<<"Warning: modifying unsolvable system not possible\n";
00176       return true;
00177     };
00178     return varstat[VarIndex[v]]==1;
00179   };
00180 
00181   int lower_bound(var v) {
00182     if(unsolvable()) {
00183       cout<<"Warning: modifying unsolvable system not possible\n";
00184       return 0;
00185     };
00186     if(one_fixed(v)) return 1;
00187     return 0;
00188   };
00189 
00190   int upper_bound(var v) {
00191     if(unsolvable()) {
00192       cout<<"Warning: modifying unsolvable system not possible\n";
00193       return 0;
00194     };
00195     if(zero_fixed(v)) return 0;
00196     return 1;
00197   };
00198 
00199   // The following function is used to determine if something has changed
00200   int number_of_fixed_variables() {
00201     return number_of_fixed;
00202   };
00203 
00204   bool unsolvable() {
00205     return number_of_fixed==-1;
00206   };
00207 
00208   // The following function can be used to modify the status of a variable
00209 
00210   // Fix a variable to a value; 
00211   // list<var> specifices the list of fixed variables used to derive the 
00212   // fixation.
00213   void fix(var v, int i, list<var>& L) {
00214     if(unsolvable()) {
00215       cout<<"Warning: modifying unsolvable system not possible\n";
00216       return;
00217     };
00218     if(fixed(v)) {
00219       cout<<"Warning: fixing fixed variable ignored\n";
00220       return;
00221     };
00222     if((i!=0) && (i!=1)) {
00223       cout<<"Error: fixing not to 0 or 1 forbitten\n";
00224       return;
00225     };
00226 
00227     varstat[VarIndex[v]]=i;
00228     number_of_fixed++;
00229     int vn=VarIndex[v];
00230     var v1;
00231     forall(v1,L) usedvars[vn].append(v1);
00232     cp_sym_constraint* c;
00233     forall(c, NV[VarIndex[v]]) {
00234       c->notify(v, *this);
00235     };
00236   };
00237 
00238 
00239   void unsolvable(list<var>& L) {
00240     if(unsolvable()) {
00241       cout<<"Warning: modifying unsolvable system not possible\n";
00242       return;
00243     };
00244     number_of_fixed=-1;
00245     var v;
00246     forall(v, L) infeasvars.append(v);
00247     return;
00248   };
00249   // Indicate that the problem is unsolvable, list<var> as above
00250 
00251   const list<var>& vars_for_infeas() {
00252     return infeasvars;
00253   };
00254 
00255   const list<var>& vars_for_var(var v) {
00256     if((unsolvable()==true) && (!fixed(v))) // it could be fixed to other value 
00257       return vars_for_infeas();
00258     return usedvars[VarIndex[v]];
00259   };
00260 };
00261 
00262 
00263 class cp_basic_constraint : public cp_sym_constraint {
00264   private:
00265   cons c;
00266   row r;
00267   cons_sense s;
00268   int rhs;
00269   int lb, ub;
00270   map<var, int> coeff;
00271   list<var> ULB, UUB;
00272   int uubn, ulbn, ubn, lbn;
00273   bool wa;
00274 
00275   public:
00276   cp_basic_constraint(cons c_) {
00277     c=c_;
00278     c.cons_pointer()->non_zero_entries(r);
00279     s=c.sense();
00280     wa=false;
00281     row_entry re;
00282     forall(re, r) coeff[re.Var]=(int) re.coeff;
00283     rhs=(int) c.rhs();
00284   };
00285 
00286   void init(infer_status& is) {
00287     // init lb and ub
00288     wa=false;
00289     ULB.clear();
00290     UUB.clear();
00291     row_entry re;
00292     ub=0; lb=0;
00293     forall(re, r) {
00294       if(re.coeff>0) {
00295         ub+=((int) re.coeff)*is.upper_bound(re.Var);
00296         lb+=((int) re.coeff)*is.lower_bound(re.Var);
00297       };
00298       if(re.coeff<0) {
00299         ub+=((int) re.coeff)*is.lower_bound(re.Var);
00300         lb+=((int) re.coeff)*is.upper_bound(re.Var);
00301       }
00302       if(is.one_fixed(re.Var)) {
00303         if(re.coeff>0) ULB.append(re.Var);
00304         else UUB.append(re.Var);
00305       };
00306       if(is.zero_fixed(re.Var)) {
00307         if(re.coeff>0) UUB.append(re.Var);
00308         else ULB.append(re.Var);
00309       }
00310       if(!is.fixed(re.Var)) is.notify(this, re.Var);
00311     };
00312     if((s!=Greater) && (lb>=rhs) && (wa==false)) { 
00313       is.make_active(ULB.size(), this);
00314       wa=true;
00315     }
00316     if((s!=Less) && (ub<=rhs) && (wa==false)) {
00317       is.make_active(UUB.size(), this);
00318       wa=true;
00319     };
00320     ulbn=ULB.size();
00321     uubn=UUB.size();
00322     ubn=ub; lbn=lb;
00323   };
00324 
00325   void reinit() {
00326     while(ULB.size()>ulbn) {
00327       ULB.pop_back();
00328     };
00329     while(UUB.size()>uubn) {
00330       UUB.pop_back();
00331     };
00332     ub=ubn; lb=lbn;
00333     wa=false;
00334   };
00335 
00336   void notify(var v, infer_status& is) {
00337     int d=coeff[v];
00338     if(is.one_fixed(v)) {
00339       if(d>0) lb+=d;
00340       else ub+=d;
00341     } else {
00342       if(d>0) ub-=d;
00343       else lb-=d;
00344     };
00345     if(is.one_fixed(v)) {
00346       if(d>0) ULB.append(v);
00347       else UUB.append(v);
00348     };
00349     if(is.zero_fixed(v)) {
00350       if(d>0) UUB.append(v);
00351       else ULB.append(v);
00352     }
00353     if((s!=Greater) && (lb>=rhs) && (wa==false)) { 
00354       is.make_active(ULB.size(), this);
00355       wa=true;
00356     }
00357     if((s!=Less) && (ub<=rhs) && (wa==false)) {
00358       is.make_active(UUB.size(), this);
00359       wa=true;
00360     };
00361   };
00362 
00363   void infer(infer_status& is) {
00364     row_entry re;
00365     int i=0;
00366     if(s!=Greater) {
00367       if(lb>rhs) {
00368         is.unsolvable(ULB);
00369         return;
00370       }
00371       if(lb==rhs) {
00372         forall(re, r) if(!is.fixed(re.Var)) {
00373           if(lb+re.coeff>rhs) is.fix(re.Var, 0, ULB);
00374           if(lb-re.coeff>rhs) is.fix(re.Var, 1, ULB);
00375         };
00376         i++;
00377       }
00378     };
00379     if(s!=Less) {
00380       if(ub<rhs) {
00381         is.unsolvable(UUB);
00382         return;
00383       };
00384       if(ub==rhs) {
00385         forall(re, r) if(!is.fixed(re.Var)) {
00386           if(ub-re.coeff<rhs) is.fix(re.Var, 1, UUB);
00387           if(ub+re.coeff<rhs) is.fix(re.Var, 0, UUB);
00388         };
00389         i++;
00390       };
00391     };
00392     if (i==0) cout<<"Error\n";
00393   };
00394 };
00395 
00396 class edge_compare : public leda_cmp_base<edge> {
00397  private:
00398   edge_array<double>& E;
00399   
00400  public:
00401   edge_compare(edge_array<double>& E_) : E(E_) {};
00402 
00403   int operator()(const edge& e1, const edge& e2) const {
00404     if(E[e1]>E[e2]) return 1;
00405     if(E[e1]<E[e2]) return -1;
00406     return 0;
00407   };
00408 };
00409 
00410 class cp_cuts : public sym_constraint {
00411 
00412  private:
00413 
00414   list<cp_sym_constraint*> CSCL;
00415   h_array<int, var> IndexVar;
00416   h_array<var, int> VarIndex;
00417 
00418  public:
00419   
00420   // add a symbolic constraint with cp-functionallity to derive cuts
00421   void add_sym_constraint(cp_sym_constraint* c) {
00422     CSCL.append(c);
00423   };
00424   
00425   void infer(infer_status& is, var v, int b) {
00426     cp_sym_constraint* c;
00427     forall(c, CSCL) c->reinit();
00428 
00429     list<var> dummy;
00430 
00431     is.fix(v, b, dummy);
00432 
00433     while(is.has_active()) {
00434       c=is.get_active();
00435       c->infer(is);
00436     };
00437 
00438   };
00439 
00440   status separate(subproblem& s) {
00441 
00442     cout<<"HIER\n";
00443 
00444     int i=0;
00445     var v;
00446     forall_variables(v,s) {
00447       IndexVar[i]=v;
00448       VarIndex[v]=i;
00449       i++;
00450     };
00451 
00452     list<cons_obj*> CL;
00453 
00454     bool st=false;
00455 
00456     bool B[s.nVar()];
00457     int number_of_fixed_variables=0;
00458 
00459     for(int i=0; i<s.nVar(); i++) {
00460       if((s.value(IndexVar[i])<EPS) || (s.value(IndexVar[i])>1-EPS)) {
00461         B[i]=true;
00462         number_of_fixed_variables++;
00463       } else B[i]=false;
00464     };
00465 
00466     graph G;
00467     int uf=s.nVar()-number_of_fixed_variables;
00468     node_array<var> GtV(G,2*uf, nil);
00469     node_array<double> W(G, 2*uf, 0);
00470     node_array<bool> Si(G, 2*uf, false);
00471     node VtG0[s.nVar()], VtG1[s.nVar()];
00472     edge e;
00473     edge_map<list<var> > EM(G);
00474     edge_map<double> EV(G);
00475 
00476     for(int j=0; j<s.nVar(); j++) if(!B[j]) {
00477       VtG0[j]=G.new_node();
00478       VtG1[j]=G.new_node();
00479       GtV[VtG0[j]]=IndexVar[j];
00480       GtV[VtG1[j]]=IndexVar[j];
00481       W[VtG0[j]]=1-s.value(IndexVar[j]);
00482       W[VtG1[j]]=s.value(IndexVar[j]);
00483       Si[VtG1[j]]=true;
00484     };
00485 
00486     list<node> IEC;
00487     list<var> dummy;
00488     infer_status is(VarIndex);
00489 
00490     cp_sym_constraint* c;
00491     is.init(s,true);
00492     forall(c, CSCL) c->init(is);
00493 
00494     for(int j=0; j<s.nVar(); j++) if(!B[j]) {
00495       bool b0=false, b1=false;
00496 
00497       if(G.degree(VtG1[j])==0) {
00498         is.reinit(s,true);
00499         b0=true;
00500         infer(is, IndexVar[j], 0);
00501         //};
00502         
00503         if((!is.unsolvable()) && (b0)) {
00504           for(int k=0; k<s.nVar(); k++) 
00505             if((k!=j) && (!B[k]) && (is.fixed(IndexVar[k]))) {
00506               if(is.zero_fixed(IndexVar[k])) { 
00507                 if(j<k) e=G.new_edge(VtG0[j], VtG1[k]);
00508                 else e=G.new_edge(VtG1[k], VtG0[j]);
00509                 search_vars(is, IndexVar[k], EM, EV, e, s, IndexVar[j]);
00510                 //if(EV[e]>0.1) G.del_edge(e);
00511               }
00512               if(is.one_fixed(IndexVar[k])) {
00513                 if(j<k) e=G.new_edge(VtG0[j], VtG0[k]);
00514                 else e=G.new_edge(VtG0[k], VtG0[j]);
00515                 search_vars(is, IndexVar[k], EM, EV, e, s, IndexVar[j]);
00516                 //if(EV[e]>0.1) G.del_edge(e);
00517               }
00518             };
00519         };
00520 
00521         if(is.unsolvable()) {
00522           if(add_infeasible_constraint(s, is, IndexVar[j],0, CL))
00523             st=true;
00524         };
00525         is.remove();
00526       };
00527 
00528       if(G.degree(VtG0[j])==0) {
00529         is.reinit(s,true);
00530         b1=true;
00531         infer(is, IndexVar[j], 1);
00532         //};
00533 
00534         if((!is.unsolvable()) && (b1)) {
00535           for(int k=0; k<s.nVar(); k++) 
00536             if((k!=j) && (!B[k]) && (is.fixed(IndexVar[k]))) {
00537               if(is.zero_fixed(IndexVar[k])) {
00538                 if(j<k) e=G.new_edge(VtG1[j], VtG1[k]);
00539                 else e=G.new_edge(VtG1[k], VtG1[j]);
00540                 search_vars(is, IndexVar[k], EM, EV, e, s, IndexVar[j]);
00541                 //if(EV[e]>0.1) G.del_edge(e);
00542               }
00543               if(is.one_fixed(IndexVar[k])) {
00544                 if(j<k) e=G.new_edge(VtG1[j], VtG0[k]);
00545                 else e=G.new_edge(VtG0[k], VtG1[j]);
00546                 search_vars(is, IndexVar[k], EM, EV, e, s, IndexVar[j]);
00547                 //if(EV[e]>0.1) G.del_edge(e);
00548               }
00549             };
00550         };
00551 
00552         if(is.unsolvable()) {
00553           if(add_infeasible_constraint(s, is, IndexVar[j],1, CL))
00554             st=true;
00555         };
00556       };
00557       is.remove();
00558     };
00559 
00560     cout<<G.number_of_nodes()<<" "<<G.number_of_edges()<<" size of G\n";
00561     Make_Simple(G);
00562     cout<<G.number_of_nodes()<<" "<<G.number_of_edges()<<" size of G\n";
00563     clique_separate(G, W, GtV, s, Si, IEC, st, EM, EV, CL);
00564  
00565     if(st==false) return no_constraint_found;
00566 
00567     cons_obj* bc;
00568     forall(bc, CL)
00569       s.add_basic_constraint(bc);
00570 
00571     VarIndex.clear();
00572     IndexVar.clear();
00573 
00574     cout<<CL.size()<<" CONSTRAINTS\n";
00575 
00576     return constraint_found;
00577   };
00578 
00579   void search_vars(infer_status& is, var v, edge_map<list<var> >& EM, 
00580                    edge_map<double>& EV, edge e, subproblem& s, var rv) {
00581     bool B[s.nVar()];
00582     for(int i=0; i<s.nVar(); i++) B[i]=false;
00583     list<var> TD;
00584     var w,x;
00585     if(is.unsolvable()) {
00586       forall(w, is.vars_for_infeas()) TD.append(w);
00587     } else {
00588       forall(w,is.vars_for_var(v)) TD.append(w);
00589     };
00590     //TD=is.vars_for_var(v);
00591     while(!TD.empty()) {
00592       w=TD.pop();
00593       if(!B[VarIndex[w]]) {
00594         B[VarIndex[w]]=true;
00595         forall(x, is.vars_for_var(w)) TD.append(x);
00596       };
00597     };
00598     EV[e]=0;
00599     for(int i=0; i<s.nVar(); i++) 
00600       if((B[i]) && (IndexVar[i]!=rv) && (is.vars_for_var(IndexVar[i]).size()==0)) { 
00601         EM[e].append(IndexVar[i]);
00602         if(s.value(IndexVar[i])>0.5) {
00603           EV[e]-=1-s.value(IndexVar[i]);
00604         } else {
00605           EV[e]-=s.value(IndexVar[i]);
00606         };
00607       };
00608   };
00609 
00610   void clique_separate(graph& G, node_array<double>& W, node_array<var>& V, 
00611                        subproblem& s, node_array<bool>& Si, list<node>& IEC,
00612                        bool& st, 
00613                        edge_map<list<var> >& EM, edge_array<double>& EV,
00614                        list<cons_obj*>& CL) {
00615     node u,v;
00616     edge e;
00617     list<node> C;
00618     list<edge> L, L1;
00619     edge_array<double> E(G);
00620     double ic=0;
00621     forall(u, IEC) ic+=W[u];
00622     node_array<bool> CA(G, false);
00623     forall_nodes(u, G) {
00624       C.clear();
00625       double w=ic;
00626       L=G.in_edges(u);
00627       L1=G.out_edges(u);
00628       L.conc(L1);
00629       forall(e, L) {
00630         E[e]=W[G.opposite(u,e)];
00631       };
00632       edge_compare ecomp(E); 
00633       L.sort(ecomp);
00634       forall(e, L) {
00635         v=G.opposite(u,e);
00636         double un=W[v]+EV[e];
00637         int c=0;
00638         forall_inout_edges(e, v) {
00639           if(CA[G.opposite(v,e)]) { 
00640             c++;
00641             un+=EV[e];
00642           };
00643         };
00644         if((c==C.size()) && (un>0.1)) {
00645           C.append(v);
00646           CA[v]=true;
00647           w+=W[v];
00648         }
00649       }
00650       if(w+W[u]>1.1) {
00651         C.append(u);
00652         CA[u]=true;
00653         cout<<w<<" clique\n";
00654         cout<<C.size()<<" size\n";
00655         row r; int rhs=1;
00656         forall(v, C) {
00657           if(Si[v]) r+=V[v];
00658           else { r-=V[v]; rhs--; }
00659         }
00660         var v;
00661         int B[s.nVar()];
00662         for(int k=0; k<s.nVar(); k++) B[k]=0;
00663         forall_edges(e, G) if((CA[source(e)]) && (CA[target(e)])) {
00664           forall(v, EM[e]) B[VarIndex[v]]++;
00665         }
00666         int m,l=0, n;
00667         for(int k=0; k<s.nVar(); k++) if(B[k]>0) {
00668           m=0, n=0;
00669           while(n<B[k]) { m++; n+=m; };
00670           //cout<<m<<" m\n";
00671           l+=m;
00672           if(s.value(IndexVar[k])>1-EPS) { r+=m*IndexVar[k]; rhs+=m; }
00673           if(s.value(IndexVar[k])<EPS) { r-=m*IndexVar[k]; }
00674         }
00675         cout<<l<<" size of cut\n";
00676         row_entry re; double rd=0;
00677         forall(re, r) rd+=re.coeff*s.value(re.Var);
00678         cout<<rd<<" <= "<<rhs<<endl;
00679         if((rd-rhs)/l>MIN_ADD_VIOLATION) {
00680           cout<<"add cliequ\n";
00681           CL.append(r<=rhs);
00682           if((rd-rhs)/l>MIN_MAX_VIOLATION) {
00683             st=true;
00684           };
00685         }
00686       }
00687       CA[u]=false;
00688       forall(v, C) CA[v]=false;
00689     }
00690   };
00691 
00692   bool add_infeasible_constraint(subproblem& s, infer_status& is0, 
00693                                  infer_status& is1, var v1, list<cons_obj*>& CL) {
00694     array<bool> UV0(s.nVar());
00695     array<bool> UV1(s.nVar());
00696     var_map<int> VM;
00697     for(int i=0; i<s.nVar(); i++) {
00698       UV0[i]=false;
00699       UV1[i]=false;
00700       VM[i]=IndexVar[i];
00701     };
00702 
00703     list<var> U;
00704     var v;
00705     list<var> L=is0.vars_for_infeas();
00706     //cout<<L.size()<<" size of used 0\n";
00707     forall(v, L) U.append(v);
00708     mark_variables(UV0, U, VM, is0);
00709     L=is1.vars_for_infeas();
00710     //cout<<L.size()<<" size of used 1\n";
00711     forall(v, L) U.append(v);
00712     mark_variables(UV1, U, VM, is1);
00713 
00714     row r;
00715     int rhs=-1;
00716     int k=0;
00717     double vo=1;
00718 
00719     for(int i=0; i<s.nVar(); i++) {
00720       if((UV0[i]) || (UV1[i])) {
00721         if(s.value(IndexVar[i])<EPS) {
00722           r-=IndexVar[i];
00723           vo-=s.value(IndexVar[i]);
00724           k++;
00725         };
00726         if (s.value(IndexVar[i])>1-EPS) {
00727           r+=IndexVar[i];
00728           rhs++;
00729           vo-=1-s.value(IndexVar[i]);
00730           k++;
00731         }
00732       }
00733     };
00734 
00735     if(strengthen_cut(s, r, rhs)) vo+=0.5;
00736 
00737     cout<<vo/k<<" violation\n";
00738     if(vo/k<MIN_ADD_VIOLATION) return false;
00739     CL.append(r<=rhs);
00740     if(vo/k<MIN_MAX_VIOLATION) return false;
00741     return true;
00742 
00743   };
00744 
00745   bool add_infeasible_constraint(subproblem& s, infer_status& is, list<cons_obj*>& CL) { 
00746     array<bool> UV(s.nVar());
00747     var_map<int> VM;
00748     for(int i=0; i<s.nVar(); i++) {
00749       UV[i]=false;
00750       VM[i]=IndexVar[i];
00751     };
00752 
00753     list<var> U;
00754     var v;
00755     list<var> L=is.vars_for_infeas();
00756     //cout<<L.size()<<" size of used\n";
00757     forall(v, L) U.append(v);
00758     mark_variables(UV, U, VM, is);
00759 
00760     row r;
00761     int rhs=-1;
00762     int k=0;
00763     double vo=1;
00764 
00765     for(int i=0; i<s.nVar(); i++) {
00766       if(UV[i]) {
00767         if(s.value(IndexVar[i])<EPS) {
00768           r-=IndexVar[i];
00769           k++;
00770           vo-=s.value(IndexVar[i]);
00771         };
00772         if (s.value(IndexVar[i])>1-EPS) {
00773           r+=IndexVar[i];
00774           rhs++;
00775           k++;
00776           vo-=1-s.value(IndexVar[i]);
00777         }
00778       }
00779     };
00780 
00781     strengthen_cut(s, r, rhs);
00782 
00783     cout<<vo/k<<" violation\n";
00784     if(vo/k<MIN_ADD_VIOLATION) return false;
00785     CL.append(r<=rhs);
00786     if(vo/k<MIN_MAX_VIOLATION) return false;
00787     return true;
00788 
00789   };
00790 
00791   bool add_infeasible_constraint(subproblem& s, infer_status& is, 
00792                                  var v0, int f, list<cons_obj*>& CL) { 
00793     array<bool> UV(s.nVar());
00794     var_map<int> VM;
00795     for(int i=0; i<s.nVar(); i++) {
00796       UV[i]=false;
00797       VM[i]=IndexVar[i];
00798     };
00799 
00800     list<var> U;
00801     var v;
00802     list<var> L=is.vars_for_infeas();
00803     //cout<<L.size()<<" size of used\n";
00804     forall(v, L) U.append(v);
00805     mark_variables(UV, U, VM, is);
00806 
00807     row r;
00808     double vo;
00809     int rhs;
00810     if(f==0) {
00811       r-=v0;
00812       rhs=-1;
00813       vo=1-s.value(v0);
00814     } else {
00815       r+=v0;
00816       rhs=0;
00817       vo=s.value(v0);
00818     };
00819     int k=1;
00820 
00821     for(int i=0; i<s.nVar(); i++) {
00822       if(UV[i]) {
00823         if(s.value(IndexVar[i])<EPS) {
00824           r-=IndexVar[i];
00825           k++;
00826           vo-=s.value(IndexVar[i]);
00827         };
00828         if (s.value(IndexVar[i])>1-EPS) {
00829           r+=IndexVar[i];
00830           rhs++;
00831           k++;
00832           vo-=1-s.value(IndexVar[i]);
00833         }
00834       }
00835     };
00836 
00837     strengthen_cut(s, r, rhs);
00838 
00839     cout<<vo/k<<" violation\n";
00840     if(vo/k<MIN_ADD_VIOLATION) return false;
00841     CL.append(r<=rhs);
00842     if(vo/k<MIN_MAX_VIOLATION) return false;
00843     return true;
00844   };
00845 
00846   void mark_variables(array<bool>& UV, list<var>& L, var_map<int>& VM, 
00847                       infer_status& is) {
00848     var v;
00849     list<var> L1;
00850     while(L.size()>0) {
00851       v=L.pop();
00852       if(!UV[VarIndex[v]]) {
00853         UV[VarIndex[v]]=true;
00854         L1=is.vars_for_var(v);
00855         if(L1.size()>0)
00856           forall(v, L1) L.append(v);
00857       }
00858     };
00859   };
00860 
00861   bool strengthen_cut(subproblem& s, row& r, int& rhs) {
00862     return false;
00863     /*
00864     row_entry re;
00865     list<var> dummy;
00866     bool B0[s.nVar()];
00867     bool B1[s.nVar()];
00868     for(int i=0; i<s.nVar(); i++) {
00869       B0[i]=true;
00870       B1[i]=true;
00871     };
00872     forall(re, r) {
00873       B0[VarIndex[re.Var]]=false;
00874       B1[VarIndex[re.Var]]=false;
00875       infer_status is(s, false);
00876       row_entry re1;
00877       forall(re1, r) {
00878         if(re1.Var!=re.Var) {
00879           if(re1.coeff==1) is.fix(re1.Var, 1, dummy);
00880           if(re1.coeff==-1) is.fix(re1.Var, 0, dummy);
00881         } else {
00882           if(re1.coeff==1) is.fix(re1.Var, 0, dummy);
00883           if(re1.coeff==-1) is.fix(re1.Var, 1, dummy);
00884         }
00885       };
00886       infer(is);
00887       if(!is.unsolvable()) {
00888         for(int i=0; i<s.nVar(); i++) {
00889           if(!is.one_fixed(IndexVar[i])) B1[i]=false;
00890           if(!is.zero_fixed(IndexVar[i])) B0[i]=false;
00891         };
00892       };
00893     };
00894     double v=-1; 
00895     int j;
00896     for(int i=0; i<s.nVar(); i++) {
00897       if((B0[i]) && (s.value(IndexVar[i])>v)) {
00898         v=s.value(IndexVar[i]);
00899         j=i;
00900       };
00901       if((B1[i]) && (1-s.value(IndexVar[i])>v)) {
00902         v=1-s.value(IndexVar[i]);
00903         j=i;
00904       }; 
00905     }; 
00906     if(v<=0.01) return false;
00907       
00908     cout<<v<<" Jaa\n";
00909     if(B0[j]) r+=IndexVar[j];
00910     else { r-=IndexVar[j]; rhs--; }
00911     return true;
00912     */
00913   };
00914 };
00915
cp_cuts.h
