Loading [MathJax]/extensions/MathZoom.js
Scippy

SCIP

Solving Constraint Integer Programs

branch_relpscost.c
Go to the documentation of this file.
1/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
2/* */
3/* This file is part of the program and library */
4/* SCIP --- Solving Constraint Integer Programs */
5/* */
6/* Copyright (c) 2002-2025 Zuse Institute Berlin (ZIB) */
7/* */
8/* Licensed under the Apache License, Version 2.0 (the "License"); */
9/* you may not use this file except in compliance with the License. */
10/* You may obtain a copy of the License at */
11/* */
12/* http://www.apache.org/licenses/LICENSE-2.0 */
13/* */
14/* Unless required by applicable law or agreed to in writing, software */
15/* distributed under the License is distributed on an "AS IS" BASIS, */
16/* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. */
17/* See the License for the specific language governing permissions and */
18/* limitations under the License. */
19/* */
20/* You should have received a copy of the Apache-2.0 license */
21/* along with SCIP; see the file LICENSE. If not visit scipopt.org. */
22/* */
23/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
24
25/**@file branch_relpscost.c
26 * @ingroup DEFPLUGINS_BRANCH
27 * @brief reliable pseudo costs branching rule
28 * @author Tobias Achterberg
29 * @author Timo Berthold
30 * @author Gerald Gamrath
31 * @author Marc Pfetsch
32 */
33
34/*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
35
38#include "scip/treemodel.h"
39#include "scip/cons_and.h"
40#include "scip/pub_branch.h"
41#include "scip/pub_cons.h"
42#include "scip/pub_message.h"
43#include "scip/pub_misc.h"
44#include "scip/pub_sol.h"
45#include "scip/pub_tree.h"
46#include "scip/pub_var.h"
47#include "scip/scip_branch.h"
48#include "scip/scip_cons.h"
49#include "scip/scip_general.h"
50#include "scip/scip_lp.h"
51#include "scip/scip_mem.h"
52#include "scip/scip_message.h"
53#include "scip/scip_nlp.h"
54#include "scip/scip_numerics.h"
55#include "scip/scip_param.h"
56#include "scip/scip_prob.h"
58#include "scip/scip_sol.h"
60#include "scip/scip_tree.h"
61#include "scip/scip_var.h"
62#include "scip/prop_symmetry.h"
63#include "scip/symmetry.h"
64#include <string.h>
65
66#define BRANCHRULE_NAME "relpscost"
67#define BRANCHRULE_DESC "reliability branching on pseudo cost values"
68#define BRANCHRULE_PRIORITY 10000
69#define BRANCHRULE_MAXDEPTH -1
70#define BRANCHRULE_MAXBOUNDDIST 1.0
71
72#define DEFAULT_CONFLICTWEIGHT 0.01 /**< weight in score calculations for conflict score */
73#define DEFAULT_CONFLENGTHWEIGHT 0.0 /**< weight in score calculations for conflict length score*/
74#define DEFAULT_INFERENCEWEIGHT 0.0001 /**< weight in score calculations for inference score */
75#define DEFAULT_CUTOFFWEIGHT 0.0001 /**< weight in score calculations for cutoff score */
76#define DEFAULT_GMIAVGEFFWEIGHT 0.0 /**< weight in score calculations of average GMI cut normed efficacies */
77#define DEFAULT_GMILASTEFFWEIGHT 0.00001 /**< weight in score calculations of last GMI cut normed efficacy */
78#define DEFAULT_PSCOSTWEIGHT 1.0 /**< weight in score calculations for pseudo cost score */
79#define DEFAULT_NLSCOREWEIGHT 0.1 /**< weight in score calculations for nlcount score */
80#define DEFAULT_MINRELIABLE 1.0 /**< minimal value for minimum pseudo cost size to regard pseudo cost value as reliable */
81#define DEFAULT_MAXRELIABLE 5.0 /**< maximal value for minimum pseudo cost size to regard pseudo cost value as reliable */
82#define DEFAULT_SBITERQUOT 0.5 /**< maximal fraction of strong branching LP iterations compared to normal iterations */
83#define DEFAULT_SBITEROFS 100000 /**< additional number of allowed strong branching LP iterations */
84#define DEFAULT_MAXLOOKAHEAD 9 /**< maximal number of further variables evaluated without better score */
85#define DEFAULT_INITCAND 100 /**< maximal number of candidates initialized with strong branching per node */
86#define DEFAULT_INITITER 0 /**< iteration limit for strong branching initialization of pseudo cost entries (0: auto) */
87#define DEFAULT_MAXBDCHGS 5 /**< maximal number of bound tightenings before the node is reevaluated (-1: unlimited) */
88#define DEFAULT_MAXPROPROUNDS -2 /**< maximum number of propagation rounds to be performed during strong branching
89 * before solving the LP (-1: no limit, -2: parameter settings) */
90#define DEFAULT_PROBINGBOUNDS TRUE /**< should valid bounds be identified in a probing-like fashion during strong
91 * branching (only with propagation)? */
92#define DEFAULT_USERELERRORFORRELIABILITY FALSE /**< should reliability be based on relative errors? */
93#define DEFAULT_LOWERRORTOL 0.05 /**< lowest tolerance beneath which relative errors are reliable */
94#define DEFAULT_HIGHERRORTOL 1.0 /**< highest tolerance beneath which relative errors are reliable */
95#define DEFAULT_USEHYPTESTFORRELIABILITY FALSE /**< should the strong branching decision be based on a hypothesis test? */
96#define DEFAULT_USEDYNAMICCONFIDENCE FALSE /**< should the confidence level be adjusted dynamically? */
97#define DEFAULT_STORESEMIINITCOSTS FALSE /**< should strong branching result be considered for pseudo costs if the other direction was infeasible? */
98#define DEFAULT_USESBLOCALINFO FALSE /**< should the scoring function use only local cutoff and inference information obtained for strong branching candidates? */
99#define DEFAULT_CONFIDENCELEVEL 2 /**< The confidence level for statistical methods, between 0 (Min) and 4 (Max). */
100#define DEFAULT_SKIPBADINITCANDS TRUE /**< should branching rule skip candidates that have a low probability to be
101 * better than the best strong-branching or pseudo-candidate? */
102#define DEFAULT_STARTRANDSEED 5 /**< start random seed for random number generation */
103#define DEFAULT_RANDINITORDER FALSE /**< should slight perturbation of scores be used to break ties in the prior scores? */
104#define DEFAULT_USESMALLWEIGHTSITLIM FALSE /**< should smaller weights be used for pseudo cost updates after hitting the LP iteration limit? */
105#define DEFAULT_DYNAMICWEIGHTS TRUE /**< should the weights of the branching rule be adjusted dynamically during solving based
106 * infeasible and objective leaf counters? */
107#define DEFAULT_DEGENERACYAWARE 1 /**< should degeneracy be taken into account to update weights and skip strong branching? (0: off, 1: after root, 2: always)*/
108
109/* symmetry handling */
110#define DEFAULT_FILTERCANDSSYM FALSE /**< Use symmetry to filter branching candidates? */
111#define DEFAULT_TRANSSYMPSCOST FALSE /**< Transfer pscost information to symmetric variables if filtering is performed? */
112
113/** branching rule data */
114struct SCIP_BranchruleData
115{
116 SCIP_Real conflictweight; /**< weight in score calculations for conflict score */
117 SCIP_Real conflengthweight; /**< weight in score calculations for conflict length score */
118 SCIP_Real inferenceweight; /**< weight in score calculations for inference score */
119 SCIP_Real cutoffweight; /**< weight in score calculations for cutoff score */
120 SCIP_Real gmiavgeffweight; /**< weight in score calculations of average GMI normed cut efficacies */
121 SCIP_Real gmilasteffweight; /**< weight in score calculations of last GMI cut normalized efficacy */
122 SCIP_Real pscostweight; /**< weight in score calculations for pseudo cost score */
123 SCIP_Real nlscoreweight; /**< weight in score calculations for nlcount score */
124 SCIP_Real minreliable; /**< minimal value for minimum pseudo cost size to regard pseudo cost value as reliable */
125 SCIP_Real maxreliable; /**< maximal value for minimum pseudo cost size to regard pseudo cost value as reliable */
126 SCIP_Real sbiterquot; /**< maximal fraction of strong branching LP iterations compared to normal iterations */
127 int sbiterofs; /**< additional number of allowed strong branching LP iterations */
128 int maxlookahead; /**< maximal number of further variables evaluated without better score */
129 int initcand; /**< maximal number of candidates initialized with strong branching per node */
130 int inititer; /**< iteration limit for strong branching initialization of pseudo cost entries (0: auto) */
131 int maxbdchgs; /**< maximal number of bound tightenings before the node is reevaluated (-1: unlimited) */
132 int maxproprounds; /**< maximum number of propagation rounds to be performed during strong branching
133 * before solving the LP (-1: no limit, -2: parameter settings) */
134 SCIP_Bool probingbounds; /**< should valid bounds be identified in a probing-like fashion during strong
135 * branching (only with propagation)? */
136 SCIP_Bool userelerrorforreliability; /**< should reliability be based on relative errors? */
137 SCIP_Real lowerrortol; /**< lowest tolerance beneath which relative errors are reliable */
138 SCIP_Real higherrortol; /**< highest tolerance beneath which relative errors are reliable */
139 SCIP_Bool usehyptestforreliability; /**< should the strong branching decision be based on a hypothesis test? */
140 SCIP_Bool usedynamicconfidence; /**< should the confidence level be adjusted dynamically? */
141 SCIP_Bool storesemiinitcosts; /**< should strong branching result be considered for pseudo costs if the
142 * other direction was infeasible? */
143 SCIP_Bool usesblocalinfo; /**< should the scoring function disregard cutoffs for variable if sb-lookahead was feasible ? */
144 SCIP_Bool skipbadinitcands; /**< should branching rule skip candidates that have a low probability to be
145 * better than the best strong-branching or pseudo-candidate? */
146 SCIP_Bool dynamicweights; /**< should the weights of the branching rule be adjusted dynamically during
147 * solving based on objective and infeasible leaf counters? */
148 int degeneracyaware; /**< should degeneracy be taken into account to update weights and skip strong branching? (0: off, 1: after root, 2: always) */
149 int confidencelevel; /**< The confidence level for statistical methods, between 0 (Min) and 4 (Max). */
150 int* nlcount; /**< array to store nonlinear count values */
151 int nlcountsize; /**< length of nlcount array */
152 int nlcountmax; /**< maximum entry in nlcount array or 1 if NULL */
153 SCIP_Bool randinitorder; /**< should slight perturbation of scores be used to break ties in the prior scores? */
154 SCIP_RANDNUMGEN* randnumgen; /**< random number generator */
155 int startrandseed; /**< start random seed for random number generation */
156 SCIP_Bool usesmallweightsitlim; /**< should smaller weights be used for pseudo cost updates after hitting the LP iteration limit? */
157 SCIP_TREEMODEL* treemodel; /**< Parameters for the Treemodel branching rules */
158
159 /* for symmetry */
160 SCIP_Bool filtercandssym; /**< Use symmetry to filter branching candidates? */
161 SCIP_Bool transsympscost; /**< Transfer pscost information to symmetric variables? */
162
163 SCIP_Bool nosymmetry; /**< No symmetry present? */
164 int* orbits; /**< array of non-trivial orbits */
165 int* orbitbegins; /**< array containing begin positions of new orbits in orbits array */
166 int norbits; /**< pointer to number of orbits currently stored in orbits */
167 int* varorbitmap; /**< array for storing indices of the containing orbit for each variable */
168 int* orbitrep; /**< representative variable of each orbit */
169 SCIP_VAR** permvars; /**< variables on which permutations act */
170 int npermvars; /**< number of variables for permutations */
171 SCIP_HASHMAP* permvarmap; /**< map of variables to indices in permvars array */
172};
173
174/*
175 * local methods
176 */
177
178/** initialize orbits */
179static
181 SCIP* scip, /**< SCIP data structure */
182 SCIP_BRANCHRULEDATA* branchruledata /**< branching rule data */
183 )
184{
185 int** permstrans = NULL;
186 int* components = NULL;
187 int* componentbegins = NULL;
188 int* vartocomponent = NULL;
189 int ncomponents = 0;
190 int nperms = -1;
191
192 assert( scip != NULL );
193 assert( branchruledata != NULL );
194 assert( branchruledata->filtercandssym );
195
196 /* exit if no symmetry or orbits already available */
197 if( branchruledata->nosymmetry || branchruledata->orbits != NULL )
198 return SCIP_OKAY;
199
200 assert( branchruledata->orbitbegins == NULL );
201 assert( branchruledata->varorbitmap == NULL );
202 assert( branchruledata->orbitrep == NULL );
203
204 /* obtain symmetry including permutations */
205 SCIP_CALL( SCIPgetSymmetry(scip, &branchruledata->npermvars, &branchruledata->permvars, &branchruledata->permvarmap,
206 &nperms, NULL, &permstrans, NULL, NULL, &components, &componentbegins, &vartocomponent, &ncomponents) );
207
208 /* turn off symmetry handling if there is no symmetry or the number of variables is not equal */
209 if( nperms <= 0 || branchruledata->npermvars != SCIPgetNVars(scip) )
210 {
211 branchruledata->nosymmetry = TRUE;
212 return SCIP_OKAY;
213 }
214 assert( branchruledata->permvars != NULL );
215 assert( branchruledata->permvarmap != NULL );
216 assert( branchruledata->npermvars > 0 );
217 assert( permstrans != NULL );
218 assert( components != NULL );
219 assert( componentbegins != NULL );
220 assert( vartocomponent != NULL );
221 assert( ncomponents > 0 );
222
223 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &branchruledata->orbits, branchruledata->npermvars) );
224 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &branchruledata->orbitbegins, branchruledata->npermvars) );
225 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &branchruledata->varorbitmap, branchruledata->npermvars) );
226 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &branchruledata->orbitrep, branchruledata->npermvars) );
227
228 /* Compute orbits on all variables, since this might help for branching and this computation is only done once. */
229 SCIP_CALL( SCIPcomputeOrbitsComponentsSym(scip, branchruledata->npermvars, permstrans, nperms,
230 components, componentbegins, vartocomponent, ncomponents,
231 branchruledata->orbits, branchruledata->orbitbegins, &branchruledata->norbits, branchruledata->varorbitmap) );
232 assert( branchruledata->norbits < branchruledata->npermvars );
233
234 return SCIP_OKAY;
235}
236
237/** filter out symmetric variables from branching variables */
238static
240 SCIP* scip, /**< SCIP data structure */
241 SCIP_BRANCHRULEDATA* branchruledata, /**< branching rule data */
242 SCIP_VAR** origbranchcands, /**< original branching candidates */
243 SCIP_Real* origbranchcandssol, /**< original solution value for the branching candidates */
244 SCIP_Real* origbranchcandsfrac,/**< original fractional part of the branching candidates */
245 int norigbranchcands, /**< original number of branching candidates */
246 SCIP_VAR** branchcands, /**< branching candidates */
247 SCIP_Real* branchcandssol, /**< solution value for the branching candidates */
248 SCIP_Real* branchcandsfrac, /**< fractional part of the branching candidates */
249 int* branchorbitidx, /**< array of indices of orbit of branching candidates */
250 int* nbranchcands /**< pointer to store number of branching candidates */
251 )
252{
253 int i;
254
255 assert( scip != NULL );
256 assert( branchruledata != NULL );
257 assert( origbranchcands != NULL );
258 assert( origbranchcandssol != NULL );
259 assert( origbranchcandsfrac != NULL );
260 assert( branchcands != NULL );
261 assert( branchcandssol != NULL );
262 assert( branchcandsfrac != NULL );
263 assert( nbranchcands != NULL );
264
265 assert( ! branchruledata->nosymmetry );
266 assert( branchruledata->orbitbegins != NULL );
267 assert( branchruledata->orbits != NULL );
268 assert( branchruledata->permvarmap != NULL );
269 assert( branchruledata->varorbitmap != NULL );
270 assert( branchruledata->orbitrep != NULL );
271 assert( branchruledata->norbits < branchruledata->npermvars );
272
273 /* init representatives (used to see whether variable is the first in an orbit) */
274 for( i = 0; i < branchruledata->norbits; ++i )
275 branchruledata->orbitrep[i] = -1;
276
277 /* loop through branching variables, determine orbit and whether they are the first ones */
278 *nbranchcands = 0;
279 for( i = 0; i < norigbranchcands; ++i )
280 {
281 int orbitidx = -1;
282 int varidx;
283
284 varidx = SCIPhashmapGetImageInt(branchruledata->permvarmap, (void*) origbranchcands[i]);
285 if( varidx != INT_MAX )
286 {
287 assert( 0 <= varidx && varidx < branchruledata->npermvars );
288 orbitidx = branchruledata->varorbitmap[varidx];
289 }
290 assert( -1 <= orbitidx && orbitidx < branchruledata->norbits );
291
292 /* Check whether the variable is not present (can happen if variable was added after computing symmetries or is in
293 * a singleton orbit). */
294 if( orbitidx == -1 )
295 {
296 branchcands[*nbranchcands] = origbranchcands[i];
297 branchcandssol[*nbranchcands] = origbranchcandssol[i];
298 branchcandsfrac[*nbranchcands] = origbranchcandsfrac[i];
299 branchorbitidx[*nbranchcands] = -1;
300 ++(*nbranchcands);
301 }
302 else if( branchruledata->orbitrep[orbitidx] == -1 )
303 {
304 /* if variable is the first in a nontrivial orbit */
305 assert( 0 <= varidx && varidx < branchruledata->npermvars );
306 branchruledata->orbitrep[orbitidx] = varidx;
307 branchcands[*nbranchcands] = origbranchcands[i];
308 branchcandssol[*nbranchcands] = origbranchcandssol[i];
309 branchcandsfrac[*nbranchcands] = origbranchcandsfrac[i];
310 branchorbitidx[*nbranchcands] = orbitidx;
311 ++(*nbranchcands);
312 }
313 }
314
315 SCIPdebugMsg(scip, "Filtered out %d variables by symmetry.\n", norigbranchcands - *nbranchcands);
316
317 return SCIP_OKAY;
318}
319
320/** updates the pseudo costs of the given variable and all its symmetric variables */
321static
323 SCIP* scip, /**< SCIP data structure */
324 SCIP_BRANCHRULEDATA* branchruledata, /**< branching rule data */
325 SCIP_VAR* branchvar, /**< branching variable candidate */
326 int* branchorbitidx, /**< array of orbit indices */
327 int branchvaridx, /**< index of variable in branchorbitidx */
328 SCIP_Real solvaldelta, /**< difference of variable's new LP value - old LP value */
329 SCIP_Real objdelta, /**< difference of new LP's objective value - old LP's objective value */
330 SCIP_Real weight /**< weight in (0,1] of this update in pseudo cost sum */
331 )
332{
333 int orbitidx;
334 int j;
335
336 assert( scip != NULL );
337 assert( branchruledata != NULL );
338
339 if( branchruledata->nosymmetry || ! branchruledata->transsympscost || branchorbitidx == NULL )
340 {
341 /* use original update function */
342 SCIP_CALL( SCIPupdateVarPseudocost(scip, branchvar, solvaldelta, objdelta, weight) );
343 return SCIP_OKAY;
344 }
345
346 assert( branchruledata->orbitbegins != NULL );
347 assert( branchruledata->orbits != NULL );
348 assert( 0 <= branchvaridx && branchvaridx < branchruledata->npermvars );
349
350 orbitidx = branchorbitidx[branchvaridx];
351 if( orbitidx < 0 )
352 {
353 /* only update given variable */
354 SCIP_CALL( SCIPupdateVarPseudocost(scip, branchvar, solvaldelta, objdelta, weight) );
355 return SCIP_OKAY;
356 }
357 assert( 0 <= orbitidx && orbitidx < branchruledata->norbits );
358
359 /* loop through orbit containing variable and update pseudo costs for all variables */
360 for( j = branchruledata->orbitbegins[orbitidx]; j < branchruledata->orbitbegins[orbitidx+1]; ++j )
361 {
362 SCIP_VAR* var;
363 int idx;
364
365 idx = branchruledata->orbits[j];
366 assert( 0 <= idx && idx < branchruledata->npermvars );
367
368 var = branchruledata->permvars[idx];
369 assert( var != NULL );
370
371 if( SCIPvarIsActive(var) )
372 {
373 SCIP_CALL( SCIPupdateVarPseudocost(scip, var, solvaldelta, objdelta, weight) );
374 }
375 }
376
377 return SCIP_OKAY;
378}
379
380/**! [SnippetCodeStyleDeclaration] */
381
382/** counts number of nonlinear constraints in which each variable appears */
383static
385 SCIP* scip, /**< SCIP data structure */
386 int* nlcount, /**< pointer to array for storing count values */
387 int nlcountsize, /**< buffer for storing length of nlcount array */
388 int* nlcountmax /**< buffer for storing maximum value in nlcount array */
389 )
390{
391 SCIP_CONSHDLR* andconshdlr;
392 SCIP_VAR** vars;
393 int nvars;
394 int i;
395
396/**! [SnippetCodeStyleDeclaration] */
397
398 assert(scip != NULL);
399 assert(nlcount != NULL);
400 assert(nlcountmax != NULL);
401
402 SCIP_CALL( SCIPgetVarsData(scip, &vars, &nvars, NULL, NULL, NULL, NULL) );
403 assert(nlcountsize >= nvars);
404
405 /* get nonlinearity for constraints in NLP */
407 {
408 assert(SCIPgetNNLPVars(scip) == nvars);
410 }
411 else
412 {
413 BMSclearMemoryArray(nlcount, nvars);
414 }
415
416 /* increase counters for and constraints */
417 andconshdlr = SCIPfindConshdlr(scip, "and");
418 if( andconshdlr != NULL )
419 {
420 int c;
421
422 for( c = 0; c < SCIPconshdlrGetNActiveConss(andconshdlr); ++c )
423 {
424 SCIP_CONS* andcons;
425 SCIP_VAR** andvars;
426 SCIP_VAR* andres;
427 int probindex;
428 int nandvars;
429 int v;
430
431 /* get constraint and variables */
432 andcons = SCIPconshdlrGetConss(andconshdlr)[c];
433 nandvars = SCIPgetNVarsAnd(scip, andcons);
434 andvars = SCIPgetVarsAnd(scip, andcons);
435 andres = SCIPgetResultantAnd(scip, andcons);
436
437 /* get active index of resultant */
438 probindex = SCIPvarGetProbindex(SCIPvarGetProbvar(andres));
439
440 /* the resultant might be deleted */
441 if( probindex >= 0 )
442 ++nlcount[probindex];
443
444 /**! [SnippetCodeStyleIfFor] */
445 for( v = 0; v < nandvars; ++v )
446 {
447 /* get active index of operator */
448 probindex = SCIPvarGetProbindex(SCIPvarGetProbvar(andvars[v]));
449
450 /* the operator might be deleted */
451 if( probindex >= 0 )
452 ++nlcount[probindex];
453 }
454 /**! [SnippetCodeStyleIfFor] */
455 }
456 }
457
458 /* compute maximum count value */
459 *nlcountmax = 1;
460 for( i = 0; i < nvars; ++i )
461 {
462 if( *nlcountmax < nlcount[i] )
463 *nlcountmax = nlcount[i];
464 }
465
466 return SCIP_OKAY;
467}
468
469static
471 SCIP* scip, /**< SCIP data structure */
472 SCIP_BRANCHRULEDATA* branchruledata /**< branching rule data */
473 )
474{
475 int nvars;
476
477 assert(scip != NULL);
478 assert(branchruledata != NULL);
479
480 nvars = SCIPgetNVars(scip);
481
482 /**@todo test whether we want to apply this as if problem has only and constraints */
483 /**@todo update changes in and constraints */
484 if( branchruledata->nlscoreweight > 0.0 ) /* && SCIPisNLPConstructed(scip) */
485 {
486 if( branchruledata->nlcount == NULL )
487 {
488 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &branchruledata->nlcount, nvars) );
489 branchruledata->nlcountsize = nvars;
490
491 SCIP_CALL( countNonlinearities(scip, branchruledata->nlcount, branchruledata->nlcountsize, &branchruledata->nlcountmax) );
492 }
493 else if( branchruledata->nlcountsize < nvars )
494 {
495 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &branchruledata->nlcount, branchruledata->nlcountsize, nvars) );
496 /**@todo should we update nlcounts for new variables? */
497 BMSclearMemoryArray(&(branchruledata->nlcount[branchruledata->nlcountsize]), nvars - branchruledata->nlcountsize); /*lint !e866*/
498 branchruledata->nlcountsize = nvars;
499 }
500 assert(branchruledata->nlcount != NULL);
501 assert(branchruledata->nlcountsize == nvars);
502 assert(branchruledata->nlcountmax >= 1);
503 }
504 else
505 {
506 SCIPfreeBlockMemoryArrayNull(scip, &branchruledata->nlcount, branchruledata->nlcountsize);
507 branchruledata->nlcountsize = 0;
508 branchruledata->nlcountmax = 1;
509 }
510
511 return SCIP_OKAY;
512}
513
514
515/** calculates nlscore value between 0 and 1 */
516static
518 SCIP* scip, /**< SCIP data structure */
519 int* nlcount, /**< array to store count values */
520 int nlcountmax, /**< maximum value in nlcount array */
521 int probindex /**< index of branching candidate */
522 )
523{
524 if( nlcountmax >= 1 && nlcount != NULL )
525 {
526 SCIP_Real nlscore;
527
528 assert(scip != NULL);
529 assert(probindex >= 0);
530 assert(probindex < SCIPgetNVars(scip));
531
532 nlscore = nlcount[probindex] / (SCIP_Real)nlcountmax;
533
534 assert(nlscore >= 0.0);
535 assert(nlscore <= 1.0);
536 return nlscore;
537 }
538 else
539 return 0.0;
540}
541
542/** calculates an overall score value for the given individual score values */
543static
545 SCIP* scip, /**< SCIP data structure */
546 SCIP_BRANCHRULEDATA* branchruledata, /**< branching rule data */
547 SCIP_Real conflictscore, /**< conflict score of current variable */
548 SCIP_Real avgconflictscore, /**< average conflict score */
549 SCIP_Real conflengthscore, /**< conflict length score of current variable */
550 SCIP_Real avgconflengthscore, /**< average conflict length score */
551 SCIP_Real inferencescore, /**< inference score of current variable */
552 SCIP_Real avginferencescore, /**< average inference score */
553 SCIP_Real cutoffscore, /**< cutoff score of current variable */
554 SCIP_Real avgcutoffscore, /**< average cutoff score */
555 SCIP_Real gmieffscore, /**< normalized-eff of avg GMI cuts from row when var was frac and basic */
556 SCIP_Real lastgmieffscore, /**< last normalized gmieffscore when var was frac and basic */
557 SCIP_Real pscostscore, /**< pscost score of current variable */
558 SCIP_Real avgpscostscore, /**< average pscost score */
559 SCIP_Real nlscore, /**< nonlinear score of current variable between 0 and 1 */
560 SCIP_Real frac, /**< fractional value of variable in current solution */
561 SCIP_Real degeneracyfactor /**< factor to apply because of degeneracy */
562 )
563{
564 SCIP_Real score;
565 SCIP_Real dynamicfactor;
566
567 assert(branchruledata != NULL);
568 assert(0.0 < frac && frac < 1.0);
569
570 if( branchruledata->dynamicweights )
571 {
572 dynamicfactor = (SCIPgetNInfeasibleLeaves(scip) + 1.0) / (SCIPgetNObjlimLeaves(scip) + 1.0);
573 }
574 else
575 dynamicfactor = 1.0;
576
577 dynamicfactor *= degeneracyfactor;
578
579 score = dynamicfactor * (branchruledata->conflictweight * (1.0 - 1.0/(1.0+conflictscore/avgconflictscore))
580 + branchruledata->conflengthweight * (1.0 - 1.0/(1.0+conflengthscore/avgconflengthscore))
581 + branchruledata->inferenceweight * (1.0 - 1.0/(1.0+inferencescore/avginferencescore))
582 + branchruledata->cutoffweight * (1.0 - 1.0/(1.0+cutoffscore/avgcutoffscore))
583 + branchruledata->gmiavgeffweight * gmieffscore + branchruledata->gmilasteffweight * lastgmieffscore)
584 + branchruledata->pscostweight / dynamicfactor * (1.0 - 1.0/(1.0+pscostscore/avgpscostscore))
585 + branchruledata->nlscoreweight * nlscore;
586
587 /* avoid close to integral variables */
588 if( MIN(frac, 1.0 - frac) < 10.0 * SCIPfeastol(scip) )
589 score *= 1e-6;
590
591 return score;
592}
593
594/** adds given index and direction to bound change arrays */
595static
597 SCIP* scip, /**< SCIP data structure */
598 int** bdchginds, /**< pointer to bound change index array */
599 SCIP_BOUNDTYPE** bdchgtypes, /**< pointer to bound change types array */
600 SCIP_Real** bdchgbounds, /**< pointer to bound change new bounds array */
601 int* nbdchgs, /**< pointer to number of bound changes */
602 int ind, /**< index to store in bound change index array */
603 SCIP_BOUNDTYPE type, /**< type of the bound change to store in bound change type array */
604 SCIP_Real bound /**< new bound to store in bound change new bounds array */
605 )
606{
607 assert(bdchginds != NULL);
608 assert(bdchgtypes != NULL);
609 assert(bdchgbounds != NULL);
610 assert(nbdchgs != NULL);
611
612 SCIP_CALL( SCIPreallocBufferArray(scip, bdchginds, (*nbdchgs) + 1) );
613 SCIP_CALL( SCIPreallocBufferArray(scip, bdchgtypes, (*nbdchgs) + 1) );
614 SCIP_CALL( SCIPreallocBufferArray(scip, bdchgbounds, (*nbdchgs) + 1) );
615 assert(*bdchginds != NULL);
616 assert(*bdchgtypes != NULL);
617 assert(*bdchgbounds != NULL);
618
619 (*bdchginds)[*nbdchgs] = ind;
620 (*bdchgtypes)[*nbdchgs] = type;
621 (*bdchgbounds)[*nbdchgs] = bound;
622 (*nbdchgs)++;
623
624 return SCIP_OKAY;
625}
626
627/**! [SnippetCodeStyleStaticAsserts] */
628/** frees bound change arrays */
629static
631 SCIP* scip, /**< SCIP data structure */
632 int** bdchginds, /**< pointer to bound change index array */
633 SCIP_BOUNDTYPE** bdchgtypes, /**< pointer to bound change types array */
634 SCIP_Real** bdchgbounds, /**< pointer to bound change new bounds array */
635 int* nbdchgs /**< pointer to number of bound changes */
636 )
637{
638 assert(scip != NULL);
639 assert(bdchginds != NULL);
640 assert(bdchgtypes != NULL);
641 assert(bdchgbounds != NULL);
642 assert(nbdchgs != NULL);
643/**! [SnippetCodeStyleStaticAsserts] */
644
645 SCIPfreeBufferArrayNull(scip, bdchgbounds);
646 SCIPfreeBufferArrayNull(scip, bdchgtypes);
647 SCIPfreeBufferArrayNull(scip, bdchginds);
648 *nbdchgs = 0;
649}
650
651/** applies bound changes stored in bound change arrays */
652static
654 SCIP* scip, /**< SCIP data structure */
655 SCIP_VAR** vars, /**< problem variables */
656 int* bdchginds, /**< bound change index array */
657 SCIP_BOUNDTYPE* bdchgtypes, /**< bound change types array */
658 SCIP_Real* bdchgbounds, /**< bound change new bound array */
659 int nbdchgs, /**< number of bound changes */
660 SCIP_RESULT* result /**< result pointer */
661 )
662{
663#ifndef NDEBUG
664 SCIP_BRANCHRULE* branchrule;
665 SCIP_BRANCHRULEDATA* branchruledata;
666#endif
667 SCIP_Bool infeasible;
668 SCIP_Bool tightened;
669 int i;
670
671 assert(vars != NULL);
672
673#ifndef NDEBUG
674 /* find branching rule */
676 assert(branchrule != NULL);
677
678 /* get branching rule data */
679 branchruledata = SCIPbranchruleGetData(branchrule);
680 assert(branchruledata != NULL);
681#endif
682
683 SCIPdebugMsg(scip, "applying %d bound changes\n", nbdchgs);
684
685 for( i = 0; i < nbdchgs; ++i )
686 {
687 int v;
688
689 v = bdchginds[i];
690
691 SCIPdebugMsg(scip, " -> <%s> [%g,%g]\n",
692 SCIPvarGetName(vars[v]), SCIPvarGetLbLocal(vars[v]), SCIPvarGetUbLocal(vars[v]));
693
694 if( bdchgtypes[i] == SCIP_BOUNDTYPE_LOWER )
695 {
696 /* change lower bound of variable to given bound */
697 SCIP_CALL( SCIPtightenVarLb(scip, vars[v], bdchgbounds[i], TRUE, &infeasible, &tightened) );
698 if( infeasible )
699 {
700 *result = SCIP_CUTOFF;
701 return SCIP_OKAY;
702 }
703
704 /* if we did propagation, the bound change might already have been added */
705 assert(tightened || (branchruledata->maxproprounds != 0));
706 }
707 else
708 {
709 assert(bdchgtypes[i] == SCIP_BOUNDTYPE_UPPER);
710
711 /* change upper bound of variable to given bound */
712 SCIP_CALL( SCIPtightenVarUb(scip, vars[v], bdchgbounds[i], TRUE, &infeasible, &tightened) );
713 if( infeasible )
714 {
715 *result = SCIP_CUTOFF;
716 return SCIP_OKAY;
717 }
718
719 /* if we did propagation, the bound change might already have been added */
720 assert(tightened || (branchruledata->maxproprounds != 0));
721 }
722 SCIPdebugMsg(scip, " -> [%g,%g]\n", SCIPvarGetLbLocal(vars[v]), SCIPvarGetUbLocal(vars[v]));
723 }
724
725 return SCIP_OKAY;
726}
727
728/** execute reliability pseudo cost branching */
729static
731 SCIP* scip, /**< SCIP data structure */
732 SCIP_BRANCHRULE* branchrule, /**< branching rule */
733 SCIP_VAR** branchcands, /**< branching candidates */
734 SCIP_Real* branchcandssol, /**< solution value for the branching candidates */
735 SCIP_Real* branchcandsfrac, /**< fractional part of the branching candidates */
736 int* branchorbitidx, /**< indices of orbit (or NULL) */
737 int nbranchcands, /**< number of branching candidates */
738 SCIP_Bool executebranch, /**< execute a branching step or run probing only */
739 SCIP_RESULT* result /**< pointer to the result of the execution */
740 )
741{ /*lint --e{715}*/
742 SCIP_BRANCHRULEDATA* branchruledata;
743 SCIP_Real lpobjval;
744 SCIP_Real bestsbdown;
745 SCIP_Real bestsbup;
746 SCIP_Real provedbound;
747 SCIP_Bool bestsbdownvalid;
748 SCIP_Bool bestsbupvalid;
749 SCIP_Bool bestsbdowncutoff;
750 SCIP_Bool bestsbupcutoff;
751 SCIP_Bool bestisstrongbranch;
752 SCIP_Bool allcolsinlp;
753 SCIP_Bool exactsolve;
754 int ninitcands;
755 int bestcand;
756
757 /* remember which variables strong branching is performed on, and the
758 * recorded lp bound changes that are observed */
759 SCIP_Real* sbdown = NULL;
760 SCIP_Real* sbup = NULL;
761 SCIP_Bool* sbdownvalid = NULL;
762 SCIP_Bool* sbupvalid = NULL;
763
764 *result = SCIP_DIDNOTRUN;
765
767
768 /* get branching rule data */
769 branchruledata = SCIPbranchruleGetData(branchrule);
770 assert(branchruledata != NULL);
771
772 /* get current LP objective bound of the local sub problem and global cutoff bound */
773 lpobjval = SCIPgetLPObjval(scip);
774
775 /* check, if we want to solve the problem exactly, meaning that strong branching information is not useful
776 * for cutting off sub problems and improving lower bounds of children
777 */
778 exactsolve = SCIPisExactSolve(scip);
779
780 /* check, if all existing columns are in LP, and thus the strong branching results give lower bounds */
781 allcolsinlp = SCIPallColsInLP(scip);
782
783 bestcand = -1;
784 bestisstrongbranch = FALSE;
785 bestsbdown = SCIP_INVALID;
786 bestsbup = SCIP_INVALID;
787 bestsbdownvalid = FALSE;
788 bestsbupvalid = FALSE;
789 bestsbdowncutoff = FALSE;
790 bestsbupcutoff = FALSE;
791 provedbound = lpobjval;
792
793 /* Allocate memory to store the lp bounds of the up and down children
794 * for those of the variables that we performed sb on
795 */
796 SCIP_CALL( SCIPallocBufferArray(scip, &sbdown, nbranchcands) );
797 SCIP_CALL( SCIPallocBufferArray(scip, &sbup, nbranchcands) );
798 SCIP_CALL( SCIPallocBufferArray(scip, &sbdownvalid, nbranchcands) );
799 SCIP_CALL( SCIPallocBufferArray(scip, &sbupvalid, nbranchcands) );
800
801 if( nbranchcands == 1 )
802 {
803 /* only one candidate: nothing has to be done */
804 bestcand = 0;
805 SCIPdebug(ninitcands = 0);
806 sbdownvalid[0] = FALSE;
807 sbupvalid[0] = FALSE;
808 }
809 else
810 {
811 SCIP_VAR** vars;
812 int* initcands;
813 SCIP_Real* initcandscores;
814 SCIP_Real* newlbs = NULL;
815 SCIP_Real* newubs = NULL;
816 SCIP_Real* mingains = NULL;
817 SCIP_Real* maxgains = NULL;
818 /* scores computed from pseudocost branching */
819 SCIP_Real* scores = NULL;
820 SCIP_Real* scoresfrompc = NULL;
821 SCIP_Real* scoresfromothers = NULL;
822 int* bdchginds;
823 SCIP_BOUNDTYPE* bdchgtypes;
824 SCIP_Real* bdchgbounds;
825 int maxninitcands;
826 int nuninitcands;
827 int nbdchgs;
828 int nbdconflicts;
829 SCIP_Real avgconflictscore;
830 SCIP_Real avgconflengthscore;
831 SCIP_Real avginferencescore;
832 SCIP_Real avgcutoffscore;
833 SCIP_Real avgpscostscore;
834 SCIP_Real bestpsscore;
835 SCIP_Real bestpsfracscore;
836 SCIP_Real bestpsdomainscore;
837 SCIP_Real bestsbscore;
838 SCIP_Real bestuninitsbscore;
839 SCIP_Real bestsbfracscore;
840 SCIP_Real bestsbdomainscore;
841 SCIP_Real prio;
842 SCIP_Real reliable;
843 SCIP_Real maxlookahead;
844 SCIP_Real lookahead;
845 SCIP_Real relerrorthreshold;
846 SCIP_Bool initstrongbranching;
847 SCIP_Bool propagate;
848 SCIP_Bool probingbounds;
849 SCIP_Longint nodenum;
850 SCIP_Longint nlpiterationsquot;
851 SCIP_Longint nsblpiterations;
852 SCIP_Longint maxnsblpiterations;
853 int bestsolidx;
854 int maxbdchgs;
855 int bestpscand;
856 int bestsbcand;
857 int bestuninitsbcand;
858 int inititer;
859 int nvars;
860 int i;
861 int c;
863 SCIP_Real degeneracyfactor = 1.0;
864
865 /* get LP degeneracy information and compute a factor to change weighting of pseudo cost score vs. other scores */
866 if( branchruledata->degeneracyaware > 0 && (SCIPgetDepth(scip) > 0 || branchruledata->degeneracyaware > 1) )
867 {
868 SCIP_Real degeneracy;
869 SCIP_Real varconsratio;
870
871 /* get LP degeneracy information */
872 SCIP_CALL( SCIPgetLPDualDegeneracy(scip, &degeneracy, &varconsratio) );
873
874 assert(degeneracy >= 0.0);
875 assert(degeneracy <= 1.0);
876 assert(varconsratio >= 1.0);
877
878 /* increase factor for a degeneracy >= 80% */
879 if( degeneracy >= 0.8 )
880 {
881 degeneracy = 10.0 * (degeneracy - 0.7);
882 degeneracyfactor = degeneracyfactor * pow(10.0,degeneracy);
883 }
884 /* increase factor for a variable-constraint ratio >= 2.0 */
885 if( varconsratio >= 2.0 )
886 {
887 degeneracyfactor *= 10.0 * varconsratio;
888 }
889 }
890
891 vars = SCIPgetVars(scip);
892 nvars = SCIPgetNVars(scip);
893
894 bestsolidx = SCIPgetBestSol(scip) == NULL ? -1 : SCIPsolGetIndex(SCIPgetBestSol(scip));
895
896 /* get average conflict, inference, and pseudocost scores */
897 avgconflictscore = SCIPgetAvgConflictScore(scip);
898 avgconflictscore = MAX(avgconflictscore, 0.1);
899 avgconflengthscore = SCIPgetAvgConflictlengthScore(scip);
900 avgconflengthscore = MAX(avgconflengthscore, 0.1);
901 avginferencescore = SCIPgetAvgInferenceScore(scip);
902 avginferencescore = MAX(avginferencescore, 0.1);
903 avgcutoffscore = SCIPgetAvgCutoffScore(scip);
904 avgcutoffscore = MAX(avgcutoffscore, 0.1);
905 avgpscostscore = SCIPgetAvgPseudocostScore(scip);
906 avgpscostscore = MAX(avgpscostscore, 0.1);
907
908 /* get nonlinear counts according to parameters */
909 SCIP_CALL( branchruledataEnsureNlcount(scip, branchruledata) );
910
911 initstrongbranching = FALSE;
912
913 /* check whether propagation should be performed */
914 propagate = (branchruledata->maxproprounds != 0);
915
916 /* check whether valid bounds should be identified in probing-like fashion */
917 probingbounds = propagate && branchruledata->probingbounds;
918
919 /* get maximal number of candidates to initialize with strong branching; if the current solutions is not basic,
920 * we cannot warmstart the simplex algorithm and therefore don't initialize any candidates
921 */
922 maxninitcands = MIN(nbranchcands, branchruledata->initcand);
923 if( !SCIPisLPSolBasic(scip) )
924 maxninitcands = 0;
925
926 /* calculate maximal number of strong branching LP iterations; if we used too many, don't apply strong branching
927 * any more; also, if degeneracy is too high, don't run strong branching at this node
928 */
929 nlpiterationsquot = (SCIP_Longint)(branchruledata->sbiterquot * SCIPgetNNodeLPIterations(scip));
930 maxnsblpiterations = nlpiterationsquot + branchruledata->sbiterofs + SCIPgetNRootStrongbranchLPIterations(scip);
931 nsblpiterations = SCIPgetNStrongbranchLPIterations(scip);
932 if( nsblpiterations > maxnsblpiterations || degeneracyfactor >= 10.0 )
933 maxninitcands = 0;
934
935 /* get buffer for storing the unreliable candidates */
936 SCIP_CALL( SCIPallocBufferArray(scip, &initcands, maxninitcands+1) ); /* allocate one additional slot for convenience */
937 SCIP_CALL( SCIPallocBufferArray(scip, &initcandscores, maxninitcands+1) );
938 ninitcands = 0;
939
940 /* Allocate memory for the down and up gains, and the computed pseudocost scores */
941 SCIP_CALL( SCIPallocBufferArray(scip, &mingains, nbranchcands) );
942 SCIP_CALL( SCIPallocBufferArray(scip, &maxgains, nbranchcands) );
943 SCIP_CALL( SCIPallocBufferArray(scip, &scores, nbranchcands) );
944 SCIP_CALL( SCIPallocBufferArray(scip, &scoresfrompc, nbranchcands) );
945 SCIP_CALL( SCIPallocBufferArray(scip, &scoresfromothers, nbranchcands) );
946
947 /* get current node number */
948 nodenum = SCIPgetNNodes(scip);
949
950 /* initialize bound change arrays */
951 bdchginds = NULL;
952 bdchgtypes = NULL;
953 bdchgbounds = NULL;
954 nbdchgs = 0;
955 nbdconflicts = 0;
956 maxbdchgs = branchruledata->maxbdchgs;
957 if( maxbdchgs == -1 )
958 maxbdchgs = INT_MAX;
959
960 /* calculate value used as reliability */
961 prio = (maxnsblpiterations - nsblpiterations)/(nsblpiterations + 1.0);
962 prio = MIN(prio, 1.0);
963 prio = MAX(prio, (nlpiterationsquot - nsblpiterations)/(nsblpiterations + 1.0));
964 reliable = (1.0-prio) * branchruledata->minreliable + prio * branchruledata->maxreliable;
965
966 /* calculate the threshold for the relative error in the same way; low tolerance is more strict than higher tolerance */
967 relerrorthreshold = (1.0 - prio) * branchruledata->higherrortol + prio * branchruledata->lowerrortol;
968
969 clevel = (SCIP_CONFIDENCELEVEL)branchruledata->confidencelevel;
970 /* determine the confidence level for hypothesis testing based on value of prio */
971 if( branchruledata->usedynamicconfidence )
972 {
973 /* with decreasing priority, use a less strict confidence level */
974 if( prio >= 0.9 )
976 else if( prio >= 0.7 )
978 else if( prio >= 0.5 )
980 else if( prio >= 0.3 )
982 else
984 }
985
986 /* search for the best pseudo cost candidate, while remembering unreliable candidates in a sorted buffer */
987 nuninitcands = 0;
988 bestpscand = -1;
989 bestpsscore = -SCIPinfinity(scip);
990 bestpsfracscore = -SCIPinfinity(scip);
991 bestpsdomainscore = -SCIPinfinity(scip);
992
993 /* search for the best candidate first */
994 if( branchruledata->usehyptestforreliability )
995 {
996 for( c = 0; c < nbranchcands; ++c )
997 {
998 SCIP_Real conflictscore;
999 SCIP_Real conflengthscore;
1000 SCIP_Real inferencescore;
1001 SCIP_Real cutoffscore;
1002 SCIP_Real gmieffscore;
1003 SCIP_Real lastgmieffscore;
1004 SCIP_Real pscostscore;
1005 SCIP_Real nlscore;
1006 SCIP_Real score;
1007
1008 conflictscore = SCIPgetVarConflictScore(scip, branchcands[c]);
1009 conflengthscore = SCIPgetVarConflictlengthScore(scip, branchcands[c]);
1010 inferencescore = SCIPgetVarAvgInferenceScore(scip, branchcands[c]);
1011 cutoffscore = SCIPgetVarAvgCutoffScore(scip, branchcands[c]);
1012 gmieffscore = SCIPgetVarAvgGMIScore(scip, branchcands[c]);
1013 lastgmieffscore = SCIPgetVarLastGMIScore(scip, branchcands[c]);
1014 nlscore = calcNlscore(scip, branchruledata->nlcount, branchruledata->nlcountmax, SCIPvarGetProbindex(branchcands[c]));
1015 pscostscore = SCIPgetVarPseudocostScore(scip, branchcands[c], branchcandssol[c]);
1016
1017 /* replace the pseudo cost score with the already calculated one;
1018 * @todo: use old data for strong branching with propagation?
1019 */
1020 if( SCIPgetVarStrongbranchNode(scip, branchcands[c]) == nodenum )
1021 {
1022 SCIP_Real down;
1023 SCIP_Real up;
1024 SCIP_Real lastlpobjval;
1025 SCIP_Real downgain;
1026 SCIP_Real upgain;
1027
1028 /* use the score of the strong branching call at the current node */
1029 SCIP_CALL( SCIPgetVarStrongbranchLast(scip, branchcands[c], &down, &up, NULL, NULL, NULL, &lastlpobjval) );
1030 downgain = MAX(down - lastlpobjval, 0.0);
1031 upgain = MAX(up - lastlpobjval, 0.0);
1032 pscostscore = SCIPgetBranchScore(scip, branchcands[c], downgain, upgain);
1033
1034 SCIPdebugMsg(scip, " -> strong branching on variable <%s> already performed (down=%g (%+g), up=%g (%+g), pscostscore=%g)\n",
1035 SCIPvarGetName(branchcands[c]), down, downgain, up, upgain, pscostscore);
1036 }
1037
1038 score = calcScore(scip, branchruledata, conflictscore, avgconflictscore, conflengthscore, avgconflengthscore,
1039 inferencescore, avginferencescore, cutoffscore, avgcutoffscore, gmieffscore, lastgmieffscore,
1040 pscostscore, avgpscostscore, nlscore, branchcandsfrac[c], degeneracyfactor);
1041
1042 /* check for better score of candidate */
1043 if( SCIPisSumGE(scip, score, bestpsscore) )
1044 {
1045 SCIP_Real fracscore;
1046 SCIP_Real domainscore;
1047
1048 fracscore = MIN(branchcandsfrac[c], 1.0 - branchcandsfrac[c]);
1049 domainscore = -(SCIPvarGetUbLocal(branchcands[c]) - SCIPvarGetLbLocal(branchcands[c]));
1050 if( SCIPisSumGT(scip, score, bestpsscore)
1051 || SCIPisSumGT(scip, fracscore, bestpsfracscore)
1052 || (SCIPisSumGE(scip, fracscore, bestpsfracscore) && domainscore > bestpsdomainscore) )
1053 {
1054 bestpscand = c;
1055 bestpsscore = score;
1056 bestpsfracscore = fracscore;
1057 bestpsdomainscore = domainscore;
1058 }
1059 }
1060 }
1061 }
1062
1063 for( c = 0; c < nbranchcands; ++c )
1064 {
1065 SCIP_Real conflictscore;
1066 SCIP_Real conflengthscore;
1067 SCIP_Real inferencescore;
1068 SCIP_Real cutoffscore;
1069 SCIP_Real gmieffscore;
1070 SCIP_Real lastgmieffscore;
1071 SCIP_Real pscostscore;
1072 SCIP_Real nlscore;
1073 SCIP_Real score;
1074 SCIP_Bool usesb;
1075 SCIP_Real downgain;
1076 SCIP_Real upgain;
1077 SCIP_Real fracpart;
1078
1079 assert(branchcands[c] != NULL);
1080 assert(!SCIPisFeasIntegral(scip, branchcandssol[c]));
1081 assert(!SCIPisFeasIntegral(scip, SCIPvarGetLPSol(branchcands[c])));
1082
1083 /* Record the variables current pseudocosts. These may be overwritten if
1084 * strong branching is performed.
1085 */
1086 sbdownvalid[c] = FALSE;
1087 sbupvalid[c] = FALSE;
1088 fracpart = SCIPfeasFrac(scip, SCIPvarGetLPSol(branchcands[c]));
1089 downgain = SCIPgetVarPseudocostVal(scip, branchcands[c], 0.0 - fracpart);
1090 upgain = SCIPgetVarPseudocostVal(scip, branchcands[c], 1.0 - fracpart);
1091 mingains[c] = MIN(downgain, upgain);
1092 maxgains[c] = MAX(downgain, upgain);
1093
1094 /* get conflict, inference, cutoff, nonlinear, and pseudo cost scores for candidate */
1095 conflictscore = SCIPgetVarConflictScore(scip, branchcands[c]);
1096 conflengthscore = SCIPgetVarConflictlengthScore(scip, branchcands[c]);
1097 inferencescore = SCIPgetVarAvgInferenceScore(scip, branchcands[c]);
1098 cutoffscore = SCIPgetVarAvgCutoffScore(scip, branchcands[c]);
1099 gmieffscore = SCIPgetVarAvgGMIScore(scip, branchcands[c]);
1100 lastgmieffscore = SCIPgetVarLastGMIScore(scip, branchcands[c]);
1101 nlscore = calcNlscore(scip, branchruledata->nlcount, branchruledata->nlcountmax, SCIPvarGetProbindex(branchcands[c]));
1102 pscostscore = SCIPgetVarPseudocostScore(scip, branchcands[c], branchcandssol[c]);
1103 usesb = FALSE;
1104
1105 /* don't use strong branching on variables that have already been initialized at the current node;
1106 * instead replace the pseudo cost score with the already calculated one;
1107 * @todo: use old data for strong branching with propagation?
1108 */
1109 if( SCIPgetVarStrongbranchNode(scip, branchcands[c]) == nodenum )
1110 {
1111 SCIP_Real down;
1112 SCIP_Real up;
1113 SCIP_Real lastlpobjval;
1114
1115 /* use the score of the strong branching call at the current node */
1116 SCIP_CALL( SCIPgetVarStrongbranchLast(scip, branchcands[c], &down, &up, NULL, NULL, NULL, &lastlpobjval) );
1117 downgain = MAX(down - lastlpobjval, 0.0);
1118 upgain = MAX(up - lastlpobjval, 0.0);
1119 pscostscore = SCIPgetBranchScore(scip, branchcands[c], downgain, upgain);
1120
1121 mingains[c] = MIN(downgain, upgain);
1122 maxgains[c] = MAX(downgain, upgain);
1123
1124 SCIPdebugMsg(scip, " -> strong branching on variable <%s> already performed (down=%g (%+g), up=%g (%+g), pscostscore=%g)\n",
1125 SCIPvarGetName(branchcands[c]), down, downgain, up, upgain, pscostscore);
1126 }
1127 else if( maxninitcands > 0 )
1128 {
1129 SCIP_Real downsize;
1130 SCIP_Real upsize;
1131 SCIP_Real size;
1132
1133 /* check, if the pseudo cost score of the variable is reliable */
1136 size = MIN(downsize, upsize);
1137
1138 /* determine if variable is considered reliable based on the current reliability setting */
1139 /* check fixed number threshold (aka original) reliability first */
1140 assert(!branchruledata->usehyptestforreliability || bestpscand >= 0);
1141 usesb = FALSE;
1142 if( size < reliable )
1143 usesb = TRUE;
1144 else if( branchruledata->userelerrorforreliability && branchruledata->usehyptestforreliability )
1145 {
1146 if( !SCIPisVarPscostRelerrorReliable(scip, branchcands[c], relerrorthreshold, clevel) &&
1147 !SCIPsignificantVarPscostDifference(scip, branchcands[bestpscand], branchcandsfrac[bestpscand],
1148 branchcands[c], branchcandsfrac[c], SCIP_BRANCHDIR_DOWNWARDS, clevel, TRUE) &&
1149 !SCIPsignificantVarPscostDifference(scip, branchcands[bestpscand], 1 - branchcandsfrac[bestpscand],
1150 branchcands[c], 1 - branchcandsfrac[c], SCIP_BRANCHDIR_UPWARDS, clevel, TRUE) )
1151 usesb = TRUE;
1152 }
1153 /* check if relative error is tolerable */
1154 else if( branchruledata->userelerrorforreliability &&
1155 !SCIPisVarPscostRelerrorReliable(scip, branchcands[c], relerrorthreshold, clevel))
1156 usesb = TRUE;
1157 /* check if best pseudo-candidate is significantly better in both directions, use strong-branching otherwise */
1158 else if( branchruledata->usehyptestforreliability &&
1159 !SCIPsignificantVarPscostDifference(scip, branchcands[bestpscand], branchcandsfrac[bestpscand],
1160 branchcands[c], branchcandsfrac[c], SCIP_BRANCHDIR_DOWNWARDS, clevel, TRUE) &&
1161 !SCIPsignificantVarPscostDifference(scip, branchcands[bestpscand], 1 - branchcandsfrac[bestpscand],
1162 branchcands[c], 1 - branchcandsfrac[c], SCIP_BRANCHDIR_UPWARDS, clevel, TRUE))
1163 usesb = TRUE;
1164
1165 /* count the number of variables that are completely uninitialized */
1166 if( size < 0.1 )
1167 nuninitcands++;
1168 }
1169
1170 /* combine the five score values */
1171 scoresfrompc[c] = calcScore(scip, branchruledata, 0.0, avgconflictscore, 0.0, avgconflengthscore,
1172 0.0, avginferencescore, 0.0, avgcutoffscore, 0.0, 0.0,
1173 pscostscore, avgpscostscore, 0.0, branchcandsfrac[c], degeneracyfactor);
1174 scoresfromothers[c] = calcScore(scip, branchruledata, conflictscore, avgconflictscore, conflengthscore, avgconflengthscore,
1175 inferencescore, avginferencescore, cutoffscore, avgcutoffscore, gmieffscore, lastgmieffscore,
1176 0.0, avgpscostscore, nlscore, branchcandsfrac[c], degeneracyfactor);
1177 score = scoresfrompc[c] + scoresfromothers[c];
1178 scores[c] = score;
1179 /*score = calcScore(scip, branchruledata, conflictscore, avgconflictscore, conflengthscore, avgconflengthscore,
1180 inferencescore, avginferencescore, cutoffscore, avgcutoffscore, gmieffscore, lastgmieffscore,
1181 pscostscore, avgpscostscore, nlscore, branchcandsfrac[c], degeneracyfactor);*/
1182 if( usesb )
1183 {
1184 int j;
1185
1186 mingains[c] = 0;
1187 maxgains[c] = 0;
1188 scoresfrompc[c] = 0;
1189 scoresfromothers[c] = 0;
1190
1191 /* assign a random score to this uninitialized candidate */
1192 if( branchruledata->randinitorder )
1193 score += SCIPrandomGetReal(branchruledata->randnumgen, 0.0, 1e-4);
1194
1195 /* pseudo cost of variable is not reliable: insert candidate in initcands buffer */
1196 for( j = ninitcands; j > 0 && score > initcandscores[j-1]; --j )
1197 {
1198 initcands[j] = initcands[j-1];
1199 initcandscores[j] = initcandscores[j-1];
1200 }
1201 initcands[j] = c;
1202 initcandscores[j] = score;
1203 ninitcands++;
1204 ninitcands = MIN(ninitcands, maxninitcands);
1205 }
1206 /* in the case of hypothesis reliability, the best pseudo candidate has been determined already */
1207 else if( !branchruledata->usehyptestforreliability )
1208 {
1209 /* variable will keep its pseudo cost value: check for better score of candidate */
1210 if( SCIPisSumGE(scip, score, bestpsscore) )
1211 {
1212 SCIP_Real fracscore;
1213 SCIP_Real domainscore;
1214
1215 fracscore = MIN(branchcandsfrac[c], 1.0 - branchcandsfrac[c]);
1216 domainscore = -(SCIPvarGetUbLocal(branchcands[c]) - SCIPvarGetLbLocal(branchcands[c]));
1217 if( SCIPisSumGT(scip, score, bestpsscore)
1218 || SCIPisSumGT(scip, fracscore, bestpsfracscore)
1219 || (SCIPisSumGE(scip, fracscore, bestpsfracscore) && domainscore > bestpsdomainscore) )
1220 {
1221 bestpscand = c;
1222 bestpsscore = score;
1223 bestpsfracscore = fracscore;
1224 bestpsdomainscore = domainscore;
1225 }
1226 }
1227 }
1228 }
1229
1230 /* in the special case that only the best pseudo candidate was selected for strong branching, skip the strong branching */
1231 if( branchruledata->usehyptestforreliability && ninitcands == 1 )
1232 {
1233 ninitcands = 0;
1234 SCIPdebugMsg(scip, "Only one single candidate for initialization-->Skipping strong branching\n");
1235 }
1236
1237 /* initialize unreliable candidates with strong branching until maxlookahead is reached,
1238 * search best strong branching candidate
1239 */
1240 maxlookahead = (SCIP_Real)branchruledata->maxlookahead * (1.0 + (SCIP_Real)nuninitcands/(SCIP_Real)nbranchcands);
1241 inititer = branchruledata->inititer;
1242 if( inititer == 0 )
1243 {
1244 SCIP_Longint nlpiterations;
1245 SCIP_Longint nlps;
1246
1247 /* iteration limit is set to twice the average number of iterations spent to resolve a dual feasible SCIP_LP;
1248 * at the first few nodes, this average is not very exact, so we better increase the iteration limit on
1249 * these very important nodes
1250 */
1251 nlpiterations = SCIPgetNDualResolveLPIterations(scip);
1253 if( nlps == 0 )
1254 {
1255 nlpiterations = SCIPgetNNodeInitLPIterations(scip);
1256 nlps = SCIPgetNNodeInitLPs(scip);
1257 if( nlps == 0 )
1258 {
1259 nlpiterations = 1000;
1260 nlps = 1;
1261 }
1262 }
1263 assert(nlps >= 1);
1264 inititer = (int)(2*nlpiterations / nlps);
1265 inititer = (int)((SCIP_Real)inititer * (1.0 + 20.0/nodenum));
1266 inititer = MAX(inititer, 10);
1267 inititer = MIN(inititer, 500);
1268 }
1269
1270 SCIPdebugMsg(scip, "strong branching (reliable=%g, %d/%d cands, %d uninit, maxcands=%d, maxlookahead=%g, maxbdchgs=%d, inititer=%d, iterations:%" SCIP_LONGINT_FORMAT "/%" SCIP_LONGINT_FORMAT ", basic:%u)\n",
1271 reliable, ninitcands, nbranchcands, nuninitcands, maxninitcands, maxlookahead, maxbdchgs, inititer,
1273
1274 bestsbcand = -1;
1275 bestsbscore = -SCIPinfinity(scip);
1276 bestsbfracscore = -SCIPinfinity(scip);
1277 bestsbdomainscore = -SCIPinfinity(scip);
1278 bestuninitsbscore = -SCIPinfinity(scip);
1279 bestuninitsbcand = -1;
1280 lookahead = 0.0;
1281 for( i = 0; i < ninitcands && lookahead < maxlookahead && nbdchgs + nbdconflicts < maxbdchgs
1282 && (i < (int) maxlookahead || SCIPgetNStrongbranchLPIterations(scip) < maxnsblpiterations); ++i )
1283 {
1284 SCIP_Real down;
1285 SCIP_Real up;
1286 SCIP_Real downgain;
1287 SCIP_Real upgain;
1288 SCIP_Bool downvalid;
1289 SCIP_Bool upvalid;
1290 SCIP_Longint ndomredsdown;
1291 SCIP_Longint ndomredsup;
1292 SCIP_Bool lperror;
1293 SCIP_Bool downinf;
1294 SCIP_Bool upinf;
1295 SCIP_Bool downconflict;
1296 SCIP_Bool upconflict;
1297
1298 /* get candidate number to initialize */
1299 c = initcands[i];
1300 assert(!SCIPisFeasIntegral(scip, branchcandssol[c]));
1301
1302 if( branchruledata->skipbadinitcands )
1303 {
1304 SCIP_Bool skipsb = FALSE;
1305 /* if the current best candidate is a candidate found by strong branching, determine if candidate pseudo-costs are
1306 * significantly smaller in at least one direction, in which case we safe the execution of strong-branching for now
1307 */
1308 if( bestsbscore > bestpsscore && bestsbscore > bestuninitsbscore && bestsbupvalid && bestsbdownvalid )
1309 {
1310 assert(bestsbcand != -1);
1311 assert(bestsbup != SCIP_INVALID && bestsbdown != SCIP_INVALID); /*lint !e777 lint doesn't like comparing floats */
1312
1313 /* test if the variable is unlikely to produce a better gain than the currently best one. Skip strong-branching
1314 * in such a case
1315 */
1316 if( SCIPpscostThresholdProbabilityTest(scip, branchcands[c], branchcandsfrac[c], bestsbdown,
1318 || SCIPpscostThresholdProbabilityTest(scip, branchcands[c], 1.0 - branchcandsfrac[c], bestsbup,
1319 SCIP_BRANCHDIR_UPWARDS, clevel) )
1320 skipsb = TRUE;
1321 }
1322 /* the currently best candidate is also a pseudo-candidate; apply significance test and skip candidate if it
1323 * is significantly worse in at least one direction
1324 */
1325 else if( bestpscand != -1 && bestpsscore > bestuninitsbscore )
1326 {
1327 if( SCIPsignificantVarPscostDifference(scip, branchcands[bestpscand], branchcandsfrac[bestpscand],
1328 branchcands[c], branchcandsfrac[c], SCIP_BRANCHDIR_DOWNWARDS, clevel, TRUE)
1329 || SCIPsignificantVarPscostDifference(scip, branchcands[bestpscand], 1.0 - branchcandsfrac[bestpscand],
1330 branchcands[c], 1.0 - branchcandsfrac[c], SCIP_BRANCHDIR_UPWARDS, clevel, TRUE) )
1331 skipsb = TRUE;
1332 }
1333 /* compare against the best init cand that has been skipped already */
1334 else if( bestuninitsbcand != -1 )
1335 {
1336 if( SCIPsignificantVarPscostDifference(scip, branchcands[bestuninitsbcand], branchcandsfrac[bestuninitsbcand],
1337 branchcands[c], branchcandsfrac[c], SCIP_BRANCHDIR_DOWNWARDS, clevel, TRUE)
1338 || SCIPsignificantVarPscostDifference(scip, branchcands[bestuninitsbcand], 1.0 - branchcandsfrac[bestuninitsbcand],
1339 branchcands[c], 1.0 - branchcandsfrac[c], SCIP_BRANCHDIR_UPWARDS, clevel, TRUE) )
1340 skipsb = TRUE;
1341 }
1342
1343 /* skip candidate, update the best score of an unitialized candidate */
1344 if( skipsb )
1345 {
1346 if( bestuninitsbcand == -1 )
1347 {
1348 bestuninitsbcand = c;
1349 bestuninitsbscore = initcandscores[i];
1350 }
1351 continue;
1352 }
1353 }
1354 SCIPdebugMsg(scip, "init pseudo cost (%g/%g) of <%s> at %g (score:%g) with strong branching (%d iterations) -- %" SCIP_LONGINT_FORMAT "/%" SCIP_LONGINT_FORMAT " iterations\n",
1357 SCIPvarGetName(branchcands[c]), branchcandssol[c], initcandscores[i],
1358 inititer, SCIPgetNStrongbranchLPIterations(scip), maxnsblpiterations);
1359
1360 /* use strong branching on candidate */
1361 if( !initstrongbranching )
1362 {
1363 initstrongbranching = TRUE;
1364
1365 SCIP_CALL( SCIPstartStrongbranch(scip, propagate) );
1366
1367 /* create arrays for probing-like bound tightening */
1368 if( probingbounds )
1369 {
1370 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &newlbs, nvars) );
1371 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &newubs, nvars) );
1372 }
1373 }
1374
1375 if( propagate )
1376 {
1377 /* apply strong branching */
1378 SCIP_CALL( SCIPgetVarStrongbranchWithPropagation(scip, branchcands[c], branchcandssol[c], lpobjval, inititer,
1379 branchruledata->maxproprounds, &down, &up, &downvalid, &upvalid, &ndomredsdown, &ndomredsup, &downinf, &upinf,
1380 &downconflict, &upconflict, &lperror, newlbs, newubs) );
1381 }
1382 else
1383 {
1384 /* apply strong branching */
1385 SCIP_CALL( SCIPgetVarStrongbranchFrac(scip, branchcands[c], inititer, FALSE,
1386 &down, &up, &downvalid, &upvalid, &downinf, &upinf, &downconflict, &upconflict, &lperror) );
1387
1388 ndomredsdown = ndomredsup = 0;
1389 }
1390
1391 /* check for an error in strong branching */
1392 if( lperror )
1393 {
1394 if( !SCIPisStopped(scip) )
1395 {
1397 "(node %" SCIP_LONGINT_FORMAT ") error in strong branching call for variable <%s> with solution %g\n",
1398 SCIPgetNNodes(scip), SCIPvarGetName(branchcands[c]), branchcandssol[c]);
1399 }
1400 break;
1401 }
1402
1403 /* Strong branching might have found a new primal solution which updated the cutoff bound. In this case, the
1404 * provedbound computed before can be higher than the cutoffbound and the current node can be cut off.
1405 * Additionally, also if the value for the current best candidate is valid and exceeds the new cutoff bound,
1406 * we want to change the domain of this variable rather than branching on it.
1407 */
1408 if( SCIPgetBestSol(scip) != NULL && SCIPsolGetIndex(SCIPgetBestSol(scip)) != bestsolidx )
1409 {
1410 bestsolidx = SCIPsolGetIndex(SCIPgetBestSol(scip));
1411
1412 SCIPdebugMsg(scip, " -> strong branching on variable <%s> lead to a new incumbent\n",
1413 SCIPvarGetName(branchcands[c]));
1414
1415 /* proved bound for current node is larger than new cutoff bound -> cut off current node */
1416 if( SCIPisGE(scip, provedbound, SCIPgetCutoffbound(scip)) )
1417 {
1418 SCIPdebugMsg(scip, " -> node can be cut off (provedbound=%g, cutoff=%g)\n", provedbound, SCIPgetCutoffbound(scip));
1419
1420 *result = SCIP_CUTOFF;
1421 break; /* terminate initialization loop, because node is infeasible */
1422 }
1423 /* proved bound for down child of best candidate is larger than cutoff bound
1424 * -> increase lower bound of best candidate
1425 * we must only do this if the LP is complete, i.e., we are not doing column generation
1426 */
1427
1428 else if( bestsbcand != -1 && allcolsinlp )
1429 {
1430 if( !bestsbdowncutoff && bestsbdownvalid && SCIPisGE(scip, bestsbdown, SCIPgetCutoffbound(scip)) )
1431 {
1432 bestsbdowncutoff = TRUE;
1433
1434 SCIPdebugMsg(scip, " -> valid dual bound for down child of best candidate <%s> is higher than new cutoff bound (valid=%u, bestsbdown=%g, cutoff=%g)\n",
1435 SCIPvarGetName(branchcands[bestsbcand]), bestsbdownvalid, bestsbdown, SCIPgetCutoffbound(scip));
1436
1437 SCIPdebugMsg(scip, " -> increase lower bound of best candidate <%s> to %g\n",
1438 SCIPvarGetName(branchcands[bestsbcand]), SCIPfeasCeil(scip, branchcandssol[bestsbcand]));
1439
1440 SCIP_CALL( addBdchg(scip, &bdchginds, &bdchgtypes, &bdchgbounds, &nbdchgs, SCIPvarGetProbindex(branchcands[bestsbcand]),
1441 SCIP_BOUNDTYPE_LOWER, SCIPfeasCeil(scip, branchcandssol[bestsbcand])) );
1442 }
1443 /* proved bound for up child of best candidate is larger than cutoff bound -> decrease upper bound of best candidate */
1444 else if( !bestsbupcutoff && bestsbupvalid && SCIPisGE(scip, bestsbup, SCIPgetCutoffbound(scip)) )
1445 {
1446 bestsbupcutoff = TRUE;
1447
1448 SCIPdebugMsg(scip, " -> valid dual bound for up child of best candidate <%s> is higher than new cutoff bound (valid=%u, bestsbup=%g, cutoff=%g)\n",
1449 SCIPvarGetName(branchcands[bestsbcand]), bestsbupvalid, bestsbup, SCIPgetCutoffbound(scip));
1450
1451 SCIPdebugMsg(scip, " -> decrease upper bound of best candidate <%s> to %g\n",
1452 SCIPvarGetName(branchcands[bestsbcand]), SCIPfeasFloor(scip, branchcandssol[bestsbcand]));
1453
1454 SCIP_CALL( addBdchg(scip, &bdchginds, &bdchgtypes, &bdchgbounds, &nbdchgs, SCIPvarGetProbindex(branchcands[bestsbcand]),
1455 SCIP_BOUNDTYPE_UPPER, SCIPfeasFloor(scip, branchcandssol[bestsbcand])) );
1456 }
1457 }
1458 }
1459
1460 /* evaluate strong branching */
1461 down = MAX(down, lpobjval);
1462 up = MAX(up, lpobjval);
1463 downgain = down - lpobjval;
1464 upgain = up - lpobjval;
1465 assert(!allcolsinlp || exactsolve || !downvalid || downinf == SCIPisGE(scip, down, SCIPgetCutoffbound(scip)));
1466 assert(!allcolsinlp || exactsolve || !upvalid || upinf == SCIPisGE(scip, up, SCIPgetCutoffbound(scip)));
1467 assert(downinf || !downconflict);
1468 assert(upinf || !upconflict);
1469
1470 /* @todo: store pseudo cost only for valid bounds?
1471 * depending on the user parameter choice of storesemiinitcosts, pseudo costs are also updated in single directions,
1472 * if the node in the other direction was infeasible or cut off
1473 */
1474 if( !downinf
1475#ifdef WITH_LPSOLSTAT
1477#endif
1478 && (!upinf || branchruledata->storesemiinitcosts) )
1479 {
1480 SCIP_Real weight;
1481
1482 /* smaller weights are given if the strong branching hit the time limit in the corresponding direction */
1483 if( branchruledata->usesmallweightsitlim )
1485 else
1486 weight = 1.0;
1487
1488 /* update pseudo cost values */
1489 SCIP_CALL( SCIPupdateVarPseudocostSymmetric(scip, branchruledata, branchcands[c], branchorbitidx, c, 0.0 - branchcandsfrac[c], downgain, weight) );
1490 }
1491 if( !upinf
1492#ifdef WITH_LPSOLSTAT
1494#endif
1495 && (!downinf || branchruledata->storesemiinitcosts) )
1496 {
1497 SCIP_Real weight;
1498
1499 /* smaller weights are given if the strong branching hit the time limit in the corresponding direction */
1500 if( branchruledata->usesmallweightsitlim )
1502 else
1503 weight = 1.0;
1504
1505 /* update pseudo cost values */
1506 SCIP_CALL( SCIPupdateVarPseudocostSymmetric(scip, branchruledata, branchcands[c], branchorbitidx, c, 1.0 - branchcandsfrac[c], upgain, weight) );
1507 }
1508
1509 /* the minimal lower bound of both children is a proved lower bound of the current subtree */
1510 if( allcolsinlp && !exactsolve && downvalid && upvalid )
1511 {
1512 SCIP_Real minbound;
1513
1514 minbound = MIN(down, up);
1515 provedbound = MAX(provedbound, minbound);
1516 assert((downinf && upinf) || SCIPisLT(scip, provedbound, SCIPgetCutoffbound(scip)));
1517
1518 /* save probing-like bounds detected during strong branching */
1519 if( probingbounds && ( !downinf || !upinf ) )
1520 {
1521 int v;
1522
1523 assert(newlbs != NULL);
1524 assert(newubs != NULL);
1525
1526 for( v = 0; v < nvars; ++v )
1527 {
1528 if( SCIPisGT(scip, newlbs[v], SCIPvarGetLbLocal(vars[v])) )
1529 {
1530 SCIPdebugMsg(scip, "better lower bound for variable <%s>: %.9g -> %.9g (by strongbranching on <%s>)\n",
1531 SCIPvarGetName(vars[v]), SCIPvarGetLbLocal(vars[v]), newlbs[v], SCIPvarGetName(branchcands[c]));
1532
1533 SCIP_CALL( addBdchg(scip, &bdchginds, &bdchgtypes, &bdchgbounds, &nbdchgs, v,
1534 SCIP_BOUNDTYPE_LOWER, newlbs[v]) );
1535 }
1536 if( SCIPisLT(scip, newubs[v], SCIPvarGetUbLocal(vars[v])) )
1537 {
1538 SCIPdebugMsg(scip, "better upper bound for variable <%s>: %.9g -> %.9g (by strongbranching on <%s>)\n",
1539 SCIPvarGetName(vars[v]), SCIPvarGetUbLocal(vars[v]), newubs[v], SCIPvarGetName(branchcands[c]));
1540
1541 SCIP_CALL( addBdchg(scip, &bdchginds, &bdchgtypes, &bdchgbounds, &nbdchgs, v,
1542 SCIP_BOUNDTYPE_UPPER, newubs[v]) );
1543 }
1544 }
1545 }
1546 }
1547
1548 /* check if there are infeasible roundings */
1549 if( downinf || upinf )
1550 {
1551 assert(allcolsinlp || propagate);
1552 assert(!exactsolve);
1553
1554 if( downinf && upinf )
1555 {
1556 /* both roundings are infeasible -> node is infeasible */
1557 SCIPdebugMsg(scip, " -> variable <%s> is infeasible in both directions (conflict: %u/%u)\n",
1558 SCIPvarGetName(branchcands[c]), downconflict, upconflict);
1559 *result = SCIP_CUTOFF;
1560 break; /* terminate initialization loop, because node is infeasible */
1561 }
1562 else
1563 {
1564 /* rounding is infeasible in one direction -> round variable in other direction */
1565 SCIPdebugMsg(scip, " -> variable <%s> is infeasible in %s branch (conflict: %u/%u)\n",
1566 SCIPvarGetName(branchcands[c]), downinf ? "downward" : "upward", downconflict, upconflict);
1567 SCIP_CALL( addBdchg(scip, &bdchginds, &bdchgtypes, &bdchgbounds, &nbdchgs, SCIPvarGetProbindex(branchcands[c]),
1569 (downinf ? SCIPfeasCeil(scip, branchcandssol[c]) : SCIPfeasFloor(scip, branchcandssol[c]))) );
1570 if( nbdchgs + nbdconflicts >= maxbdchgs )
1571 break; /* terminate initialization loop, because enough roundings are performed */
1572 }
1573 }
1574 else
1575 {
1576 SCIP_Real conflictscore;
1577 SCIP_Real conflengthscore;
1578 SCIP_Real inferencescore;
1579 SCIP_Real cutoffscore;
1580 SCIP_Real gmieffscore;
1581 SCIP_Real lastgmieffscore;
1582 SCIP_Real pscostscore;
1583 SCIP_Real nlscore;
1584 SCIP_Real score;
1585
1586 mingains[c] = MIN(downgain, upgain);
1587 maxgains[c] = MAX(downgain, upgain);
1588
1589 sbdown[c] = down;
1590 sbup[c] = up;
1591 sbdownvalid[c] = downvalid;
1592 sbupvalid[c] = upvalid;
1593
1594 /* check for a better score */
1595 conflictscore = SCIPgetVarConflictScore(scip, branchcands[c]);
1596 conflengthscore = SCIPgetVarConflictlengthScore(scip, branchcands[c]);
1597 nlscore = calcNlscore(scip, branchruledata->nlcount, branchruledata->nlcountmax, SCIPvarGetProbindex(branchcands[c]));
1598
1599 /* optionally, use only local information obtained via strong branching for this candidate, i.e., local
1600 * domain reductions and no cutoff score
1601 */
1602 inferencescore = branchruledata->usesblocalinfo ? SCIPgetBranchScore(scip, branchcands[c], (SCIP_Real)ndomredsdown, (SCIP_Real)ndomredsup)
1603 : SCIPgetVarAvgInferenceScore(scip, branchcands[c]);
1604 cutoffscore = branchruledata->usesblocalinfo ? 0.0 : SCIPgetVarAvgCutoffScore(scip, branchcands[c]);
1605 gmieffscore = branchruledata->usesblocalinfo ? 0.0 : SCIPgetVarAvgGMIScore(scip, branchcands[c]);
1606 lastgmieffscore = branchruledata->usesblocalinfo ? 0.0 : SCIPgetVarLastGMIScore(scip, branchcands[c]);
1607 pscostscore = SCIPgetBranchScore(scip, branchcands[c], downgain, upgain);
1608
1609 scoresfrompc[c] = calcScore(scip, branchruledata, 0.0, avgconflictscore, 0.0, avgconflengthscore,
1610 0.0, avginferencescore, 0.0, avgcutoffscore, 0.0, 0.0, pscostscore,
1611 avgpscostscore, 0.0, branchcandsfrac[c], degeneracyfactor);
1612 scoresfromothers[c] = calcScore(scip, branchruledata, conflictscore, avgconflictscore, conflengthscore, avgconflengthscore,
1613 inferencescore, avginferencescore, cutoffscore, avgcutoffscore, gmieffscore,
1614 lastgmieffscore, 0.0, avgpscostscore, nlscore, branchcandsfrac[c],
1615 degeneracyfactor);
1616 score = scoresfrompc[c] + scoresfromothers[c];
1617 scores[c] = score;
1618 /*score = calcScore(scip, branchruledata, conflictscore, avgconflictscore, conflengthscore, avgconflengthscore,
1619 inferencescore, avginferencescore, cutoffscore, avgcutoffscore, gmieffscore, lastgmieffscore,
1620 pscostscore, avgpscostscore, nlscore, branchcandsfrac[c], degeneracyfactor);*/
1621
1622 if( SCIPisSumGE(scip, score, bestsbscore) )
1623 {
1624 SCIP_Real fracscore;
1625 SCIP_Real domainscore;
1626
1627 fracscore = MIN(branchcandsfrac[c], 1.0 - branchcandsfrac[c]);
1628 domainscore = -(SCIPvarGetUbLocal(branchcands[c]) - SCIPvarGetLbLocal(branchcands[c]));
1629 if( SCIPisSumGT(scip, score, bestsbscore)
1630 || SCIPisSumGT(scip, fracscore, bestsbfracscore)
1631 || (SCIPisSumGE(scip, fracscore, bestsbfracscore) && domainscore > bestsbdomainscore) )
1632 {
1633 bestsbcand = c;
1634 bestsbscore = score;
1635 bestsbdown = down;
1636 bestsbup = up;
1637 bestsbdownvalid = downvalid;
1638 bestsbupvalid = upvalid;
1639 bestsbdowncutoff = FALSE;
1640 bestsbupcutoff = FALSE;
1641 bestsbfracscore = fracscore;
1642 bestsbdomainscore = domainscore;
1643 lookahead = 0.0;
1644 }
1645 else
1646 lookahead += 0.5;
1647 }
1648 else
1649 lookahead += 1.0;
1650
1651 SCIPdebugMsg(scip, " -> variable <%s> (solval=%g, down=%g (%+g,valid=%u), up=%g (%+g,valid=%u), score=%g/ %g/%g %g/%g %g -> %g)\n",
1652 SCIPvarGetName(branchcands[c]), branchcandssol[c], down, downgain, downvalid, up, upgain, upvalid,
1653 pscostscore, conflictscore, conflengthscore, inferencescore, cutoffscore, gmieffscore, score);
1654 }
1655 }
1656#ifdef SCIP_DEBUG
1657 if( bestsbcand >= 0 )
1658 {
1659 SCIPdebugMsg(scip, " -> best: <%s> (%g / %g / %g), lookahead=%g/%g\n",
1660 SCIPvarGetName(branchcands[bestsbcand]), bestsbscore, bestsbfracscore, bestsbdomainscore,
1661 lookahead, maxlookahead);
1662 }
1663#endif
1664
1665 if( initstrongbranching )
1666 {
1667 if( probingbounds )
1668 {
1669 assert(newlbs != NULL);
1670 assert(newubs != NULL);
1671
1672 SCIPfreeBlockMemoryArray(scip, &newubs, nvars);
1673 SCIPfreeBlockMemoryArray(scip, &newlbs, nvars);
1674 }
1675
1677
1679 {
1680 assert(SCIPhasCurrentNodeLP(scip));
1681 *result = SCIP_CUTOFF;
1682 }
1683 }
1684
1685 if( *result != SCIP_CUTOFF )
1686 {
1687 /* get the score of the best uninitialized strong branching candidate */
1688 if( i < ninitcands && bestuninitsbcand == -1 )
1689 bestuninitsbscore = initcandscores[i];
1690
1691 /* if the best pseudo cost candidate is better than the best uninitialized strong branching candidate,
1692 * compare it to the best initialized strong branching candidate
1693 */
1694 if( bestpsscore > bestuninitsbscore && SCIPisSumGT(scip, bestpsscore, bestsbscore) )
1695 {
1696 bestcand = bestpscand;
1697 bestisstrongbranch = FALSE;
1698 }
1699 else if( bestsbcand >= 0 )
1700 {
1701 bestcand = bestsbcand;
1702 bestisstrongbranch = TRUE;
1703 }
1704 else
1705 {
1706 /* no candidate was initialized, and the best score is the one of the first candidate in the initialization
1707 * queue
1708 */
1709 assert(ninitcands >= 1);
1710 bestcand = initcands[0];
1711 bestisstrongbranch = FALSE;
1712 }
1713
1714 /* update lower bound of current node */
1715 if( allcolsinlp && !exactsolve )
1716 {
1717 assert(SCIPisLT(scip, provedbound, SCIPgetCutoffbound(scip)));
1718 SCIP_CALL( SCIPupdateLocalLowerbound(scip, provedbound) );
1719 assert(SCIPisGE(scip, SCIPgetLocalLowerbound(scip), provedbound));
1720 }
1721 }
1722
1723 /* apply domain reductions */
1724 if( nbdchgs > 0 )
1725 {
1726 if( *result != SCIP_CUTOFF )
1727 {
1728 SCIP_CALL( applyBdchgs(scip, vars, bdchginds, bdchgtypes, bdchgbounds, nbdchgs, result) );
1729 if( *result != SCIP_CUTOFF )
1730 *result = SCIP_REDUCEDDOM;
1731 }
1732 freeBdchgs(scip, &bdchginds, &bdchgtypes, &bdchgbounds, &nbdchgs);
1733 }
1734
1735 /* Apply the Treemodel branching rule to potentially select a better branching candidate than the current one. */
1736 if( *result != SCIP_CUTOFF && *result != SCIP_REDUCEDDOM && *result != SCIP_CONSADDED && SCIPtreemodelIsEnabled(scip, branchruledata->treemodel) )
1737 {
1738 SCIP_Real smallpscost;
1739 SCIP_Bool usetreemodel;
1740
1741 usetreemodel = TRUE;
1742
1743 /* If the pseudocosts are zero, use SCIPs best variable since the Treemodel is not applicable */
1744 if( SCIPisZero(scip, maxgains[bestcand]))
1745 {
1746 usetreemodel = FALSE;
1747 }
1748
1749 /* If SCIPs best candidate was selected due to hybrid branching scores
1750 * rather than because of pseudocosts, then we keep it.
1751 */
1752 SCIP_CALL( SCIPgetRealParam(scip, "branching/treemodel/smallpscost", &smallpscost) );
1753 if( usetreemodel == TRUE && avgpscostscore <= smallpscost )
1754 {
1755 int cand;
1756 for( cand = 0; cand < nbranchcands; ++cand )
1757 {
1758 if( scoresfrompc[cand] > scoresfrompc[bestcand] )
1759 {
1760 usetreemodel = FALSE;
1761 break;
1762 }
1763 }
1764 }
1765
1766 if( usetreemodel == TRUE )
1767 {
1769 scip, /* SCIP data structure */
1770 branchruledata->treemodel, /* branching rule */
1771 branchcands, /* branching candidate storage */
1772 mingains, /* minimum gain of rounding downwards or upwards */
1773 maxgains, /* maximum gain of rounding downwards or upwards */
1774 scoresfromothers, /* scores from other branching methods */
1775 nbranchcands, /* the number of branching candidates */
1776 &bestcand /* the best branching candidate found by SCIP */
1777 ) );
1778 }
1779 }
1780
1781 /* free buffer for the lp gains and pseudocost scores */
1782 SCIPfreeBufferArray(scip, &scoresfromothers);
1783 SCIPfreeBufferArray(scip, &scoresfrompc);
1784 SCIPfreeBufferArray(scip, &scores);
1785 SCIPfreeBufferArray(scip, &maxgains);
1786 SCIPfreeBufferArray(scip, &mingains);
1787
1788 /* free buffer for the unreliable candidates */
1789 SCIPfreeBufferArray(scip, &initcandscores);
1790 SCIPfreeBufferArray(scip, &initcands);
1791 }
1792
1793 /* if no domain could be reduced, create the branching */
1794 if( *result != SCIP_CUTOFF && *result != SCIP_REDUCEDDOM && *result != SCIP_CONSADDED && executebranch )
1795 {
1796 SCIP_NODE* downchild;
1797 SCIP_NODE* upchild;
1798 SCIP_VAR* var;
1799 SCIP_Real val;
1800
1801 assert(*result == SCIP_DIDNOTRUN);
1802 assert(0 <= bestcand && bestcand < nbranchcands);
1803 assert(!SCIPisFeasIntegral(scip, branchcandssol[bestcand]));
1804 assert(!allcolsinlp || SCIPisLT(scip, provedbound, SCIPgetCutoffbound(scip)));
1805 assert(!bestsbdowncutoff && !bestsbupcutoff);
1806
1807 var = branchcands[bestcand];
1808 val = branchcandssol[bestcand];
1809
1810 /* perform the branching */
1811 SCIPdebugMsg(scip, " -> %d (%d) cands, sel cand %d: var <%s> (sol=%g, down=%g (%+g), up=%g (%+g), sb=%u, psc=%g/%g [%g])\n",
1812 nbranchcands, ninitcands, bestcand, SCIPvarGetName(var), branchcandssol[bestcand],
1813 bestsbdown, bestsbdown - lpobjval, bestsbup, bestsbup - lpobjval, bestisstrongbranch,
1816 SCIPgetVarPseudocostScoreCurrentRun(scip, var, branchcandssol[bestcand]));
1817 SCIP_UNUSED(bestisstrongbranch);
1818 SCIP_CALL( SCIPbranchVarVal(scip, var, val, &downchild, NULL, &upchild) );
1819 assert(downchild != NULL);
1820 assert(upchild != NULL);
1821
1822 /* update the lower bounds in the children */
1823 if( allcolsinlp && !exactsolve )
1824 {
1825 if( sbdownvalid[bestcand] )
1826 {
1827 assert(SCIPisLT(scip, sbdown[bestcand], SCIPgetCutoffbound(scip)));
1828 SCIP_CALL( SCIPupdateNodeLowerbound(scip, downchild, sbdown[bestcand]) );
1829 assert(SCIPisGE(scip, SCIPgetNodeLowerbound(scip, downchild), provedbound));
1830 }
1831 if( sbupvalid[bestcand] )
1832 {
1833 assert(SCIPisLT(scip, sbup[bestcand], SCIPgetCutoffbound(scip)));
1834 SCIP_CALL( SCIPupdateNodeLowerbound(scip, upchild, sbup[bestcand]) );
1835 assert(SCIPisGE(scip, SCIPgetNodeLowerbound(scip, upchild), provedbound));
1836 }
1837 }
1838
1839 SCIPdebugMsg(scip, " -> down child's lowerbound: %g\n", SCIPnodeGetLowerbound(downchild));
1840 SCIPdebugMsg(scip, " -> up child's lowerbound : %g\n", SCIPnodeGetLowerbound(upchild));
1841
1843
1844 *result = SCIP_BRANCHED;
1845 }
1846
1847 /* free buffer for the strong branching lp gains */
1848 SCIPfreeBufferArray(scip, &sbupvalid);
1849 SCIPfreeBufferArray(scip, &sbdownvalid);
1850 SCIPfreeBufferArray(scip, &sbup);
1851 SCIPfreeBufferArray(scip, &sbdown);
1852
1853 return SCIP_OKAY;
1854}
1855
1856
1857/*
1858 * Callback methods
1859 */
1860
1861/** copy method for branchrule plugins (called when SCIP copies plugins) */
1862static
1863SCIP_DECL_BRANCHCOPY(branchCopyRelpscost)
1864{ /*lint --e{715}*/
1865 assert(scip != NULL);
1866 assert(branchrule != NULL);
1867 assert(strcmp(SCIPbranchruleGetName(branchrule), BRANCHRULE_NAME) == 0);
1868
1869 /* call inclusion method of branchrule */
1871
1872 return SCIP_OKAY;
1873}
1874
1875/** destructor of branching rule to free user data (called when SCIP is exiting) */
1876static
1877SCIP_DECL_BRANCHFREE(branchFreeRelpscost)
1878{ /*lint --e{715}*/
1879 SCIP_BRANCHRULEDATA* branchruledata;
1880 branchruledata = SCIPbranchruleGetData(branchrule);
1881
1882 /* free Treemodel parameter data structure */
1883 SCIP_CALL( SCIPtreemodelFree(scip, &branchruledata->treemodel) );
1884
1885 /* free branching rule data */
1886 SCIPfreeBlockMemory(scip, &branchruledata);
1887 SCIPbranchruleSetData(branchrule, NULL);
1888
1889 return SCIP_OKAY;
1890}
1891
1892
1893/** solving process initialization method of branching rule (called when branch and bound process is about to begin) */
1894static
1895SCIP_DECL_BRANCHINITSOL(branchInitsolRelpscost)
1896{ /*lint --e{715}*/
1897 SCIP_BRANCHRULEDATA* branchruledata;
1898
1899 /* initialize branching rule data */
1900 branchruledata = SCIPbranchruleGetData(branchrule);
1901 branchruledata->nlcount = NULL;
1902 branchruledata->nlcountsize = 0;
1903 branchruledata->nlcountmax = 1;
1904 assert(branchruledata->startrandseed >= 0);
1905
1906 /* create a random number generator */
1907 SCIP_CALL( SCIPcreateRandom(scip, &branchruledata->randnumgen,
1908 (unsigned int)branchruledata->startrandseed, TRUE) );
1909
1910 return SCIP_OKAY;
1911}
1912
1913
1914/** solving process deinitialization method of branching rule (called before branch and bound process data is freed) */
1915static
1916SCIP_DECL_BRANCHEXITSOL(branchExitsolRelpscost)
1917{ /*lint --e{715}*/
1918 SCIP_BRANCHRULEDATA* branchruledata;
1919
1920 /* free memory in branching rule data */
1921 branchruledata = SCIPbranchruleGetData(branchrule);
1922 SCIPfreeBlockMemoryArrayNull(scip, &branchruledata->nlcount, branchruledata->nlcountsize);
1923
1924 /* free random number generator */
1925 SCIPfreeRandom(scip, &branchruledata->randnumgen);
1926
1927 SCIPfreeBlockMemoryArrayNull(scip, &branchruledata->orbitrep, branchruledata->npermvars);
1928 SCIPfreeBlockMemoryArrayNull(scip, &branchruledata->varorbitmap, branchruledata->npermvars);
1929 SCIPfreeBlockMemoryArrayNull(scip, &branchruledata->orbitbegins, branchruledata->npermvars);
1930 SCIPfreeBlockMemoryArrayNull(scip, &branchruledata->orbits, branchruledata->npermvars);
1931 branchruledata->nosymmetry = FALSE;
1932 branchruledata->norbits = 0;
1933 branchruledata->permvars = NULL;
1934 branchruledata->permvarmap = NULL;
1935 branchruledata->npermvars = 0;
1936
1937 return SCIP_OKAY;
1938}
1939
1940
1941/** branching execution method for fractional LP solutions */
1942static
1943SCIP_DECL_BRANCHEXECLP(branchExeclpRelpscost)
1944{ /*lint --e{715}*/
1945 SCIP_BRANCHRULEDATA* branchruledata;
1946 SCIP_VAR** lpcands;
1947 SCIP_Real* lpcandssol;
1948 SCIP_Real* lpcandsfrac;
1949 SCIP_VAR** filteredlpcands;
1950 SCIP_Real* filteredlpcandssol;
1951 SCIP_Real* filteredlpcandsfrac;
1952 SCIP_Bool runfiltering;
1953 int* filteredlpcandsorbitidx = NULL;
1954 int nfilteredlpcands;
1955 int nlpcands;
1956
1957 assert(branchrule != NULL);
1958 assert(strcmp(SCIPbranchruleGetName(branchrule), BRANCHRULE_NAME) == 0);
1959 assert(scip != NULL);
1960 assert(result != NULL);
1961
1962 SCIPdebugMsg(scip, "Execlp method of relpscost branching in node %" SCIP_LONGINT_FORMAT "\n", SCIPnodeGetNumber(SCIPgetCurrentNode(scip)));
1963
1965 {
1966 *result = SCIP_DIDNOTRUN;
1967 SCIPdebugMsg(scip, "Could not apply relpscost branching, as the current LP was not solved to optimality.\n");
1968
1969 return SCIP_OKAY;
1970 }
1971
1972 /* get branching candidates */
1973 SCIP_CALL( SCIPgetLPBranchCands(scip, &lpcands, &lpcandssol, &lpcandsfrac, NULL, &nlpcands, NULL) );
1974 assert(nlpcands > 0);
1975
1976 branchruledata = SCIPbranchruleGetData(branchrule);
1977 assert(branchruledata != NULL);
1978
1979 /* determine whether we should run filtering */
1980 runfiltering = ! branchruledata->nosymmetry && branchruledata->filtercandssym && SCIPgetSubscipDepth(scip) == 0 && ! SCIPinDive(scip) && ! SCIPinProbing(scip);
1981
1982 /* init orbits if necessary */
1983 if( runfiltering )
1984 {
1985 SCIP_CALL( initOrbits(scip, branchruledata) );
1986 }
1987
1988 /* determine fractional variables (possibly filter by using symmetries) */
1989 if( runfiltering && branchruledata->norbits != 0 )
1990 {
1991 SCIP_CALL( SCIPallocBufferArray(scip, &filteredlpcands, nlpcands) );
1992 SCIP_CALL( SCIPallocBufferArray(scip, &filteredlpcandssol, nlpcands) );
1993 SCIP_CALL( SCIPallocBufferArray(scip, &filteredlpcandsfrac, nlpcands) );
1994 SCIP_CALL( SCIPallocBufferArray(scip, &filteredlpcandsorbitidx, nlpcands) );
1995
1996 /* determine filtered fractional variables */
1997 SCIP_CALL( filterSymmetricVariables(scip, branchruledata, lpcands, lpcandssol, lpcandsfrac, nlpcands,
1998 filteredlpcands, filteredlpcandssol, filteredlpcandsfrac, filteredlpcandsorbitidx, &nfilteredlpcands) );
1999 }
2000 else
2001 {
2002 /* No orbits available. Copy all (unfiltered) branching candidates, because they will be updated w.r.t. the strong branching LP solution */
2003 SCIP_CALL( SCIPduplicateBufferArray(scip, &filteredlpcands, lpcands, nlpcands) );
2004 SCIP_CALL( SCIPduplicateBufferArray(scip, &filteredlpcandssol, lpcandssol, nlpcands) );
2005 SCIP_CALL( SCIPduplicateBufferArray(scip, &filteredlpcandsfrac, lpcandsfrac, nlpcands) );
2006 nfilteredlpcands = nlpcands;
2007 }
2008
2009 /* execute branching rule */
2010 SCIP_CALL( execRelpscost(scip, branchrule, filteredlpcands, filteredlpcandssol, filteredlpcandsfrac, filteredlpcandsorbitidx, nfilteredlpcands, TRUE, result) );
2011
2012 SCIPfreeBufferArrayNull(scip, &filteredlpcandsorbitidx);
2013 SCIPfreeBufferArray(scip, &filteredlpcandsfrac);
2014 SCIPfreeBufferArray(scip, &filteredlpcandssol);
2015 SCIPfreeBufferArray(scip, &filteredlpcands);
2016
2017 return SCIP_OKAY;
2018}
2019
2020
2021/*
2022 * branching specific interface methods
2023 */
2024
2025/** creates the reliable pseudo cost branching rule and includes it in SCIP */
2027 SCIP* scip /**< SCIP data structure */
2028 )
2029{
2030 SCIP_BRANCHRULEDATA* branchruledata;
2031 SCIP_BRANCHRULE* branchrule;
2032
2033 /* create relpscost branching rule data */
2034 SCIP_CALL( SCIPallocBlockMemory(scip, &branchruledata) );
2035
2036 branchruledata->nosymmetry = FALSE;
2037 branchruledata->orbits = NULL;
2038 branchruledata->orbitbegins = NULL;
2039 branchruledata->orbitrep = NULL;
2040 branchruledata->varorbitmap = NULL;
2041 branchruledata->norbits = 0;
2042 branchruledata->permvars = NULL;
2043 branchruledata->npermvars = 0;
2044 branchruledata->permvarmap = NULL;
2045
2046 /* include branching rule */
2048 BRANCHRULE_MAXDEPTH, BRANCHRULE_MAXBOUNDDIST, branchruledata) );
2049
2050 assert(branchrule != NULL);
2051
2052 /* set non-fundamental callbacks via specific setter functions*/
2053 SCIP_CALL( SCIPsetBranchruleCopy(scip, branchrule, branchCopyRelpscost) );
2054 SCIP_CALL( SCIPsetBranchruleFree(scip, branchrule, branchFreeRelpscost) );
2055 SCIP_CALL( SCIPsetBranchruleInitsol(scip, branchrule, branchInitsolRelpscost) );
2056 SCIP_CALL( SCIPsetBranchruleExitsol(scip, branchrule, branchExitsolRelpscost) );
2057 SCIP_CALL( SCIPsetBranchruleExecLp(scip, branchrule, branchExeclpRelpscost) );
2058
2059 /* relpscost branching rule parameters */
2061 "branching/relpscost/conflictweight",
2062 "weight in score calculations for conflict score",
2063 &branchruledata->conflictweight, TRUE, DEFAULT_CONFLICTWEIGHT, SCIP_REAL_MIN, SCIP_REAL_MAX, NULL, NULL) );
2065 "branching/relpscost/conflictlengthweight",
2066 "weight in score calculations for conflict length score",
2067 &branchruledata->conflengthweight, TRUE, DEFAULT_CONFLENGTHWEIGHT, SCIP_REAL_MIN, SCIP_REAL_MAX, NULL, NULL) );
2069 "branching/relpscost/inferenceweight",
2070 "weight in score calculations for inference score",
2071 &branchruledata->inferenceweight, TRUE, DEFAULT_INFERENCEWEIGHT, SCIP_REAL_MIN, SCIP_REAL_MAX, NULL, NULL) );
2073 "branching/relpscost/cutoffweight",
2074 "weight in score calculations for cutoff score",
2075 &branchruledata->cutoffweight, TRUE, DEFAULT_CUTOFFWEIGHT, SCIP_REAL_MIN, SCIP_REAL_MAX, NULL, NULL) );
2077 "branching/relpscost/gmiavgeffweight",
2078 "weight in score calculations for average GMI cuts normalized efficacy",
2079 &branchruledata->gmiavgeffweight, TRUE, DEFAULT_GMIAVGEFFWEIGHT, SCIP_REAL_MIN, SCIP_REAL_MAX, NULL, NULL) );
2081 "branching/relpscost/gmilasteffweight",
2082 "weight in score calculations for last GMI cuts normalized efficacy",
2083 &branchruledata->gmilasteffweight, TRUE, DEFAULT_GMILASTEFFWEIGHT, SCIP_REAL_MIN, SCIP_REAL_MAX, NULL, NULL) );
2085 "branching/relpscost/pscostweight",
2086 "weight in score calculations for pseudo cost score",
2087 &branchruledata->pscostweight, TRUE, DEFAULT_PSCOSTWEIGHT, SCIP_REAL_MIN, SCIP_REAL_MAX, NULL, NULL) );
2089 "branching/relpscost/nlscoreweight",
2090 "weight in score calculations for nlcount score",
2091 &branchruledata->nlscoreweight, TRUE, DEFAULT_NLSCOREWEIGHT, SCIP_REAL_MIN, SCIP_REAL_MAX, NULL, NULL) );
2093 "branching/relpscost/minreliable",
2094 "minimal value for minimum pseudo cost size to regard pseudo cost value as reliable",
2095 &branchruledata->minreliable, TRUE, DEFAULT_MINRELIABLE, 0.0, SCIP_REAL_MAX, NULL, NULL) );
2097 "branching/relpscost/maxreliable",
2098 "maximal value for minimum pseudo cost size to regard pseudo cost value as reliable",
2099 &branchruledata->maxreliable, TRUE, DEFAULT_MAXRELIABLE, 0.0, SCIP_REAL_MAX, NULL, NULL) );
2101 "branching/relpscost/sbiterquot",
2102 "maximal fraction of strong branching LP iterations compared to node relaxation LP iterations",
2103 &branchruledata->sbiterquot, FALSE, DEFAULT_SBITERQUOT, 0.0, SCIP_REAL_MAX, NULL, NULL) );
2105 "branching/relpscost/sbiterofs",
2106 "additional number of allowed strong branching LP iterations",
2107 &branchruledata->sbiterofs, FALSE, DEFAULT_SBITEROFS, 0, INT_MAX, NULL, NULL) );
2109 "branching/relpscost/maxlookahead",
2110 "maximal number of further variables evaluated without better score",
2111 &branchruledata->maxlookahead, TRUE, DEFAULT_MAXLOOKAHEAD, 1, INT_MAX, NULL, NULL) );
2113 "branching/relpscost/initcand",
2114 "maximal number of candidates initialized with strong branching per node",
2115 &branchruledata->initcand, FALSE, DEFAULT_INITCAND, 0, INT_MAX, NULL, NULL) );
2117 "branching/relpscost/inititer",
2118 "iteration limit for strong branching initializations of pseudo cost entries (0: auto)",
2119 &branchruledata->inititer, FALSE, DEFAULT_INITITER, 0, INT_MAX, NULL, NULL) );
2121 "branching/relpscost/maxbdchgs",
2122 "maximal number of bound tightenings before the node is reevaluated (-1: unlimited)",
2123 &branchruledata->maxbdchgs, TRUE, DEFAULT_MAXBDCHGS, -1, INT_MAX, NULL, NULL) );
2125 "branching/relpscost/maxproprounds",
2126 "maximum number of propagation rounds to be performed during strong branching before solving the LP (-1: no limit, -2: parameter settings)",
2127 &branchruledata->maxproprounds, TRUE, DEFAULT_MAXPROPROUNDS, -2, INT_MAX, NULL, NULL) );
2129 "branching/relpscost/probingbounds",
2130 "should valid bounds be identified in a probing-like fashion during strong branching (only with propagation)?",
2131 &branchruledata->probingbounds, TRUE, DEFAULT_PROBINGBOUNDS, NULL, NULL) );
2132
2133 SCIP_CALL( SCIPaddBoolParam(scip, "branching/relpscost/userelerrorreliability",
2134 "should reliability be based on relative errors?", &branchruledata->userelerrorforreliability, TRUE, DEFAULT_USERELERRORFORRELIABILITY,
2135 NULL, NULL) );
2136
2137 SCIP_CALL( SCIPaddRealParam(scip, "branching/relpscost/lowerrortol", "low relative error tolerance for reliability",
2138 &branchruledata->lowerrortol, TRUE, DEFAULT_LOWERRORTOL, 0.0, SCIP_REAL_MAX, NULL, NULL) );
2139
2140 SCIP_CALL( SCIPaddRealParam(scip, "branching/relpscost/higherrortol", "high relative error tolerance for reliability",
2141 &branchruledata->higherrortol, TRUE, DEFAULT_HIGHERRORTOL, 0.0, SCIP_REAL_MAX, NULL, NULL) );
2142
2143/**! [SnippetCodeStyleParenIndent] */
2144
2145 SCIP_CALL( SCIPaddBoolParam(scip, "branching/relpscost/storesemiinitcosts",
2146 "should strong branching result be considered for pseudo costs if the other direction was infeasible?",
2147 &branchruledata->storesemiinitcosts, TRUE, DEFAULT_STORESEMIINITCOSTS,
2148 NULL, NULL) );
2149
2150/**! [SnippetCodeStyleParenIndent] */
2151
2152 SCIP_CALL( SCIPaddBoolParam(scip, "branching/relpscost/usesblocalinfo",
2153 "should the scoring function use only local cutoff and inference information obtained for strong branching candidates?",
2154 &branchruledata->usesblocalinfo, TRUE, DEFAULT_USESBLOCALINFO,
2155 NULL, NULL) );
2156
2157 SCIP_CALL( SCIPaddBoolParam(scip, "branching/relpscost/usehyptestforreliability",
2158 "should the strong branching decision be based on a hypothesis test?",
2159 &branchruledata->usehyptestforreliability, TRUE, DEFAULT_USEHYPTESTFORRELIABILITY,
2160 NULL, NULL) );
2161
2162 SCIP_CALL( SCIPaddBoolParam(scip, "branching/relpscost/usedynamicconfidence",
2163 "should the confidence level be adjusted dynamically?",
2164 &branchruledata->usedynamicconfidence, TRUE, DEFAULT_USEDYNAMICCONFIDENCE,
2165 NULL, NULL) );
2166 SCIP_CALL( SCIPaddBoolParam(scip, "branching/relpscost/skipbadinitcands",
2167 "should branching rule skip candidates that have a low probability to "
2168 "be better than the best strong-branching or pseudo-candidate?",
2169 &branchruledata->skipbadinitcands, TRUE, DEFAULT_SKIPBADINITCANDS,
2170 NULL, NULL) );
2172 "branching/relpscost/confidencelevel",
2173 "the confidence level for statistical methods, between 0 (Min) and 4 (Max).",
2174 &branchruledata->confidencelevel, TRUE, DEFAULT_CONFIDENCELEVEL, 0, 4, NULL, NULL) );
2175
2176 SCIP_CALL( SCIPaddBoolParam(scip, "branching/relpscost/randinitorder",
2177 "should candidates be initialized in randomized order?",
2178 &branchruledata->randinitorder, TRUE, DEFAULT_RANDINITORDER,
2179 NULL, NULL) );
2180
2181 SCIP_CALL( SCIPaddBoolParam(scip, "branching/relpscost/usesmallweightsitlim",
2182 "should smaller weights be used for pseudo cost updates after hitting the LP iteration limit?",
2183 &branchruledata->usesmallweightsitlim, TRUE, DEFAULT_USESMALLWEIGHTSITLIM,
2184 NULL, NULL) );
2185
2186 SCIP_CALL( SCIPaddBoolParam(scip, "branching/relpscost/dynamicweights",
2187 "should the weights of the branching rule be adjusted dynamically during solving based on objective and infeasible leaf counters?",
2188 &branchruledata->dynamicweights, TRUE, DEFAULT_DYNAMICWEIGHTS,
2189 NULL, NULL) );
2190 SCIP_CALL( SCIPaddIntParam(scip, "branching/relpscost/degeneracyaware",
2191 "should degeneracy be taken into account to update weights and skip strong branching? (0: off, 1: after root, 2: always)",
2192 &branchruledata->degeneracyaware, TRUE, DEFAULT_DEGENERACYAWARE, 0, 2,
2193 NULL, NULL) );
2194 SCIP_CALL( SCIPaddIntParam(scip, "branching/relpscost/startrandseed", "start seed for random number generation",
2195 &branchruledata->startrandseed, TRUE, DEFAULT_STARTRANDSEED, 0, INT_MAX, NULL, NULL) );
2196
2197 SCIP_CALL( SCIPaddBoolParam(scip, "branching/relpscost/filtercandssym",
2198 "Use symmetry to filter branching candidates?",
2199 &branchruledata->filtercandssym, TRUE, DEFAULT_FILTERCANDSSYM, NULL, NULL) );
2200
2201 SCIP_CALL( SCIPaddBoolParam(scip, "branching/relpscost/transsympscost",
2202 "Transfer pscost information to symmetric variables?",
2203 &branchruledata->transsympscost, TRUE, DEFAULT_TRANSSYMPSCOST, NULL, NULL) );
2204
2205 /* initialise the Treemodel parameters */
2206 SCIP_CALL( SCIPtreemodelInit(scip, &branchruledata->treemodel) );
2207
2208 return SCIP_OKAY;
2209}
2210
2211/** execution reliability pseudo cost branching with the given branching candidates */
2213 SCIP* scip, /**< SCIP data structure */
2214 SCIP_VAR** branchcands, /**< branching candidates */
2215 SCIP_Real* branchcandssol, /**< solution value for the branching candidates */
2216 SCIP_Real* branchcandsfrac, /**< fractional part of the branching candidates */
2217 int nbranchcands, /**< number of branching candidates */
2218 SCIP_Bool executebranching, /**< perform a branching step after probing */
2219 SCIP_RESULT* result /**< pointer to the result of the execution */
2220 )
2221{
2222 SCIP_BRANCHRULE* branchrule;
2223
2224 assert(scip != NULL);
2225 assert(result != NULL);
2226
2227 /* find branching rule */
2229 assert(branchrule != NULL);
2230
2231 /* execute branching rule */
2232 SCIP_CALL( execRelpscost(scip, branchrule, branchcands, branchcandssol, branchcandsfrac, NULL, nbranchcands, executebranching, result) );
2233
2234 return SCIP_OKAY;
2235}
static long bound
#define BRANCHRULE_DESC
#define DEFAULT_DEGENERACYAWARE
#define DEFAULT_USESBLOCALINFO
#define DEFAULT_SKIPBADINITCANDS
#define DEFAULT_SBITERQUOT
static SCIP_Real calcNlscore(SCIP *scip, int *nlcount, int nlcountmax, int probindex)
#define DEFAULT_HIGHERRORTOL
#define BRANCHRULE_PRIORITY
#define DEFAULT_PROBINGBOUNDS
#define DEFAULT_INITITER
#define DEFAULT_USESMALLWEIGHTSITLIM
static SCIP_RETCODE applyBdchgs(SCIP *scip, SCIP_VAR **vars, int *bdchginds, SCIP_BOUNDTYPE *bdchgtypes, SCIP_Real *bdchgbounds, int nbdchgs, SCIP_RESULT *result)
#define DEFAULT_STARTRANDSEED
#define DEFAULT_FILTERCANDSSYM
#define DEFAULT_USEDYNAMICCONFIDENCE
#define DEFAULT_DYNAMICWEIGHTS
#define DEFAULT_MAXBDCHGS
#define DEFAULT_USEHYPTESTFORRELIABILITY
#define DEFAULT_GMIAVGEFFWEIGHT
static SCIP_RETCODE SCIPupdateVarPseudocostSymmetric(SCIP *scip, SCIP_BRANCHRULEDATA *branchruledata, SCIP_VAR *branchvar, int *branchorbitidx, int branchvaridx, SCIP_Real solvaldelta, SCIP_Real objdelta, SCIP_Real weight)
#define DEFAULT_INFERENCEWEIGHT
#define BRANCHRULE_NAME
#define DEFAULT_RANDINITORDER
static SCIP_RETCODE initOrbits(SCIP *scip, SCIP_BRANCHRULEDATA *branchruledata)
#define DEFAULT_MAXLOOKAHEAD
#define DEFAULT_SBITEROFS
static SCIP_DECL_BRANCHFREE(branchFreeRelpscost)
#define DEFAULT_NLSCOREWEIGHT
#define DEFAULT_CONFLICTWEIGHT
static SCIP_RETCODE filterSymmetricVariables(SCIP *scip, SCIP_BRANCHRULEDATA *branchruledata, SCIP_VAR **origbranchcands, SCIP_Real *origbranchcandssol, SCIP_Real *origbranchcandsfrac, int norigbranchcands, SCIP_VAR **branchcands, SCIP_Real *branchcandssol, SCIP_Real *branchcandsfrac, int *branchorbitidx, int *nbranchcands)
#define DEFAULT_GMILASTEFFWEIGHT
#define DEFAULT_USERELERRORFORRELIABILITY
#define DEFAULT_CUTOFFWEIGHT
static SCIP_DECL_BRANCHCOPY(branchCopyRelpscost)
static SCIP_RETCODE countNonlinearities(SCIP *scip, int *nlcount, int nlcountsize, int *nlcountmax)
static SCIP_DECL_BRANCHINITSOL(branchInitsolRelpscost)
#define DEFAULT_TRANSSYMPSCOST
static SCIP_RETCODE addBdchg(SCIP *scip, int **bdchginds, SCIP_BOUNDTYPE **bdchgtypes, SCIP_Real **bdchgbounds, int *nbdchgs, int ind, SCIP_BOUNDTYPE type, SCIP_Real bound)
#define DEFAULT_LOWERRORTOL
static SCIP_RETCODE execRelpscost(SCIP *scip, SCIP_BRANCHRULE *branchrule, SCIP_VAR **branchcands, SCIP_Real *branchcandssol, SCIP_Real *branchcandsfrac, int *branchorbitidx, int nbranchcands, SCIP_Bool executebranch, SCIP_RESULT *result)
static SCIP_RETCODE branchruledataEnsureNlcount(SCIP *scip, SCIP_BRANCHRULEDATA *branchruledata)
#define DEFAULT_CONFLENGTHWEIGHT
#define DEFAULT_CONFIDENCELEVEL
#define DEFAULT_MINRELIABLE
static SCIP_Real calcScore(SCIP *scip, SCIP_BRANCHRULEDATA *branchruledata, SCIP_Real conflictscore, SCIP_Real avgconflictscore, SCIP_Real conflengthscore, SCIP_Real avgconflengthscore, SCIP_Real inferencescore, SCIP_Real avginferencescore, SCIP_Real cutoffscore, SCIP_Real avgcutoffscore, SCIP_Real gmieffscore, SCIP_Real lastgmieffscore, SCIP_Real pscostscore, SCIP_Real avgpscostscore, SCIP_Real nlscore, SCIP_Real frac, SCIP_Real degeneracyfactor)
#define DEFAULT_STORESEMIINITCOSTS
#define DEFAULT_MAXPROPROUNDS
static SCIP_DECL_BRANCHEXITSOL(branchExitsolRelpscost)
#define DEFAULT_INITCAND
#define BRANCHRULE_MAXDEPTH
#define DEFAULT_MAXRELIABLE
#define DEFAULT_PSCOSTWEIGHT
#define BRANCHRULE_MAXBOUNDDIST
static void freeBdchgs(SCIP *scip, int **bdchginds, SCIP_BOUNDTYPE **bdchgtypes, SCIP_Real **bdchgbounds, int *nbdchgs)
static SCIP_DECL_BRANCHEXECLP(branchExeclpRelpscost)
reliable pseudo costs branching rule
Constraint handler for AND constraints, .
#define NULL
Definition: def.h:266
#define SCIP_Longint
Definition: def.h:157
#define SCIP_UNUSED(x)
Definition: def.h:427
#define SCIP_REAL_MAX
Definition: def.h:173
#define SCIP_INVALID
Definition: def.h:192
#define SCIP_Bool
Definition: def.h:91
#define MIN(x, y)
Definition: def.h:242
#define SCIP_Real
Definition: def.h:172
#define TRUE
Definition: def.h:93
#define FALSE
Definition: def.h:94
#define MAX(x, y)
Definition: def.h:238
#define SCIP_LONGINT_FORMAT
Definition: def.h:164
#define SCIP_REAL_MIN
Definition: def.h:174
#define SCIP_CALL(x)
Definition: def.h:373
SCIP_RETCODE SCIPexecRelpscostBranching(SCIP *scip, SCIP_VAR **branchcands, SCIP_Real *branchcandssol, SCIP_Real *branchcandsfrac, int nbranchcands, SCIP_Bool executebranching, SCIP_RESULT *result)
SCIP_RETCODE SCIPincludeBranchruleRelpscost(SCIP *scip)
SCIP_VAR * SCIPgetResultantAnd(SCIP *scip, SCIP_CONS *cons)
Definition: cons_and.c:5260
int SCIPgetNVarsAnd(SCIP *scip, SCIP_CONS *cons)
Definition: cons_and.c:5211
SCIP_VAR ** SCIPgetVarsAnd(SCIP *scip, SCIP_CONS *cons)
Definition: cons_and.c:5235
int SCIPgetSubscipDepth(SCIP *scip)
Definition: scip_copy.c:2605
SCIP_Bool SCIPisStopped(SCIP *scip)
Definition: scip_general.c:734
SCIP_Bool SCIPisExactSolve(SCIP *scip)
Definition: scip_general.c:621
SCIP_RETCODE SCIPgetVarsData(SCIP *scip, SCIP_VAR ***vars, int *nvars, int *nbinvars, int *nintvars, int *nimplvars, int *ncontvars)
Definition: scip_prob.c:1866
int SCIPgetNVars(SCIP *scip)
Definition: scip_prob.c:1992
SCIP_VAR ** SCIPgetVars(SCIP *scip)
Definition: scip_prob.c:1947
int SCIPhashmapGetImageInt(SCIP_HASHMAP *hashmap, void *origin)
Definition: misc.c:3284
SCIP_RETCODE SCIPupdateNodeLowerbound(SCIP *scip, SCIP_NODE *node, SCIP_Real newbound)
Definition: scip_prob.c:3762
SCIP_RETCODE SCIPupdateLocalLowerbound(SCIP *scip, SCIP_Real newbound)
Definition: scip_prob.c:3697
SCIP_Real SCIPgetLocalLowerbound(SCIP *scip)
Definition: scip_prob.c:3586
SCIP_Real SCIPgetNodeLowerbound(SCIP *scip, SCIP_NODE *node)
Definition: scip_prob.c:3623
void SCIPverbMessage(SCIP *scip, SCIP_VERBLEVEL msgverblevel, FILE *file, const char *formatstr,...)
Definition: scip_message.c:225
#define SCIPdebugMsg
Definition: scip_message.h:78
SCIP_RETCODE SCIPaddIntParam(SCIP *scip, const char *name, const char *desc, int *valueptr, SCIP_Bool isadvanced, int defaultvalue, int minvalue, int maxvalue, SCIP_DECL_PARAMCHGD((*paramchgd)), SCIP_PARAMDATA *paramdata)
Definition: scip_param.c:83
SCIP_RETCODE SCIPaddRealParam(SCIP *scip, const char *name, const char *desc, SCIP_Real *valueptr, SCIP_Bool isadvanced, SCIP_Real defaultvalue, SCIP_Real minvalue, SCIP_Real maxvalue, SCIP_DECL_PARAMCHGD((*paramchgd)), SCIP_PARAMDATA *paramdata)
Definition: scip_param.c:139
SCIP_RETCODE SCIPgetRealParam(SCIP *scip, const char *name, SCIP_Real *value)
Definition: scip_param.c:307
SCIP_RETCODE SCIPaddBoolParam(SCIP *scip, const char *name, const char *desc, SCIP_Bool *valueptr, SCIP_Bool isadvanced, SCIP_Bool defaultvalue, SCIP_DECL_PARAMCHGD((*paramchgd)), SCIP_PARAMDATA *paramdata)
Definition: scip_param.c:57
SCIP_RETCODE SCIPsetBranchruleExecLp(SCIP *scip, SCIP_BRANCHRULE *branchrule, SCIP_DECL_BRANCHEXECLP((*branchexeclp)))
Definition: scip_branch.c:256
SCIP_BRANCHRULE * SCIPfindBranchrule(SCIP *scip, const char *name)
Definition: scip_branch.c:304
SCIP_RETCODE SCIPsetBranchruleCopy(SCIP *scip, SCIP_BRANCHRULE *branchrule, SCIP_DECL_BRANCHCOPY((*branchcopy)))
Definition: scip_branch.c:160
SCIP_RETCODE SCIPincludeBranchruleBasic(SCIP *scip, SCIP_BRANCHRULE **branchruleptr, const char *name, const char *desc, int priority, int maxdepth, SCIP_Real maxbounddist, SCIP_BRANCHRULEDATA *branchruledata)
Definition: scip_branch.c:123
const char * SCIPbranchruleGetName(SCIP_BRANCHRULE *branchrule)
Definition: branch.c:1971
SCIP_BRANCHRULEDATA * SCIPbranchruleGetData(SCIP_BRANCHRULE *branchrule)
Definition: branch.c:1849
SCIP_RETCODE SCIPsetBranchruleFree(SCIP *scip, SCIP_BRANCHRULE *branchrule, SCIP_DECL_BRANCHFREE((*branchfree)))
Definition: scip_branch.c:176
SCIP_RETCODE SCIPsetBranchruleExitsol(SCIP *scip, SCIP_BRANCHRULE *branchrule, SCIP_DECL_BRANCHEXITSOL((*branchexitsol)))
Definition: scip_branch.c:240
void SCIPbranchruleSetData(SCIP_BRANCHRULE *branchrule, SCIP_BRANCHRULEDATA *branchruledata)
Definition: branch.c:1859
SCIP_RETCODE SCIPsetBranchruleInitsol(SCIP *scip, SCIP_BRANCHRULE *branchrule, SCIP_DECL_BRANCHINITSOL((*branchinitsol)))
Definition: scip_branch.c:224
SCIP_RETCODE SCIPbranchVarVal(SCIP *scip, SCIP_VAR *var, SCIP_Real val, SCIP_NODE **downchild, SCIP_NODE **eqchild, SCIP_NODE **upchild)
Definition: scip_branch.c:1133
SCIP_RETCODE SCIPgetLPBranchCands(SCIP *scip, SCIP_VAR ***lpcands, SCIP_Real **lpcandssol, SCIP_Real **lpcandsfrac, int *nlpcands, int *npriolpcands, int *nfracimplvars)
Definition: scip_branch.c:402
SCIP_Real SCIPgetBranchScore(SCIP *scip, SCIP_VAR *var, SCIP_Real downgain, SCIP_Real upgain)
Definition: scip_branch.c:856
SCIP_CONSHDLR * SCIPfindConshdlr(SCIP *scip, const char *name)
Definition: scip_cons.c:941
int SCIPconshdlrGetNActiveConss(SCIP_CONSHDLR *conshdlr)
Definition: cons.c:4670
SCIP_CONS ** SCIPconshdlrGetConss(SCIP_CONSHDLR *conshdlr)
Definition: cons.c:4593
SCIP_RETCODE SCIPgetLPDualDegeneracy(SCIP *scip, SCIP_Real *degeneracy, SCIP_Real *varconsratio)
Definition: scip_lp.c:2792
SCIP_Bool SCIPinDive(SCIP *scip)
Definition: scip_lp.c:2775
SCIP_Bool SCIPhasCurrentNodeLP(SCIP *scip)
Definition: scip_lp.c:83
SCIP_LPSOLSTAT SCIPgetLPSolstat(SCIP *scip)
Definition: scip_lp.c:168
SCIP_Bool SCIPallColsInLP(SCIP *scip)
Definition: scip_lp.c:649
SCIP_Real SCIPgetLPObjval(SCIP *scip)
Definition: scip_lp.c:247
SCIP_Bool SCIPisLPSolBasic(SCIP *scip)
Definition: scip_lp.c:667
#define SCIPfreeBlockMemoryArray(scip, ptr, num)
Definition: scip_mem.h:110
#define SCIPallocBufferArray(scip, ptr, num)
Definition: scip_mem.h:124
#define SCIPreallocBufferArray(scip, ptr, num)
Definition: scip_mem.h:128
#define SCIPfreeBufferArray(scip, ptr)
Definition: scip_mem.h:136
#define SCIPduplicateBufferArray(scip, ptr, source, num)
Definition: scip_mem.h:132
#define SCIPallocBlockMemoryArray(scip, ptr, num)
Definition: scip_mem.h:93
#define SCIPreallocBlockMemoryArray(scip, ptr, oldnum, newnum)
Definition: scip_mem.h:99
#define SCIPfreeBlockMemory(scip, ptr)
Definition: scip_mem.h:108
#define SCIPfreeBlockMemoryArrayNull(scip, ptr, num)
Definition: scip_mem.h:111
#define SCIPfreeBufferArrayNull(scip, ptr)
Definition: scip_mem.h:137
#define SCIPallocBlockMemory(scip, ptr)
Definition: scip_mem.h:89
SCIP_Bool SCIPisNLPConstructed(SCIP *scip)
Definition: scip_nlp.c:110
int SCIPgetNNLPVars(SCIP *scip)
Definition: scip_nlp.c:201
SCIP_RETCODE SCIPgetNLPVarsNonlinearity(SCIP *scip, int *nlcount)
Definition: scip_nlp.c:223
SCIP_Real SCIPnodeGetLowerbound(SCIP_NODE *node)
Definition: tree.c:7531
SCIP_Longint SCIPnodeGetNumber(SCIP_NODE *node)
Definition: tree.c:7511
SCIP_Bool SCIPinProbing(SCIP *scip)
Definition: scip_probing.c:97
SCIP_SOL * SCIPgetBestSol(SCIP *scip)
Definition: scip_sol.c:2165
int SCIPsolGetIndex(SCIP_SOL *sol)
Definition: sol.c:2835
SCIP_Real SCIPgetAvgInferenceScore(SCIP *scip)
SCIP_Real SCIPgetAvgConflictlengthScore(SCIP *scip)
SCIP_Longint SCIPgetNInfeasibleLeaves(SCIP *scip)
SCIP_Longint SCIPgetNNodes(SCIP *scip)
SCIP_Longint SCIPgetNNodeInitLPs(SCIP *scip)
SCIP_Longint SCIPgetNStrongbranchLPIterations(SCIP *scip)
SCIP_Longint SCIPgetNNodeLPIterations(SCIP *scip)
SCIP_Real SCIPgetAvgConflictScore(SCIP *scip)
SCIP_Longint SCIPgetNDualResolveLPIterations(SCIP *scip)
SCIP_Real SCIPgetAvgPseudocostScore(SCIP *scip)
SCIP_Real SCIPgetCutoffbound(SCIP *scip)
SCIP_Longint SCIPgetNObjlimLeaves(SCIP *scip)
SCIP_Longint SCIPgetNNodeInitLPIterations(SCIP *scip)
SCIP_Longint SCIPgetNDualResolveLPs(SCIP *scip)
SCIP_Longint SCIPgetNRootStrongbranchLPIterations(SCIP *scip)
SCIP_Real SCIPgetAvgCutoffScore(SCIP *scip)
SCIP_RETCODE SCIPcomputeOrbitsComponentsSym(SCIP *scip, int npermvars, int **permstrans, int nperms, int *components, int *componentbegins, int *vartocomponent, int ncomponents, int *orbits, int *orbitbegins, int *norbits, int *varorbitmap)
Definition: symmetry.c:420
SCIP_Real SCIPinfinity(SCIP *scip)
SCIP_Bool SCIPisGE(SCIP *scip, SCIP_Real val1, SCIP_Real val2)
SCIP_Real SCIPfeasFrac(SCIP *scip, SCIP_Real val)
SCIP_Real SCIPfeasCeil(SCIP *scip, SCIP_Real val)
SCIP_Real SCIPfeasFloor(SCIP *scip, SCIP_Real val)
SCIP_Bool SCIPisFeasIntegral(SCIP *scip, SCIP_Real val)
SCIP_Real SCIPfeastol(SCIP *scip)
SCIP_Bool SCIPisSumGE(SCIP *scip, SCIP_Real val1, SCIP_Real val2)
SCIP_Bool SCIPisGT(SCIP *scip, SCIP_Real val1, SCIP_Real val2)
SCIP_Bool SCIPisZero(SCIP *scip, SCIP_Real val)
SCIP_Bool SCIPisLT(SCIP *scip, SCIP_Real val1, SCIP_Real val2)
SCIP_Bool SCIPisSumGT(SCIP *scip, SCIP_Real val1, SCIP_Real val2)
int SCIPgetDepth(SCIP *scip)
Definition: scip_tree.c:672
SCIP_NODE * SCIPgetCurrentNode(SCIP *scip)
Definition: scip_tree.c:91
SCIP_RETCODE SCIPgetVarStrongbranchFrac(SCIP *scip, SCIP_VAR *var, int itlim, SCIP_Bool idempotent, SCIP_Real *down, SCIP_Real *up, SCIP_Bool *downvalid, SCIP_Bool *upvalid, SCIP_Bool *downinf, SCIP_Bool *upinf, SCIP_Bool *downconflict, SCIP_Bool *upconflict, SCIP_Bool *lperror)
Definition: scip_var.c:2919
SCIP_RETCODE SCIPtightenVarLb(SCIP *scip, SCIP_VAR *var, SCIP_Real newbound, SCIP_Bool force, SCIP_Bool *infeasible, SCIP_Bool *tightened)
Definition: scip_var.c:5326
SCIP_Bool SCIPpscostThresholdProbabilityTest(SCIP *scip, SCIP_VAR *var, SCIP_Real frac, SCIP_Real threshold, SCIP_BRANCHDIR dir, SCIP_CONFIDENCELEVEL clevel)
Definition: scip_var.c:9170
SCIP_Bool SCIPvarIsActive(SCIP_VAR *var)
Definition: var.c:17775
SCIP_Real SCIPgetVarPseudocostCountCurrentRun(SCIP *scip, SCIP_VAR *var, SCIP_BRANCHDIR dir)
Definition: scip_var.c:9061
SCIP_Real SCIPgetVarAvgInferenceScore(SCIP *scip, SCIP_VAR *var)
Definition: scip_var.c:9584
SCIP_RETCODE SCIPendStrongbranch(SCIP *scip)
Definition: scip_var.c:2744
SCIP_Real SCIPgetVarPseudocostCurrentRun(SCIP *scip, SCIP_VAR *var, SCIP_BRANCHDIR dir)
Definition: scip_var.c:9007
SCIP_Real SCIPvarGetUbLocal(SCIP_VAR *var)
Definition: var.c:18171
SCIP_Bool SCIPisVarPscostRelerrorReliable(SCIP *scip, SCIP_VAR *var, SCIP_Real threshold, SCIP_CONFIDENCELEVEL clevel)
Definition: scip_var.c:9189
SCIP_VAR * SCIPvarGetProbvar(SCIP_VAR *var)
Definition: var.c:12245
SCIP_RETCODE SCIPtightenVarUb(SCIP *scip, SCIP_VAR *var, SCIP_Real newbound, SCIP_Bool force, SCIP_Bool *infeasible, SCIP_Bool *tightened)
Definition: scip_var.c:5443
int SCIPvarGetProbindex(SCIP_VAR *var)
Definition: var.c:17795
const char * SCIPvarGetName(SCIP_VAR *var)
Definition: var.c:17446
SCIP_Longint SCIPgetVarStrongbranchNode(SCIP *scip, SCIP_VAR *var)
Definition: scip_var.c:4283
SCIP_LPSOLSTAT SCIPgetLastStrongbranchLPSolStat(SCIP *scip, SCIP_BRANCHDIR branchdir)
Definition: scip_var.c:4111
SCIP_Real SCIPgetVarPseudocostVal(SCIP *scip, SCIP_VAR *var, SCIP_Real solvaldelta)
Definition: scip_var.c:8925
SCIP_Real SCIPgetVarConflictlengthScore(SCIP *scip, SCIP_VAR *var)
Definition: scip_var.c:9414
SCIP_Bool SCIPsignificantVarPscostDifference(SCIP *scip, SCIP_VAR *varx, SCIP_Real fracx, SCIP_VAR *vary, SCIP_Real fracy, SCIP_BRANCHDIR dir, SCIP_CONFIDENCELEVEL clevel, SCIP_Bool onesided)
Definition: scip_var.c:9140
SCIP_RETCODE SCIPgetVarStrongbranchWithPropagation(SCIP *scip, SCIP_VAR *var, SCIP_Real solval, SCIP_Real lpobjval, int itlim, int maxproprounds, SCIP_Real *down, SCIP_Real *up, SCIP_Bool *downvalid, SCIP_Bool *upvalid, SCIP_Longint *ndomredsdown, SCIP_Longint *ndomredsup, SCIP_Bool *downinf, SCIP_Bool *upinf, SCIP_Bool *downconflict, SCIP_Bool *upconflict, SCIP_Bool *lperror, SCIP_Real *newlbs, SCIP_Real *newubs)
Definition: scip_var.c:3396
SCIP_Real SCIPvarGetLPSol(SCIP_VAR *var)
Definition: var.c:18479
SCIP_Real SCIPgetVarPseudocostScoreCurrentRun(SCIP *scip, SCIP_VAR *var, SCIP_Real solval)
Definition: scip_var.c:9252
SCIP_Real SCIPvarGetLbLocal(SCIP_VAR *var)
Definition: var.c:18161
SCIP_RETCODE SCIPupdateVarPseudocost(SCIP *scip, SCIP_VAR *var, SCIP_Real solvaldelta, SCIP_Real objdelta, SCIP_Real weight)
Definition: scip_var.c:8891
SCIP_Real SCIPgetVarAvgGMIScore(SCIP *scip, SCIP_VAR *var)
Definition: scip_var.c:9988
SCIP_Real SCIPgetVarLastGMIScore(SCIP *scip, SCIP_VAR *var)
Definition: scip_var.c:10043
SCIP_Real SCIPgetVarConflictScore(SCIP *scip, SCIP_VAR *var)
Definition: scip_var.c:9352
SCIP_Real SCIPgetVarAvgCutoffScore(SCIP *scip, SCIP_VAR *var)
Definition: scip_var.c:9838
SCIP_RETCODE SCIPgetVarStrongbranchLast(SCIP *scip, SCIP_VAR *var, SCIP_Real *down, SCIP_Real *up, SCIP_Bool *downvalid, SCIP_Bool *upvalid, SCIP_Real *solval, SCIP_Real *lpobjval)
Definition: scip_var.c:4133
SCIP_Real SCIPgetVarPseudocostScore(SCIP *scip, SCIP_VAR *var, SCIP_Real solval)
Definition: scip_var.c:9214
SCIP_RETCODE SCIPstartStrongbranch(SCIP *scip, SCIP_Bool enablepropagation)
Definition: scip_var.c:2686
void SCIPfreeRandom(SCIP *scip, SCIP_RANDNUMGEN **randnumgen)
SCIP_Real SCIPrandomGetReal(SCIP_RANDNUMGEN *randnumgen, SCIP_Real minrandval, SCIP_Real maxrandval)
Definition: misc.c:10133
SCIP_RETCODE SCIPcreateRandom(SCIP *scip, SCIP_RANDNUMGEN **randnumgen, unsigned int initialseed, SCIP_Bool useglobalseed)
memory allocation routines
#define BMSclearMemoryArray(ptr, num)
Definition: memory.h:130
SCIP_RETCODE SCIPgetSymmetry(SCIP *scip, int *npermvars, SCIP_VAR ***permvars, SCIP_HASHMAP **permvarmap, int *nperms, int ***perms, int ***permstrans, SCIP_Real *log10groupsize, SCIP_Bool *binvaraffected, int **components, int **componentbegins, int **vartocomponent, int *ncomponents)
propagator for symmetry handling
public methods for branching rules
public methods for managing constraints
public methods for message output
#define SCIPdebug(x)
Definition: pub_message.h:93
public data structures and miscellaneous methods
public methods for primal CIP solutions
public methods for branch and bound tree
public methods for problem variables
public methods for branching rule plugins and branching
public methods for constraint handler plugins and constraints
general public methods
public methods for the LP relaxation, rows and columns
public methods for memory management
public methods for message handling
public methods for nonlinear relaxation
public methods for numerical tolerances
public methods for SCIP parameter handling
public methods for global and local (sub)problems
public methods for random numbers
public methods for solutions
public methods for querying solving statistics
public methods for the branch-and-bound tree
public methods for SCIP variables
methods for handling symmetries
SCIP_RETCODE SCIPtreemodelSelectCandidate(SCIP *scip, SCIP_TREEMODEL *treemodel, SCIP_VAR **branchcands, SCIP_Real *mingains, SCIP_Real *maxgains, SCIP_Real *tiebreakerscore, int nbranchcands, int *bestcand)
Definition: treemodel.c:912
SCIP_RETCODE SCIPtreemodelInit(SCIP *scip, SCIP_TREEMODEL **treemodel)
Definition: treemodel.c:826
SCIP_RETCODE SCIPtreemodelFree(SCIP *scip, SCIP_TREEMODEL **treemodel)
Definition: treemodel.c:884
SCIP_Bool SCIPtreemodelIsEnabled(SCIP *scip, SCIP_TREEMODEL *treemodel)
Definition: treemodel.c:900
Branching rules based on the Single-Variable-Branching (SVB) model.
struct SCIP_BranchruleData SCIP_BRANCHRULEDATA
Definition: type_branch.h:57
@ SCIP_BRANCHDIR_DOWNWARDS
Definition: type_history.h:43
@ SCIP_BRANCHDIR_UPWARDS
Definition: type_history.h:44
@ SCIP_BOUNDTYPE_UPPER
Definition: type_lp.h:57
@ SCIP_BOUNDTYPE_LOWER
Definition: type_lp.h:56
enum SCIP_BoundType SCIP_BOUNDTYPE
Definition: type_lp.h:59
@ SCIP_LPSOLSTAT_OPTIMAL
Definition: type_lp.h:43
@ SCIP_LPSOLSTAT_INFEASIBLE
Definition: type_lp.h:44
@ SCIP_LPSOLSTAT_OBJLIMIT
Definition: type_lp.h:46
@ SCIP_LPSOLSTAT_ITERLIMIT
Definition: type_lp.h:47
@ SCIP_VERBLEVEL_HIGH
Definition: type_message.h:56
@ SCIP_CONFIDENCELEVEL_MAX
Definition: type_misc.h:51
@ SCIP_CONFIDENCELEVEL_MEDIUM
Definition: type_misc.h:49
@ SCIP_CONFIDENCELEVEL_HIGH
Definition: type_misc.h:50
@ SCIP_CONFIDENCELEVEL_MIN
Definition: type_misc.h:47
@ SCIP_CONFIDENCELEVEL_LOW
Definition: type_misc.h:48
enum SCIP_Confidencelevel SCIP_CONFIDENCELEVEL
Definition: type_misc.h:53
@ SCIP_DIDNOTRUN
Definition: type_result.h:42
@ SCIP_CUTOFF
Definition: type_result.h:48
@ SCIP_REDUCEDDOM
Definition: type_result.h:51
@ SCIP_CONSADDED
Definition: type_result.h:52
@ SCIP_BRANCHED
Definition: type_result.h:54
enum SCIP_Result SCIP_RESULT
Definition: type_result.h:61
@ SCIP_OKAY
Definition: type_retcode.h:42
enum SCIP_Retcode SCIP_RETCODE
Definition: type_retcode.h:63