(!******************************************************
Mosel Example Problems
======================
file f5tour3.mos
````````````````
Tour planning - traveling salesman problem (TSP)
(See "Applications of optimization with Xpress-MP",
Section 9.6: Planning a flight tour)
Combination with CP as start solution does not
work with the iterative MIP algorithm in 'f5tour.mos',
the MIP model used here therefore contains the complete
problem definition.
ALG 1+2: use CP with local-search style enumeration
to generate a start solution for the MIP search
For small instances (up to approx. 50 cities) pure MIP tends
to beat the combined approach. Beyond approx. 120 cities
the problem grows too large for the formulations used here.
(c) 2008 Fair Isaac Corporation
author: S. Heipcke, Dec. 2007, rev. Nov. 2017
*******************************************************!)
model "F-5 Tour planning (MIP + CP)"
uses "mmxprs", "mmsystem", "kalis"
parameters
DATAFILE="Data/gr96.dat"
ALG=2 ! 0: pure MIP (no start solution)
! 1: CP heuristic for initial solution
! 2: LP as initial solution
end-parameters
forward procedure print_CP_sol
forward procedure print_MIP_sol
declarations
starttime: real
CITIES: range ! Set of cities
end-declarations
starttime:=gettime
initializations from DATAFILE
CITIES
end-initializations
finalize(CITIES)
CITIES1:= CITIES-{min(i in CITIES) i}
! **** MIP model ****
declarations
NCITIES=getsize(CITIES)
DIST: array(CITIES,CITIES) of integer ! Distance between cities
NEXTC: array(CITIES) of integer ! Next city after i in the solution
fly: array(CITIES,CITIES) of mpvar ! 1 if flight from i to j
y: array(CITIES) of mpvar ! Variables for excluding subtours
end-declarations
initializations from DATAFILE
DIST
end-initializations
! Objective: total distance
TotalDist:= sum(i,j in CITIES | i<>j) DIST(i,j)*fly(i,j)
! Visit every city once
forall(i in CITIES) sum(j in CITIES | i<>j) fly(i,j) = 1
forall(j in CITIES) sum(i in CITIES | i<>j) fly(i,j) = 1
forall(i,j in CITIES | i<>j) fly(i,j) is_binary
! Exclude subtours
forall(i in CITIES, j in CITIES1 | i<>j) y(j) >= y(i) + 1 - NCITIES * (1 - fly(i,j))
! **** CP model ****
setparam("KALIS_DEFAULT_LB", 0)
setparam("KALIS_DEFAULT_UB", NCITIES-1)
declarations
D: array(CITIES) of integer ! Auxiliary array
flyc: array(CITIES) of cpvar ! Successor of a city
end-declarations
MAXD:= sum(i in CITIES) max(j in CITIES) DIST(i,j)
setparam("KALIS_DEFAULT_UB", MAXD)
declarations
totaldist: cpvar ! Total travel distance (objective)
end-declarations
MAXD:= max(i,j in CITIES) DIST(i,j)
setparam("KALIS_DEFAULT_UB", MAXD)
declarations
dist: array(CITIES) of cpvar ! Distance to successor of city i
fsol: array(CITIES) of integer ! CP solution values
end-declarations
! Assign values to 'flyc' variables as to obtain a single cycle,
! 'totaldist' is the distance travelled
cycle(flyc, totaldist, DIST)
! Element constraints, defining 'dist' variables used in search
! (redundant constraints)
forall(i in CITIES) do
forall(j in CITIES) D(j):= DIST(i,j)
dist(i) = element(D, flyc(i))
end-do
! **** Solving: CP start solution ****
declarations
CORDER: array(1..NCITIES) of integer ! Ordered city indices
dsol: array(CITIES) of real ! LP solution values
allsol: array(mpvar) of real ! MIP start solution
sdist: integer ! CP objective function value
end-declarations
cp_set_solution_callback("print_CP_sol")
setparam("KALIS_MAX_COMPUTATION_TIME", 60)
! Find a first solution
case ALG of
1: do ! Initial solution generated by CP:
cp_set_branching(assign_var(KALIS_LARGEST_MAX, KALIS_MIN_TO_MAX, dist))
if not (cp_find_next_sol) then
writeln("No initial solution found")
cpsol:=false
end-if
forall(i in CITIES) settarget(flyc(i), getsol(flyc(i)))
forall(i in CITIES) settarget(dist(i), getsol(dist(i)))
forall(i in CITIES) fsol(i):=flyc(i).sol
cpsol:=true
sdist:= getsol(totaldist)
cp_reset_search
totaldist <= sdist - maxlist(sdist * 0.001,1)
cp_set_branching(assign_var(KALIS_LARGEST_MAX, KALIS_NEAREST_VALUE, dist))
setparam("KALIS_MAX_COMPUTATION_TIME", 10)
end-do
2: do ! Using LP results as target values for CP:
minimize(XPRS_LIN,TotalDist) ! Solve the LP-relaxation
forall(i in CITIES) do
minj:=i; minv:=0.0
forall(j in CITIES | i<>j)
if getsol(fly(i,j))>minv then
minj:=j; minv:=getsol(fly(i,j))
end-if
dsol(i):= minv
settarget(flyc(i), minj)
settarget(dist(i),
round(sum(j in CITIES | i<>j) DIST(i,j)*getsol(fly(i,j))) )
end-do
qsort(SYS_DOWN, dsol, CORDER)
forall(i in CITIES)
Strategy(i):= assign_var(KALIS_LARGEST_MAX, KALIS_NEAREST_VALUE,
{dist(CORDER(i+1))})
cp_set_branching(Strategy)
end-do
end-case
if ALG<>0 then
! Solve the CP problem: minimize total distance using nearest-neighbor search
! on distance variables
cpstart:=gettime
while (gettime-cpstart<5 and cp_find_next_sol) do
cpsol:=true
sdist:= getsol(totaldist)
forall(i in CITIES) settarget(dist(i),getsol(dist(i)))
forall(i in CITIES) fsol(i):=flyc(i).sol
cp_reset_search
totaldist <= sdist - maxlist(sdist * 0.001,1)
setparam("KALIS_MAX_COMPUTATION_TIME", 10)
end-do
if cpsol then
loadprob(TotalDist) ! Load best CP solution as MIP start sol.
forall(i in CITIES) allsol(fly(i,fsol(i))):=1
res:=loadmipsol(allsol)
else
writeln("(", gettime-starttime, "s) No CP start solution found.")
end-if
end-if
! **** Solving: MIP ****
setparam("XPRS_VERBOSE", true)
setparam("XPRS_MIPLOG", -500)
setparam("XPRS_GOMCUTS",0)
setparam("XPRS_HEURSEARCHEFFORT", 2) ! Increased local search effort
setparam("XPRS_HEURSEARCHROOTSELECT", 7)
setparam("XPRS_HEURSEARCHTREESELECT", 3)
if ALG<>0 then
setparam("XPRS_MAXNODE", 1000)
else
setparam("XPRS_MAXNODE", 1500)
end-if
minimize(TotalDist) ! Solve the MIP problem
if getparam("XPRS_MIPSOLS") > 0 then
print_MIP_sol
else
writeln("(", gettime-starttime, "s) No MIP solutions.")
end-if
!-----------------------------------------------------------------
! Print the current CP solution
procedure print_CP_sol
declarations
ALLCITIES: set of integer
end-declarations
writeln("(", gettime-starttime, "s) Total distance: ", getsol(totaldist))
ALLCITIES:={}
forall(i in CITIES) do
if(i not in ALLCITIES) then
write(i)
first:=i
repeat
ALLCITIES+={first}
write(" - ", flyc(first).sol)
first:= flyc(first).sol
until first=i
writeln
end-if
end-do
end-procedure
! Print the current MIP solution
procedure print_MIP_sol
declarations
ALLCITIES: set of integer
end-declarations
writeln("(", gettime-starttime, "s) Total distance: ", getobjval)
forall(i in CITIES)
forall(j in CITIES)
if(getsol(fly(i,j))>0.99) then
NEXTC(i):=j
break
end-if
ALLCITIES:={}
forall(i in CITIES) do
if(i not in ALLCITIES) then
write(i)
first:=i
repeat
ALLCITIES+={first}
write(" - ", NEXTC(first))
first:=NEXTC(first)
until first=i
writeln
end-if
end-do
end-procedure
end-model
|
(!******************************************************
Mosel Example Problems
======================
file f5tour4.mos
````````````````
Tour planning - traveling salesman problem (TSP)
(See "Applications of optimization with Xpress-MP",
Section 9.6: Planning a flight tour)
Generate and load heuristic CP solutions during MIP search
using the 'optimized node' and 'integer solution' callbacks.
For small instances (up to approx. 50 cities) pure MIP tends
to beat the combined approach. Beyond approx. 120 cities
the problem grows too large for the formulations used here.
(c) 2008 Fair Isaac Corporation
author: S. Heipcke, Dec. 2007, rev. Nov. 2017
*******************************************************!)
model "F-5 Tour planning (MIP + CP)"
uses "mmxprs", "mmsystem", "kalis"
parameters
DATAFILE="Data/gr96.dat"
end-parameters
forward public procedure update_CP
forward public function solve_CP: boolean
forward procedure print_CP_sol
forward procedure print_MIP_sol
declarations
starttime: real
CITIES: range ! Set of cities
end-declarations
starttime:=gettime
initializations from DATAFILE
CITIES
end-initializations
finalize(CITIES)
CITIES1:= CITIES-{min(i in CITIES) i}
! **** MIP model ****
declarations
NCITIES=getsize(CITIES)
DIST: array(CITIES,CITIES) of integer ! Distance between cities
NEXTC: array(CITIES) of integer ! Next city after i in the solution
fly: array(CITIES,CITIES) of mpvar ! 1 if flight from i to j
y: array(CITIES) of mpvar ! Variables for excluding subtours
end-declarations
initializations from DATAFILE
DIST
end-initializations
! Objective: total distance
TotalDist:= sum(i,j in CITIES | i<>j) DIST(i,j)*fly(i,j)
! Visit every city once
forall(i in CITIES) sum(j in CITIES | i<>j) fly(i,j) = 1
forall(j in CITIES) sum(i in CITIES | i<>j) fly(i,j) = 1
forall(i,j in CITIES | i<>j) fly(i,j) is_binary
! Exclude subtours
forall(i in CITIES, j in CITIES1 | i<>j) y(j) >= y(i) + 1 - NCITIES * (1 - fly(i,j))
! **** CP model ****
setparam("KALIS_DEFAULT_LB", 0)
setparam("KALIS_DEFAULT_UB", NCITIES-1)
declarations
D: array(CITIES) of integer ! Auxiliary array
flyc: array(CITIES) of cpvar ! Successor of a city
end-declarations
MAXD:= sum(i in CITIES) max(j in CITIES) DIST(i,j)
setparam("KALIS_DEFAULT_UB", MAXD)
declarations
totaldist: cpvar ! Total travel distance (objective)
end-declarations
MAXD:= max(i,j in CITIES) DIST(i,j)
setparam("KALIS_DEFAULT_UB", MAXD)
declarations
dist: array(CITIES) of cpvar ! Distance to successor of city i
fsol: array(CITIES) of integer ! CP solution values
end-declarations
! Assign values to 'flyc' variables as to obtain a single cycle,
! 'totaldist' is the distance travelled
cycle(flyc, totaldist, DIST)
! Element constraints, defining 'dist' variables used in search
! (redundant constraints)
forall(i in CITIES) do
forall(j in CITIES) D(j):= DIST(i,j)
dist(i) = element(D, flyc(i))
end-do
! **** Solving ****
declarations
CORDER: array(1..NCITIES) of integer ! Ordered city indices
dsol: array(CITIES) of real ! LP solution values
allsol: array(mpvar) of real ! MIP start solution
sdist: integer ! CP objective function value
end-declarations
! cp_set_solution_callback("print_CP_sol")
setparam("XPRS_VERBOSE", true)
setparam("XPRS_MIPLOG", -25)
setparam("XPRS_GOMCUTS",0)
setparam("XPRS_HEURSEARCHEFFORT", 2) ! Increased local search effort
setparam("XPRS_HEURSEARCHROOTSELECT", 0)
setparam("XPRS_HEURSEARCHTREESELECT", 7)
setparam("XPRS_MAXNODE", 300)
setcallback(XPRS_CB_OPTNODE, "solve_CP")
setcallback(XPRS_CB_INTSOL, "update_CP")
!setparam("XPRS_HEURSTRATEGY", 0)
setparam("XPRS_NODESELECTION", 2) ! Breadth-first search
setparam("XPRS_MIPTHREADS", 1) ! Desactivate parallel MIP
minimize(TotalDist) ! Solve the MIP problem
if getparam("XPRS_MIPSOLS") > 0 then
print_MIP_sol
else
writeln("(", gettime-starttime, "s) No MIP solutions.")
end-if
!-----------------------------------------------------------------
! Update CP with solutions found by MIP
public procedure update_CP
totaldist <= getparam("XPRS_LPOBJVAL") - 1
end-procedure
! Generate and load CP solutions into MIP
public function solve_CP:boolean
if getparam("XPRS_NODES")=1 or getparam("XPRS_NODEDEPTH") MOD 4 = 0 then
writeln("Solving CP")
forall(i in CITIES) do
minj:=i; minv:=0.0
forall(j in CITIES | i<>j)
if getsol(fly(i,j))>minv then
minj:=j; minv:=getsol(fly(i,j))
end-if
dsol(i):= minv
settarget(flyc(i), minj)
settarget(dist(i),
round(sum(j in CITIES | i<>j) DIST(i,j)*getsol(fly(i,j))) )
end-do
qsort(SYS_DOWN, dsol, CORDER)
forall(i in CITIES)
Strategy(i):= assign_var(KALIS_LARGEST_MAX, KALIS_NEAREST_VALUE,
{dist(CORDER(i+1))})
cp_set_branching(Strategy)
setparam("KALIS_MAX_COMPUTATION_TIME", ceil(NCITIES/10))
setparam("KALIS_MAX_SOLUTIONS", 1)
if cp_find_next_sol then
writeln("(", gettime-starttime, "s) CP solution: ", getsol(totaldist))
sdist:= getsol(totaldist)
forall(i in CITIES) fsol(i):=flyc(i).sol
forall(i in CITIES) allsol(fly(i,fsol(i))):=1
cp_reset_search
res:=loadmipsol(allsol)
forall(i in CITIES) allsol(fly(i,fsol(i))):=0
totaldist <= sdist - maxlist(sdist * 0.0001,1)
end-if
end-if
end-function
!-----------------------------------------------------------------
! Print the current CP solution
procedure print_CP_sol
declarations
ALLCITIES: set of integer
end-declarations
writeln("(", gettime-starttime, "s) Total distance: ", getsol(totaldist))
ALLCITIES:={}
forall(i in CITIES) do
if(i not in ALLCITIES) then
write(i)
first:=i
repeat
ALLCITIES+={first}
write(" - ", flyc(first).sol)
first:= flyc(first).sol
until first=i
writeln
end-if
end-do
end-procedure
! Print the current MIP solution
procedure print_MIP_sol
declarations
ALLCITIES: set of integer
end-declarations
writeln("(", gettime-starttime, "s) Total distance: ", getobjval)
forall(i in CITIES)
forall(j in CITIES)
if(getsol(fly(i,j))>0.99) then
NEXTC(i):=j
break
end-if
ALLCITIES:={}
forall(i in CITIES) do
if(i not in ALLCITIES) then
write(i)
first:=i
repeat
ALLCITIES+={first}
write(" - ", NEXTC(first))
first:=NEXTC(first)
until first=i
writeln
end-if
end-do
end-procedure
end-model
|
