This paper introduces a new variable depth search method for the Quadratic Assignment Problem. The new method considers the cost of edges assignment as the criterion to decide which vertices to exchange during local search moves. It also presents the results of an extensive experimental study that compares the performance of local search and variable depth search algorithms for the Quadratic Assignment Problem. The investigation presented here contributes to a better understanding of the potential of these techniques, which are widely used as intensification tools in more sophisticated heuristic methods, such as evolutionary algorithms. Different algorithms presented in the literature were implemented and compared to the proposed methods. The results of a computational experiment with 161 benchmark instances are reported. Different statistical tests are applied in order to analyze the results provided by the experiments.

Local search algorithms are popular methods to find approximate solutions to hard optimization problems. At present many metaheuristics are extensions or generalizations of local search procedures, such as Tabu Search, Simulated Annealing, GRASP and Variable Neighborhood Search, among others (Glover & Kochenberger, 2003). Local search is also often embedded in evolutionary algorithms to deal with intensification tasks, where the objective is to combine the power of evolutionary operators in determining interesting regions of the search space with the power of local search procedures to conduct a thorough search in a specified region (Grosan & Abraham, 2007). In this sense, investigations of efficient local search methods contribute to the development of better heuristics to hard optimization problems, since the development of more efficient search strategies will result in more powerful heuristics. In its basic form, the local search algorithm starts from an initial solution and searches among the neighboring solutions one that improves the starting solution. If such a neighbor is found, the process is re-started with the new solution replacing the current one as the starting solution. The iterations continue until no solution with better quality is found in the neighborhood of the current solution. The best solution found by the algorithm is said to be a local optimum with respect to the neighborhood structure used in the search. Usually, the neighborhood structure depends on the problem under consideration and is a determinant factor for the performance of local search algorithms. Variable Depth Search (VDS) is an important variation of local search. It was firstly presented for graph partitioning (Kernighan & Lin, 1970) and for the Traveling Salesman Problem (Lin & Kernighan, 1973). Lin and Kernighan's heuristic for the Traveling Salesman Problem is considered to be one of the most effective methods to symmetric instances (Johnson & McGeoch, 2002). The main feature of VDS is to perform subsequences of local search moves such that large portions of the space of solutions are explored in reasonable processing times. Usually, the search is based on simple neighborhood structures and a sequence of solutions is obtained with the application of simple local search moves. The number of solutions generated in each sequence varies during the search. In most cases, the effectiveness of variable depth searches is highly dependent on the choices of the exchanged elements at each iteration step. Applications of VDS were proposed for the generalized assignment problem (Yagiura et al., 1998) and for bi-partitioning signal flow graphs (de Kock et al., 1995).

Given the trade-off between processing time and quality of solution, a question arises on whether it is advantageous or not to use VDS instead of simple local search approaches. This is the first point investigated in this paper that focuses on the Quadratic Assignment Problem (QAP) which is a strong NP-hard combinatorial optimization problem (Sahni & Gonzalez, 1976) with extensive practical application and an important test ground for most heuristics. A wide range of algorithms based on metaheuristic strategies have been proposed for this problem, many of them have local search methods as an important component (Loiola et al., 2007). Although a huge number of papers are dedicated to the QAP, as far as the authors' knowledge concerns, no paper addressed the investigation of the benefits of using variable depth search methods rather than simple local search methods. This is the first point examined in this paper that compares strategies based on local search and on VDS for the QAP. Two VDS algorithms were proposed for the investigated problem (Li et al., 1994; Rego et al., 2006) which use the 2-exchange neighborhood structure. Those algorithms are compared to the VLSN (Very Large Scale Neighborhood) (Ahuja et al., 2007) and to an iterated local search algorithm that uses the classical 2-exchange neighborhood and was embedded on an efficient Evolution Strategies algorithm (Stützle, 2006). The results of a computational experiment with 92 benchmark instances show that both VDS approaches outperform the local search methods and therefore are a better option as intensification tools for other heuristics.

The VDS algorithms previously proposed in the literature exchange vertices systematically and do not consider any other information besides the value of the objective function to determine which vertices are exchanged. Since the elements that contribute in the computation of the objective function of the QAP are the costs of the assignment of edges resultant from the assignment of vertices, a new variable depth search heuristic that aims at taking advantage of the information of the cost of edge assignments to choose which vertices are exchanged during the search is proposed here. A computational experiment with 161 benchmark instances shows that the consideration of the cost of the edge assignments is beneficial for several families of instances.

This paper is organized in other five sections besides this one. An overview of the Quadratic Assignment Problem is presented in Section 2. Section 3 presents the methodology used in the computational experiments reported in subsequent sections. Existing VDS and local search algorithms are compared in Section 4. Section 5 introduces the new VDS algorithm and reports the results of computational experiments. Concluding remarks are presented in Section 6.

2 THE QUADRATIC ASSIGNMENT PROBLEM

The Quadratic Assignment Problem was first introduced by Koopmans & Beckman (1957) in the context of facility-location problems. Given square matrices of order n, F = (fij) andD = (dij), the problem consists in finding a permutation ρ of the set N ={1,..., n} that minimizes the cost c(ρ) given in equation 1.

In the location theory, the elements of matrix F represent flows of materials between facilities and the elements of matrix D represent distances between locations. If matrices F and D are symmetric then the QAP is said to be symmetric, otherwise it is asymmetric. In terms of Graph Theory, the QAP can be thought as an assignment between the vertices of two complete graphs of order n, whose weighted adjacency matrices correspond to F and D. The assignment of vertices yields an assignment of the edges of the correspondent graphs. The cost of the assignment of edge (ρ(i), ρ( j)) to edge (i, j) is represented by fρ(i)ρ( j)dijin equation (1). The objective is to find an assignment of vertices which minimizes the cost of the assignment of edges given by the sum of the costs of each edge assignment.

Since, usually, real world problems modeled by the QAP are larger than the ones that exact algorithms are able to solve, a wide variety of heuristic approaches were proposed for handling near optimal solutions for this problem (Loiola et al., 2004, 2007). Heuristic algorithms developed according to a myriad of metaheuristic techniques were proposed for the QAP. A review of such algorithms is presented by Loiola et al. (2004, 2007). Some recent works include Tabu Search (James et al., 2009a, 2009b; Fescioglu-Unver & Kokar, 2011), Differential Evolution (Davendra et al., 2009), Iterated Local Search (Ramkumar et al., 2009), Ant Colony (Puris et al., 2010) and Lagrangian Smoothing Algorithm (Xia, 2010). Paul (2010) presents a comparison between Simulated Annealing and Tabu Search algorithms for the QAP.

Most successful heuristics presented for the QAP are based on, extend or embed local search methods with the 2-exchange neighborhood structure (Drezner, 2005; Drezner, 2008; James et al., 2009; Fescioglu-Unver & Kokar, 2011). Therefore, an investigation on the potential of local search procedures and variations is very useful in the context of deciding which method to use to compose more sophisticated heuristics.

3 MATERIALS AND METHODS

The next two sections present computational experiments that consider two sets of test cases: QAPLIB (Burkard et al., 1991) (http://www.opt.math.tu-graz.ac.at/qaplib/inst.html), Taixxe and Dre instances (Drezner et al., 2005) (http://mistic.heig-vd.ch/taillard/problemes.dir/qap.dir/qap.html). TheQAPLIB is a benchmark library of QAP instances with different structures. Those instances can be divided in four classes (Taillard, 1995): instances randomly generated from a uniform distribution, instances whose distance matrix is based on the Manhattan distance on a grid and the elements of the flow matrix are randomly generated, not necessarily from a uniform distribution, instances that arise from practical applications and instances that have the same characteristics of real-life problems, that is, they are randomly generated, but not from a uniform distribution and the distance matrix corresponds to the Euclidean distance between n points in the plane. Taixxe and Dre instances are introduced by Drezner et al. (2005). Those instances were designed to be difficult to solve by algorithms that depend on transposition neighborhoods.

The VDS algorithms utilized in the computational experiments reported in this paper, including the algorithms presented by Li et al. (1994) and Rego et al. (2006), were implemented under an Iterated Local Search (ILS) framework. ILS is an effective way of using repeated runs of local search. This method has been proposed, independently, by several researchers who gave it different names, such as iterated descent, large-step Markov chains, and iterated local search among others (Lourenço et al., 2002). This search strategy leads to a randomized walk in the space of local optima and, in many cases, is more effective than the random re-starts. In iterated local search algorithms a starting solution, ρ0, is generated and a local search is applied to it. The resulting solution, ρ1, is disturbed and the local search procedure is applied to ρ1, resulting in ρ2. An acceptance criterion is applied to decide which solution will be the input for a new application of the local search. These steps are repeated until a given condition is met.

The platform used to implement the iterated local search algorithms is an Intel Core 2, 2.4 GHz, 2 Gb RAM, running Linux. Twenty five independent executions of each algorithm were performed for each instance. The stopping criterion adopted in the experiments is a pre-specified maximum processing time.

In Section 4, existing VDS like algorithms (Li et al., 1994; Rego et al., 2006) are compared to an effective ILS algorithm that uses the classical 2-exchange neighborhood (Stützle, 2006) and to a multi-start algorithm that implements k-exchange neighborhoods (Ahuja et al., 2007). Except for the latter, all algorithms were implemented by the authors. The results reported in the experiments for the multi-start algorithm were presented by Ahuja et al. (2007). The experiments were performed on 92 symmetric and asymmetric QAPLIB test cases, including all families of instances.

The computational experiment concerning the VDS algorithms proposed in this paper were performed on 161 symmetric and asymmetric QAP instances with n ranging from 20 to 150: 92 QAPLIB instances, 60 Taixxe instances and 9 Dre instances.

In general, the mean and the median costs of the solutions found by each algorithm on the QAPLIB instances were close. It did not occur with Taixxe instances. Therefore, we opted to present the mean when dealing to QAPLIB instances and the median to Taixxe and Dre instances. The choice to present the median for the latter classes of instances was made since this statistic presents an upper limit to the minimal values found in at least 50% of the independent executions. The results are presented in terms of percent difference from the best known solution calculated with equation 2, where val denotes the value of the mean, median or best solution accordingly and best denotes the value of the optimal solution or the best known value reported in the literature when the optimum is not known.

Statistical analyses were done with the Kruskal-Wallis' test, the Mann-Whitney U test (Conover, 1999) and the test for proportions comparison proposed by Taillard et al. (2008) (two-tailed test). In the latter case the proportions refer to number of successes obtained by each heuristic method considering a pre-specified computational effort. Success, here, means to achieve the lowest value when average or median solutions are compared and to achieve the highest value when dealing with number of best solutions. Algorithms are compared in pairs. The input data for proportions comparison tests are the number of instances considered in a given experiment and the number of instances where each algorithm is successful in relation to the other.

To illustrate the potential of the proposed algorithms, their results are compared to the ones produced with an exact algorithm proposed recently by Zhang et al. (2010). The comparison is performed on 17 sparse instances with n up to 32.

4 LOCAL SEARCH AND VARIATIONS

This section reviews some approaches of local search and variations applied to the QAP. The neighborhood structures used in these algorithms are based on the exchange of elements of a current solution. The exchanges are done systematically, without any information, except for the value of the objective function calculated with equation 1. The performances of these algorithms are compared and conclusions are drawn. This section is divided into two parts. The different heuristics are presented in Subsection 4.1 and the computational experiment is presented inSubsection 4.2.

4.1 Local Search for the QAP

Exchange neighborhoods are very popular for problems that can be represented as permutations of n elements, such as the QAP. The local search procedures used in the most successful heuristics proposed to tackle the investigated problem rely on the 2-exchange neighborhood structure. Given a solution ρ, the neighborhood א(ρ) obtained through the 2-exchange structure, is defined by the set of permutations which can be reached by exchanging two elements of ρ, as defined in expression 3.

There are n(n  1)/2 possible combinations of two locations of a permutation with n elements and, given two elements, there is only one way of exchanging them. An advantage of this neighborhood for the QAP is that the difference between the costs of two neighbor solutions can be computed in O(n) (Taillard, 1991). This neighborhood is the most widely used in several heuristic algorithms proposed for the examined problem.

The simplest way to use a local search algorithm is to perform several re-starts of the procedure from different points of the search space. Multi-start applications with the 2-exchange neighborhood are reported to the QAP in a number of papers (Fleurent & Glover, 1999; Misevicius, 1997; Misevicius & Riskus, 1999). Resende et al. (1996) present a greedy randomized adaptive search (GRASP) for the QAP. An extension is presented by Oliveira et al. (2004) who improve the performance of the GRASP algorithm with the inclusion of a path-relinking procedure. They present experimental results for 92 instances of the QAPLIB (Burkard et al., 1991) with n up to 64. Lourenço et al. (2003) observe that as the instance size increases, the probability of random re-starts methods to find solutions with significantly lower costs decreases. The Iterated Local Search method, ILS, is an alternative to prevent that effect. The method is thought to perform a biased sampling of the search space by performing small perturbations in the solution generated by the local search procedure and using this disturbed solution as the starting solution of the next iteration (Lourenço et al., 2003). An ILS approach with the 2-exchange neighborhood is proposed for the QAP by Stützle (2006). The general ILS framework and the ILS algorithm proposed for the QAP are reviewed in the following, since the ILS algorithm proposed by Stützle (2006) is tested in the computational experiments presented in this paper and also the ILS architecture is used in the implementation of the proposed methods. A general framework of the ILS algorithm is presented in Figure 1.

An initial solution is generated and a local search is applied to it in steps 1 and 2, respectively. The ILS algorithm may use or not a set to store the history of the search, called history in Figure 1. This set may contain, for example, interesting solutions generated during the search. The history of the search may influence the decisions made in steps 4 and 6. Although history can be considered, in many ILS implementations, the decisions regarding the perturbation of the current solution (step 4) and acceptance criterion (step 6) do not depend on it. The solution generated in the first two steps is disturbed and a local search is applied to the new solution in steps 4 and 5, respectively. In order to re-start the search, an acceptance criterion chooses from which solution the search will be continued. At least, two new decisions, beyond those of the local search method, have to be made in order to use an iterated local search algorithm: the method to disturb solutions and the criteria to accept a solution.

Four ILS algorithmic versions were investigated with the major difference between them regarding the criterion to choose the starting solution of the next iteration, called acceptance criterion. Those versions are called: Better, Restart, RandomWalk and LSMC. In Better, the best solution among ρ0 and ρ2 is the input for the next iteration. The same acceptance criterion is used in Restart, but a new random solution is returned by the acceptance criterion if an improved solution is not found for a fixed number of iterations. In RandomWalk, ρ0 is always replaced by ρ2 irrespective of the cost of the latter. Finally, in LSMC a simulated annealing based criterion is used. In this version, ρ2 is accepted with a probability p that depends on the difference between the costs of ρ0 and ρ2 and on a temperature. The latter version obtains the overall best solutions when compared to the other three versions and was implemented in our computational experiments for the comparison with the other algorithms. The implementation followed the directions given in the paper it was presented. The local search procedure uses the first pivoting rule and "don't look bits". Stützle (2006) presents experimental results of the iterated local search algorithms for 38 instances with n ranging from 20 to 100. The ILS algorithm is then embedded in an evolutionary algorithm that presents high quality results.

Multi-exchange neighborhoods are natural extensions of the 2-exchange. Ahuja et al. (2007) investigate multi-exchange neighborhoods for the QAP. They defined ρ', a k-exchange neighbor of ρ, as the permutation obtained from ρ by a cyclic sequence of k elements. Given a sequence of locations i1,..., ik, the multi-exchange move corresponds to assign facility ρ(ij) to location ij + 1, for all j≠k and to assign facility ρ(ik) to location i1. The local search algorithm iteratively examines all cyclic sequences of increasing values for k, where the maximum is aspecified parameter. The size of the neighborhood grows exponentially with k, thus this neighborhood structure is classified as a Very Large Scale Neighborhood, for which an exhaustive examination of the whole neighborhood is impracticable in computational terms. With the use of an improvement graph, Ahuja et al. (2007) show that, in average, they can compute the cost of a k-exchange move in O(k). Due to the high processing times needed to the full enumeration scheme the authors established a maximum value for k equals 4. They test multi-start implementations of the proposed neighborhood structure against a multi-start implementation of the 2-exchange neighborhood in 132 test cases of the QAPLIB, with fixed processing times. They show that their neighborhood obtains consistently better results than the multi-start local search with the 2-exchange neighborhood.

In general, for fixed values, as k grows, the computational effort to deal with k-exchange neighborhoods rises rapidly and, usually, the improvement in the quality of solutions does not compensate a greater effort. An effective alternative for the fixed size neighborhood structures is presented by Kernighan & Lin (1970) and Lin & Kernighan (1973) for the Partitioning Problem and the Traveling Salesman Problem, respectively. The idea is to perform sequences of moves, such that each move is likely to lead to a better solution instead of testing all moves of a given neighborhood. This method is called Variable Depth Search (VDS). The VDS algorithm examines iteratively, for growing values of k, whether exchanging k elements may result in a better neighboring solution. A gain function is computed and if it is likely that the previous solution can be improved, then the algorithm tests if k can be extended to k + 1. When no more gain can be made, the algorithm performs the exchange of the k elements. The general framework of a VDS algorithm is shown in Figure 2.

Similar to other local search algorithms, an initial solution ρ0 is generated in step 1. The current best solution is stored in ρbestinitially set to ρ0. The variable gain stores the result of the gain function. The variable bestgain maintains the best value achieved by the gain function during the iterations of the inner loop. Whenever bestgain is updated, variable k is set to the value of variable i, the number of the current iteration. The value of k indicates the length of the best sequence of exchanges the algorithm has to perform after the inner loop is executed. Lists Lout and Lin contain elements that were withdrawn and included in the current solution, respectively. The choice of elements xi(that is withdrawn) and yi(that is included in the current solution) aims at maximizing the improvement of the current solution within the sequence of moves. Variable ρ1 stores the neighboring solutions found during the VDS iterations. A test is done in step 14 to find out if, given conditions established in advance, the value in gain still indicates that an improvement exists with the current sequence of moves. The inner loop also finishes if no elements exist to be exchanged due to the restrictions imposed by the lists of prohibited elements what is verified in procedure element exist( ). In step 15 a new solution ρ2 is generated with the replacement of elements x1,..., xkby y1,..., ykin solution ρ0. The main loop continues iterating until a given stop criterion is met.

A VDS algorithm for the QAP was presented by Li et al. (1994) who investigated the complexity of local search for the QAP with a 2-exchange neighborhood which they proved to be PLS-complete. Their algorithm starts with a permutation generated at random, with uniform probability. For a current permutation ρ0, a sequence of permutations, ρ1,...,ρs, is built, each of them obtained from the previous one with a 2-exchange move with cost lower than the current permutation. They used a cumulative gain function, G(i), for permutation ρi,where G(i) = c(ρi)  c(ρ0). At each re-start of the local search, ρ0 is replaced by the best ρi ,the one with the best G(i). The stopping condition was the inexistence of a sequence for the current permutation. Another VDS like method for the QAP is presented by Rego et al. (2006) who propose an ejection chain algorithm. Their method builds sequences of exchanges until no promising permutations are likely to exist or until a maximum of n moves. Initially, the best 2-exchange move for each facility ρ(i), i = 1,..., n, is identified. Let j be the best location for facility ρ(i), then this facility is assigned to location j. The facility that previously occupied position j is said to be ejected and has to be assigned to other location. The method prohibits elements to be moved twice. Thus the sequence can grow until all n facilities are moved. Each element of the move sequence corresponds to one local search step in a QAP instance of q elements, where q varies from n to 1.

4.2 Experiments Comparing Local Search Algorithms

The algorithms by Li et al. (1994) and Rego et al. (2006) were implemented under an iterated local search framework and are named IT-LPR and IT-EC, respectively. The solution generated by the algorithm presented by Rego et al. (2006) is not guaranteed to be a local optimum for the 2-exchange neighborhood. Therefore, a final 2-exchange local search step was included at the end of the iterated local search iterations. Preliminary experiments showed that this modification improved significantly the performance of this algorithm.

The VDS algorithms are compared to the LSMC (Stützle, 2006) and to the VLSN (Ahuja et al., 2007). For the experiments performed in this paper, we tested two methods to update the "don't look bits" in the LSMC. The first method followed the directions given in the reference paper by Stützle (2006) where only the "don't look bits" of exchanged elements are reset to 0. In the second method all "don't look bits" are reset to 0 after a solution is disturbed. The results of the experimentation made to compare these methods showed that the first method produces worse solutions than the second method for the majority of the tested instances. Thus, for comparison purposes, the best result found for each instance by one of these two versions is reported for LSMC.

The results listed in columns related to the VLSN were reported by Ahuja et al. (2007) who used an IBM RS6000 platform (333 MHz). Their stopping criterion was 1 hour for problems with n > 40 and 2 hours for problems greater than 40. The stop condition for LSMC, IT-LPR and IT-EC is the maximum processing time (in seconds) shown in the last column of Tables 1 and 2. These tables present the results for symmetric and asymmetric QAPLIB instances, respectively. The name of the instance is presented in the first column followed by the best known solution in the second column. Columns Av and %best present, respectively, the average percent deviation from the best known solution and the percentage of best known solutions obtained in 25 independent executions. The results obtained for instances Esc32x, Sko100x and Bur26x are grouped. In these cases the best known solution is not presented in the correspondent column. The best result obtained for each instance by one of the compared algorithms is bold.

VLSN produces the best average percent deviation on instances Esc32e, Esc32g and Tai100a. The LSMC presents the best results on instances Esc32c, Esc32e, Esc32g and Esc64a. The results in Table 1 show that the two VDS algorithms present, in general, the best performance for the 58 symmetric instances of the computational experiment. The IT-LPR presents the best average results on 18 instances. The IT-EC presents the overall best performance regarding quality of solution. It presents 46 best average results including instances Esc32c, Esc32d, Esc32e, Esc32g and the whole Sko100 set. The number of best known solutions produced by the two VDS algorithms is similar with some advantage for the IT-EC. The average concerning the 58 symmetric instances and the 25 independent executions are 4.02, 4.90, 23.65 and 28.39 for the VLSN, LSMC, IT-LPR and IT-EC, respectively.

The Kruskal-Wallis' statistical test was applied to determine if significant differences exist between the results presented by the three ILS algorithms. The VLSN could not be tested, once the results for each instance are not available. The results of the statistical test considering significance level 0.05 showed that the null hypothesis (samples come from identical populations) cannot be rejected only for instances Esc32c, Esc32e, Esc32g, Esc64a and Esc128, indicating that for the other 53 instances the results produced by at least one algorithm were significantly different from the other two. In all cases the test indicated that the LSMC performed significantly worse than the other algorithms. The Mann-Whitney U-test was applied to the results of the two VDS algorithms. With significance level 0.05, the test showed that the IT-EC produced significantly superior results on instances Chr22a, Chr22b, Chr25a, Kra30b, Nug27, Nug28, Tai25a, Tai35a, Sko42, Sko49, Sko56, Sko64, Sko72, Sko81, Sko90, Sko100a-f, Ste36a-c, Tho30, Tho40, Tho150, Wil50 and Wil100. Significant better results were presented by the IT-LPR on instances: Rou20, Tai50a, Tai60a, Tai80a and Tai100a.

In order to apply the test for proportions comparison (Taillard et al., 2008) it is necessary to define "success" In this paper, one algorithm is considered successful on one instance when compared to other, if the former produces a better average value than the latter. With significance level 0.01, the statistical test showed that the IT-EC was the best approach when compared to each of the other three algorithms. The VLSN performed conclusively worse than the others. The comparison between the IT-LPR and to the LSMC showed that that the former outperformed the latter.

Similar results were observed on the asymmetric instances presented in Table 2. The best results for these instances were produced by the IT-LPR and the IT-EC. The IT-LPR presented the best average results on 15 instances and the IT-EC presented 26 best average results including Bur26 instances. The average of best known solutions found by each algorithm on the 34 asymmetric instances are 1.55, 3.12, 25.19 and 41.98 for the VLSN, LSMC, IT-LPR and IT-EC, respectively. With significance level 0.05, the Kruskal-Wallis' test showed that the results of, at least, one algorithm is significantly different from the others. The statistical test showed that the two VDS algorithms performed significantly better than the LSMC. Concerning the comparison of the two VDS algorithms, significant best results were obtained by the IT-EC on instances Bur26a-h, Tai25b, Tai30b, Tai35b, Tai40b and Tai100b. IT-LPR presented significantbest results on instances Lipa80a, Lipa80b and Lipa90b.

The test for proportions comparison with significance level 0.01 showed that the two VDS algorithms outperform the VLSN and the LSMC. The null hypothesis cannot be rejected (even with significance level 0.05) in the comparison between VLSN and LSMC. The test also pointed out a significant better result of the IT-EC over the IT-LPR.

5 GUIDING VDS WITH THE COSTS OF EDGE ASSIGNMENTS

In terms of Graph Theory, each parcel of equation 1 corresponds to the cost of an assignment of two edges of complete graphs of order n. From this viewpoint, the objective function seeks at minimizing the sum of products of edge weights. The problem of assigning edges is solved in polynomial time, but only few assignments of edges, regarding the total number of such assignments, lead to feasible solutions concerning the assignment of vertices when n > 3 (Rangel & Abreu, 2007). Nevertheless, the costs of the edges assignments can be used to guide the search indicating good candidates to be exchanged. In this sense, the edges are thought as the exchanging elements, taking due care to maintain feasibility of the solutions produced during the search. Given a solution, ρ, the flow edge efwith terminal vertices ρ(i) and ρ( j) and the location edge edwith terminal vertices i and j, the assignment of edges is represented by (ef, ed).

Consider that edge efis assigned to ed, the assignment of ef to a new location e'd implies in freeing the edge e'f, previously assigned to e'd. Therefore, the basic idea consists in, iteratively, canceling an edge assignment and to re-assign the freed flow edge to another distance edge.

In this paper, two alternatives are considered for the new assignment. In the first alternative, the assignment of one terminal vertex of efis undone and the assignment of the other terminal vertex remains fixed. In the other alternative both assignments of the terminal vertices of edge efare undone. In both cases, the freed vertex(vertices) is(are) assigned to a new location(s). This new assignment induces a new solution, ρ1. The difference between the costs of the original solution, ρ, and the solution induced by the basic move, ρ1, is calculated and the gain function is updated. The first and the second alternatives lead to one and two algorithmic versions, respectively. The general framework of the variable-depth search based on the first alternative is presented in algorithm VDS QAP.

After the initial solution, ρ0, is generated, the first assignment of edges to be undone is chosen in procedure edge_assignment( ). A list of prohibited moves, LProhib, stores the pair of edges (ef, ed) of the assignment that is undone. This list is initialized with (ef, ed). Its maximum size is limited to #sizelistelements. Inside the loop a new distance edge, e d, is chosen. The flow edge efis assigned to e'd in procedure new_edge( ). The edge e'f previously assigned to e'd is set free. A new solution, ρ2, is generated by exchanging edges ef and e'f . Clearly, the edge assignments depend on the terminal vertices of the corresponding edges. The distinct algorithmic versions proposed in this paper differ in this step, where different strategies to exchange the edges and to assign the terminal vertices are adopted. The variation between the costs of the two solutions is calculated and incremented in variable gain. The best value of this variable is stored in variable bestgain and the position in the sequence of exchanges is stored in variable k. The list of prohibited assignments is updated with (e'f , e'd). If there are less than #sizelistelements in LProhib, then the new pair is added to it. Otherwise, the new element replaces the "oldest" element of LProhib. The main loop iterates while positive gains exist (variable gain is greater than 0). Alternatively a stopping criterion can also consider a maximum number of iterations of the main loop. Finally, the initial solution, ρ0, is updated with the best sequence of edges exchange in procedure exchange_sequence( ).

Two important decisions to be made when implementing VDS QAP are: to define how many and which assignments are canceled at each iteration step and to define how the new locations for the freed edges are chosen. A number of alternatives exist for these and other implementation issues. Depending on the decisions made by the developer, algorithms with wide varying behaviors are expected to be created. In the next two sections, three algorithmic versions are presented where the details concerning those decisions are explained.

In the algorithms presented in the next sections, the choice of the first edge assignment that is undone is based on the cost of the assignments of each pair of edges in the initial solution. Given a solution, ρ, the weight of the assignment (ef, ed) of edge ef with terminal vertices ρ(i) and ρ( j) to edge edwith terminal vertices i and j is given by w(ef, ed) = fρ(i)ρ( j)dij. In order to select the first pair, (ef, ed), a restricted list of candidates is built with the lsize% pairs with thegreatest weights. One element of that list is randomly chosen with a probability that is directly proportional to the weight of the correspondent assignment.

In order to determine the best value for lsize, a preliminary experiment was done with 42 QAPLIB instances with n between 20 and 150. Four values were investigated for lsize: 10, 30, 50 and 70.

It was observed that the best results were obtained with 10 or 30, with the latter value yielding slightly better results than the former.

The sequence of movements is terminated if gain < 0. In the implementations reported here an additional stopping criterion was added to establish a maximum number of iterations being twice the size of LProhib without the improvement of the best solution found up to the current iteration. The size of LProhib was set to 150 after preliminary experiments.

5.1 One Terminal Vertex is set Free  VDS1

In the first variant of the proposed algorithm, denoted as VDS1, only one terminal vertex of edge efis set free. Thus, it is necessary to decide which terminal vertex assignment remains fixed and which is undone. It corresponds to choose which location of edge edbecomes unoccupied. This decision is made simultaneously with the choice of the new location for the free facility. Both terminal vertices of edge efare tested.

In a first experiment to examine the potential of the variant VDS1 in relation to the 2-exchange neighborhood, 100 random re-starts of each method were applied to all QAPLIB instances with n ranging from 20 to 30. In this experiment, each released vertex was tested in all locations not prohibited for it. The minimal, average and median solution reached by each method were computed, as well as the average processing time. The results are shown in Table 3 where the values of the minimal, average and median solutions are reported in columns Min, Av and Med, respectively. The average processing time in seconds is exhibited in column T .

The results presented in Table 3 show that the minimal values obtained by VDS1 are better than the ones produced by the 2-exchange on 32 of the 38 instances, there are 4 ties and the latter method finds 2 minimal values better than the former. The proposed method also finds the best average solutions on 36 instances. All median values produced by the VDS1 are lower than the median values obtained with the other method. A summary of the improvements in average solutions obtained by VDS1 over the other algorithm per instance class is shown in Table 4. The values exhibited in Table 4 were calculated on classes with at least three instances. Columns Class, Best and Improvement show, respectively, the class name, which algorithm produced the best average value and the improvement of the average solution over the other algorithm.

Table 4 shows that the 2-exchange local search exhibits better results on Taixxb class, where an improvement of 4.01% is obtained. The proposed approach outperforms the 2-exchange algorithm on the remaining instances with the lowest improvement being 16.67% on class Lipa.

Once normality cannot be assumed, the Mann-Whitney test was applied to the results produced by the investigated methods to test whether or not significant differences exist. With significance level 0.01, the null hypothesis cannot be rejected only for four instances: the three Taixxb instances and the Chr20c. Statistical significant differences were found on the remaining instances. Table 3 also shows that smallest processing times are, in general, obtained with the 2-exchange local search. Nevertheless, when the same statistical investigation is applied to the processing times, with significance level 0.05, significant differences exist in favor of the 2-exchange on 10 of the 38 instances only. It means that on the majority of instances, 28 of 38, the statistical test did not indicate significant differences regarding computational times.

In the ILS version of VDS1, two matrices are used to establish a preference order among the localities to which a free facility can be assigned. In a pre-processing phase, two matrices n × n  1 of ranks, RFand RD, are built from F and D, respectively. In order to generate these matrices, a list is created for each facility(location) that contains the remaining facilities(locations) in non-decreasing(non-increasing) order of flow(distance). Element RF[p][q] is the position of facility q in the list built for facility p. Element RD[p][q] is the locality that is in the q-th position of the list correspondent to location p.

Let ρ(i) and ρ(j) be the free and the fixed vertices, respectively, of edge ef. A list of possible locations for ρ(i), Lposs, is built with basis on RDand element RF[ρ(j)][ρ(i)]. Since j is the occupied location, the first candidate as the new location of facility ρ(i) to be included in Lpossis the element in position RF[ρ(j)][ρ(i)] of the list of locations corresponding to j, that is, element RD[ j][RF[ρ(j)][ρ(i)]]. The other elements of Lpossare RD[ j][RF[ρ(j)][ρ(i) + 1], RD[ j][RF[ρ(j)][ρ(i)  1], ..., RD[ j][RF[ρ(j)][ρ(i) + k], RD[ j][RF[ρ(j)][ρ(i)  k]. If it is the case that ρ(i)  k = 0 and ρ(i) + k < n  1 or ρ(i)] k > 0 and ρ(i) + k = n  1, Lpossis completed with elements RD[ j][RF[ρ( j)][ρ(i) + k + 1], ..., RD[ j][RF[ρ(j)][n  1], or RD[ j][RF[ρ(j)][ρ(i)  k], ..., RD[ j][RF[ρ(j)][1], respectively. Therefore, Lpossis mounted with the locations with best chance of being a good location for ρ(i) according to the values of the corresponding edges. The algorithm seeks the new location for ρ(i), analyzing the exchange of ρ(i) and ρ(h), h being the current location examined of list Lposs. In the algorithm described in this paper, a first improvement rule was used. Let bestloc be the best location found for ρ(i). Then a new solution ρ1 with cost c(ρ1) is induced by the interchange of vertices ρ(i) and ρ(bestloc). The same analysis is done for vertex ρ(j) resulting in solution ρ2 with cost c(ρ2). The assignment that induces the solution with the lowest cost among ρ1 and ρ2 is assumed as the exchange move.

5.2 Two Terminal Vertexes are Set Free  VDS2

In this second variant, denoted as VDS2, both terminal vertices of edge efare re-assigned to other locations. The strategy adopted in VDS1 to define the new location for the free vertex, tends to be time consuming in this case, once two vertices need to be re-located. Thus, in the VDS2 versions, the choice of the distance edge e'd to which efis assigned has a random element. Two versions of VDS2 are investigated. These variants do not use the list of forbidden movements. This occurs because the probability of undoing previous movements is very small due to the randomness on the choice of the new assignments.

The first version, named VDS2A, in order to define e'd, one location is selected at random. Let us consider that edge edhas no facilities assigned to its terminal vertices, since the current facilities will be removed. One facility, fac1, is chosen at random with equal probability to be assigned to one terminal vertex of ed. Let i and j be the terminal vertices of edge ed and suppose that fac1 is assigned to location i. If j is the k-th element in the list of vertex i, regarding the distances between i and the other locations, then the candidate facilities to be assigned to the second free location are close to the k-th element of the corresponding list of fac1. In the algorithm implemented in this paper, the facilities in positions k ± a in the list of fac1 are the candidates to be the facility that will be assigned to location j. In this version, one of these facilities is selected at random with uniform probability. Preliminary tests evaluated values 0, 1 and 2 for variable a. The experiments showed that the results with a = 0 were inferior to both a = 1 and a = 2. A small advantage was observed for a = 2 with no significant difference in processing time. Once the two new facilities are chosen, the algorithm determines the best assignment among the two pair of facilities, the ones that left the locations and the new ones.

In the second VDS2 version, VDS2R, e'd is selected at random among the possible distance edges.

The absence of the best improvement criterion employed in VDS1 reflects on the quality of solutions, but a gain regarding the processing times is obtained. Like in the IT-EC, to guarantee that the solution produced by each algorithm is a local optimum for the 2-exchange neighborhood, a final 2-exchange local search step was included at the end of each ILS iteration.

Tables 5 and 6 present the results of the computational experiments regarding symmetric and asymmetric QAPLIB instances for the iterated local search versions of the three proposed VDS algorithms. The Kruskal-Wallis' statistical test with significance level 0.05 was applied to determine if differences exist between the results presented by the three new algorithms.

The results showed that significant favorable results occur for IT-VDS1 against the other two algorithms on all Bur26 instances. This is also the case for instances Chr20b, Esc128, Sko100a, Tai20a, Tai30a, Tai80a and Tai100b concerning the IT-VDS2A and for instances Lipa80a, Tai40a and Tho150 regarding the IT-VDS2R. The algorithm IT-VDS2A presented significant better results than the algorithm IT-VDS2R on all Bur26 instances and on instances Esc32a, Tai40a and Lipa60a. Better results than the IT-VDS1 were presented by the IT-VDS2A on instances Nug21, Nug24 and Tai100a. Superior results for IT-VDS2R against IT-VDS1 were pointed out on instances Nug24, Sko100d and Lipa30a, and against IT-VDS2A on instances Chr20b, Esc128, Scr20, Sko56, Sko72, Sko100a, Sko100d, Ste36a, Tai80a, Tho40 and Lipa40b.

The test for proportions comparison with significance level 0.05 shows that, except for the pair IT-VDS1 and IT-VDS2R, significant differences do exist among the other pairs of algorithms. Considering the averages presented in Tables 1 and 5, that test also shows that the IT-LPR presents significant worse results than the three proposed VDS algorithms in the symmetric QAPLIB instances. IT-VDS2A is outperformed by IT-VDS1, IT-VDS2R and IT-EC. The algorithms ITVDS1 and ITVDS2R present better results than the IT-EC. The algorithm IT-VDS1 shows a particular good performance on instances Taixxa where the algorithm presents the best results on six of the nine instances. The algorithms IT-VDS2R and IT-VDS2A present 22 and 13 best average results and 13 and 11 best percentages of the best known solution, respectively. Considering the 25 independent executions of each algorithm for each instance of the experiment, the algorithm ITVDS2R obtains, in average, 31.12 % of the best known solutions, followed by IT-VDS1, 29.92% and IT-VDS2A, 28.90%. These three results are better than the ones presented by algorithms IT-LPR and IT-EC.

Considering the average values presented in Table 6, the test for proportions comparison with significance level 0.05 shows that, the algorithm IT-VDS2R presents significantly better results than the other VDS algorithms on the set of asymmetric instances. The test also points out the best performance of IT-VDS2R over IT-LPR and IT-EC and that the latter outperforms the IT-VDS1. A significant difference is also found in favor of IT-VDS1 against IT-VDS2A. With the specified significance level, it is not possible to point out differences between the other results.

The results concerning the number of best known solutions of the asymmetric instances obtained by the proposed algorithms show that, in average, the algorithm ITVDS2R obtains the best performance with 42.22% of the best known solutions, followed by IT-VDS1 with 42.20%, IT-EC with 41.98%, IT-VDS2A with 38.77% and IT-LPR with 25.19%.

5.3 Comparison of the proposed VDS heuristics with an exact algorithm

To illustrate the difference between the results obtained with exact algorithms and those obtained with the methods proposed in this paper, a comparison is presented in Table 7. The comparison is made with an exact algorithm recently proposed by Zhang et al. (2010) to solve sparse QAP instances. Columns Min present the percent difference from the optimal solution and the best solution found by the corresponding algorithm. Column T (s) of the exact algorithm shows the processing time to obtain the optimal solution or the maximum processing time given as stop criterion for the algorithm. An asterisk means that the exact algorithm reached the maximum processing time without solving the corresponding instance. The exact algorithm was executed by its authors on a laptop Intel Pentium M-1.70GHz processor, 1.23Gb RAM. Columns T (s) of the VDS algorithms show the average processing time to find the optimal solution. An asterisk in those columns means that the corresponding heuristic algorithm reached the maximum processing time without finding the optimal solution. Values 0.00 in columns T (s) mean that the average processing time was less than 0.01s. Columns %best show the percentage of best known solutions obtained in the 25 independent executions.

The two ITVDS2 algorithms did not find the optimal solution of the Chr20b instance and no algorithm found the optimal solution of Esc32a instance. The results produced by the exact algorithm on Chr instances are very effective, as shown in Table 7. However on the remaining instances, the heuristic algorithms are always preferable to the exact algorithm. Even on the Esc32e instance, where the exact algorithm takes only 5.9s to solve it, the heuristic algorithms take, in average, less than 0.01s to find the optimal solution and find the optimum in all executions. The exact algorithm cannot find the optimal solution of eight instances. From those instances, the heuristic algorithm misses the optimum only on instance Esc32a, yet the difference to the optimal value is much smaller as well as the processing time. Another point to consider is that the heuristic algorithms find optimal solutions on instances larger than those considered by the exact algorithms currently proposed for the QAP.

5.4 Testing with Drexx instances

Instances with n between 15 and 72 were tested. Table 8 shows averages and percentage of optimal solutions obtained by each algorithm. The maximum processing times in seconds are presented in column T (s).

As shown in Table 8, IT-LPR behaves poorly on this set of instances presenting results significantly inferior to the other algorithms. The Kruskal-Wallis' test performed on the results of the other four algorithms show that, with significance level 0.05, no significant differences were found between their performances on instances Dre from 15 to 30. Significant differences in favor of IT-EC were found on instances Dre42, 56 and 72 in comparison to the other three algorithms. The test for proportions comparison confirms the best performance of IT-EC on Dre instances in comparison to the IT-VDS2R. Evidences of significant differences are not pointed out by the proportions comparison test in pairwise comparisons of the other algorithms.

5.5 Testing with Taixxe instances

Instances with n = 27, 45 and 75 are tested. The maximum processing times given to the algorithms are 10 seconds (n = 27), 30 seconds (n = 45) and 200 seconds (n = 75). Tables 9, 10 and 11 show the results for classes n = 27, 45 and 75, respectively. The first column contains the identification of each test case. The second column shows the optimum or the best known solution of each instance. Columns Md show the percent deviation of the median from the optimum or best known solution. Columns %best show the percentage of optimal or best known solutions found by each algorithm.

With significance level 0.05, the Kruskal-Wallis' test pointed out significant differences on the performance of the algorithms for all instances of class 27, where the inferior performance of IT-LPR was identified. Due to this fact, other comparisons were carried on the results obtained only with the other four algorithms. The success for this set of instances is established on the number of best solutions found by each algorithm. With significance level 0.05, the test for proportion comparison showed that IT-EC and IT-VDS1 are outperformed by IT-VDS2R. In general, IT-VDS1 is outperformed by IT-EC and IT-VDS2A, although the hypothesis test does not point out significant differences. IT-VDS2A and IT-EC performs similarly.

Significant differences between the results presented by the IT-LPR and the other algorithms in class Tai45e instances were also pointed out by the Kruskal-Wallis' test. The Mann-Whitney test with significance level 0.05 indicated that better results were found by algorithm IT-VDS2A against IT-EC on instances 7 and 11, by algorithm IT-EC against IT-VDS2A on instance 9, by IT-VDS1 against IT-EC on instance 9, and by algorithm IT-VDS2R against IT-VDS2A on instances 6 and 14. In average, the best performance regarding the number of solutions with cost equals to BKS is obtained by the IT-VDS1, followed by the IT-VDS2A, IT-EC and IT-VDS2R as shown in Table 10. The test for proportions comparison indicated that IT-VDS1 performs better than IT-VDS2R in the Tai45e instances with confidence level 93.38%. The test did not detect significant differences among the performances of the remaining algorithms.

The Kruskal-Wallis' test pointed out significant differences between the results presented by at least one algorithm on class Tai75e instances. The Mann-Whitney U-test showed that the IT-LPR performed worse than the IT-VDS1 and IT-VDS2R on 13 instances, the IT-EC on 14 instances and the IT-VDS2A on 17 instances. The test for proportions comparison considering the medians also showed significant differences in favor of IT-VDS1, IT-VDS2A, IT-VDS2R and IT-EC against IT-LPR. At significance level 0.05, the test for proportions comparison showed that IT-EC, IT-VDS2R and IT-VDS2A are outperformed by IT-VDS1. The test does not indicate significant differences among the other three algorithms.

In general, it was observed that the relative performance of IT-VDS1 improves for growing sizes of Taixxe instances, being better than IT-LPR on all sets and better than the IT-EC on the set n = 75. Moreover, the relative performance of the IT-VDS2R worsens as the sizes of the instances grow. This fact shows that the systematic search produces better results on these classes for increasing sizes of instances. All VDS algorithms performed better than the IT-LPR.

6 CONCLUDING REMARKS

This paper presented an extensive experimental investigation that compared local search and variable depth search algorithms for the Quadratic Assignment Problem. Moreover, variable depth search algorithms that use information concerning the cost of the edges assignments were presented. This study contributes to a better understanding of the potential of these techniques, which are widely used as intensification tools in more sophisticated heuristic methods, such as evolutionary algorithms. The experiments were performed with 161 benchmark QAP instances belonging to classes designed with different structures. Nonparametric statistical tests were applied to the data generated by the computational experiments to support conclusion about the performance of the algorithms.

Existing VDS-like algorithms for the QAP (Li et al., 1994; Rego et al., 2006) were implemented under an iterated local search framework and submitted to a computational experiment more extensive than the ones reported in their original papers. These algorithms were compared with recent approaches that use the classical 2-exchange neighborhood (Stützle, 2006) and k-exchange neighborhoods Ahuja et al. (2007). They were also compared to the VDS algorithms proposed in this paper.

The experiments showed that, given the same computational effort, variable depth search is always preferable to pure local search approaches. The test showed that both previously proposed variable depth search methods outperformed LSMC (Stützle, 2006) and VLSN (Ahuja et al.,2007). Among these four methods, the best performance was obtained by the iterated local search version of the algorithm presented by Rego et al. (2006).

The paper introduced VDS like algorithms that consider the cost of edge assignments as a criterion to decide which elements are exchanged while searching the space of solutions of QAP instances. The experiments showed that the proposed criterion was beneficial for several classes of benchmark instances, indicating that the approaches introduced here are preferable to the existing ones on those classes. Statistical tests showed that the proposed algorithms are always preferable to the approach presented by Li et al. (1994). Those tests also confirmed that the IT-VDS1 version is better than an ILS version of the algorithm presented by Rego et al. (2006) on QAPLIB classes Ste, Chr, Taixxa, Nug and Tai75e being outperformed by the latter on Lipaxxa class.

Future works will investigate statistical indicators in the choice of edges to be exchanged in variable depth search algorithms.

ACKNOWLEDGEMENTS

We thank anonymous reviewers for their suggestions and useful remarks which contributed to improve the paper. We also acknowledge the support of CNPq to this research under grants 303538/200-2 and 300778/2010-4.