1 
2 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
3 /*                                                                           */
4 /*                  This file is part of the program and library             */
5 /*         SCIP --- Solving Constraint Integer Programs                      */
6 /*                                                                           */
7 /*    Copyright (C) 2002-2021 Konrad-Zuse-Zentrum                            */
8 /*                            fuer Informationstechnik Berlin                */
9 /*                                                                           */
10 /*  SCIP is distributed under the terms of the ZIB Academic License.         */
11 /*                                                                           */
12 /*  You should have received a copy of the ZIB Academic License              */
13 /*  along with SCIP; see the file COPYING. If not visit scipopt.org.         */
14 /*                                                                           */
15 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
16 
17 /**@file   pricer_stp.c
18  * @brief  stp variable pricer
19  * @author Daniel Rehfeldt
20  */
21 
22 #include <assert.h>
23 #include <string.h>
24 #include "scip/cons_linear.h"
25 #include <stdio.h>
26 #include <stdlib.h>
27 #include "pricer_stp.h"
28 #include "grph.h"
29 
30 
31 #define PRICER_NAME            "stp"
32 #define PRICER_DESC            "pricer for stp"
33 #define PRICER_PRIORITY        1
34 #define PRICER_DELAY           TRUE    /* only call pricer if all problem variables have non-negative reduced costs */
35 
36 /*
37  * Data structures
38  */
39 
40 
41 /** variable pricer data */
42 struct SCIP_PricerData
43 {
44    SCIP*                 scip;               /**< SCIP data structure */
45    SCIP_CONS**           edgecons;           /**< array containing all edge constraints */
46    SCIP_CONS**           pathcons;           /**< array containing all path constraints */
47    SCIP_Real*            pi;                 /**< array of the dual solutions to the edge constraints*/
48    SCIP_Real*            mi;	             /**< array of the dual solutions to the path constraint*/
49    SCIP_Real             lowerbound;         /**< lower bound computed by the pricer */
50    int* 	         realterms;          /**< array of all terminals except the root */
51    int*                  ncreatedvars;       /**< array (1,...,realnterms); in ncreatedvars[t] the number of
52                                                 created variables (during pricing) pertaining to terminal t is stored */
53    int                   root;               /**< root node */
54    int                   nedges;             /**< number of edges */
55    int                   realnterms;         /**< number of terminals except the root */
56    SCIP_Bool             bigt;               /**< stores whether the 'T' model is being used */
57 };
58 
59 
60 /*
61  * Local methods
62  */
63 
64 
65 /*
66  * Callback methods of variable pricer
67  */
68 
69 /** copy method for pricer plugins (called when SCIP copies plugins) */
70 static
SCIP_DECL_PRICERCOPY(pricerCopyStp)71 SCIP_DECL_PRICERCOPY(pricerCopyStp)
72 { /*lint --e{715}*/
73    SCIPdebugPrintf("pricerCopy \n");
74    assert(pricer != NULL);
75    assert(strcmp(SCIPpricerGetName(pricer), PRICER_NAME) == 0);
76 
77    return SCIP_OKAY;
78 }
79 
80 /** destructor of variable pricer to free user data (called when SCIP is exiting) */
81 /**! [SnippetPricerFreeSTP] */
82 static
SCIP_DECL_PRICERFREE(pricerFreeStp)83 SCIP_DECL_PRICERFREE(pricerFreeStp)
84 {
85    SCIP_PRICERDATA* pricerdata;
86    SCIPdebugPrintf("pricerfree \n");
87    assert(scip != NULL);
88 
89    /* get pricerdata */
90    pricerdata = SCIPpricerGetData(pricer);
91 
92    /* free memory for pricerdata*/
93    if ( pricerdata != NULL )
94    {
95       SCIPfreeMemory(scip, &pricerdata);
96    }
97 
98    SCIPpricerSetData(pricer, NULL);
99 
100    return SCIP_OKAY;
101 }
102 /**! [SnippetPricerFreeSTP] */
103 
104 /** initialization method of variable pricer (called after problem was transformed) */
105 static
SCIP_DECL_PRICERINIT(pricerInitStp)106 SCIP_DECL_PRICERINIT(pricerInitStp)
107 {
108    SCIP_PRICERDATA* pricerdata;
109    SCIPdebugPrintf("pricer Init \n");
110    assert(scip != NULL);
111    assert(pricer != NULL);
112 
113    pricerdata = SCIPpricerGetData(pricer);
114    pricerdata->edgecons = SCIPprobdataGetEdgeConstraints(scip);
115    pricerdata->pathcons = SCIPprobdataGetPathConstraints(scip);
116    pricerdata->root = SCIPprobdataGetRoot(scip);
117    pricerdata->nedges = SCIPprobdataGetNEdges(scip);
118    pricerdata->realnterms = SCIPprobdataGetRNTerms(scip);
119    pricerdata->realterms = SCIPprobdataGetRTerms(scip);
120    pricerdata->bigt = SCIPprobdataIsBigt(scip);
121    pricerdata->lowerbound = 0;
122    return SCIP_OKAY;
123 }
124 
125 /** solving process initialization method of variable pricer (called when branch and bound process is about to begin) */
126 static
SCIP_DECL_PRICERINITSOL(pricerInitsolStp)127 SCIP_DECL_PRICERINITSOL(pricerInitsolStp)
128 {
129    SCIP_PRICERDATA* pricerdata;
130    SCIPdebugPrintf("pricerinitsol \n");
131    assert(scip != NULL);
132    assert(pricer != NULL);
133 
134    pricerdata = SCIPpricerGetData(pricer);
135    assert(pricerdata != NULL);
136 
137    /* allocate memory */
138    if( !pricerdata->bigt )
139    {
140       SCIP_CALL( SCIPallocMemoryArray(scip, &(pricerdata->pi), SCIPprobdataGetNEdges(scip) * SCIPprobdataGetRNTerms(scip)) );
141    }
142    else
143    {
144       SCIP_CALL( SCIPallocMemoryArray(scip, &(pricerdata->pi), SCIPprobdataGetNEdges(scip)) );
145    }
146 
147    SCIP_CALL( SCIPallocMemoryArray(scip, &(pricerdata->mi), SCIPprobdataGetRNTerms(scip)) );
148    SCIP_CALL( SCIPallocMemoryArray(scip, &pricerdata->ncreatedvars, SCIPprobdataGetRNTerms(scip)) );
149    BMSclearMemoryArray(pricerdata->ncreatedvars, SCIPprobdataGetRNTerms(scip));
150 
151    return SCIP_OKAY;
152 }
153 
154 /** solving process deinitialization method of variable pricer (called before branch and bound process data is freed) */
155 static
SCIP_DECL_PRICEREXITSOL(pricerExitsolStp)156 SCIP_DECL_PRICEREXITSOL(pricerExitsolStp)
157 {
158    SCIP_PRICERDATA* pricerdata;
159 
160    assert(scip != NULL);
161    assert(pricer != NULL);
162 
163    pricerdata = SCIPpricerGetData(pricer);
164    assert(pricerdata != NULL);
165 
166    /* free memory */
167    SCIPfreeMemoryArray(scip, &(pricerdata->mi));
168    SCIPfreeMemoryArray(scip, &(pricerdata->pi));
169    SCIPfreeMemoryArray(scip, &pricerdata->ncreatedvars);
170 
171    return SCIP_OKAY;
172 }
173 
174 /** method for either Farkas or Redcost pricing */
175 static
pricing(SCIP * scip,SCIP_PRICER * pricer,SCIP_Real * lowerbound,SCIP_Bool farkas)176 SCIP_RETCODE pricing(
177    SCIP*                 scip,               /**< SCIP data structure */
178    SCIP_PRICER*          pricer,             /**< pricer */
179    SCIP_Real*            lowerbound,         /**< lowerbound pointer */
180    SCIP_Bool             farkas              /**< TRUE: Farkas pricing; FALSE: Redcost pricing */
181    )
182 {
183    SCIP_PRICERDATA* pricerdata; /* the data of the pricer */
184    SCIP_PROBDATA* probdata;
185    GRAPH* graph;
186    SCIP_VAR* var;
187    PATH* path;
188    SCIP_Real* edgecosts;  /* edgecosts of the current subproblem */
189    char varname[SCIP_MAXSTRLEN];
190    SCIP_Real newlowerbound = -SCIPinfinity(scip);
191    SCIP_Real redcost;   /* reduced cost */
192    int tail;
193    int e;
194    int t;
195    int i;
196 
197    assert(scip != NULL);
198    assert(pricer != NULL);
199 
200    /* get pricer data */
201    pricerdata = SCIPpricerGetData(pricer);
202    assert(pricerdata != NULL);
203 
204    /* get problem data */
205    probdata = SCIPgetProbData(scip);
206    assert(probdata != NULL);
207 
208    SCIPdebugMessage("solstat=%d\n", SCIPgetLPSolstat(scip));
209 
210    if( !farkas && SCIPgetLPSolstat(scip) == SCIP_LPSOLSTAT_OPTIMAL )
211       newlowerbound = SCIPgetSolTransObj(scip, NULL);
212 
213    SCIPdebug( SCIP_CALL( SCIPprintSol(scip, NULL, NULL, FALSE) ) );
214 
215 # if 0
216    if ( pricerdata->lowerbound <= 4 )
217    {
218       char label[SCIP_MAXSTRLEN];
219       (void)SCIPsnprintf(label, SCIP_MAXSTRLEN, "X%g.gml", pricerdata->lowerbound);
220       SCIP_CALL( SCIPprobdataPrintGraph(scip, label , NULL, TRUE) );
221       pricerdata->lowerbound++;
222    }
223 #endif
224    /* get the graph*/
225    graph = SCIPprobdataGetGraph(probdata);
226 
227    /* get dual solutions and save them in mi and pi */
228    for( t = 0; t < pricerdata->realnterms; ++t )
229    {
230       if( farkas )
231       {
232 	 pricerdata->mi[t] = SCIPgetDualfarkasLinear(scip, pricerdata->pathcons[t]);
233       }
234       else
235       {
236          pricerdata->mi[t] = SCIPgetDualsolLinear(scip, pricerdata->pathcons[t]);
237          assert(!SCIPisNegative(scip, pricerdata->mi[t]));
238       }
239    }
240 
241    for( e = 0; e < pricerdata->nedges; ++e )
242    {
243       if( !pricerdata->bigt )
244       {
245          for( t = 0; t < pricerdata->realnterms; ++t )
246          {
247             if( farkas )
248 	    {
249                pricerdata->pi[t * pricerdata->nedges + e] = SCIPgetDualfarkasLinear(
250                   scip, pricerdata->edgecons[t * pricerdata->nedges + e]);
251 	    }
252             else
253 	    {
254                pricerdata->pi[t * pricerdata->nedges + e] = SCIPgetDualsolLinear(
255                   scip, pricerdata->edgecons[t * pricerdata->nedges + e]);
256 	    }
257          }
258       }
259       else
260       {
261          if( farkas )
262 	 {
263 	    pricerdata->pi[e] = SCIPgetDualfarkasLinear(
264                scip, pricerdata->edgecons[e]);
265 	 }
266 	 else
267 	 {
268 	    pricerdata->pi[e] = SCIPgetDualsolLinear(
269                scip, pricerdata->edgecons[e]);
270 	 }
271       }
272    }
273 
274    SCIP_CALL( SCIPallocMemoryArray(scip, &path, graph->knots) );
275    SCIP_CALL( SCIPallocMemoryArray(scip, &edgecosts, pricerdata->nedges) );
276 
277    if( pricerdata->bigt )
278    {
279       for( e = 0; e < pricerdata->nedges; ++e )
280       {
281          edgecosts[e] = (-pricerdata->pi[e]);
282       }
283    }
284    /* find shortest r-t (r root, t terminal) paths and create corresponding variables iff reduced cost < 0 */
285    for( t = 0; t < pricerdata->realnterms; ++t )
286    {
287       for( e = 0; e < pricerdata->nedges; ++e )
288       {
289 	 if( !pricerdata->bigt )
290 	 {
291             edgecosts[e] = (-pricerdata->pi[t * pricerdata->nedges + e]);
292 	 }
293 
294          assert(!SCIPisNegative(scip, edgecosts[e]));
295       }
296 
297       for( i = 0; i < graph->knots; i++ )
298          graph->mark[i] = 1;
299 
300       graph_path_exec(scip, graph, FSP_MODE, pricerdata->root, edgecosts, path);
301 
302       /* compute reduced cost of shortest path to terminal t */
303       redcost = 0.0;
304       tail = pricerdata->realterms[t];
305       while( tail != pricerdata->root )
306       {
307          redcost += edgecosts[path[tail].edge];
308 	 tail = graph->tail[path[tail].edge];
309       }
310       redcost -= pricerdata->mi[t];
311 
312       if( !farkas && SCIPgetLPSolstat(scip) == SCIP_LPSOLSTAT_OPTIMAL )
313       {
314          newlowerbound += redcost;
315       }
316       /* check if reduced cost < 0 */
317       if( SCIPisNegative(scip, redcost) )
318       {
319 	 /* create variable to the shortest path (having reduced cost < 0) */
320          var = NULL;
321 	 sprintf(varname, "PathVar%d_%d", t, pricerdata->ncreatedvars[t]);
322          ++(pricerdata->ncreatedvars[t]);
323 
324          SCIP_CALL( SCIPcreateVarBasic(scip, &var, varname, 0.0, SCIPinfinity(scip), 0.0, SCIP_VARTYPE_CONTINUOUS) );
325          SCIP_CALL( SCIPaddPricedVar(scip, var, -redcost) );
326          tail = pricerdata->realterms[t];
327          while( tail != pricerdata->root )
328          {
329             /* add variable to constraints */
330 	    if( !pricerdata->bigt )
331 	    {
332 	       SCIP_CALL( SCIPaddCoefLinear(scip, pricerdata->edgecons[t * pricerdata->nedges + path[tail].edge], var, 1.0) );
333 	    }
334 	    else
335 	    {
336 	       SCIP_CALL( SCIPaddCoefLinear(scip, pricerdata->edgecons[path[tail].edge], var, 1.0) );
337 	    }
338 
339 	    tail = graph->tail[path[tail].edge];
340          }
341          SCIP_CALL( SCIPaddCoefLinear(scip, pricerdata->pathcons[t], var, 1.0) );
342       }
343    }
344 
345    if( !farkas && SCIPgetLPSolstat(scip) == SCIP_LPSOLSTAT_OPTIMAL )
346       *lowerbound = newlowerbound;
347 
348    SCIPfreeMemoryArray(scip, &edgecosts);
349    SCIPfreeMemoryArray(scip, &path);
350 
351    return SCIP_OKAY;
352 }
353 
354 /** reduced cost pricing method of variable pricer for feasible LPs */
355 static
SCIP_DECL_PRICERREDCOST(pricerRedcostStp)356 SCIP_DECL_PRICERREDCOST(pricerRedcostStp)
357 { /*lint --e{715}*/
358    SCIP_CALL( pricing(scip, pricer, lowerbound, FALSE) );
359 
360    /* set result pointer */
361    *result = SCIP_SUCCESS;
362    return SCIP_OKAY;
363 }
364 
365 /** Farkas pricing method of variable pricer for infeasible LPs */
366 static
SCIP_DECL_PRICERFARKAS(pricerFarkasStp)367 SCIP_DECL_PRICERFARKAS(pricerFarkasStp)
368 { /*lint --e{715}*/
369    SCIP_CALL( pricing(scip, pricer,  NULL, TRUE) );
370 
371    return SCIP_OKAY;
372 }
373 
374 /*
375  * variable pricer specific interface methods
376  */
377 
378 /** creates the stp variable pricer and includes it in SCIP */
SCIPincludePricerStp(SCIP * scip)379 SCIP_RETCODE SCIPincludePricerStp(
380    SCIP*                 scip                /**< SCIP data structure */
381    )
382 {
383 
384    SCIP_PRICERDATA* pricerdata;
385    SCIP_PRICER* pricer;
386    SCIPdebugPrintf("include Pricer \n");
387 
388    /* create stp variable pricer data */
389    SCIP_CALL( SCIPallocMemory(scip, &pricerdata) );
390    pricerdata->scip = scip;
391 
392    /* include variable pricer */
393    pricer = NULL;
394    SCIP_CALL( SCIPincludePricerBasic(scip, &pricer, PRICER_NAME, PRICER_DESC, PRICER_PRIORITY, PRICER_DELAY,
395          pricerRedcostStp, pricerFarkasStp, pricerdata) );
396    assert(pricer != NULL);
397 
398    /* set non fundamental callbacks via setter functions */
399    SCIP_CALL( SCIPsetPricerCopy(scip, pricer, pricerCopyStp) );
400    SCIP_CALL( SCIPsetPricerFree(scip, pricer, pricerFreeStp) );
401    SCIP_CALL( SCIPsetPricerInit(scip, pricer, pricerInitStp) );
402    SCIP_CALL( SCIPsetPricerInitsol(scip, pricer, pricerInitsolStp) );
403    SCIP_CALL( SCIPsetPricerExitsol(scip, pricer, pricerExitsolStp) );
404 
405    return SCIP_OKAY;
406 }
407