This issue investigates the idea suggested in CoPhy paper that employs the well known technique of Lagrangian Relaxation to transform the formulated BIP to an equivalent BIP that is easier to solve. The trick is to move the constraint of \sum{x{qkia}} = y_{qk} into the objective function. However, it is not straightforward to implement this trick and we need to check the correctness of other variants of BIP.
I encountered this issue when handling node failure constraint, which is generally harder to solve (in compare to the query imbalance and node imbalance constraint). For instance with n = 5, B = 0.5x, CPLEX solver cannot find a solution within 600 seconds. In addition, if I allow CPLEX to run in at most 900 seconds, the machine automatically shuts down (due to the heat issue).