diff --git a/extern/zstr/strict_fstream.hpp b/extern/zstr/strict_fstream.hpp index 6c3ea85d90..953311a907 100644 --- a/extern/zstr/strict_fstream.hpp +++ b/extern/zstr/strict_fstream.hpp @@ -17,7 +17,7 @@ namespace strict_fstream { -// Help people out a bit, it seems like this is a common recommenation since +// Help people out a bit, it seems like this is a common recommendation since // musl breaks all over the place. #if defined(__NEED_size_t) && !defined(__MUSL__) #warning "It seems to be recommended to patch in a define for __MUSL__ if you use musl globally: https://www.openwall.com/lists/musl/2013/02/10/5" diff --git a/src/ipm/ipx/crossover.h b/src/ipm/ipx/crossover.h index 3587014353..7e16b00eb1 100644 --- a/src/ipm/ipx/crossover.h +++ b/src/ipm/ipx/crossover.h @@ -44,7 +44,7 @@ class Crossover { // as long as the Crossover object is used. Crossover(const Control& control); - // First runs the dual push phase; if this was succesful, then runs the + // First runs the dual push phase; if this was successful, then runs the // primal push phase. // // weights: Must either be NULL or an array of size n+m. diff --git a/src/ipm/ipx/forrest_tomlin.h b/src/ipm/ipx/forrest_tomlin.h index 32c373080b..dc60bfd6ab 100644 --- a/src/ipm/ipx/forrest_tomlin.h +++ b/src/ipm/ipx/forrest_tomlin.h @@ -11,7 +11,7 @@ namespace ipx { // Generic implementation of the Forrest-Tomlin update [1] that can be used with // any LU factorization. The implementation does not exploit hypersparsity, -// which exludes its use for such problems. For non-hypersparse problems the +// which excludes its use for such problems. For non-hypersparse problems the // implementation is better suited than BASICLU, however, because it stores L // and U in compressed form with permuted indices; hence solving triangular // systems with a dense rhs/lhs accesses memory contiguously. BASICLU could not diff --git a/src/ipm/ipx/indexed_vector.h b/src/ipm/ipx/indexed_vector.h index 6ad46af292..07a0864be1 100644 --- a/src/ipm/ipx/indexed_vector.h +++ b/src/ipm/ipx/indexed_vector.h @@ -24,7 +24,7 @@ namespace ipx { // otherwise. // // When modifying the vector changes its pattern (e.g. by writing to v[i] for an -// arbitray index i), you have to invalidate the pattern or provide the new one. +// arbitrary index i), you have to invalidate the pattern or provide the new one. class IndexedVector { public: diff --git a/src/ipm/ipx/iterate.cc b/src/ipm/ipx/iterate.cc index 755e143f8f..c407eb3de4 100644 --- a/src/ipm/ipx/iterate.cc +++ b/src/ipm/ipx/iterate.cc @@ -279,7 +279,7 @@ void Iterate::Postprocess() { // For fixed variables compute xl[j] and xu[j] from x[j]. If the lower and // upper bound are equal, set zl[j] or zu[j] such that the variable is dual - // feasibile. Otherwise leave them zero. + // feasible. Otherwise leave them zero. for (Int j = 0; j < n+m; j++) { if (StateOf(j) == State::fixed) { xl_[j] = x_[j] - lb[j]; diff --git a/src/ipm/ipx/iterate.h b/src/ipm/ipx/iterate.h index 53b05cce37..0f1ce529d4 100644 --- a/src/ipm/ipx/iterate.h +++ b/src/ipm/ipx/iterate.h @@ -77,7 +77,7 @@ class Iterate { double zl(Int j) const { return zl_[j]; } double zu(Int j) const { return zu_[j]; } - // Returns const rerefences to the residual vectors + // Returns const references to the residual vectors // rb = b-AI*x, // rl = lb-x+xl, // ru = ub-x-xu, @@ -155,7 +155,7 @@ class Iterate { double presidual() const; double dresidual() const; - // copmlementarity() returns the sum of the pairwise complementarity + // complementarity() returns the sum of the pairwise complementarity // products xl[j]*zl[j] and xu[j]*zu[j] from all barrier terms. mu() // returns the average, mu_min() the minimum and mu_max() the maximum. double complementarity() const; diff --git a/src/ipm/ipx/lu_factorization.h b/src/ipm/ipx/lu_factorization.h index 6effb8b094..5bc4340678 100644 --- a/src/ipm/ipx/lu_factorization.h +++ b/src/ipm/ipx/lu_factorization.h @@ -38,7 +38,7 @@ class LuFactorization { // kLuDependencyTol as absolute pivot tolerance and to // remove columns from the active submatrix // immediately when all entries became smaller than - // the abolute pivot tolerance. Need not be supported + // the absolute pivot tolerance. Need not be supported // by the implementation. // @L, @U: return the matrix factors with sorted indices. The objects are // resized as necessary. diff --git a/src/ipm/ipx/lu_update.h b/src/ipm/ipx/lu_update.h index 2ef13e0228..299ec4d8c0 100644 --- a/src/ipm/ipx/lu_update.h +++ b/src/ipm/ipx/lu_update.h @@ -25,7 +25,7 @@ class LuUpdate { // kLuDependencyTol as absolute pivot tolerance and to // remove columns from the active submatrix // immediately when all entries became smaller than - // the abolute pivot tolerance. Need not be supported + // the absolute pivot tolerance. Need not be supported // by the implementation. // // Factorize() cannot fail other than for out of memory, in which case diff --git a/src/ipm/ipx/sparse_utils.h b/src/ipm/ipx/sparse_utils.h index 54bdbc782c..db45600ef7 100644 --- a/src/ipm/ipx/sparse_utils.h +++ b/src/ipm/ipx/sparse_utils.h @@ -16,7 +16,7 @@ namespace ipx { // have previously been unmarked; they are marked now and newtop // is returned. // @marked, @marker Node i is "marked" iff @marked[i] == @marker. -// @work worksapce of size # rows of A. +// @work workspace of size # rows of A. // // The code has been copied and adapted from cs_dfs.c, included in the CSPARSE // package [1]. diff --git a/src/lp_data/HighsInterface.cpp b/src/lp_data/HighsInterface.cpp index 7f7a201f99..d5ffd342a0 100644 --- a/src/lp_data/HighsInterface.cpp +++ b/src/lp_data/HighsInterface.cpp @@ -2198,7 +2198,7 @@ void Highs::formIllConditioningLp0(HighsLp& ill_conditioning_lp, ill_conditioning_matrix.value_.push_back(1.0); ill_conditioning_matrix.start_.push_back( HighsInt(ill_conditioning_matrix.index_.size())); - // Subracting x_- with cost 1 + // Subtracting x_- with cost 1 ill_conditioning_lp.col_cost_.push_back(1); ill_conditioning_lp.col_lower_.push_back(0); ill_conditioning_lp.col_upper_.push_back(kHighsInf); @@ -2312,7 +2312,7 @@ void Highs::formIllConditioningLp1(HighsLp& ill_conditioning_lp, } assert(ill_conditioning_lp.num_col_ == incumbent_num_row); if (constraint) { - // Add the identiy matrix for constraint y - u + w = 0 + // Add the identity matrix for constraint y - u + w = 0 for (HighsInt iRow = 0; iRow < incumbent_num_row; iRow++) { ill_conditioning_matrix.index_.push_back(iRow); ill_conditioning_matrix.value_.push_back(1.0); diff --git a/src/lp_data/HighsSolve.cpp b/src/lp_data/HighsSolve.cpp index 29377e7fec..d839bfb88d 100644 --- a/src/lp_data/HighsSolve.cpp +++ b/src/lp_data/HighsSolve.cpp @@ -254,7 +254,7 @@ HighsStatus solveUnconstrainedLp(const HighsOptions& options, const HighsLp& lp, if (lp.num_row_ > 0) { // Assign primal, dual and basis status for rows, checking for - // infeasiblility + // infeasibility for (HighsInt iRow = 0; iRow < lp.num_row_; iRow++) { double primal_infeasibility = 0; double lower = lp.row_lower_[iRow]; diff --git a/src/mip/HighsDomain.cpp b/src/mip/HighsDomain.cpp index 8f084fe1f0..a00e8784c0 100644 --- a/src/mip/HighsDomain.cpp +++ b/src/mip/HighsDomain.cpp @@ -46,6 +46,29 @@ static double activityContributionMax(double coef, const double& lb, } } +static double computeDelta(double val, double oldbound, double newbound, + double inf, HighsInt& numinfs) { + // if bounds are huge, HighsCDouble should be used when computing bound + // differences. todo: qualify usage of HighsCDouble in this function. + if (oldbound == inf) { + --numinfs; + return newbound * val; + } else if (newbound == inf) { + ++numinfs; + return -oldbound * val; + } else { + return (newbound - oldbound) * val; + } +} + +static inline double boundRange(double upper_bound, double lower_bound, + double tolerance, HighsVarType var_type) { + double range = upper_bound - lower_bound; + return range - (var_type == HighsVarType::kContinuous + ? std::max(0.3 * range, 1000.0 * tolerance) + : tolerance); +} + HighsDomain::HighsDomain(HighsMipSolver& mipsolver) : mipsolver(&mipsolver) { col_lower_ = mipsolver.model_->col_lower_; col_upper_ = mipsolver.model_->col_upper_; @@ -369,14 +392,11 @@ void HighsDomain::CutpoolPropagation::recomputeCapacityThreshold(HighsInt cut) { if (domain->col_upper_[arindex[i]] == domain->col_lower_[arindex[i]]) continue; - double boundRange = - domain->col_upper_[arindex[i]] - domain->col_lower_[arindex[i]]; - - boundRange -= domain->variableType(arindex[i]) == HighsVarType::kContinuous - ? std::max(0.3 * boundRange, 1000.0 * domain->feastol()) - : domain->feastol(); - - double threshold = std::fabs(arvalue[i]) * boundRange; + double threshold = + std::fabs(arvalue[i]) * boundRange(domain->col_upper_[arindex[i]], + domain->col_lower_[arindex[i]], + domain->feastol(), + domain->variableType(arindex[i])); capacityThreshold_[cut] = std::max({capacityThreshold_[cut], threshold, domain->feastol()}); @@ -461,17 +481,8 @@ void HighsDomain::CutpoolPropagation::updateActivityLbChange(HighsInt col, cutpool->getMatrix().forEachPositiveColumnEntry( col, [&](HighsInt row, double val) { assert(val > 0); - double deltamin; - - if (oldbound == -kHighsInf) { - --activitycutsinf_[row]; - deltamin = newbound * val; - } else if (newbound == -kHighsInf) { - ++activitycutsinf_[row]; - deltamin = -oldbound * val; - } else { - deltamin = (newbound - oldbound) * val; - } + double deltamin = computeDelta(val, oldbound, newbound, -kHighsInf, + activitycutsinf_[row]); activitycuts_[row] += deltamin; if (deltamin <= 0) { @@ -504,18 +515,8 @@ void HighsDomain::CutpoolPropagation::updateActivityLbChange(HighsInt col, cutpool->getMatrix().forEachPositiveColumnEntry( col, [&](HighsInt row, double val) { assert(val > 0); - double deltamin; - - if (oldbound == -kHighsInf) { - --activitycutsinf_[row]; - deltamin = newbound * val; - } else if (newbound == -kHighsInf) { - ++activitycutsinf_[row]; - deltamin = -oldbound * val; - } else { - deltamin = (newbound - oldbound) * val; - } - activitycuts_[row] += deltamin; + activitycuts_[row] += computeDelta(val, oldbound, newbound, + -kHighsInf, activitycutsinf_[row]); if (domain->infeasible_reason.index == row) return false; @@ -541,17 +542,8 @@ void HighsDomain::CutpoolPropagation::updateActivityUbChange(HighsInt col, cutpool->getMatrix().forEachNegativeColumnEntry( col, [&](HighsInt row, double val) { assert(val < 0); - double deltamin; - - if (oldbound == kHighsInf) { - --activitycutsinf_[row]; - deltamin = newbound * val; - } else if (newbound == kHighsInf) { - ++activitycutsinf_[row]; - deltamin = -oldbound * val; - } else { - deltamin = (newbound - oldbound) * val; - } + double deltamin = computeDelta(val, oldbound, newbound, kHighsInf, + activitycutsinf_[row]); activitycuts_[row] += deltamin; if (deltamin <= 0) { @@ -582,18 +574,8 @@ void HighsDomain::CutpoolPropagation::updateActivityUbChange(HighsInt col, cutpool->getMatrix().forEachNegativeColumnEntry( col, [&](HighsInt row, double val) { assert(val < 0); - double deltamin; - - if (oldbound == kHighsInf) { - --activitycutsinf_[row]; - deltamin = newbound * val; - } else if (newbound == kHighsInf) { - ++activitycutsinf_[row]; - deltamin = -oldbound * val; - } else { - deltamin = (newbound - oldbound) * val; - } - activitycuts_[row] += deltamin; + activitycuts_[row] += computeDelta(val, oldbound, newbound, kHighsInf, + activitycutsinf_[row]); if (domain->infeasible_reason.index == row) return false; @@ -780,12 +762,11 @@ void HighsDomain::ObjectivePropagation::recomputeCapacityThreshold() { for (HighsInt i = partitionStarts[numPartitions]; i < numObjNzs; ++i) { HighsInt col = objNonzeros[i]; - double boundRange = (domain->col_upper_[col] - domain->col_lower_[col]); - boundRange -= domain->variableType(col) == HighsVarType::kContinuous - ? std::max(0.3 * boundRange, 1000.0 * domain->feastol()) - : domain->feastol(); - capacityThreshold = - std::max(capacityThreshold, std::fabs(cost[col]) * boundRange); + capacityThreshold = std::max( + capacityThreshold, + std::fabs(cost[col]) * + boundRange(domain->col_upper_[col], domain->col_lower_[col], + domain->feastol(), domain->variableType(col))); } } @@ -793,11 +774,11 @@ void HighsDomain::ObjectivePropagation::updateActivityLbChange( HighsInt col, double oldbound, double newbound) { if (cost[col] <= 0.0) { if (cost[col] != 0.0 && newbound < oldbound) { - double boundRange = domain->col_upper_[col] - newbound; - boundRange -= domain->variableType(col) == HighsVarType::kContinuous - ? std::max(0.3 * boundRange, 1000.0 * domain->feastol()) - : domain->feastol(); - capacityThreshold = std::max(capacityThreshold, -cost[col] * boundRange); + capacityThreshold = + std::max(capacityThreshold, + -cost[col] * boundRange(domain->col_upper_[col], newbound, + domain->feastol(), + domain->variableType(col))); isPropagated = false; } debugCheckObjectiveLower(); @@ -821,11 +802,11 @@ void HighsDomain::ObjectivePropagation::updateActivityLbChange( debugCheckObjectiveLower(); if (newbound < oldbound) { - double boundRange = (domain->col_upper_[col] - domain->col_lower_[col]); - boundRange -= domain->variableType(col) == HighsVarType::kContinuous - ? std::max(0.3 * boundRange, 1000.0 * domain->feastol()) - : domain->feastol(); - capacityThreshold = std::max(capacityThreshold, cost[col] * boundRange); + capacityThreshold = std::max( + capacityThreshold, + cost[col] * boundRange(domain->col_upper_[col], + domain->col_lower_[col], domain->feastol(), + domain->variableType(col))); } else if (numInfObjLower == 0 && objectiveLower > domain->mipsolver->mipdata_->upper_limit) { domain->infeasible_ = true; @@ -914,11 +895,10 @@ void HighsDomain::ObjectivePropagation::updateActivityUbChange( HighsInt col, double oldbound, double newbound) { if (cost[col] >= 0.0) { if (cost[col] != 0.0 && newbound > oldbound) { - double boundRange = newbound - domain->col_lower_[col]; - boundRange -= domain->variableType(col) == HighsVarType::kContinuous - ? std::max(0.3 * boundRange, 1000.0 * domain->feastol()) - : domain->feastol(); - capacityThreshold = std::max(capacityThreshold, cost[col] * boundRange); + capacityThreshold = std::max( + capacityThreshold, + cost[col] * boundRange(newbound, domain->col_lower_[col], + domain->feastol(), domain->variableType(col))); isPropagated = false; } debugCheckObjectiveLower(); @@ -942,11 +922,11 @@ void HighsDomain::ObjectivePropagation::updateActivityUbChange( debugCheckObjectiveLower(); if (newbound > oldbound) { - double boundRange = (domain->col_upper_[col] - domain->col_lower_[col]); - boundRange -= domain->variableType(col) == HighsVarType::kContinuous - ? std::max(0.3 * boundRange, 1000.0 * domain->feastol()) - : domain->feastol(); - capacityThreshold = std::max(capacityThreshold, -cost[col] * boundRange); + capacityThreshold = std::max( + capacityThreshold, + -cost[col] * boundRange(domain->col_upper_[col], + domain->col_lower_[col], domain->feastol(), + domain->variableType(col))); } else if (numInfObjLower == 0 && objectiveLower > domain->mipsolver->mipdata_->upper_limit) { domain->infeasible_ = true; @@ -1287,7 +1267,6 @@ void HighsDomain::computeMinActivity(HighsInt start, HighsInt end, activitymin += contributionmin; } - activitymin.renormalize(); } else { activitymin = 0.0; ninfmin = 0; @@ -1305,9 +1284,8 @@ void HighsDomain::computeMinActivity(HighsInt start, HighsInt end, else activitymin += contributionmin; } - - activitymin.renormalize(); } + activitymin.renormalize(); } void HighsDomain::computeMaxActivity(HighsInt start, HighsInt end, @@ -1333,8 +1311,6 @@ void HighsDomain::computeMaxActivity(HighsInt start, HighsInt end, else activitymax += contributionmin; } - - activitymax.renormalize(); } else { activitymax = 0.0; ninfmax = 0; @@ -1352,9 +1328,8 @@ void HighsDomain::computeMaxActivity(HighsInt start, HighsInt end, else activitymax += contributionmin; } - - activitymax.renormalize(); } + activitymax.renormalize(); } double HighsDomain::adjustedUb(HighsInt col, HighsCDouble boundVal, @@ -1362,13 +1337,10 @@ double HighsDomain::adjustedUb(HighsInt col, HighsCDouble boundVal, double bound; if (mipsolver->variableType(col) != HighsVarType::kContinuous) { - bound = std::floor(double(boundVal + mipsolver->mipdata_->feastol)); - if (bound < col_upper_[col] && - col_upper_[col] - bound > - 1000.0 * mipsolver->mipdata_->feastol * std::fabs(bound)) - accept = true; - else - accept = false; + bound = static_cast(floor(boundVal + mipsolver->mipdata_->feastol)); + accept = bound < col_upper_[col] && + col_upper_[col] - bound > + 1000.0 * mipsolver->mipdata_->feastol * std::fabs(bound); } else { if (std::fabs(double(boundVal) - col_lower_[col]) <= mipsolver->mipdata_->epsilon) @@ -1397,13 +1369,10 @@ double HighsDomain::adjustedLb(HighsInt col, HighsCDouble boundVal, double bound; if (mipsolver->variableType(col) != HighsVarType::kContinuous) { - bound = std::ceil(double(boundVal - mipsolver->mipdata_->feastol)); - if (bound > col_lower_[col] && - bound - col_lower_[col] > - 1000.0 * mipsolver->mipdata_->feastol * std::fabs(bound)) - accept = true; - else - accept = false; + bound = static_cast(ceil(boundVal - mipsolver->mipdata_->feastol)); + accept = bound > col_lower_[col] && + bound - col_lower_[col] > + 1000.0 * mipsolver->mipdata_->feastol * std::fabs(bound); } else { if (std::fabs(col_upper_[col] - double(boundVal)) <= mipsolver->mipdata_->epsilon) @@ -1517,14 +1486,10 @@ HighsInt HighsDomain::propagateRowLower(const HighsInt* Rindex, void HighsDomain::updateThresholdLbChange(HighsInt col, double newbound, double val, double& threshold) { if (newbound != col_upper_[col]) { - double boundRange = (col_upper_[col] - newbound); - - boundRange -= - variableType(col) == HighsVarType::kContinuous - ? std::max(0.3 * boundRange, 1000.0 * mipsolver->mipdata_->feastol) - : mipsolver->mipdata_->feastol; - - double thresholdNew = std::fabs(val) * boundRange; + double thresholdNew = + std::fabs(val) * boundRange(col_upper_[col], newbound, + mipsolver->mipdata_->feastol, + variableType(col)); // the new threshold is now the maximum of the new threshold and the current // one @@ -1536,14 +1501,10 @@ void HighsDomain::updateThresholdLbChange(HighsInt col, double newbound, void HighsDomain::updateThresholdUbChange(HighsInt col, double newbound, double val, double& threshold) { if (newbound != col_lower_[col]) { - double boundRange = (newbound - col_lower_[col]); - - boundRange -= - variableType(col) == HighsVarType::kContinuous - ? std::max(0.3 * boundRange, 1000.0 * mipsolver->mipdata_->feastol) - : mipsolver->mipdata_->feastol; - - double thresholdNew = std::fabs(val) * boundRange; + double thresholdNew = + std::fabs(val) * boundRange(newbound, col_lower_[col], + mipsolver->mipdata_->feastol, + variableType(col)); // the new threshold is now the maximum of the new threshold and the current // one @@ -1567,16 +1528,9 @@ void HighsDomain::updateActivityLbChange(HighsInt col, double oldbound, for (HighsInt i = start; i != end; ++i) { if (mip->a_matrix_.value_[i] > 0) { - double deltamin; - if (oldbound == -kHighsInf) { - --activitymininf_[mip->a_matrix_.index_[i]]; - deltamin = newbound * mip->a_matrix_.value_[i]; - } else if (newbound == -kHighsInf) { - ++activitymininf_[mip->a_matrix_.index_[i]]; - deltamin = -oldbound * mip->a_matrix_.value_[i]; - } else { - deltamin = (newbound - oldbound) * mip->a_matrix_.value_[i]; - } + double deltamin = + computeDelta(mip->a_matrix_.value_[i], oldbound, newbound, -kHighsInf, + activitymininf_[mip->a_matrix_.index_[i]]); activitymin_[mip->a_matrix_.index_[i]] += deltamin; #ifndef NDEBUG @@ -1618,16 +1572,9 @@ void HighsDomain::updateActivityLbChange(HighsInt col, double oldbound, mip->row_upper_[mip->a_matrix_.index_[i]] != kHighsInf) markPropagate(mip->a_matrix_.index_[i]); } else { - double deltamax; - if (oldbound == -kHighsInf) { - --activitymaxinf_[mip->a_matrix_.index_[i]]; - deltamax = newbound * mip->a_matrix_.value_[i]; - } else if (newbound == -kHighsInf) { - ++activitymaxinf_[mip->a_matrix_.index_[i]]; - deltamax = -oldbound * mip->a_matrix_.value_[i]; - } else { - deltamax = (newbound - oldbound) * mip->a_matrix_.value_[i]; - } + double deltamax = + computeDelta(mip->a_matrix_.value_[i], oldbound, newbound, -kHighsInf, + activitymaxinf_[mip->a_matrix_.index_[i]]); activitymax_[mip->a_matrix_.index_[i]] += deltamax; #ifndef NDEBUG @@ -1684,29 +1631,13 @@ void HighsDomain::updateActivityLbChange(HighsInt col, double oldbound, std::swap(oldbound, newbound); for (HighsInt i = start; i != end; ++i) { if (mip->a_matrix_.value_[i] > 0) { - double deltamin; - if (oldbound == -kHighsInf) { - --activitymininf_[mip->a_matrix_.index_[i]]; - deltamin = newbound * mip->a_matrix_.value_[i]; - } else if (newbound == -kHighsInf) { - ++activitymininf_[mip->a_matrix_.index_[i]]; - deltamin = -oldbound * mip->a_matrix_.value_[i]; - } else { - deltamin = (newbound - oldbound) * mip->a_matrix_.value_[i]; - } - activitymin_[mip->a_matrix_.index_[i]] += deltamin; + activitymin_[mip->a_matrix_.index_[i]] += + computeDelta(mip->a_matrix_.value_[i], oldbound, newbound, + -kHighsInf, activitymininf_[mip->a_matrix_.index_[i]]); } else { - double deltamax; - if (oldbound == -kHighsInf) { - --activitymaxinf_[mip->a_matrix_.index_[i]]; - deltamax = newbound * mip->a_matrix_.value_[i]; - } else if (newbound == -kHighsInf) { - ++activitymaxinf_[mip->a_matrix_.index_[i]]; - deltamax = -oldbound * mip->a_matrix_.value_[i]; - } else { - deltamax = (newbound - oldbound) * mip->a_matrix_.value_[i]; - } - activitymax_[mip->a_matrix_.index_[i]] += deltamax; + activitymax_[mip->a_matrix_.index_[i]] += + computeDelta(mip->a_matrix_.value_[i], oldbound, newbound, + -kHighsInf, activitymaxinf_[mip->a_matrix_.index_[i]]); } } @@ -1736,16 +1667,9 @@ void HighsDomain::updateActivityUbChange(HighsInt col, double oldbound, for (HighsInt i = start; i != end; ++i) { if (mip->a_matrix_.value_[i] > 0) { - double deltamax; - if (oldbound == kHighsInf) { - --activitymaxinf_[mip->a_matrix_.index_[i]]; - deltamax = newbound * mip->a_matrix_.value_[i]; - } else if (newbound == kHighsInf) { - ++activitymaxinf_[mip->a_matrix_.index_[i]]; - deltamax = -oldbound * mip->a_matrix_.value_[i]; - } else { - deltamax = (newbound - oldbound) * mip->a_matrix_.value_[i]; - } + double deltamax = + computeDelta(mip->a_matrix_.value_[i], oldbound, newbound, kHighsInf, + activitymaxinf_[mip->a_matrix_.index_[i]]); activitymax_[mip->a_matrix_.index_[i]] += deltamax; #ifndef NDEBUG @@ -1790,17 +1714,9 @@ void HighsDomain::updateActivityUbChange(HighsInt col, double oldbound, // propagateinds_.push_back(mip->a_matrix_.index_[i]); } } else { - double deltamin; - if (oldbound == kHighsInf) { - --activitymininf_[mip->a_matrix_.index_[i]]; - deltamin = newbound * mip->a_matrix_.value_[i]; - } else if (newbound == kHighsInf) { - ++activitymininf_[mip->a_matrix_.index_[i]]; - deltamin = -oldbound * mip->a_matrix_.value_[i]; - } else { - deltamin = (newbound - oldbound) * mip->a_matrix_.value_[i]; - } - + double deltamin = + computeDelta(mip->a_matrix_.value_[i], oldbound, newbound, kHighsInf, + activitymininf_[mip->a_matrix_.index_[i]]); activitymin_[mip->a_matrix_.index_[i]] += deltamin; #ifndef NDEBUG @@ -1860,30 +1776,13 @@ void HighsDomain::updateActivityUbChange(HighsInt col, double oldbound, std::swap(oldbound, newbound); for (HighsInt i = start; i != end; ++i) { if (mip->a_matrix_.value_[i] > 0) { - double deltamax; - if (oldbound == kHighsInf) { - --activitymaxinf_[mip->a_matrix_.index_[i]]; - deltamax = newbound * mip->a_matrix_.value_[i]; - } else if (newbound == kHighsInf) { - ++activitymaxinf_[mip->a_matrix_.index_[i]]; - deltamax = -oldbound * mip->a_matrix_.value_[i]; - } else { - deltamax = (newbound - oldbound) * mip->a_matrix_.value_[i]; - } - activitymax_[mip->a_matrix_.index_[i]] += deltamax; + activitymax_[mip->a_matrix_.index_[i]] += + computeDelta(mip->a_matrix_.value_[i], oldbound, newbound, + kHighsInf, activitymaxinf_[mip->a_matrix_.index_[i]]); } else { - double deltamin; - if (oldbound == kHighsInf) { - --activitymininf_[mip->a_matrix_.index_[i]]; - deltamin = newbound * mip->a_matrix_.value_[i]; - } else if (newbound == kHighsInf) { - ++activitymininf_[mip->a_matrix_.index_[i]]; - deltamin = -oldbound * mip->a_matrix_.value_[i]; - } else { - deltamin = (newbound - oldbound) * mip->a_matrix_.value_[i]; - } - - activitymin_[mip->a_matrix_.index_[i]] += deltamin; + activitymin_[mip->a_matrix_.index_[i]] += + computeDelta(mip->a_matrix_.value_[i], oldbound, newbound, + kHighsInf, activitymininf_[mip->a_matrix_.index_[i]]); } } @@ -1908,13 +1807,9 @@ void HighsDomain::recomputeCapacityThreshold(HighsInt row) { if (col_upper_[col] == col_lower_[col]) continue; - double boundRange = col_upper_[col] - col_lower_[col]; - - boundRange -= variableType(col) == HighsVarType::kContinuous - ? std::max(0.3 * boundRange, 1000.0 * feastol()) - : feastol(); - - double threshold = std::fabs(mipsolver->mipdata_->ARvalue_[i]) * boundRange; + double threshold = std::fabs(mipsolver->mipdata_->ARvalue_[i]) * + boundRange(col_upper_[col], col_lower_[col], feastol(), + variableType(col)); capacityThreshold_[row] = std::max({capacityThreshold_[row], threshold, feastol()}); @@ -2180,13 +2075,13 @@ void HighsDomain::setDomainChangeStack( if (k == stacksize) return; // For redundant branching bound changes we need to be more careful due to - // symmetry handling. If these boundchanges are redundant simply because the - // corresponding subtree was enumerated and hence the global bound updated, - // then we still need to keep their status as branching variables for - // computing correct stabilizers. - // They can, however, be safely dropped if they are either strictly - // redundant in the global domain, or if there is already a local bound - // change that makes the branching change redundant. + // symmetry handling. If these bound changes are redundant simply because + // the corresponding subtree was enumerated and hence the global bound + // updated, then we still need to keep their status as branching variables + // for computing correct stabilizers. They can, however, be safely dropped + // if they are either strictly redundant in the global domain, or if there + // is already a local bound change that makes the branching change + // redundant. if (domchgstack[k].boundtype == HighsBoundType::kLower) { if (domchgstack[k].boundval <= col_lower_[domchgstack[k].column]) { if (domchgstack[k].boundval < col_lower_[domchgstack[k].column]) @@ -3037,6 +2932,7 @@ bool HighsDomain::ConflictSet::resolveLinearGeq(HighsCDouble M, double Mupper, relaxUb = std::floor(relaxUb); if (relaxUb - ub <= localdom.feastol()) continue; + locdomchg.domchg.boundval = relaxUb; if (relaxUb - gub >= -localdom.mipsolver->mipdata_->epsilon) { diff --git a/src/mip/HighsModkSeparator.h b/src/mip/HighsModkSeparator.h index e7c2c9005f..8e6eb88bac 100644 --- a/src/mip/HighsModkSeparator.h +++ b/src/mip/HighsModkSeparator.h @@ -20,7 +20,7 @@ * cut. * * If a row contains continuous variables that sit at zero after bound - * substitution, then those rows are included in the congurence system, as the + * substitution, then those rows are included in the congruence system, as the * presence of such variables does not reduce the cuts violation when applying * the MIR procedure. In order to handle their presence the row must simply be * scaled, such that all integer variables that have a non-zero solution value diff --git a/src/pdlp/cupdlp/cupdlp_defs.h b/src/pdlp/cupdlp/cupdlp_defs.h index bdc4257a73..068c150935 100644 --- a/src/pdlp/cupdlp/cupdlp_defs.h +++ b/src/pdlp/cupdlp/cupdlp_defs.h @@ -154,7 +154,7 @@ struct CUPDLP_CSC_MATRIX { cupdlp_int *colMatIdx; cupdlp_float *colMatElem; - // Used to aviod implementing NormInf on cuda + // Used to avoid implementing NormInf on cuda cupdlp_float MatElemNormInf; #ifndef CUPDLP_CPU // Pointers to GPU vectors diff --git a/src/presolve/HPresolve.cpp b/src/presolve/HPresolve.cpp index b25b87af00..f2d549e865 100644 --- a/src/presolve/HPresolve.cpp +++ b/src/presolve/HPresolve.cpp @@ -1456,77 +1456,75 @@ HPresolve::Result HPresolve::runProbing(HighsPostsolveStack& postsolve_stack) { for (const auto& binvar : binaries) { HighsInt i = std::get<3>(binvar); - if (cliquetable.getSubstitution(i) != nullptr) continue; - - if (domain.isBinary(i)) { - // when a large percentage of columns have been deleted, stop this round - // of probing - // if (numDel > std::max(model->num_col_ * 0.2, 1000.)) break; - if (numDel > - std::max(1000., (model->num_row_ + model->num_col_) * 0.05)) { - probingEarlyAbort = true; - break; - } + if (cliquetable.getSubstitution(i) != nullptr || !domain.isBinary(i)) + continue; - // break in case of too many new implications to not spent ages in - // probing - if (cliquetable.isFull() || - cliquetable.numCliques() - numCliquesStart > - std::max(HighsInt{1000000}, 2 * numNonzeros()) || - implications.getNumImplications() - numImplicsStart > - std::max(HighsInt{1000000}, 2 * numNonzeros())) - break; + // when a large percentage of columns have been deleted, stop this round + // of probing + // if (numDel > std::max(model->num_col_ * 0.2, 1000.)) break; + probingEarlyAbort = + numDel > + std::max(HighsInt{1000}, (model->num_row_ + model->num_col_) / 20); + if (probingEarlyAbort) break; + + // break in case of too many new implications to not spent ages in + // probing + if (cliquetable.isFull() || + cliquetable.numCliques() - numCliquesStart > + std::max(HighsInt{1000000}, 2 * numNonzeros()) || + implications.getNumImplications() - numImplicsStart > + std::max(HighsInt{1000000}, 2 * numNonzeros())) + break; - // if (numProbed % 10 == 0) - // printf( - // "numprobed=%d numDel=%d newcliques=%d " - // "numNeighbourhoodQueries=%ld " - // "splayContingent=%ld\n", - // numProbed, numDel, cliquetable.numCliques() - numCliquesStart, - // cliquetable.numNeighbourhoodQueries, splayContingent); - if (cliquetable.numNeighbourhoodQueries > splayContingent) break; - - if (probingContingent - numProbed < 0) break; - - HighsInt numBoundChgs = 0; - HighsInt numNewCliques = -cliquetable.numCliques(); - if (!implications.runProbing(i, numBoundChgs)) continue; - probingContingent += numBoundChgs; - numNewCliques += cliquetable.numCliques(); - numNewCliques = std::max(numNewCliques, HighsInt{0}); - while (domain.getChangedCols().size() != numChangedCols) { - if (domain.isFixed(domain.getChangedCols()[numChangedCols++])) - ++probingNumDelCol; - } - HighsInt newNumDel = probingNumDelCol - numDelStart + - implications.substitutions.size() + - cliquetable.getSubstitutions().size(); - - if (newNumDel > numDel) { - probingContingent += numDel; - if (!mipsolver->submip) { - splayContingent += 100 * (newNumDel + numDelStart); - splayContingent += 1000 * numNewCliques; - } - numDel = newNumDel; - numFail = 0; - } else if (mipsolver->submip || numNewCliques == 0) { - splayContingent -= 100 * numFail; - ++numFail; - } else { + // if (numProbed % 10 == 0) + // printf( + // "numprobed=%d numDel=%d newcliques=%d " + // "numNeighbourhoodQueries=%ld " + // "splayContingent=%ld\n", + // numProbed, numDel, cliquetable.numCliques() - numCliquesStart, + // cliquetable.numNeighbourhoodQueries, splayContingent); + if (cliquetable.numNeighbourhoodQueries > splayContingent) break; + + if (probingContingent - numProbed < 0) break; + + HighsInt numBoundChgs = 0; + HighsInt numNewCliques = -cliquetable.numCliques(); + if (!implications.runProbing(i, numBoundChgs)) continue; + probingContingent += numBoundChgs; + numNewCliques += cliquetable.numCliques(); + numNewCliques = std::max(numNewCliques, HighsInt{0}); + while (domain.getChangedCols().size() != numChangedCols) { + if (domain.isFixed(domain.getChangedCols()[numChangedCols++])) + ++probingNumDelCol; + } + HighsInt newNumDel = probingNumDelCol - numDelStart + + implications.substitutions.size() + + cliquetable.getSubstitutions().size(); + + if (newNumDel > numDel) { + probingContingent += numDel; + if (!mipsolver->submip) { + splayContingent += 100 * (newNumDel + numDelStart); splayContingent += 1000 * numNewCliques; - numFail = 0; } + numDel = newNumDel; + numFail = 0; + } else if (mipsolver->submip || numNewCliques == 0) { + splayContingent -= 100 * numFail; + ++numFail; + } else { + splayContingent += 1000 * numNewCliques; + numFail = 0; + } - ++numProbed; - numProbes[i] += 1; + ++numProbed; + numProbes[i] += 1; - // printf("nprobed: %" HIGHSINT_FORMAT ", numCliques: %" HIGHSINT_FORMAT - // "\n", nprobed, - // cliquetable.numCliques()); - if (domain.infeasible()) { - return Result::kPrimalInfeasible; - } + // printf("nprobed: %" HIGHSINT_FORMAT ", numCliques: %" HIGHSINT_FORMAT + // "\n", nprobed, + // cliquetable.numCliques()); + if (domain.infeasible()) { + return Result::kPrimalInfeasible; } } @@ -5195,8 +5193,7 @@ HighsInt HPresolve::strengthenInequalities() { // do not run on very dense rows as this could get expensive if (rowsize[row] > - std::max(HighsInt{1000}, - HighsInt(0.05 * (model->num_col_ - numDeletedCols)))) + std::max(HighsInt{1000}, (model->num_col_ - numDeletedCols) / 20)) continue; // printf("strengthening knapsack of %" HIGHSINT_FORMAT " vars\n", diff --git a/src/simplex/SimplexStruct.h b/src/simplex/SimplexStruct.h index e29ce4e1ee..eb49c8ca96 100644 --- a/src/simplex/SimplexStruct.h +++ b/src/simplex/SimplexStruct.h @@ -25,7 +25,7 @@ struct SimplexBasis { // The basis for the simplex method consists of basicIndex, // nonbasicFlag and nonbasicMove. If HighsSimplexStatus has_basis // is true then it is assumed that basicIndex_ and nonbasicFlag_ are - // self-consistent and correpond to the dimensions of an associated + // self-consistent and correspond to the dimensions of an associated // HighsLp, but the basis matrix B is not necessarily nonsingular. std::vector basicIndex_; std::vector nonbasicFlag_; diff --git a/src/util/HighsCDouble.h b/src/util/HighsCDouble.h index 11cc8460d7..b02622d162 100644 --- a/src/util/HighsCDouble.h +++ b/src/util/HighsCDouble.h @@ -32,7 +32,7 @@ class HighsCDouble { // Proceedings of. 2005. /// performs an exact transformation such that x + y = a + b - /// and x = double(a + b). The operation uses 6 flops (addition/substraction). + /// and x = double(a + b). The operation uses 6 flops (addition/subtraction). static void two_sum(double& x, double& y, double a, double b) { x = a + b; double z = x - a; @@ -50,7 +50,7 @@ class HighsCDouble { /// performs an exact transformation such that x + y = a * b /// and x = double(a * b). The operation uses 10 flops for - /// addition/substraction and 7 flops for multiplication. + /// addition/subtraction and 7 flops for multiplication. static void two_product(double& x, double& y, double a, double b) { x = a * b; double a1, a2, b1, b2;