1% BEGIN LICENSE BLOCK
2% Version: CMPL 1.1
3%
4% The contents of this file are subject to the Cisco-style Mozilla Public
5% License Version 1.1 (the "License"); you may not use this file except
6% in compliance with the License.  You may obtain a copy of the License
7% at www.eclipse-clp.org/license.
8% 
9% Software distributed under the License is distributed on an "AS IS"
10% basis, WITHOUT WARRANTY OF ANY KIND, either express or implied.  See
11% the License for the specific language governing rights and limitations
12% under the License. 
13% 
14% The Original Code is  The ECLiPSe Constraint Logic Programming System. 
15% The Initial Developer of the Original Code is  Cisco Systems, Inc. 
16% Portions created by the Initial Developer are
17% Copyright (C) 1995 - 2006 Cisco Systems, Inc.  All Rights Reserved.
18% 
19% Contributor(s): Joachim Schimpf, IC-Parc
20% 
21% END LICENSE BLOCK
22% ----------------------------------------------------------------------
23% System:	ECLiPSe Constraint Logic Programming System
24% Version:	$Id: mip.ecl,v 1.1 2012/07/31 02:17:06 jschimpf Exp $
25%
26% MIP branch and bound in ECLiPSe using lib(eplex)
27%
28% J.Schimpf, IC-Parc, 1/96
29% ----------------------------------------------------------------------
30
31% ----------------------------------------------------------------------
32:- module(mip).
33% ----------------------------------------------------------------------
34
35:- comment(categories, ["Constraints"]).
36:- comment(summary, "An example implementing MIP-style branch-and-bound").
37:- comment(author, "Joachim Schimpf, IC-Parc").
38:- comment(date, "$Date: 2012/07/31 02:17:06 $").
39:- comment(copyright, "Cisco Systems, Inc.").
40
41:- use_module(library(eplex)).
42:- lib(branch_and_bound).
43
44:- export bb_inf/3.
45:- export bb_inf/4.
46
47
48% ----------------------------------------------------------------------
49% MIP branch-and-bound
50% ----------------------------------------------------------------------
51
52bb_inf(IntVars, Expr, MinIpCost) :-
53	bb_inf(IntVars, Expr, MinIpCost, eplex).
54
55bb_inf(IntVars, Expr, MinIpCost, Pool) :-
56
57        Pool: (Expr $= IpCost),
58	% setup a simplex demon which wakes on bound changes
59	Pool:eplex_solver_setup(min(Expr), IpCost, [], 9, [deviating_bounds
60%		,pre(writeln(woken)),post(writeln(lpsol:IntVars))
61	]),
62	Pool:eplex_get(handle, ProbHandle),
63
64	% print some statistics
65	lp_get(ProbHandle, vars, VArr),
66	functor(VArr, _, NVars),
67	length(IntVars, NIntVars),
68	writeln(log_output, variables:NVars/intvars:NIntVars),
69
70	% call force_integers/1 within a branch-and-bound framework
71	int_tolerance(Delta),
72	bb_min((
73		force_integers(IntVars, ProbHandle, 0),
74		lp_get(ProbHandle, typed_solution, SolutionArr),
75		lp_get(ProbHandle, cost, IpCost)
76	    ),
77	    IpCost, SolutionArr, OptSolutionArr, MinIpCost,
78	    bb_options with [strategy:continue,delta:Delta]
79	),
80
81	% print some statistics
82	lp_get(ProbHandle, statistics, [SS,SF|_]),
83	writeln(log_output, (solver_succs=SS, solver_fails=SF)),
84	lp_cleanup(ProbHandle),
85
86	% instantiate solutions
87	( VArr = OptSolutionArr, true -> true
88	; writeln(error, "Instantiating solution lead to failure") ).
89
90
91% Keep branching as long as there are integer-variables whose lp-solutions
92% are not integral. The lp-demon will keep waking up during the process.
93
94force_integers(IntVars, ProbHandle, D0) :-
95%	( D0 > 20 -> abort ; true ),
96	( integer_violation(IntVars, BranchingVar, SplitValue, ProbHandle) ->
97%	    writeln(D0:pick(BranchingVar,IntVars)),
98	    branch(BranchingVar, SplitValue, D0),
99	    D1 is D0+1,
100	    force_integers(IntVars, ProbHandle, D1)
101	;
102	    true    
103	).
104
105
106    branch(BranchingVar, SplitValue, Depth) :-
107	S is round(SplitValue),
108	( S > SplitValue ->
109	    (
110		% try the upper sub-range first
111%		writeln(Depth:left->lwb(S,BranchingVar)),
112		set_var_bounds(BranchingVar, S, 1.0Inf)
113	    ;
114		S1 is S-1.0,
115%		writeln(Depth:right->upb(S1,BranchingVar)),
116		set_var_bounds(BranchingVar, -1.0Inf, S1)
117	    )
118	;
119	    (
120		% try the lower sub-range first
121%		writeln(Depth:left->upb(S,BranchingVar)),
122		set_var_bounds(BranchingVar, -1.0Inf, S)
123	    ;
124		S1 is S+1.0,
125%		writeln(Depth:right->lwb(S1,BranchingVar)),
126		set_var_bounds(BranchingVar, S1, 1.0Inf)
127	    )
128	).
129
130
131% ----------------------------------------------------------------------
132% Variable selection
133% ----------------------------------------------------------------------
134
135integer_violation([X|Xs], FractVar, FractSol, ProbHandle) :-
136	lp_var_solution(ProbHandle, X, Val),
137	( abs(Val - round(Val)) >= int_tolerance ->
138	    FractVar = X, FractSol = Val
139	;
140	    integer_violation(Xs, FractVar, FractSol, ProbHandle)
141	).
142
143
144/*
145integer_violation(Xs, FractVar, FractSol) :-
146	integer_violation(Xs, _, 1.0, FractVar),
147%	integer_violation(Xs, _, 0.0, FractVar),
148	lp_var_solution(FractVar, FractSol),
149%	call(get_var_index(FractVar, Idx))@eplex,
150%	Idx1 is Idx+1,
151%	writeln(x(Idx1):FractSol),
152	true.
153
154integer_violation([], BestX, BestDiff, BestX) :-
155	BestDiff >= int_tolerance,	% Did we actually find one?
156	BestDiff < 1.0.
157integer_violation([X|Xs], BestX, BestDiff, Res) :-
158	lp_var_solution(X, Val),
159	Diff is abs(Val - round(Val)),
160%	( var(X), Val \== 0.0 ->
161%	    call(get_var_index(X, Idx))@eplex,
162%	    Idx1 is Idx+1,
163%	    write('    '),
164%	    writeln(x(Idx1)=Val)
165%	;
166%	    true
167%	),
168	(
169	    Diff >= int_tolerance,
170%	    better(BestX, BestDiff, X, Diff)
171%	    Diff < BestDiff
172	    Diff =< BestDiff	% use 1.0 initially
173%	    Diff > BestDiff	% use 0.0 initially
174	->
175	    integer_violation(Xs, X, Diff, Res)
176	;
177	    integer_violation(Xs, BestX, BestDiff, Res)
178	).
179*/
180
181% prefer general integers to binaries
182better(OldX, _OldDiff, _NewX, _NewDiff) :-
183	free(OldX), !.
184better(OldX, OldDiff, NewX, NewDiff) :-
185	( get_var_bounds(OldX, _, 1.0) ->
186	    ( get_var_bounds(NewX, _, 1.0) ->
187		NewDiff =< OldDiff
188	    ;
189		true
190	    )
191	;
192	    ( get_var_bounds(NewX, _, 1.0) ->
193		fail
194	    ;
195		NewDiff =< OldDiff
196	    )
197	).
198
199
200% ----------------------------------------------------------------------
201% end_module(mip).
202% ----------------------------------------------------------------------
203
204