Changeset 994


Ignore:
Timestamp:
Jun 25, 2008 1:40:44 PM (11 years ago)
Author:
ladanyi
Message:

...

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/dynamicbranching/dynamicbranching.cpp

    r993 r994  
    109109  yb_plus_rl_minus_su = 0;
    110110 
    111   if (!isProvenOptimal &&
    112       !isDualObjectiveLimitReached &&
    113       !isIterationLimitReached &&
    114       !isAbandoned &&
    115       !isProvenDualInfeasible &&
    116       !isPrimalObjectiveLimitReached) {
    117     assert(isProvenPrimalInfeasible);
     111  if (isIterationLimitReached ||
     112      isAbandoned ||
     113      isProvenDualInfeasible ||
     114      isPrimalObjectiveLimitReached) {
     115    return; // None of these should happen
     116  }
     117  // If isDualObjectiveLimitReached is true then isProvenOptimal and
     118  // isProvenPrimalInfeasible are correct (search for mustResolve).
     119  // If isDualObjectiveLimitReached false, then the others were correct even
     120  // without resolve. So... we need some data only if we are infeasible, but
     121  // that flag is correct, so we can test it.
     122  if (isProvenPrimalInfeasible) {
    118123    getObjCoefficients = new double[model.getNumCols()];
    119124    CoinDisjointCopyN(model.getObjCoefficients(), model.getNumCols(), getObjCoefficients);
     
    191196private:
    192197  void gutsOfCopy(const DBNodeSimple& rhs);
    193   void gutsOfConstructor (OsiSolverInterface &model,
     198  void gutsOfConstructor (int index, OsiSolverInterface &model,
    194199                          int numberIntegers, int * integer,
    195200                          CoinWarmStart * basis);
     
    203208  // Constructor from current state (and list of integers)
    204209  // Also chooses branching variable (if none set to -1)
    205   DBNodeSimple (OsiSolverInterface &model,
     210  DBNodeSimple (int index, OsiSolverInterface &model,
    206211                int numberIntegers, int * integer,
    207212                CoinWarmStart * basis);
     
    257262 
    258263  // Public data
     264  int index_;
    259265  // The id of the node
    260266  int node_id_;
     
    293299
    294300DBNodeSimple::DBNodeSimple() :
     301  index_(-1),
    295302  node_id_(-1),
    296303  basis_(NULL),
     
    311318{
    312319}
    313 DBNodeSimple::DBNodeSimple(OsiSolverInterface & model,
     320DBNodeSimple::DBNodeSimple(int index, OsiSolverInterface & model,
    314321                           int numberIntegers, int * integer,CoinWarmStart * basis)
    315322{
    316   gutsOfConstructor(model,numberIntegers,integer,basis);
     323  gutsOfConstructor(index, model,numberIntegers,integer,basis);
    317324}
    318325void
    319 DBNodeSimple::gutsOfConstructor(OsiSolverInterface & model,
     326DBNodeSimple::gutsOfConstructor(int index, OsiSolverInterface & model,
    320327                                int numberIntegers, int * integer,CoinWarmStart * basis)
    321328{
     329  index_ = index;
    322330  node_id_ = -1;
    323331  basis_ = basis;
     
    582590DBNodeSimple::gutsOfCopy(const DBNodeSimple & rhs)
    583591{
     592  index_=rhs.index_;
    584593  node_id_=rhs.node_id_;
    585594  basis_ = rhs.basis_ ? rhs.basis_->clone() : NULL;
     
    962971  }
    963972   
    964   if (lpres.isProvenOptimal) {
     973  if (lpres.isProvenOptimal && !lpres.isDualObjectiveLimitReached) {
    965974    // THINK: should we do anything? like:
    966975#if 0
     
    977986   
    978987  if (lpres.isDualObjectiveLimitReached) {
    979     // Dual feasible, and in this case we don't care how we have
    980     // stopped (iteration limit, obj val limit, time limit, optimal solution,
    981     // etc.), we can just look at the reduced costs to figure out if the
    982     // grandparent is irrelevant. Remember, we are using dual simplex!
    983     double djValue = lpres.getReducedCost[GP_brvar_fullid]*direction;
    984     if (djValue > 1.0e-6) {
    985       // wants to go down
    986       if (parent_is_down_child) {
     988    // Dual feasible, Because of reswolving (search for mustResolve) the
     989    // problem here is either optimal or infeasible. If infeasible then we
     990    // just need to fall out of this, if optimal then we test the reduced
     991    // cost.
     992    if (lpres.isProvenOptimal) {
     993      double djValue = lpres.getReducedCost[GP_brvar_fullid]*direction;
     994      if (djValue > 1.0e-6) {
     995        // wants to go down
     996        if (parent_is_down_child) {
     997          return true;
     998        }
     999        if (lpres.getColLower[GP_brvar_fullid] > std::ceil(grandparent.value_)) {
     1000          return true;
     1001        }
     1002      } else if (djValue < -1.0e-6) {
     1003        // wants to go up
     1004        if (! parent_is_down_child) {
     1005          return true;
     1006        }
     1007        if (lpres.getColUpper[GP_brvar_fullid] < std::floor(grandparent.value_)) {
     1008          return true;
     1009        }
     1010      }
     1011      return false;
     1012    }
     1013  }
     1014
     1015  // Problem is really infeasible and has a dual ray.
     1016  assert(lpres.isProvenPrimalInfeasible);
     1017  const int greatgrandparent_id = grandparent.parent_;
     1018  const int x = GP_brvar_fullid; // for easier reference... we'll use s_x
     1019
     1020  /*
     1021    Now we are going to check that if we relax the GP's branching decision
     1022    then the child's LP relaxation is still infeasible. The test is done by
     1023    checking whether the column (or its negative) of the GP's branching
     1024    variable cuts off the dual ray proving the infeasibility.
     1025  */
     1026   
     1027  const double* dj = lpres.getReducedCost;
     1028  const double* obj = lpres.getObjCoefficients;
     1029  const double yA_x = dj[x] - obj[x];
     1030
     1031  if (direction > 0) { // minimization problem
     1032    if (parent_is_down_child) {
     1033      const double s_x = CoinMax(yA_x, -1e-8);
     1034      if (s_x < 1e-6) {
    9871035        return true;
    9881036      }
    989       if (lpres.getColLower[GP_brvar_fullid] > std::ceil(grandparent.value_)) {
     1037      if (lpres.yb_plus_rl_minus_su < 1e-8) {
     1038        printf("WARNING: lpres.yb_plus_rl_minus_su is not positive!\n");
     1039        return false;
     1040      }
     1041      const double max_u_x = lpres.yb_plus_rl_minus_su / s_x + lpres.getColUpper[x];
     1042      const double u_x_without_GP = greatgrandparent_id >= 0 ?
     1043        branchingTree.nodes_[greatgrandparent_id].upper_[GP_brvar] :
     1044        original_upper[GP_brvar];
     1045      return max_u_x > u_x_without_GP - 1e-8;
     1046    } else {
     1047      const double r_x = CoinMax(yA_x, -1e-8);
     1048      if (r_x < 1e-6) {
    9901049        return true;
    9911050      }
    992     } else if (djValue < -1.0e-6) {
    993       // wants to go up
    994       if (! parent_is_down_child) {
     1051      if (lpres.yb_plus_rl_minus_su < 1e-8) {
     1052        printf("WARNING: lpres.yb_plus_rl_minus_su is not positive!\n");
     1053        return false;
     1054      }
     1055      const double min_l_x = - lpres.yb_plus_rl_minus_su / r_x + lpres.getColLower[x];
     1056      const double l_x_without_GP = greatgrandparent_id >= 0 ?
     1057        branchingTree.nodes_[greatgrandparent_id].lower_[GP_brvar] :
     1058        original_lower[GP_brvar];
     1059      return min_l_x < l_x_without_GP + 1e-8;
     1060    }
     1061  } else { // maximization problem
     1062    if (parent_is_down_child) {
     1063      const double s_x = CoinMin(yA_x, 1e-8);
     1064      if (s_x > -1e-6) {
    9951065        return true;
    9961066      }
    997       if (lpres.getColUpper[GP_brvar_fullid] < std::floor(grandparent.value_)) {
     1067      if (lpres.yb_plus_rl_minus_su > -1e-8) {
     1068        printf("WARNING: lpres.yb_plus_rl_minus_su is not negative!\n");
     1069        return false;
     1070      }
     1071      const double max_u_x = lpres.yb_plus_rl_minus_su / s_x + lpres.getColUpper[x];
     1072      const double u_x_without_GP = greatgrandparent_id >= 0 ?
     1073        branchingTree.nodes_[greatgrandparent_id].upper_[GP_brvar] :
     1074        original_upper[GP_brvar];
     1075      return max_u_x > u_x_without_GP - 1e-8;
     1076    } else {
     1077      const double r_x = CoinMin(yA_x, 1e-8);
     1078      if (r_x < -1e-6) {
    9981079        return true;
    9991080      }
    1000     }
    1001     return false;
    1002   } else {
    1003     assert(lpres.isProvenPrimalInfeasible);
    1004     return false;
    1005     const int greatgrandparent_id = grandparent.parent_;
    1006     const int x = GP_brvar_fullid; // for easier reference... we'll use s_x
    1007 
    1008     /*
    1009       Now we are going to check that if we relax the GP's branching decision
    1010       then the child's LP relaxation is still infeasible. The test is done by
    1011       checking whether the column (or its negative) of the GP's branching
    1012       variable cuts off the dual ray proving the infeasibility.
    1013     */
    1014    
    1015     const double* dj = lpres.getReducedCost;
    1016     const double* obj = lpres.getObjCoefficients;
    1017     const double yA_x = dj[x] - obj[x];
    1018 
    1019     if (direction > 0) { // minimization problem
    1020       if (parent_is_down_child) {
    1021         const double s_x = CoinMax(yA_x, -1e-8);
    1022         if (s_x < 1e-6) {
    1023           return true;
    1024         }
    1025         if (lpres.yb_plus_rl_minus_su < 1e-8) {
    1026           printf("WARNING: lpres.yb_plus_rl_minus_su is not positive!\n");
    1027           return false;
    1028         }
    1029         const double max_u_x = lpres.yb_plus_rl_minus_su / s_x + lpres.getColUpper[x];
    1030         const double u_x_without_GP = greatgrandparent_id >= 0 ?
    1031           branchingTree.nodes_[greatgrandparent_id].upper_[GP_brvar] :
    1032           original_upper[GP_brvar];
    1033         return max_u_x > u_x_without_GP - 1e-8;
    1034       } else {
    1035         const double r_x = CoinMax(yA_x, -1e-8);
    1036         if (r_x < 1e-6) {
    1037           return true;
    1038         }
    1039         if (lpres.yb_plus_rl_minus_su < 1e-8) {
    1040           printf("WARNING: lpres.yb_plus_rl_minus_su is not positive!\n");
    1041           return false;
    1042         }
    1043         const double min_l_x = - lpres.yb_plus_rl_minus_su / r_x + lpres.getColLower[x];
    1044         const double l_x_without_GP = greatgrandparent_id >= 0 ?
    1045           branchingTree.nodes_[greatgrandparent_id].lower_[GP_brvar] :
    1046           original_lower[GP_brvar];
    1047         return min_l_x < l_x_without_GP + 1e-8;
    1048       }
    1049     } else { // maximization problem
    1050       if (parent_is_down_child) {
    1051         const double s_x = CoinMin(yA_x, 1e-8);
    1052         if (s_x > -1e-6) {
    1053           return true;
    1054         }
    1055         if (lpres.yb_plus_rl_minus_su > -1e-8) {
    1056           printf("WARNING: lpres.yb_plus_rl_minus_su is not negative!\n");
    1057           return false;
    1058         }
    1059         const double max_u_x = lpres.yb_plus_rl_minus_su / s_x + lpres.getColUpper[x];
    1060         const double u_x_without_GP = greatgrandparent_id >= 0 ?
    1061           branchingTree.nodes_[greatgrandparent_id].upper_[GP_brvar] :
    1062           original_upper[GP_brvar];
    1063         return max_u_x > u_x_without_GP - 1e-8;
    1064       } else {
    1065         const double r_x = CoinMin(yA_x, 1e-8);
    1066         if (r_x < -1e-6) {
    1067           return true;
    1068         }
    1069         if (lpres.yb_plus_rl_minus_su > -1e-8) {
    1070           printf("WARNING: lpres.yb_plus_rl_minus_su is not negative!\n");
    1071           return false;
    1072         }
    1073         const double min_l_x = - lpres.yb_plus_rl_minus_su / r_x + lpres.getColLower[x];
    1074         const double l_x_without_GP = greatgrandparent_id >= 0 ?
    1075           branchingTree.nodes_[greatgrandparent_id].lower_[GP_brvar] :
    1076           original_lower[GP_brvar];
    1077         return min_l_x < l_x_without_GP + 1e-8;
    1078       }
     1081      if (lpres.yb_plus_rl_minus_su > -1e-8) {
     1082        printf("WARNING: lpres.yb_plus_rl_minus_su is not negative!\n");
     1083        return false;
     1084      }
     1085      const double min_l_x = - lpres.yb_plus_rl_minus_su / r_x + lpres.getColLower[x];
     1086      const double l_x_without_GP = greatgrandparent_id >= 0 ?
     1087        branchingTree.nodes_[greatgrandparent_id].lower_[GP_brvar] :
     1088        original_lower[GP_brvar];
     1089      return min_l_x < l_x_without_GP + 1e-8;
    10791090    }
    10801091  }
     
    12541265        parent.way_ = DBNodeSimple::WAY_DOWN_UP__DOWN_DONE;
    12551266        sizeDeferred_--;
     1267        return; // No bound changes are needed on the GO side as GP is
     1268                // deleted...
    12561269      } else {
    12571270        grandparent.way_ = parent_is_down_child ?
     
    12781291        parent.way_ = DBNodeSimple::WAY_UP_DOWN__UP_DONE;
    12791292        sizeDeferred_--;
     1293        return;
    12801294      } else {
    12811295        grandparent.way_ = parent_is_down_child ?
     
    14501464   
    14511465    // Add continuous to it;
    1452     DBNodeSimple rootNode(model,numberIntegers,which,model.getWarmStart());
     1466    DBNodeSimple rootNode(0,model,numberIntegers,which,model.getWarmStart());
    14531467    // something extra may have been fixed by strong branching
    14541468    // if so go round again
    14551469    while (rootNode.variable_==numberIntegers) {
    14561470      model.resolve();
    1457       rootNode = DBNodeSimple(model,numberIntegers,which,model.getWarmStart());
     1471      rootNode = DBNodeSimple(0,model,numberIntegers,which,model.getWarmStart());
    14581472    }
    14591473    if (rootNode.objectiveValue_<1.0e100) {
     
    15981612        LPresult lpres(model);
    15991613        if (mustResolve) {
     1614          lpres.isDualObjectiveLimitReached = true;
    16001615          model.setDblParam(OsiDualObjectiveLimit, oldlimit);
    16011616        }
     
    16361651        }
    16371652        if (cnt > 0) {
     1653          printf("Before switch: opt: %c   inf: %c   dual_bd: %c\n",
     1654                 lpres.isProvenOptimal ? 'T' : 'F',
     1655                 lpres.isProvenPrimalInfeasible ? 'T' : 'F',
     1656                 lpres.isDualObjectiveLimitReached ? 'T' : 'F');
    16381657          model.resolve();
    16391658          // This is horribly looking... Get rid of it when properly
    16401659          // debugged...
    16411660          LPresult lpres1(model);
     1661          printf("After resolve: opt: %c   inf: %c   dual_bd: %c\n",
     1662                 lpres1.isProvenOptimal ? 'T' : 'F',
     1663                 lpres1.isProvenPrimalInfeasible ? 'T' : 'F',
     1664                 lpres1.isDualObjectiveLimitReached ? 'T' : 'F');
    16421665          assert(lpres.isAbandoned == model.isAbandoned());
    16431666          assert(lpres.isDualObjectiveLimitReached == model.isDualObjectiveLimitReached());
     
    16741697          exit(77);
    16751698        }
    1676         DBNodeSimple newNode(model,numberIntegers,which,ws);
     1699        DBNodeSimple newNode(numberNodes, model,numberIntegers,which,ws);
    16771700        // something extra may have been fixed by strong branching
    16781701        // if so go round again
    16791702        while (newNode.variable_==numberIntegers) {
    16801703          model.resolve();
    1681           newNode = DBNodeSimple(model,numberIntegers,which,model.getWarmStart());
     1704          newNode = DBNodeSimple(numberNodes, model,numberIntegers,which,model.getWarmStart());
    16821705          newNode.strong_branching_fixed_vars_ = true;
    16831706        }
Note: See TracChangeset for help on using the changeset viewer.