Scippy

UG

Ubiquity Generator framework

scipParaObjSelfSplitNodesel.cpp
Go to the documentation of this file.
1/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
2/* */
3/* This file is part of the program and software framework */
4/* UG --- Ubquity Generator Framework */
5/* */
6/* Copyright Written by Yuji Shinano <shinano@zib.de>, */
7/* Copyright (C) 2021-2024 by Zuse Institute Berlin, */
8/* licensed under LGPL version 3 or later. */
9/* Commercial licenses are available through <licenses@zib.de> */
10/* */
11/* This code is free software; you can redistribute it and/or */
12/* modify it under the terms of the GNU Lesser General Public License */
13/* as published by the Free Software Foundation; either version 3 */
14/* of the License, or (at your option) any later version. */
15/* */
16/* This program is distributed in the hope that it will be useful, */
17/* but WITHOUT ANY WARRANTY; without even the implied warranty of */
18/* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the */
19/* GNU Lesser General Public License for more details. */
20/* */
21/* You should have received a copy of the GNU Lesser General Public License */
22/* along with this program. If not, see <http://www.gnu.org/licenses/>. */
23/* */
24/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
25
26/**@file scipParaObjSelfSplitNodesel.cpp
27 * @brief node selector for self-split ramp-up
28 * @author Yuji Shinano
29 */
30
31/*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
32
33#include <cassert>
34#include "scip/scip.h"
35#include "scip/debug.h"
36#include "ug_bb/bbParaComm.h"
37#include "scipParaObjSelfSplitNodesel.h"
38
39#if SCIP_APIVERSION >= 101
40
41using namespace ParaSCIP;
42
43/*
44 * Callback methods of node selector
45 */
46
47/** node selection method of node selector */
48SCIP_DECL_NODESELSELECT(ScipParaObjSelfSplitNodesel::scip_select)
49{ /*lint --e{715}*/
50 assert(nodesel != NULL);
51 assert(strcmp(SCIPnodeselGetName(nodesel), "ScipParaObjSelfSplitNodeSel") == 0);
52 assert(scip != NULL);
53 assert(selnode != NULL);
54
55 *selnode = NULL;
56
57 if( sampling )
58 {
59 SCIP_NODE* node;
60 node = SCIPgetBestNode(scip);
61 if( node == NULL )
62 {
63 *selnode = NULL;
64 return SCIP_OKAY;
65 }
66 /* if( SCIPnodeGetDepth(SCIPgetBestNode(scip)) > nodeseldata->depthlimit ) */
67 if( SCIPnodeGetDepth(node) > depthlimit )
68 {
69 SCIP_NODE** leaves;
70 SCIP_NODE** children;
71 SCIP_NODE** siblings;
72 SCIP_NODE** nodes;
73 int* depths;
74 int i;
75 int offset;
76
77 int nleaves;
78 int nsiblings;
79 int nchildren;
80 int nnodes;
81
82 /* collect leaves, children and siblings data */
83 SCIP_CALL( SCIPgetOpenNodesData(scip, &leaves, &children, &siblings, &nleaves, &nchildren, &nsiblings) );
84 nnodes = nleaves + nchildren + nsiblings;
85 assert(nnodes == SCIPgetNNodesLeft(scip));
86
87 // printf(">>>>>>>>>>>>>>Ending sampling phase: %d open nodes\n", nnodes);
88
89 SCIP_CALL( SCIPallocMemoryArray(scip, &nodes, nnodes) );
90 SCIP_CALL( SCIPallocMemoryArray(scip, &depths, nnodes) );
91
92 for(i = 0; i < nchildren; i++)
93 {
94 nodes[i] = children[i];
95 depths[i] = SCIPnodeGetDepth(children[i]);
96 }
97
98 for(i = nchildren; i < nchildren+nsiblings; i++)
99 {
100 nodes[i] = siblings[i-nchildren];
101 depths[i] = SCIPnodeGetDepth(siblings[i-nchildren]);
102 }
103
104 offset = nchildren+nsiblings;
105 for(i = offset; i < nnodes; i++)
106 {
107 nodes[i] = leaves[i-offset];
108 depths[i] = SCIPnodeGetDepth(leaves[i-offset]);
109 }
110
111 SCIPsortIntPtr(depths, (void**)nodes, nnodes);
112
113 for(i = 0; i < nnodes; i++)
114 {
115 if( (i % selfsplitsize) == selfsplitrank )
116 {
117 keepParaNode( scip, depths[i], nodes[i] );
118 }
119 SCIP_CALL( SCIPcutoffNode(scip,nodes[i]));
120 }
121
122 SCIP_CALL( SCIPpruneTree(scip) );
123 assert( SCIPgetNNodesLeft(scip) == 0 );
124
125 SCIPfreeMemoryArray(scip, &depths);
126 SCIPfreeMemoryArray(scip, &nodes);
127
128 sampling = false;
129
130 return SCIP_OKAY;
131 }
132 else
133 {
134 /* siblings come before leaves at the same level. Sometimes it can occur that no leaves are left except for children */
135 *selnode = SCIPgetBestSibling(scip);
136 if( *selnode == NULL )
137 {
138 *selnode = SCIPgetBestLeaf(scip);
139 if( *selnode == NULL )
140 *selnode=SCIPgetBestChild(scip);
141 }
142 if( *selnode != NULL )
143 {
144 SCIPdebugMessage("Selecting next node number %" SCIP_LONGINT_FORMAT " at depth %d\n", SCIPnodeGetNumber(*selnode), SCIPnodeGetDepth(*selnode));
145 }
146 return SCIP_OKAY;
147 }
148 }
149
150 if( *selnode == NULL )
151 {
152 *selnode = SCIPgetBestNode(scip);
153 }
154 // SCIP_CALL( SCIPselectNodeEstimate(scip, nodesel_estimate, selnode) );
155
156 return SCIP_OKAY;
157}
158
159
160/** node comparison method of node selector */
161SCIP_DECL_NODESELCOMP(ScipParaObjSelfSplitNodesel::scip_comp)
162{ /*lint --e{715}*/
163
164 /* @todo implement and call interface functions of preferred node selectors in SCIP */
165 if( sampling )
166 {
167 int depth1;
168 int depth2;
169
170 depth1 = SCIPnodeGetDepth(node1);
171 depth2 = SCIPnodeGetDepth(node2);
172
173 /* if depths differ, prefer node with smaller depth */
174 if( depth1 < depth2 )
175 return -1;
176 else if( depth1 > depth2 )
177 return +1;
178 else
179 {
180 /* depths are equal; prefer node with smaller number */
181 SCIP_Longint number1;
182 SCIP_Longint number2;
183
184 number1 = SCIPnodeGetNumber(node1);
185 number2 = SCIPnodeGetNumber(node2);
186 assert(number1 != number2);
187
188 if( number1 < number2 )
189 return -1;
190 else
191 return +1;
192 }
193 }
194 else
195 {
196 SCIP_Real estimate1;
197 SCIP_Real estimate2;
198
199
200 estimate1 = SCIPnodeGetEstimate(node1);
201 estimate2 = SCIPnodeGetEstimate(node2);
202 if( (SCIPisInfinity(scip, estimate1) && SCIPisInfinity(scip, estimate2)) ||
203 (SCIPisInfinity(scip, -estimate1) && SCIPisInfinity(scip, -estimate2)) ||
204 SCIPisEQ(scip, estimate1, estimate2) )
205 {
206 SCIP_Real lowerbound1;
207 SCIP_Real lowerbound2;
208
209 lowerbound1 = SCIPnodeGetLowerbound(node1);
210 lowerbound2 = SCIPnodeGetLowerbound(node2);
211 if( SCIPisLT(scip, lowerbound1, lowerbound2) )
212 return -1;
213 else if( SCIPisGT(scip, lowerbound1, lowerbound2) )
214 return +1;
215 else
216 {
217 SCIP_NODETYPE nodetype1;
218 SCIP_NODETYPE nodetype2;
219
220 nodetype1 = SCIPnodeGetType(node1);
221 nodetype2 = SCIPnodeGetType(node2);
222 if( nodetype1 == SCIP_NODETYPE_CHILD && nodetype2 != SCIP_NODETYPE_CHILD )
223 return -1;
224 else if( nodetype1 != SCIP_NODETYPE_CHILD && nodetype2 == SCIP_NODETYPE_CHILD )
225 return +1;
226 else if( nodetype1 == SCIP_NODETYPE_SIBLING && nodetype2 != SCIP_NODETYPE_SIBLING )
227 return -1;
228 else if( nodetype1 != SCIP_NODETYPE_SIBLING && nodetype2 == SCIP_NODETYPE_SIBLING )
229 return +1;
230 else
231 {
232 int depth1;
233 int depth2;
234 depth1 = SCIPnodeGetDepth(node1);
235 depth2 = SCIPnodeGetDepth(node2);
236 if( depth1 < depth2 )
237 return -1;
238 else if( depth1 > depth2 )
239 return +1;
240 else
241 return 0;
242 }
243 }
244 }
245
246 if( SCIPisLT(scip, estimate1, estimate2) )
247 return -1;
248
249 assert(SCIPisGT(scip, estimate1, estimate2));
250 return +1;
251 }
252}
253
254void
255ScipParaObjSelfSplitNodesel::keepParaNode(
256 SCIP *scip,
257 int depth,
258 SCIP_NODE *node
259 )
260{
261 assert( depth == SCIPnodeGetDepth( node ) );
262 SCIP_VAR **branchVars = new SCIP_VAR*[depth];
263 SCIP_Real *branchBounds = new SCIP_Real[depth];
264 SCIP_BOUNDTYPE *boundTypes = new SCIP_BOUNDTYPE[depth];
265 int nBranchVars;
266 SCIPnodeGetAncestorBranchings( node, branchVars, branchBounds, boundTypes, &nBranchVars, depth );
267 if( nBranchVars > depth ) // did not have enough memory, then reallocate
268 {
269 delete [] branchVars;
270 delete [] branchBounds;
271 delete [] boundTypes;
272 branchVars = new SCIP_VAR*[nBranchVars];
273 branchBounds = new SCIP_Real[nBranchVars];
274 boundTypes = new SCIP_BOUNDTYPE[nBranchVars];
275 SCIPnodeGetAncestorBranchings( node, branchVars, branchBounds, boundTypes, &nBranchVars, nBranchVars );
276 }
277
278 int nVars = SCIPgetNVars(scip);
279 // SCIP_VAR **vars = SCIPgetVars(scip);
280 int *iBranchVars = new int[nBranchVars];
281 // create the variable mapping hash map
282 SCIP_HASHMAP* varmapLb;
283 SCIP_HASHMAP* varmapUb;
284 SCIP_CALL_ABORT( SCIPhashmapCreate(&varmapLb, SCIPblkmem(scip), nVars) );
285 SCIP_CALL_ABORT( SCIPhashmapCreate(&varmapUb, SCIPblkmem(scip), nVars) );
286 for( int i = 0; i < nBranchVars; i++ )
287 {
288 iBranchVars[i] = i;
289 if( boundTypes[i] == SCIP_BOUNDTYPE_LOWER )
290 {
291 if( !SCIPhashmapGetImage(varmapLb, branchVars[i]) )
292 {
293 SCIP_CALL_ABORT( SCIPhashmapInsert(varmapLb, branchVars[i], &iBranchVars[i] ) );
294 }
295 }
296 else
297 {
298 if( !SCIPhashmapGetImage(varmapUb, branchVars[i]) )
299 {
300 SCIP_CALL_ABORT( SCIPhashmapInsert(varmapUb, branchVars[i], &iBranchVars[i] ) );
301 }
302 }
303 }
304
305 /*
306 SCIP_VAR **preBranchVars = branchVars;
307 SCIP_Real *preBranchBounds = branchBounds;
308 SCIP_BOUNDTYPE *preBboundTypes = boundTypes;
309 branchVars = new SCIP_VAR*[nBranchVars+nVars*2];
310 branchBounds = new SCIP_Real[nBranchVars+nVars*2];
311 boundTypes = new SCIP_BOUNDTYPE[nBranchVars+nVars*2];
312 for( int i = 0; i < nBranchVars; i++ )
313 {
314 branchVars[i] = preBranchVars[i];
315 branchBounds[i] = preBranchBounds[i];
316 boundTypes[i] = preBboundTypes[i];
317 }
318 int *iBranchVar = NULL;
319 for( int i = 0; i < nVars; i++ )
320 {
321 if( scipParaSolver->isCopyIncreasedVariables() &&
322 scipParaSolver->getOriginalIndex(i) >= scipParaSolver->getNOrgVars() )
323 {
324 continue;
325 }
326 iBranchVar = (int *)SCIPhashmapGetImage(varmapLb, vars[i]);
327 if( iBranchVar )
328 {
329 if( EPSLT(preBranchBounds[*iBranchVar], SCIPvarGetLbGlobal(vars[i]) ,DEFAULT_NUM_EPSILON ) )
330 {
331 branchBounds[*iBranchVar] = SCIPvarGetLbGlobal(vars[i]); // node is current node
332 if ( EPSGT(branchBounds[*iBranchVar], SCIPvarGetUbGlobal(vars[i]), DEFAULT_NUM_EPSILON) ) abort();
333 }
334 }
335 else
336 {
337 if( EPSGT( SCIPvarGetLbGlobal(vars[i]), scipParaSolver->getOrgVarLb(i), MINEPSILON ) )
338 {
339 branchVars[nBranchVars] = vars[i];
340 branchBounds[nBranchVars] = SCIPvarGetLbGlobal(vars[i]);
341 boundTypes[nBranchVars] = SCIP_BOUNDTYPE_LOWER;
342 if ( EPSGT(branchBounds[nBranchVars], SCIPvarGetUbGlobal(vars[i]), DEFAULT_NUM_EPSILON) ) abort();
343 if ( EPSGT(branchBounds[nBranchVars], scipParaSolver->getOrgVarUb(i), DEFAULT_NUM_EPSILON) ) abort();
344 nBranchVars++;
345 // std::cout << "********* GLOBAL is updated!*************" << std::endl;
346 }
347 }
348 iBranchVar = (int *)SCIPhashmapGetImage(varmapUb, vars[i]);
349 if( iBranchVar )
350 {
351 if( EPSGT(preBranchBounds[*iBranchVar], SCIPvarGetUbGlobal(vars[i]) ,DEFAULT_NUM_EPSILON ) )
352 {
353 branchBounds[*iBranchVar] = SCIPvarGetUbGlobal(vars[i]);
354 if ( EPSLT(branchBounds[*iBranchVar], SCIPvarGetLbGlobal(vars[i]),DEFAULT_NUM_EPSILON) ) abort();
355 }
356 }
357 else
358 {
359 if( EPSLT( SCIPvarGetUbGlobal(vars[i]), scipParaSolver->getOrgVarUb(i), MINEPSILON ) )
360 {
361 branchVars[nBranchVars] = vars[i];
362 branchBounds[nBranchVars] = SCIPvarGetUbGlobal(vars[i]);
363 boundTypes[nBranchVars] = SCIP_BOUNDTYPE_UPPER;
364 if ( EPSLT(branchBounds[nBranchVars], SCIPvarGetLbGlobal(vars[i]),DEFAULT_NUM_EPSILON) ) abort();
365 if ( EPSLT(branchBounds[nBranchVars], scipParaSolver->getOrgVarLb(i),DEFAULT_NUM_EPSILON) ) abort();
366 nBranchVars++;
367 // std::cout << "********* GLOBAL is updated!*************" << std::endl;
368 }
369 }
370 }
371 */
372 SCIPhashmapFree(&varmapLb);
373 SCIPhashmapFree(&varmapUb);
374 /*
375 delete [] preBranchVars;
376 delete [] preBranchBounds;
377 delete [] preBboundTypes;
378 */
379 delete [] iBranchVars;
380
381 if( scipParaSolver->isCopyIncreasedVariables() ) // this may not need, but only for this error occurred so far.
382 {
383 if( !ifFeasibleInOriginalProblem(scip, nBranchVars, branchVars, branchBounds) )
384 {
385 delete [] branchVars;
386 delete [] branchBounds;
387 delete [] boundTypes;
388 return;
389 }
390 }
391
392 SCIP_CONS** addedcons = 0;
393 int addedconsssize = SCIPnodeGetNAddedConss(node);
394 int naddedconss = 0;
395 if( addedconsssize > 0 )
396 {
397 SCIP_CALL_ABORT( SCIPallocBufferArray(scip, &addedcons, addedconsssize) );
398 SCIPnodeGetAddedConss(node, addedcons, &naddedconss, addedconsssize);
399 }
400
401 assert( scipParaSolver->getParentDiffSubproblem() == 0 );
402
403 DEF_SCIP_PARA_COMM( scipParaComm, paraComm);
404 ScipParaDiffSubproblem *diffSubproblem = scipParaComm->createScipParaDiffSubproblem(
405 scip,
406 scipParaSolver,
407 nBranchVars,
408 branchVars,
409 branchBounds,
410 boundTypes,
411 naddedconss,
412 addedcons
413 );
414
415 if( naddedconss )
416 {
417 SCIPfreeBufferArray(scip, &addedcons);
418 }
419
420
421 long long n = SCIPnodeGetNumber( node );
422 // Cutoff looks mess-up this assert(SCIPisFeasGE(scip, SCIPretransformObj(scip,SCIPnodeGetLowerbound(node)), SCIPgetDualbound(scip)));
423 double dualBound = std::min(SCIPretransformObj(scip, SCIPnodeGetLowerbound( node )), SCIPgetDualbound(scip));
424 // std::cout << "SCIPretransformObj(scip, SCIPnodeGetLowerbound( node )) = " << SCIPretransformObj(scip, SCIPnodeGetLowerbound( node )) << ", SCIPgetDualbound(scip) = " << SCIPgetDualbound(scip) << std::endl;
425 // if( SCIPisObjIntegral(scip) )
426 // {
427 // dualBound = ceil(dualBound);
428 // }
429 // assert(SCIPisFeasGE(scip, SCIPnodeGetLowerbound(node) , SCIPgetLowerbound(scip)));
430 double estimateValue = SCIPnodeGetEstimate( node );
431 // assert(SCIPisFeasGE(scip, estimateValue, SCIPnodeGetLowerbound(node) ));
432#ifdef UG_DEBUG_SOLUTION
433 SCIP_Bool valid = 0;
434 SCIP_CALL_ABORT( SCIPdebugSolIsValidInSubtree(scip, &valid) );
435 diffSubproblem->setOptimalSolIndicator(valid);
436 std::cout << "* R." << scipParaSolver->getRank() << ", debug = " << SCIPdebugSolIsEnabled(scip) << ", valid = " << valid << std::endl;
437#endif
438 scipParaSolver->keepParaNode(n, depth, dualBound, estimateValue, diffSubproblem);
439
440 /** remove the node sent from SCIP environment */
441#ifdef UG_DEBUG_SOLUTION
442 if( valid )
443 {
444 SCIPdebugSolDisable(scip);
445 std::cout << "R." << paraComm->getRank() << ": disable debug, node which contains optmal solution is sent." << std::endl;
446 }
447#endif
448
449 delete [] branchVars;
450 delete [] branchBounds;
451 delete [] boundTypes;
452
453}
454
455bool
456ScipParaObjSelfSplitNodesel::ifFeasibleInOriginalProblem(
457 SCIP *scip,
458 int nBranchVars,
459 SCIP_VAR **branchVars,
460 SCIP_Real *inBranchBounds)
461{
462
463 bool feasible = true;
464 SCIP_Real *branchBounds = new SCIP_Real[nBranchVars];
465 for( int v = nBranchVars -1 ; v >= 0; --v )
466 {
467 SCIP_VAR *transformVar = branchVars[v];
468 SCIP_Real scalar = 1.0;
469 SCIP_Real constant = 0.0;
470 SCIP_CALL_ABORT( SCIPvarGetOrigvarSum(&transformVar, &scalar, &constant ) );
471 if( transformVar == NULL ) continue;
472 // assert( scalar == 1.0 && constant == 0.0 );
473 branchBounds[v] = ( inBranchBounds[v] - constant ) / scalar;
474 if( SCIPvarGetType(transformVar) != SCIP_VARTYPE_CONTINUOUS
475 && SCIPvarGetProbindex(transformVar) < scipParaSolver->getNOrgVars() )
476 {
477 if( !(SCIPisLE(scip,scipParaSolver->getOrgVarLb(SCIPvarGetProbindex(transformVar)), branchBounds[v]) &&
478 SCIPisGE(scip,scipParaSolver->getOrgVarUb(SCIPvarGetProbindex(transformVar)), branchBounds[v])) )
479 {
480 feasible = false;
481 break;
482 }
483 }
484 }
485 delete [] branchBounds;
486
487 return feasible;
488}
489
490#endif
491
bbParaComm.h
Base class of communicator for UG Framework.
DEF_SCIP_PARA_COMM
#define DEF_SCIP_PARA_COMM(scip_para_comm, comm)
Definition: scipParaCommMpi.h:201
SCIP_DECL_NODESELSELECT
SCIP_DECL_NODESELSELECT(ScipParaObjNodesel::scip_select)
Definition: scipParaObjNodesel.cpp:39
SCIP_DECL_NODESELCOMP
SCIP_DECL_NODESELCOMP(ScipParaObjNodesel::scip_comp)
Definition: scipParaObjNodesel.cpp:105
scipParaObjSelfSplitNodesel.h
node selector for self-split ramp-up