proc One_City_2_Opt_ND(tour: var Tour_Array;
basePos: Tour_Index;
neighbor: Neighbor_Matrix;
position: var City_Position_In_Tour;
numberOfNeigbors: int;
DontLook: var DLB_Array): bool =
var
improved: bool
pos_Y2, i, j: Tour_Index
X1, X2, Y1, Y2: City_Number
improved = false
block two_loops:
for direction in [forward, backward]:
if direction == forward:
i = basePos
X1 = tour[i]
X2 = tour[(i+1) mod N]
else:
i = (N + basePos - 1) mod N
X2 = tour[i]
X1 = tour[(i+1) mod N] # == tour[basePos]
for neighbor_number in 0 .. numberOfNeigbors-1:
Y1 = neighbor[X1][neighbor_number]
if direction == forward:
j = position[Y1]
Y2 = tour[(j+1) mod N]
else:
j = (N + position[Y1] - 1) mod N # pos[Y1] == j+1
Y2 = tour[j]
if (X2 == Y1) or (Y2 == X1):
continue
if Gain_From_2_Opt(X1, X2, Y1, Y2) > 0:
Set_DLB_off(DontLook, [X1, X2, Y1, Y2])
Make_2_Opt_Move(tour, i, j, position)
improved = true
break two_loops
result = improved

proc LS_2_Opt_ND(tour: var Tour_Array;
neighbor: Neighbor_Matrix;
neighborListLen: int = N-1) =
## Optimizes the given tour using 2-opt with neighbor lists
## and dDon't Look Bits (DLB)
var
locallyOptimal: bool = false
improved: bool
DontLook: DLB_Array
baseCity: City_Number
position: City_Position_In_Tour
let numberOfNeigbors = min(neighborListLen, len(neighbor[0]))
position = Create_Position_In_Tour(tour)
Set_DLB_off(DontLook, tour)
while not locallyOptimal:
locallyOptimal = true
for basePos in 0 .. N-1:
baseCity = tour[basePos]
if isDLB_on(DontLook, baseCity):
continue
improved = One_City_2_Opt_ND(tour, basePos, neighbor, position,
numberOfNeigbors, DontLook)
if not improved:
Set_DLB_on(DontLook, baseCity)
else:
locallyOptimal = false

proc One_City_2_Opt_NR(tour: var Tour_Array;
basePos: Tour_Index;
neighbor: Neighbor_Matrix;
position: var City_Position_In_Tour;
numberOfNeigbors: int): bool =
var
improved: bool
pos_Y2, i, j: Tour_Index
X1, X2, Y1, Y2: City_Number
radius: Length
improved = false
block two_loops:
for direction in [forward, backward]:
if direction == forward:
i = basePos
X1 = tour[i]
X2 = tour[(i+1) mod N]
else:
i = (N + basePos - 1) mod N
X2 = tour[i]
X1 = tour[(i+1) mod N] # == tour[basePos]
radius = distance(X1, X2)
for neighbor_number in 0 .. numberOfNeigbors-1:
Y1 = neighbor[X1][neighbor_number]
if direction == forward:
j = position[Y1]
Y2 = tour[(j+1) mod N]
else:
j = (N + position[Y1] - 1) mod N # pos[Y1] == j+1
Y2 = tour[j]
if (X2 == Y1) or (Y2 == X1):
continue
if distance(X1, Y1) > radius:
break
if Gain_From_2_Opt(X1, X2, Y1, Y2) > 0:
Make_2_Opt_Move(tour, i, j, position)
improved = true
break two_loops
result = improved

proc LS_2_Opt_NR(tour: var Tour_Array;
neighbor: Neighbor_Matrix;
neighborListLen: int = N-1) =
## Optimizes tour using 2-opt with neighbor lists and fixed radius
var
locallyOptimal: bool = false
improved: bool
position: City_Position_In_Tour
let numberOfNeigbors = min(neighborListLen, len(neighbor[0]))
position = Create_Position_In_Tour(tour)
while not locallyOptimal:
locallyOptimal = true
for basePos in 0 .. N-1:
improved = One_City_2_Opt_NR(tour, basePos, neighbor,
position, numberOfNeigbors)
if improved:
locallyOptimal = false

proc One_City_2_Opt_N(tour: var Tour_Array;
basePos: Tour_Index;
neighbor: Neighbor_Matrix;
position: var City_Position_In_Tour;
numberOfNeigbors: int): bool =
## Optimizes the given tour using 2-opt with neighours lists
# For each edge looks for improvement considering neigbours
# in ascending order of their distance
var
improved: bool
pos_Y2, i, j: Tour_Index
X1, X2, Y1, Y2: City_Number
improved = false
block two_loops:
for direction in [forward, backward]:
if direction == forward:
i = basePos
X1 = tour[i]
X2 = tour[(i+1) mod N]
else:
i = (N + basePos - 1) mod N
X2 = tour[i]
X1 = tour[(i+1) mod N] # == tour[basePos]
for neighbor_number in 0 .. numberOfNeigbors-1:
# after 2-opt move new tour would contain direct link X1-Y1,
# so we look for city Y1 among cities that are close to city X1
Y1 = neighbor[X1][neighbor_number]
if direction == forward:
j = position[Y1]
Y2 = tour[(j+1) mod N]
else:
j = (N + position[Y1] - 1) mod N # pos[Y1] == j+1
Y2 = tour[j]
if (X2 == Y1) or (Y2 == X1):
# 2-opt is for non-adjacent links only
# by constuction X1 != Y1
continue
if Gain_From_2_Opt(X1, X2, Y1, Y2) > 0:
Make_2_Opt_Move(tour, i, j, position)
improved = true
break two_loops
result = improved

proc LS_2_Opt_N(tour: var Tour_Array;
neighbor: Neighbor_Matrix;
neighborListLen: int = N-1) =
## Optimizes the given tour using 2-opt with neighbor lists
var
locallyOptimal: bool = false
improved: bool
position: City_Position_In_Tour
let numberOfNeigbors = min(neighborListLen, len(neighbor[0]))
position = Create_Position_In_Tour(tour)
while not locallyOptimal:
locallyOptimal = true
for basePos in 0 .. N-1:
improved = One_City_2_Opt_N(tour, basePos, neighbor,
position, numberOfNeigbors)
if improved:
locallyOptimal = false

Some simple statistics

Relative time, best length, average length and worst length for all tours created by NN heuristic and improved by 2-opt and 2-opt with neighbor lists (number of neighbors in parentheses). Random planar problems, N=100.

Using neighbor lists forces us to use some convenient and quick method to determine for given city (a neighbor) what is its position in current tour. We use an array indexed by city numbers and storing corresponding positions.

type
City_Position_In_Tour = array[City_Number, Tour_Index]

proc Create_Position_In_Tour(tour: Tour_Array): City_Position_In_Tour =
## Creates array with position of each city in tour
var
position: City_Position_In_Tour
for pos in 0 .. N-1:
position[tour[pos]] = pos
result = position

Note that array of positions should be updated during any action that changes the tour, be it reversing a segment:

proc Reverse_Segment(tour: var Tour_Array;
startIndex, endIndex: Tour_Index;
position: var City_Position_In_Tour) =
## Reverses order of elements in segment [startIndex, endIndex]
## of tour, and updates positon[] array used with neighbour lists
var
left, right: Tour_Index
let inversionSize = ((endIndex - startIndex + 1 + N) mod N) div 2
left = startIndex
right = endIndex
for counter in 1 .. inversionSize:
swap(tour[left], tour[right])
position[tour[left]] = left
position[tour[right]] = right
left = (left + 1) mod N
right = (N + right - 1) mod N

proc Shift_Segment(tour: var Tour_Array; i, j, k: Tour_Index;
position: var City_Position_In_Tour) =
## Shifts the segment of tour:
## cities from t[i+1] to t[j] from their current position to position
## after current city t[k], that is between cities t[k] and t[k+1].
## Accordingly updates positon[] array.
# Assumes: k, k+1 are not within the segment [i+1..j]
let
segmentSize = (j - i + N) mod N
shiftSize = ((k - i + N) - segmentSize + N) mod N
offset = i + 1 + shiftSize
var
pos, posNew: Tour_Index
segment: seq[City_Number] = newSeq[City_Number](segmentSize)
# make a copy of the segment before shift
for pos in 0 .. segmentSize-1:
segment[pos] = tour[(pos + i+1) mod N]
# shift to the left by segmentSize all cities between old position
# of right end of the segment and new position of its left end
pos = (i + 1) mod N
for counter in 1 .. shiftSize:
tour[pos] = tour[(pos + segmentSize) mod N]
position[tour[pos]] = pos
pos = (pos + 1) mod N
#end_for counter loop
# put the copy of the segment into its new place in the tour
for pos in 0 .. segmentSize-1:
posNew = (pos + offset) mod N
tour[posNew] = segment[pos]
position[tour[posNew]] = posNew

Note that neighbor lists are used to speed up the process of optimization, and using them may cost slightly worse quality of resulting tours, even if it not always occurs.

type
Neighbor_Matrix = array [City_Number, seq[City_Number]]

proc Build_Neighbors_Matrix(No_Of_Neigbors: int = N-1): Neighbor_Matrix =
## Builds array of No_Of_Neigbors nearest neighbors of each city
# neighbor[i][j] is j-th nearest neighbour of city number i.
var
currentCity, otherCity: City_Number
idx: int
neighbors_of_currentCity: seq[City_Number]
neighbor: Neighbor_Matrix
newSeq(neighbors_of_currentCity, N-1)
for currentCity in 0 .. N-1:
# create sequence of all neighbors of current city:
idx = 0
for otherCity in 0 .. N-1:
if otherCity != currentCity:
neighbors_of_currentCity[idx] = otherCity
idx = idx + 1
# Sort elements in neighbors_of_currentCity by ascending distance
# to currentCity (use any convenient method of sorting), then take
# first No_Of_Neigbors cities from the result.
sort(neighbors_of_currentCity) do (x,y: City_Number) -> int:
result = cmp(distance(currentCity, x), distance(currentCity, y))
neighbor[currentCity] = @[]
setlen(neighbor[currentCity], No_Of_Neigbors)
for idx in 0 .. No_Of_Neigbors-1:
neighbor[currentCity][idx] = neighbors_of_currentCity[idx]
result = neighbor

holds if at least one of the two conditions is met:

or

The inequity (a) means that for given (X1, X2) we can search for Y2 among those cities that are closer to X2 than is X1.

But alternatively the inequity between the two sums may be considered as a pair of inequities:

or

The inequity (c) means that for given (X1, X2) we can search for Y1 among those cities that are closer to X1 than is X2.

In fixed radius search we can use either (a) or (c), the radius is the same. When fixed radius searching is used with DLB, then probably searching around X1 for Y1 would be better – because we would rather prefer to not miss a good move for base city (X1) before activating its don’t look bit. When we want to continue searching and moving in given direction, then probably searching for Y2 around X2 could be better idea.

proc One_City_2_Opt_DLBR2(tour: var Tour_Array,
basePos: Tour_Index,
DontLook: var DLB_Array): bool =
var
improved: bool
i, j: Tour_Index
X1, X2, Y1, Y2: City_Number
radius: Length
improved = false
block two_loops:
for direction in [forward, backward]:
if direction == forward:
i = basePos
X1 = tour[i]
X2 = tour[(i+1) mod N]
else:
i = (N + basePos - 1) mod N
X2 = tour[i]
X1 = tour[basePos] # == tour[(i+1) mod N]
radius = distance(X1, X2)
for counter_2 in 0 .. N-1:
j = counter_2
Y1 = tour[j]
Y2 = tour[(j+1) mod N]
if direction == backward:
swap(Y1, Y2)
if (X1 == Y1) or (X2 == Y1) or (Y2 == X1):
continue
if distance(X1, Y1) > radius: # NOTE
continue
if Gain_From_2_Opt(X1, X2, Y1, Y2) > 0:
Set_DLB_off(DontLook, [X1, X2, Y1, Y2])
Make_2_Opt_Move(tour, i, j)
improved = true
break two_loops
result = improved

This example of code does not contain any new ideas and can be skipped by a reader.

proc One_City_2_Opt_DLBR(tour: var Tour_Array,
basePos: Tour_Index,
DontLook: var DLB_Array): bool =
var
improved: bool
i, j: Tour_Index
X1, X2, Y1, Y2: City_Number
radius: Length
improved = false
block two_loops:
for direction in [forward, backward]:
if direction == forward:
i = basePos
X1 = tour[i]
X2 = tour[(i+1) mod N]
else:
i = (N + basePos - 1) mod N
X2 = tour[i]
X1 = tour[basePos] # == tour[(i+1) mod N]
radius = distance(X1, X2)
for counter_2 in 0 .. N-1:
j = counter_2
Y1 = tour[j]
Y2 = tour[(j+1) mod N]
if direction == backward:
swap(Y1, Y2)
if (X1 == Y1) or (X2 == Y1) or (Y2 == X1):
continue
if distance(X2, Y2) > radius:
continue
if Gain_From_2_Opt(X1, X2, Y1, Y2) > 0:
Set_DLB_off(DontLook, [X1, X2, Y1, Y2])
Make_2_Opt_Move(tour, i, j)
improved = true
break two_loops
result = improved

proc LS_2_Opt_DLBR(tour: var Tour_Array) =
## Optimizes the given tour using 2-opt with Don't Look Bits (DLB)
## and fixed radius search
var
locallyOptimal: bool = false
improved: bool
DontLook: DLB_Array
baseCity: City_Number
Set_DLB_off(DontLook, tour)
while not locallyOptimal:
locallyOptimal = true
for basePos in 0 .. N-1:
baseCity = tour[basePos]
if isDLB_on(DontLook, baseCity):
continue
improved = One_City_2_Opt_DLBR(tour, basePos, DontLook)
if not improved:
Set_DLB_on(DontLook, baseCity)
else:
locallyOptimal = false

Fixed radius search uses the property of the inequity which a 2-opt move must meet to improve a tour.
We can benefit from this even more, but only when we consider links in the two orientations and the two symmetric forms.

forward (ABCD)

backward (DCBA)

implementation

look forward

look backward

Using this scheme we cannot miss an improving move even if we restrict possible choices to one inequity.

Example implementation:

type
Orientation = enum forward, backward

proc One_City_2_Opt_R(tour: var Tour_Array,
basePos: Tour_Index): bool =
var
improved: bool
i, j: Tour_Index
X1, X2, Y1, Y2: City_Number
radius: Length
improved = false
block two_loops:
for direction in [forward, backward]:
if direction == forward:
i = basePos
X1 = tour[i]
X2 = tour[(i+1) mod N]
else:
i = (N + basePos - 1) mod N
X2 = tour[i]
X1 = tour[basePos] # == tour[(i+1) mod N]
radius = distance(X1, X2)
for counter_2 in 0 .. N-1:
j = counter_2
Y1 = tour[j]
Y2 = tour[(j+1) mod N]
if direction == backward:
swap(Y1, Y2)
if (X1 == Y1) or (X2 == Y1) or (Y2 == X1):
continue
if distance(X2, Y2) > radius:
continue
if Gain_From_2_Opt(X1, X2, Y1, Y2) > 0:
Make_2_Opt_Move(tour, i, j)
improved = true
break two_loops
result = improved

proc LS_2_Opt_R(tour: var Tour_Array) =
## Optimizes the given tour using 2-opt with fixed radius search
var
locallyOptimal: bool = false
improved: bool
while not locallyOptimal:
locallyOptimal = true
for basePos in 0 .. N-1:
improved = One_City_2_Opt_R(tour, basePos)
if improved:
locallyOptimal = false

holds if at least one of the two conditions is met:

or

Simple radius search with 2-opt

Note that for simplicity of the code below we do not change orientation and instead use both conditions.

proc LS_2_Opt_simpleRadius(tour: var Tour_Array) =
var
locallyOptimal: bool = false
i, j: Tour_Index
X1, X2, Y1, Y2: City_Number
radius: Length
while not locallyOptimal:
locallyOptimal = true
block two_loops:
for counter_1 in 0 .. N-1:
i = counter_1
X1 = tour[i]
X2 = tour[(i+1) mod N]
radius = distance(X1, X2)
for counter_2 in 0 .. N-1:
j = counter_2
Y1 = tour[j]
Y2 = tour[(j+1) mod N]
if (X1 == Y1) or (X2 == Y1) or (Y2 == X1):
continue
if distance(X2, Y2) > radius:
if distance(X1, Y1) > radius:
continue
if Gain_From_2_Opt(X1, X2, Y1, Y2) > 0:
Make_2_Opt_Move(tour, i, j)
locallyOptimal = false
break two_loops

Some simple statistics

Relative time, best length, average length and worst length for all tours created by NN heuristic and improved by 2-opt and 2-opt with simple radius search algorithms. Random planar problems, N=100.

Procedures to calculate the gain from 3-opt move and to make 3-opt move assume that their parameters (cities and their indices in tour) are given in the order of ascending indices of tour array (modulo N). Therefore we need a function to recognize whether three indices are in this order, or: whether in this orientation the middle index is between the two boundary indices.

proc Between(a, x, b: Tour_Index): bool =
## Returns true if `x` is between `a` and `b` in cyclic sequence
## with established orientation of nondescending indices
# That is: when one begins a forward traversal of the tour
# at position `a', then position `x` is reached before position `b'.
# Returns true if and only if:
# a < x < b or b < a < x or x < b < a
if (b > a):
result = (x > a) and (x < b)
elif (b < a):
result = (x > a) or (x < b)
else:
result = false

Note the additional speed-up: the link between base city X1 and its tour neighbor X2 is not considered when at least one of these cities has its don’t-look bit set to on.

proc One_City_3_Opt_DLB(tour: var Tour_Array,
basePos: Tour_Index,
DontLook: var DLB_Array): bool =
type
Move3 = object
i, j, k: Tour_Index
optCase: Reconnection_3_optCase
var
improved: bool
i, j, k: Tour_Index
X1, X2, Y1, Y2, Z1, Z2: City_Number
gainExpected: Length_Gain
optCase: Reconnection_3_optCase
goodMove: Move3
improved = false
block three_loops:
for direction in [forward, backward]:
if direction == forward:
i = basePos
else:
i = (N + basePos - 1) mod N
X1 = tour[i]
X2 = tour[(i+1) mod N]
if isDLB_on(DontLook,X1) or isDLB_on(DontLook,X2): # NOTE
continue
for counter_2 in 0 .. N-1:
j = counter_2 # second cut after j
Y1 = tour[j]
Y2 = tour[(j+1) mod N]
if (X1 == Y1) or (X2 == Y1) or (Y2 == X1):
continue
# examine 2-opt move (see: opt3_case_1 .. opt3_case_3)
gainExpected = Gain_From_2_Opt(X1, X2, Y1, Y2)
if gainExpected > 0:
goodMove = Move3(i: i, j: j, k: j,
optCase: opt3_case_1)
improved = true
break three_loops
for counter_3 in 0 .. N-1:
k = counter_3 # third cut after k
Z1 = tour[k]
Z2 = tour[(k+1) mod N]
if (X1 == Z1) or (Y1 == Z1):
continue
if not Between(i, j, k):
continue
# examine pure 3-opt cases
for optCase in [opt3_case_6, opt3_case_7]:
gainExpected = Gain_From_3_Opt(X1, X2, Y1, Y2, Z1, Z2,
optCase)
if gainExpected > 0:
goodMove = Move3(i: i, j: j, k: k,
optCase: optCase)
improved = true
break three_loops
# end_block three_loops
if improved:
Set_DLB_off(DontLook, [X1, X2, Y1, Y2])
if goodMove.optCase > opt3_case_3:
Set_DLB_off(DontLook, [Z1, Z2])
Make_3_Opt_Move(tour,
goodMove.i, goodMove.j, goodMove.k,
goodMove.optCase)
result = improved

proc LS_3_Opt_DLB(tour: var Tour_Array) =
## Optimizes the given tour using 3-opt with Don't Look Bits (DLB)
var
locallyOptimal: bool = false
improved: bool
DontLook: DLB_Array
baseCity: City_Number
Set_DLB_off(DontLook, tour)
while not locallyOptimal:
locallyOptimal = true
for basePos in 0 .. N-1:
baseCity = tour[basePos]
if isDLB_on(DontLook, baseCity):
continue
improved = One_City_3_Opt_DLB(tour, basePos, DontLook)
if not improved:
Set_DLB_on(DontLook, baseCity)
else:
locallyOptimal = false

Some simple statistics

Relative time, best length, average length and worst length for all tours created by NN heuristic and improved by 3-opt and 3-opt with DLB algorithms. Random planar problems, N=100.

proc One_City_25_Opt_DLB(tour: var Tour_Array,
basePos: Tour_Index,
DontLook: var DLB_Array): bool =
var
improved: bool
i, j: Tour_Index
X1, X2, Y1, Y2, V0: City_Number
improved = false
block two_loops:
for direction in [forward, backward]: # NOTE: both tour neighbors
# consider link between t[i] and t[i+1],
# then between t[i-1] and t[i]
if direction == forward:
i = basePos
else:
i = (N + basePos - 1) mod N
X1 = tour[i]
X2 = tour[(i+1) mod N]
for counter_2 in 0 .. N-1:
j = counter_2
Y1 = tour[j]
Y2 = tour[(j+1) mod N]
if (X1 == Y1) or (X2 == Y1) or (Y2 == X1):
continue
if Gain_From_2_Opt(X1, X2, Y1, Y2) > 0:
Make_2_Opt_Move(tour, i, j)
improved = true
else:
V0 = tour[(i+2) mod N]
if V0 != Y1:
if Gain_From_Node_Shift(X1, X2, V0, Y1, Y2) > 0:
Set_DLB_off(DontLook, V0)
Make_Node_Shift_Move(tour, (i+1) mod N, j)
improved = true
else:
V0 = tour[(N+j-1) mod N]
if V0 != X2:
if Gain_From_Node_Shift(V0, Y1, Y2, X1, X2) > 0:
Set_DLB_off(DontLook, V0)
Make_Node_Shift_Move(tour, j, i)
improved = true
#end_if Gain_From_2_Opt...
if improved:
Set_DLB_off(DontLook, [X1, X2, Y1, Y2])
break two_loops
result = improved

proc LS_25_Opt_DLB(tour: var Tour_Array) =
## Optimizes the given tour using 2.5-opt with Don't Look Bits (DLB)
var
locallyOptimal: bool = false
improved: bool
DontLook: DLB_Array
baseCity: City_Number
Set_DLB_off(DontLook, tour)
while not locallyOptimal:
locallyOptimal = true
for basePos in 0 .. N-1:
baseCity = tour[basePos]
if isDLB_on(DontLook, baseCity):
continue
improved = One_City_25_Opt_DLB(tour, basePos, DontLook)
if not improved:
Set_DLB_on(DontLook, baseCity)
else:
locallyOptimal = false

proc isDLB_on(DontLook: DLB_Array,
city: City_Number): bool =
return DontLook[city]
proc Set_DLB_on(DontLook: var DLB_Array,
city: City_Number) =
DontLook[city] = true
proc Set_DLB_off(DontLook: var DLB_Array,
city: City_Number) =
DontLook[city] = false
proc Set_DLB_off(DontLook: var DLB_Array,
cities: openArray[City_Number]) =
# turns off DLB for given cities
for city in 0 .. len(cities)-1:
DontLook[cities[city]] = false

proc One_City_2_Opt_DLB(tour: var Tour_Array,
basePos: Tour_Index,
DontLook: var DLB_Array): bool =
var
improved: bool
i, j: Tour_Index
X1, X2, Y1, Y2: City_Number
improved = false
block two_loops:
for direction in [forward, backward]: # NOTE: both tour neighbors
# consider link between t[i] and t[i+1],
# then between t[i-1] and t[i]
if direction == forward:
i = basePos
else:
i = (N + basePos - 1) mod N
X1 = tour[i]
X2 = tour[(i+1) mod N]
for counter_2 in 0 .. N-1:
j = counter_2
Y1 = tour[j]
Y2 = tour[(j+1) mod N]
if (X1 == Y1) or (X2 == Y1) or (Y2 == X1):
continue
if Gain_From_2_Opt(X1, X2, Y1, Y2) > 0:
Set_DLB_off(DontLook, [X1, X2, Y1, Y2])
Make_2_Opt_Move(tour, i, j)
improved = true
break two_loops
result = improved

proc LS_2_Opt_DLB(tour: var Tour_Array) =
## Optimizes the given tour using 2-opt with Don't Look Bits (DLB)
var
locallyOptimal: bool = false
improved: bool
DontLook: DLB_Array
baseCity: City_Number
Set_DLB_off(DontLook, tour)
while not locallyOptimal:
locallyOptimal = true
for basePos in 0 .. N-1:
baseCity = tour[basePos]
if isDLB_on(DontLook, baseCity):
continue
improved = One_City_2_Opt_DLB(tour, basePos, DontLook)
if not improved:
Set_DLB_on(DontLook, baseCity)
else:
locallyOptimal = false

A local improvement algorithm can be made faster when we notice that the running time includes the time to needed perform the move, but also the time needed to examine possible good moves and find an improving move (if there is one).

If for a given city X1 we previously did not find any improving move and it still has the same neighbors (none of the links between this city and its neighbors has been changed), then chances that we will find an improving move for this city in current iteration are small (although non-zero, see below).

The concept based on this observation is known as don’t look bits (DLB). We use special flag for each of the cities. Initially all the flags are turned off, which means that we allow searching for improving moves starting from given city. When a search for improving move starting from city C fails, then the bit for this city is turned on. When a move is performed in which this city is involved as one of the endpoints of the exchanged links, then the bit for this city is turned off. When we consider candidates for X1, the first city of an improving move, we skip all cities which have their don’t-look bits set to on.

Comments

Consider also what happens with DLB when we examine moves “in one direction” only (see 2-opt basic algorithn). In such case we consider only one of the two links of X1, that is the link to its successor, the link between tour[i] and tour[i+1]. If we failed to found any improving move, then we turn the don’t-look bit for X1 on. If then in some subsequent iteration a move is performed in which X1 is one of the endpoints of the exchanged links, then we turn the don’t-look bit off. It works as expected. But what if a move changes the orientation of the segment in which X1 is inside, but is not at one of its ends? The don’t-look bit remains unchanged, although now the previous predecessor of X1 is its current successor, and we examined X1 only with its succcessor to verify that there are no improvements for X1. So now there could be some improving move for X1 in “forward” orientation, but it will not be examined, because the don’t-look bit for this city is still on.

That is why in the code below we consider both tour neighbors of a given city.

For the same reason of not missing improving moves we do not limit searching to the cases where j>i, but cosider all possible pairs of (i, j).

Aside from the above, it may seem impossible that new, improving move for given city could appear when both neighbors of this city remain unchanged, but let us take a closer look:

However we should notice that in this case cities Y1 and Z1 have their don’t-look bits set to off after making the move. Therefore the new, possibly improving move: from X0 and involving the link between Y1 and Z1, will be examined in one of the subsequent iterations, that is as a move from Y1 or as a move from Z1.

Using Dont' Looks Bits with a simple 2-opt

type
DLB_Array = array[City_Number, bool]
Orientation = enum forward, backward

proc LS_2_Opt_DLB_General_Idea(tour: var Tour_Array) =
## Optimizes the given tour using 2-opt with Don't Look Bits (DLB)
var
locallyOptimal: bool = false
i, j: Tour_Index
pos_2_Limit: Tour_Index
X1, X2, Y1, Y2: City_Number
DontLook: DLB_Array
for city in 0 .. N-1:
DontLook[city] = false
while not locallyOptimal:
locallyOptimal = true
block three_loops:
for counter_1 in 0 .. N-1:
if DontLook[tour[counter_1]]:
continue
for direction in [forward, backward]: # NOTE: both tour neighbors
# consider link between t[i] and t[i+1],
# then between t[i-1] and t[i]
if direction == forward:
i = counter_1
else:
i = (N + counter_1 - 1) mod N
X1 = tour[i]
X2 = tour[(i+1) mod N]
for counter_2 in 0 .. N-1:
j = counter_2 # 2nd cut between t[j] and t[j+1]
Y1 = tour[j]
Y2 = tour[(j+1) mod N]
if (X1 == Y1) or (X2 == Y1) or (Y2 == X1):
# the same city or adjacent cities
continue
if Gain_From_2_Opt(X1, X2, Y1, Y2) > 0:
DontLook[X1] = false
DontLook[X2] = false
DontLook[Y1] = false
DontLook[Y2] = false
Make_2_Opt_Move(tour, i, j)
locallyOptimal = false
break three_loops
# end_for backward loop
DontLook[tour[counter_1]] = true
# end_for counter_1 loop

It is expected that the resulting tours will be different from the ones from 2-opt First Take procedure.

Relative time, best length, average length and worst length for all tours created by NN heuristic and improved by 2-opt and 2-opt with DLBs algorithms. Random planar problems, N=100.

2.5-opt algorithm, also known as 2h-opt, is an extension of 2-opt, that takes into account two simple alternative moves of 3-opt type. In its most internal loop it examines not only a 2-opt move, but also two Node Shifts:

before

after

variant A

variant B

variant C

2.5-optimization

proc LS_25_Opt_Take_First(tour: var Tour_Array) =
## Optimizes the given tour using 2.5-opt
# Shortens tour by repeating 2.5-opt moves until no improvement
# can by done: in every iteration immediatelly makes permanent
# the first move found that gives any length gain.
var
locallyOptimal: bool = false
i, j: Tour_Index
X1, X2, Y1, Y2, V0: City_Number
while not locallyOptimal:
locallyOptimal = true
block two_loops:
for counter_1 in 0 .. N-3:
i = counter_1
X1 = tour[i]
X2 = tour[(i+1) mod N]
for counter_2 in (counter_1+2) .. N-1:
j = counter_2
Y1 = tour[j]
Y2 = tour[(j+1) mod N]
if Gain_From_2_Opt(X1, X2, Y1, Y2) > 0:
Make_2_Opt_Move(tour, i, j)
locallyOptimal = false
break two_loops
else:
V0 = tour[(i+2) mod N]
if V0 != Y1:
if Gain_From_Node_Shift(X1, X2, V0, Y1, Y2) > 0:
Make_Node_Shift_Move(tour, (i+1) mod N, j)
locallyOptimal = false
break two_loops
else:
V0 = tour[(N+j-1) mod N]
if V0 != X2:
if Gain_From_Node_Shift(V0, Y1, Y2, X1, X2) > 0:
Make_Node_Shift_Move(tour, j, i)
locallyOptimal = false
break two_loops

Note that this code could be improved to be more efficient:

var
del_2_Length: Length
dX1Y1, dX2Y2, dX2Y1: Length
...
for counter_2 in (counter_1+2) .. N-1:
...
del_2_Length = distance(X1, X2) + distance(Y1, Y2)
dX1Y1 = distance(X1, Y1)
dX2Y2 = distance(X2, Y2)
if del_2_Length - (dX1Y1 + dX2Y2) > 0:
Make_2_Opt_Move(tour, i, j)
locallyOptimal = false
break two_loops
else:
dX2Y1 = distance(X2, Y1)
V0 = tour[(i+2) mod N]
if V0 != Y1:
if (del_2_Length + distance(X2, V0)) -
(dX2Y2 + dX2Y1 + distance(X1, V0)) > 0:
Make_Node_Shift_Move(tour, (i+1) mod N, j)
locallyOptimal = false
break two_loops
else:
V0 = tour[(N+j-1) mod N]
if V0 != X2:
if (del_2_Length + distance(Y1, V0)) -
(dX1Y1 + dX2Y1 + distance(Y2, V0)) > 0:
Make_Node_Shift_Move(tour, j, i)
locallyOptimal = false
break two_loops

This reduces running times by about 25%. But since for other algorithms we used statistics for unimproved versions, the data below are also taken from runs of the basic version.

Some simple statistics

Relative time, best length, average length and worst length for all tours created by NN heuristic and improved by 2-opt and 2.5-opt algorithms. Random planar problems, N=100.

Or-opt is an algorithm proposed by Ilhan Or in 1976: a segment consisting of 3, 2 or 1 consecutive cities is shifted to other position in tour. This means that Or-opt is a restricted version of 3-opt and therefore it should be somewhat faster, but produce worse tours than full 3-optimization.

Or-opt move

Each move used by Or-opt is just Segment Shift:

proc Make_Segment_Shift_Move(tour: var Tour_Array;
i, j, k: Tour_Index) =
## Performs given Segment Shift move
Shift_Segment(tour, i, j, k)

Gain from Or-opt move

proc Gain_From_Segment_Shift(X1, X2, Y1, Y2, Z1, Z2: City_Number):
Length_Gain =
## Gain of tour length which can be obtained by performing Segment
## Shift
# Cities from X2 to Y1 would be moved from its current position,
# between X1 and Y2, to position between cities Z1 and Z2.
# Assumes: X1!=Z1
# X2==successor(X1); Y2==successor(Y1); Z2==successor(Z1)
let del_Length = distance(X1, X2) + distance(Y1, Y2) + distance(Z1, Z2)
let add_Length = distance(X1, Y2) + distance(Z1, X2) + distance(Y1, Z2)
result = del_Length - add_Length

Or-opt using the first improving move

proc LS_Or_opt_Take_First(tour: var Tour_Array) =
## Optimizes the given tour using Or-opt
# Shortens the tour by repeating Segment Shift moves for segment
# length equal 3, 2, 1 until no improvement can by done: in every
# iteration immediatelly makes permanent the first move found that
# gives any length gain.
var
locallyOptimal: bool = false
i, j, k: Tour_Index
X1, X2, Y1, Y2, Z1, Z2: City_Number
while not locallyOptimal:
locallyOptimal = true
for segmentLen in countdown(3, 1):
block two_loops:
for pos in 0 .. N-1:
i = pos
X1 = tour[i]
X2 = tour[(i + 1) mod N]
j = (i + segmentLen) mod N
Y1 = tour[j]
Y2 = tour[(j + 1) mod N]
for shift in segmentLen+1 .. N-1:
k = (i + shift) mod N
Z1 = tour[k]
Z2 = tour[(k + 1) mod N]
if Gain_From_Segment_Shift(X1, X2, Y1, Y2, Z1, Z2) > 0:
Make_Segment_Shift_Move(tour, i, j, k)
locallyOptimal = false
break two_loops

Some simple statistics

Relative time, best length, average length and worst length for all tours created by NN heuristic and improved by 2-opt and Or-opt algorithms. Random planar problems, N=100.

Shifting a segment is a generalization of the concept of shifting a node: instead of shifting a single city, we shift a whole segment (chain of cities) to other position in tour.

Node Shift is a special case of Segment Shift, when a segment has length of 1 (j==i+1).

proc Shift_Segment(tour: var Tour_Array; i, j, k: Tour_Index) =
## Shifts the segment of tour:
# cities from t[i+1] to t[j] from their current position to position
# after current city t[k], that is between cities t[k] and t[k+1].
# Assumes: k, k+1 are not within the segment [i+1..j]
let
segmentSize = (j - i + N) mod N
shiftSize = ((k - i + N) - segmentSize + N) mod N
offset = i + 1 + shiftSize
var
pos: Tour_Index
segment: seq[City_Number] = newSeq[City_Number](segmentSize)
# make a copy of the segment before shift
for counter in 0 .. segmentSize-1:
segment[pos] = tour[(pos + i+1) mod N]
# shift to the left by segmentSize all cities between old position
# of right end of the segment and new position of its left end
pos = (i + 1) mod N
for counter in 1 .. shiftSize:
tour[pos] = tour[(pos + segmentSize) mod N]
pos = (pos + 1) mod N
# put the copy of the segment into its new place in the tour
for pos in 0 .. segmentSize-1:
tour[(pos + offset) mod N] = segment[pos]

Node Shift is a simple tour optimization heuristics: it consists in shifting a city to other position in tour.

The idea:

Diagrams:

Node Shift move is obtained by removing three links and adding three new links, and therefore it is a special case of 3-opt move:

Gain from Node Shift move

proc Gain_From_Node_Shift(X0_pred, X0, X0_succ,
Y1, Y2: City_Number): Length_Gain =
## Gain of tour length which can be obtained by performing Node Shift
# City X0 would be moved from its current position, between X0_pred
# and X0_succ, to position between cities Y1 and Y2.
# Assumes: X0_pred!=Y1 (the same as: X0!=Y2);
# X0==successor(X0_pred); X0_succ==successor(X0); Y2==successor(Y1)
let del_Length = distance(X0_pred, X0) + distance(X0, X0_succ) +
distance(Y1, Y2)
let add_Length = distance(X0_pred, X0_succ) +
distance(Y1, X0) + distance(X0, Y2)
result = del_Length - add_Length

Node Shift move implementation

proc Make_Node_Shift_Move(tour: var Tour_Array;
i, j: Tour_Index) =
## Performs given Node Shift move
# Shifts the city t[i] from its current position to position
# after current city t[j], i.e. between cities t[j] and t[j+1].
# This is a special, simple case of segment shift, thus a special
# case of 3-opt move.
var
left, right: Tour_Index
X0: City_Number
let shiftSize = (j - i + 1 + N) mod N
X0 = tour[i]
left = i
for counter in 1 .. shiftSize:
right = (left + 1) mod N
tour[left] = tour[right]
left = right
# end_for counter loop
tour[j] = X0

Iterative Node Shift optimization

proc LS_Node_Shift(tour: var Tour_Array) =
## Iteratively optimizes given tour using Node_Shift moves.
# Shortens tour by repeating Node Shift moves until no improvement
# can by done: in every iteration immediatelly makes permanent
# the first move found that gives any length gain.
var
locallyOptimal: bool = false
i, j: Tour_Index
X0_pred, X0, X0_succ, Y1, Y2: City_Number
while not locallyOptimal:
locallyOptimal = true
for counter_1 in 0 .. N-1:
i = counter_1 # current position of city to shift
X0_pred = tour[(i + N - 1) mod N]
X0 = tour[i]
X0_succ = tour[(i + 1) mod N]
for counter_2 in 1 .. N-2:
j = (i + counter_2) mod N # new position of the city
Y1 = tour[j]
Y2 = tour[(j + 1) mod N]
if Gain_From_Node_Shift(X0_pred, X0, X0_succ, Y1, Y2) > 0:
Make_Node_Shift_Move(tour, i, j)
locallyOptimal = false
break

Some simple statistics

Relative time, best length, average length and worst length for all tours created by NN heuristic and improved by 2-opt and 2-opt plus Node Shift algorithms. Random planar problems, N=100.

Node Swap is a simple tour optimization heuristics: it consists in swapping two cities in tour.

The idea:

Diagrams:

Node Swap move is a special case of two subsequent 2-opt moves: the first including both cities and the second without them. It involves removing four links and adding four new links, therefore is a specific type of 4-opt move.

Gain from Node Swap move

proc Gain_From_Node_Swap(X0_pred, X0, X0_succ,
Y0_pred, Y0, Y0_succ: City_Number): Length_Gain =
## Gain of tour length which can be obtained by performing Node Swap
## City X0 would be moved from its current position, between X0_pred
## and X0_succ, to the position of city Y0 and vice versa
# Assumes: X0!=Y0
# X0=successor(X0_pred); X0_succ==successor(X0);
# Y0=successor(Y0_pred); Y0_succ==successor(Y0);
var
del_Length, add_Length: Length_Gain
if (Y0==X0_succ) or (X0_succ==Y0_pred):
del_Length = distance(X0_pred, X0) + distance(Y0, Y0_succ)
add_Length = distance(X0_pred, Y0) + distance(X0, Y0_succ)
elif (X0==Y0_succ) or (Y0_succ==X0_pred):
del_Length = distance(Y0_pred, Y0) + distance(X0, X0_succ)
add_Length = distance(Y0_pred, X0) + distance(Y0, X0_succ)
else:
del_Length = distance(X0_pred, X0) + distance(X0, X0_succ) +
distance(Y0_pred, Y0) + distance(Y0, Y0_succ)
add_Length = distance(X0_pred, Y0) + distance(Y0, X0_succ) +
distance(Y0_pred, X0) + distance(X0, Y0_succ)
#end_if
result = del_Length - add_Length

Node Swap move implementation

proc Make_Node_Swap_Move(tour: var Tour_Array;
i, j: Tour_Index) =
## Performs given Node Swap move: swaps city t[i] with city city t[j]
# This is a special, simple case of two subsequent 2-opt moves:
# the first including both cities and the second without them
# -a-t[i]-b-t[j]-a-
# -a-(t[i]-b-t[j])'-a- == -a-t[j]-b'-t[i]-a-
# -a-t[j]-(b')'-t[i]-a- == -a-t[j]-b-t[i]-a-
swap(tour[i], tour[j])

Iterative Node Swap optimization

proc LS_Node_Swap(tour: var Tour_Array) =
## Optimizes the given tour using Node_Swap.
# Shortens tour by repeating Node Swap moves until no improvement
# can by done: in every iteration immediatelly makes permanent
# the first move found that gives any length gain.
var
locallyOptimal: bool = false
i, j: Tour_Index
X0_pred, X0, X0_succ, Y0_pred, Y0, Y0_succ: City_Number
while not locallyOptimal:
locallyOptimal = true
for counter_1 in 0 .. N-1:
i = counter_1 # position of the first city
X0_pred = tour[(N + i - 1) mod N]
X0 = tour[i]
X0_succ = tour[(i + 1) mod N]
for counter_2 in counter_1+1 .. N-1:
j = counter_2 # position of the second city
Y0_pred = tour[(N + j - 1) mod N]
Y0 = tour[j]
Y0_succ = tour[(j + 1) mod N]
if Gain_From_Node_Swap(X0_pred, X0, X0_succ,
Y0_pred, Y0, Y0_succ) > 0:
Make_Node_Swap_Move(tour, i, j)
locallyOptimal = false
break

Some simple statistics

Relative time, best length, average length and worst length for all tours created by NN heuristic and improved by 2-opt and Node Swap algorithms. Random planar problems, N=100.