1 /*===========================================================================*
2  * This file is part of the Discrete Conic Optimization (DisCO) Solver.      *
3  *                                                                           *
4  * DisCO is distributed under the Eclipse Public License as part of the      *
5  * COIN-OR repository (http://www.coin-or.org).                              *
6  *                                                                           *
7  * Authors:                                                                  *
8  *          Aykut Bulut, Lehigh University                                   *
9  *          Yan Xu, Lehigh University                                        *
10  *          Ted Ralphs, Lehigh University                                    *
11  *                                                                           *
12  * Copyright (C) 2001-2018, Lehigh University, Aykut Bulut, Yan Xu, and      *
13  *                          Ted Ralphs.                                      *
14  * All Rights Reserved.                                                      *
15  *===========================================================================*/
16 
17 
18 #include "DcoBranchStrategyStrong.hpp"
19 #include "DcoModel.hpp"
20 #include "DcoMessage.hpp"
21 #include "DcoTreeNode.hpp"
22 
23 
DcoBranchStrategyStrong(DcoModel * model)24 DcoBranchStrategyStrong::DcoBranchStrategyStrong(DcoModel * model)
25   : BcpsBranchStrategy(model) {
26   setType(static_cast<int>(DcoBranchingStrategyStrong));
27 }
28 
29 // Assumes problem is not unbounded.
updateScore(BcpsBranchObject * bobject,double orig_lb,double orig_ub,double orig_obj) const30 void DcoBranchStrategyStrong::updateScore(BcpsBranchObject * bobject,
31                                             double orig_lb,
32                                             double orig_ub,
33                                             double orig_obj) const {
34   // get dco model and message stuff
35   DcoModel * dco_model = dynamic_cast<DcoModel*>(model());
36   // CoinMessageHandler * message_handler = dco_model->dcoMessageHandler_;
37   // CoinMessages * messages = dco_model->dcoMessages_;
38   DcoBranchObject * dco_bobject = dynamic_cast<DcoBranchObject*>(bobject);
39   // solve subproblem for the down branch
40   dco_model->solver()->setColUpper(bobject->index(), dco_bobject->ubDownBranch());
41   dco_model->solver()->solveFromHotStart();
42   //double down_obj = 0.5*ALPS_INFINITY;
43   double down_obj = 1.0;
44   // std::cout << "abandoned " << dco_model->solver()->isAbandoned()<< std::endl;
45   // std::cout << "optimal " << dco_model->solver()->isProvenOptimal()<< std::endl;
46   // std::cout << "primal inf" << dco_model->solver()->isProvenPrimalInfeasible()<< std::endl;
47   // std::cout << "dual inf" << dco_model->solver()->isProvenDualInfeasible()<< std::endl;
48   // std::cout << "primal obj limit" << dco_model->solver()->isPrimalObjectiveLimitReached()<< std::endl;
49   // std::cout << "dual obj limit" << dco_model->solver()->isDualObjectiveLimitReached()<< std::endl;
50   // std::cout << "iter limit" << dco_model->solver()->isIterationLimitReached()<< std::endl;
51   if (dco_model->solver()->isProvenOptimal()
52       or dco_model->solver()->isIterationLimitReached()
53       or dco_model->solver()->isDualObjectiveLimitReached()) {
54     down_obj = dco_model->solver()->getObjValue();
55   }
56   else {
57     std::cout << "prob not opt." << std::endl;
58   }
59   // restore bound
60   dco_model->solver()->setColUpper(bobject->index(), orig_ub);
61   // solve subproblem for the up branch
62   dco_model->solver()->setColLower(bobject->index(), dco_bobject->lbUpBranch());
63   dco_model->solver()->solveFromHotStart();
64   //double up_obj = 0.5*ALPS_INFINITY;
65   double up_obj = 1.0;
66   // std::cout << "abandoned " << dco_model->solver()->isAbandoned()<< std::endl;
67   // std::cout << "optimal " << dco_model->solver()->isProvenOptimal()<< std::endl;
68   // std::cout << "primal inf" << dco_model->solver()->isProvenPrimalInfeasible()<< std::endl;
69   // std::cout << "dual inf" << dco_model->solver()->isProvenDualInfeasible()<< std::endl;
70   // std::cout << "primal obj limit" << dco_model->solver()->isPrimalObjectiveLimitReached()<< std::endl;
71   // std::cout << "dual obj limit" << dco_model->solver()->isDualObjectiveLimitReached()<< std::endl;
72   // std::cout << "iter limit" << dco_model->solver()->isIterationLimitReached()<< std::endl;
73   if (dco_model->solver()->isProvenOptimal()
74       or dco_model->solver()->isIterationLimitReached()
75       or dco_model->solver()->isDualObjectiveLimitReached()) {
76     up_obj = dco_model->solver()->getObjValue();
77   }
78   else {
79     std::cout << "prob not opt." << std::endl;
80   }
81   // restore bound
82   dco_model->solver()->setColLower(bobject->index(), orig_lb);
83   // set score
84   double down_diff = fabs(orig_obj-down_obj);
85   double up_diff = fabs(orig_obj-up_obj);
86   double score = down_diff>up_diff ? down_diff : up_diff;
87   bobject->setScore(score);
88 }
89 
90 
infeas(double value) const91 double DcoBranchStrategyStrong::infeas(double value) const {
92   // get dco model and message stuff
93   DcoModel * dco_model = dynamic_cast<DcoModel*>(model());
94   // get integer tolerance parameter
95   double tolerance = dco_model->dcoPar()->entry(DcoParams::integerTol);
96   double dist_to_upper = ceil(value) - value;
97   double dist_to_lower = value - floor(value);
98   // return the minimum of distance to upper or lower
99   double res;
100   if (dist_to_upper>dist_to_lower) {
101     res = dist_to_lower;
102   }
103   else {
104     res = dist_to_upper;
105   }
106   if (res<tolerance) {
107     res = 0.0;
108   }
109   return res;
110 }
111 
createCandBranchObjects(BcpsTreeNode * node)112 int DcoBranchStrategyStrong::createCandBranchObjects(BcpsTreeNode * node) {
113   // notes(aykut)
114   // Considers all set of integers for now?
115   // What happens in case of IPM solvers?
116   // consider time limit
117   // iter limit
118 
119   // get node
120   DcoTreeNode * dco_node = dynamic_cast<DcoTreeNode*>(node);
121   // get dco model and message stuff
122   DcoModel * dco_model = dynamic_cast<DcoModel*>(model());
123   CoinMessageHandler * message_handler = dco_model->dcoMessageHandler_;
124   CoinMessages * messages = dco_model->dcoMessages_;
125   // get number of relaxed columns
126   // we assume all relaxed columns are integer variables.
127   int num_relaxed = dco_model->numRelaxedCols();
128   // get indices of relaxed object
129   int const * relaxed = dco_model->relaxedCols();
130   // current solution
131   double * sol = new double[dco_model->solver()->getNumCols()];
132   std::copy(dco_model->solver()->getColSolution(),
133             dco_model->solver()->getColSolution()
134             +dco_model->solver()->getNumCols(),
135             sol);
136 
137   // create data to keep branching objects
138   int cand_cap = dco_model->dcoPar()->entry(DcoParams::strongCandSize);
139   cand_cap = CoinMax(CoinMin(cand_cap, num_relaxed), 1);
140   BcpsBranchObject ** bobjects = new BcpsBranchObject*[cand_cap];
141   int num_bobjects = 0;
142   // pos is the index of the minimum score branch object
143   int min_pos = -1;
144   double min_score = ALPS_INFINITY;
145 
146   dco_model->solver()->markHotStart();
147   dco_model->solver()->setIntParam(OsiMaxNumIterationHotStart, 50);
148 
149   double const obj_val = dco_model->solver()->getObjValue();
150   //double const * collb = dco_model->colLB();
151   //double const * colub = dco_model->colUB();
152   double const * collb = dco_model->solver()->getColLower();
153   double const * colub = dco_model->solver()->getColUpper();
154 
155   // iterate over integer cols, create branch objects, solve corresponding
156   // relaxed problems and set scores.
157   for (int i=0; i<num_relaxed; ++i) {
158     // get corresponding variable
159     int var_index = relaxed[i];
160     // go to next if not infeasible
161     if (!infeas(sol[var_index])) {
162       continue;
163     }
164     // score is 0.0 for now, we will set it later
165     BcpsBranchObject * curr_branch_object =
166       new DcoBranchObject(var_index, 0.0, sol[var_index]);
167     updateScore(curr_branch_object, collb[var_index],
168                 colub[var_index], obj_val);
169     double curr_score = curr_branch_object->score();
170     dco_model->solver()->setColSolution(sol);
171 
172     // if we have capacity add branch object
173     // else check whether current performs better than the worst candidate
174     // if it is add it to candidates.
175     if (num_bobjects<cand_cap) {
176       bobjects[num_bobjects] = curr_branch_object;
177       if (curr_score<min_score) {
178         min_score = curr_score;
179         min_pos = num_bobjects;
180       }
181       num_bobjects++;
182     }
183     else if (curr_score>min_score) {
184       delete bobjects[min_pos];
185       bobjects[min_pos] = curr_branch_object;
186       // find new minimum score candidate
187       min_score = ALPS_INFINITY;
188       for (int k=0; k<cand_cap; ++k) {
189         if (bobjects[k]->score()<min_score) {
190           min_score = bobjects[k]->score();
191           min_pos = k;
192         }
193       }
194     }
195     else {
196       // score is not enough to be a candidate
197     }
198   }
199   delete[] sol;
200   if (num_bobjects==0) {
201     std::cout << "All columns are feasible." << std::endl;
202     throw std::exception();
203   }
204   dco_model->solver()->unmarkHotStart();
205 
206   // debug stuff
207   for (int i=0; i<num_bobjects; ++i) {
208     message_handler->message(DISCO_STRONG_REPORT, *messages)
209       << dco_model->broker()->getProcRank()
210       << bobjects[i]->index()
211       << bobjects[i]->score()
212       << CoinMessageEol;
213   }
214 
215   // add branch objects to branchObjects_
216   setBranchObjects(num_bobjects, bobjects);
217   // bobjects are now owned by BcpsBranchStrategy, do not free them.
218   //bobjects = NULL;
219   // set the branch object member of the node
220   dco_node->setBranchObject(new DcoBranchObject(bestBranchObject()));
221   // compare branch objects and keep the best one at bestBranchObject_
222   return 0;
223 }
224 
225 
226 
227 int
betterBranchObject(BcpsBranchObject const * current,BcpsBranchObject const * other)228 DcoBranchStrategyStrong::betterBranchObject(BcpsBranchObject const * current,
229                                             BcpsBranchObject const * other) {
230   // get model
231   // DcoModel * dco_model = dynamic_cast<DcoModel*>(model());
232   // CoinMessageHandler * message_handler = dco_model->dcoMessageHandler_;
233   // CoinMessages * messages = dco_model->dcoMessages_;
234   int res;
235   if (current->score()>other->score()) {
236     res = 1;
237   }
238   else {
239     res = 0;
240   }
241   return res;
242 }
243