Nouredine Melab

T-V Luong

Karima Boufaras

E-G Talbi

The-Van Luong

Towards ParadisEO-MO-GPU: a Framework for GPU-based Local Search Metaheuristics

Keywords: Software Framework, Local Search Meta-heuristics, Parallel Computing, GPU Computing, Neighborhood Exploration

This paper is a major step towards a pioneering software framework for the reusable design and implementation of parallel metaheuristics on Graphics Processing Units (GPU). The objective is to revisit the ParadisEO framework to allow its utilization on GPU accelerators. The focus is on local search metaheuristics and the parallel exploration of their neighborhood. The challenge is to make the GPU as transparent as possible for the user. The first release of the new GPU-based ParadisEO framework has been experimented on the Quadratic Assignment Problem (QAP). The preliminary results are convincing, both in terms of flexibility and easiness of reuse at implementation, and in terms of efficiency at execution on GPU.

Introduction

Nowadays, parallel metaheuristics have grown to be a highly useful paradigm to solve large-scale CPU time-intensive and complex combinatorial problems. Metaheuristics are either single-solution namely S-Metaheuristics (local search metaheuristics) or population-based namely P-Metaheuristics (e.g. evolutionary algorithms). The focus in this paper is on S-Metaheuristics. Recently, GPU accelerators have emerged as a new powerful support for massively parallel computing. Last year, we came up with the pioneering work on GPU-based S-Metaheuristics [START_REF] Luong | Local Search Algorithms on Graphics Processing Units. A Case Study: the Permutation Perceptron Problem[END_REF]. Such experience has shown that parallel combinatorial optimization on GPU is not straightforward, and requires a huge effort at design as well as at implementation level.

Indeed, the design of GPU-aware S-Metaheuristics often involves the cost of a sometimes painful apprenticeship of parallelization techniques and GPU computing technologies. In order to free from such burden those who are unfamiliar with those advanced features, optimization frameworks must integrate the up-to-date parallelization techniques and allow their transparent exploitation and deployment on GPU accelerators. To the best of our knowledge, there does not exist any software framework for GPU-based metaheuristics. In [START_REF] Cahon | ParadisEO: A Framework for the Reusable Design of Parallel and Distributed Metaheuristics[END_REF], we have proposed a framework called ParadisEO dedicated to the reusable design of parallel and distributed metaheuristics for only dedicated parallel hardware platforms. Later, we have extended the framework in [START_REF] Melab | Grid computing for parallel bioinspired algorithms[END_REF] to dynamic and heterogeneous large-scale environments using Condor-MW middleware and in [START_REF] Tantar | Building a Virtual Globus Grid in a Reconfigurable Environment -A case study: Grid[END_REF] to computational grids using Globus.

In this paper, we extend ParadisEO-MO (ParadisEO for S-Metaheuristics) to deal with GPU accelerators. The challenges and contributions consist in [START_REF] Luong | Local Search Algorithms on Graphics Processing Units. A Case Study: the Permutation Perceptron Problem[END_REF] rethinking the parallel models provided into the framework to manage efficiently the hierarchical organization of the memories (different latencies and sizes) of the GPU device as well as the interaction of this latter with the CPU ; (2) making the GPU as transparent as possible for the user minimizing his or her involvement in its management. In this paper, we propose solutions to this challenge as an extension of the ParadisEO framework. The focus is on the iteration-level parallel model of S-Metaheuristics which consists in exploring in parallel the neighborhood of a problem solution. The first release of the new GPU-based ParadisEO framework has been implemented using C++ and CUDA [START_REF] Nickolls | Scalable Parallel Programming with CUDA[END_REF] and then experimented on the QAP.

The remainder of the paper is organized as follows. Section 2 highlights the principles of parallel iteration-level S-Metaheuristics and their challenges when using GPU computing. In Section 3, we describe the major design features and architecture of ParadisEO. Section 4 presents the design and implementation of ParadisEO-MO on top of GPU called ParadisEO-MO-GPU. Section 5 shows and comments some experimental results obtained with ParadisEO-MO-GPU on the QAP. In Section 6, we conclude the paper and draw some perspectives of the presented work.

Parallel GPU-based S-Metaheuristics

An S-Metaheuristic is an iterative procedure which explores the neighborhood of a solution in order to improve its quality. The associated algorithm generates an initial (current) solution to the problem to be solved. This current solution is evaluated and its neighborhood is generated and evaluated. Based on the evaluation of the neighborhood, the best solution is selected to become the current solution. The process is repeated until a stopping criterion is found.

For large-scale combinatorial optimization problems, the neighborhood of a solution is often extremely large. Therefore, massively parallel computing is required to generate and evaluate it. The parallel generation and evaluation of the neighborhood is a master-worker and problem independent regular data-parallel application. GPU computing is very well-suited for this kind of parallel application. In the GPU (CUDA-based) model, the master is the CPU and the workers are threads executed by the processing cores of the GPU. Using GPU computing is not straightforward especially for non-experts in parallel computing. Indeed, a GPU accelerator provides a hierarchy of memories with different sizes and access latencies.

The challenge is to re-think the design of the parallel exploration and evaluation of the neighborhood taking into account the characteristics of GPU. Different issues have to be dealt with: (1) defining an efficient cooperation between CPU and GPU, which requires to share the work and to optimize the data transfer between the two components; (2) GPU computing is based on hyperthreading (massively parallel multi-threading) and the order in which the threads are executed is not known. Therefore, an efficient mapping has to be defined between each neighboring candidate solution and a thread designated by a unique identifier assigned by the GPU runtime; (3) the neighborhood has to be placed efficiently on the different memories taking into account their sizes and access latencies. From an implementation point of view, the challenge is to provide solutions to these issues in ParadisEO as transparent as possible way for the user.

3 ParadisEO-MO-based Parallel S-Metaheuristics

Parallel iteration-level model on GPU

The parallel iteration-level model is designed according to the data-parallel SPMD model of CUDA. In this model, a function code called the kernel is sent to the GPU to be executed by a large number of threads grouped into blocks. The task partitioning is such that the CPU hosts executes the whole serial part of the local search method. The GPU is in charge of the evaluation of the neighborhood of the current solution at each iteration. In order to minimize the cost of the data transfer from the CPU to GPU, the neighboring solutions are generated on GPU rather than on CPU. Indeed, only the current solution is sent to the GPU and each thread executes the same kernel. This is highly efficient for large neighborhoods. The kernel consists in generating and evaluating a neighbor. A defined mapping function allows each thread to find its corresponding neighboring solution. Once all the neighboring solutions are generated and evaluated they are sent back to the CPU where the best solution is selected. The process is iterated until a stopping criterion is satisfied. LGPL extensible open source framework based on a clear conceptual separation of the metaheuristics. ParadisEO is one of the rare frameworks that provide the most common parallel and distributed models. These models are portable on distributed-memory machines and shared-memory multi-processors as they are implemented using standard libraries such as MPI, PVM and Pthreads. The user level indicates the different problem-dependent components which must be defined: input data, the evaluation function, neighbor and neighborhood representations. The second level presents the ParadisEO-MO framework including optimization solvers embedding S-metaheuristics. The interaction is done with the ParadisEO-GPU module which automatically pilots the CUDA programming interface. The hardware level supplies the different transparent tools provided by ParadisEO-GPU such as the allocation and copy of data or the parallel generation and evaluation of the considered neighborhood. In addition to this, the platform proposes predefined neighborhood and mapping wrappers adapted to hardware constraints to deal with binary and permutation problems.

GPU-enabled ParadisEO

A case study: parallel evaluation of a neighborhood

ParadisEO-MO-GPU is illustrated in Fig. 2 through an UML sequence diagram. The scenario shows the design and implementation of parallel neighborhood evaluation of a local search on GPU. The different steps of the parallel evaluation process on GPU for each iteration are the following: 1. The neighborhood moCudaNeighborhood prepares all the steps for the parallel generation of the neighborhood on GPU. The initialization consists in setting a mapping table between GPU threads and neighbors. Then, the associated data are sent only once to the GPU since the mapping structure does not change during the execution process of local searches. The last step invokes the parallel evaluation and will be informed on its completion to retrieve the precomputed fitnesses structure. 2. Before proceeding to the parallel evaluation, the object moCudaEval configures a kernel with N threads such that each thread is associated exactly with one neighbor evaluation (N designates the neighborhood size). During the first iteration, the object allocates the neighborhood fitnesses structure where the result of the evaluated neighbors will be stored. Otherwise, in any case, it only sends to the GPU device the candidate solution which generates the neighborhood. 3. The object moCudaKernelEval modelizes the main body which will be executed by N concurrent threads on different input data. A first step consists in getting the thread identifier then the set of its associated data. This mechanism is done thanks to the mapping table previously mentioned. The second step performs the evaluation computation of the corresponding neighbor. Finally, the resulting fitness is stored in the corresponding index of the fitnesses structure. 4. The worker moCudaEvalFunc is the specific object with computes on the GPU device the corresponding evaluation neighbor and returns the provided result.

Once the entire neighborhood has been performed in parallel on GPU, the precalculated fitness structure is copied back to the CPU and given as input to the ParadisEO-MO module. This way, the local search continues the neighborhood exploration (iteration) on the CPU side. Instead of evaluating again each neighbor in a sequential manner, the corresponding fitness value will only need to be retrieved from the precomputed fitnesses structure. Hence, this mechanism has the advantage of allowing both the deployment of any S-metaheuristics and the use of toolboxes provided in ParadisEO-MO (e.g. or fitnesses landscape analysis, checkpoint monitors, etc.)

Experimentation

The QAP arises in many applications such as facility location or data analysis. Let A = (a ij ) and B = (b ij ) be n × n matrices of positive integers. Finding a solution of the QAP is equivalent to finding a permutation π = (1, 2, . . . , n) that minimizes the objective function:

z(π) = n i=1 n j=1 a ij b π(i)π(j)
As the iteration-level parallel model does not change the semantics of the sequential algorithm, the effectiveness in terms of quality of solutions is not addressed here. The objective is to assess the impact in terms of efficiency of an implementation done with ParadisEO-MO-GPU compared with an optimized version made outside the platform. To achieve this, a tabu search has been implemented in 4 different versions: 1) a ParadisEO-MO implementation on CPU and its counterpart on GPU; 2) an optimized CPU implementation and its associated GPU version. ParadisEO versions are pure object-based implementations whereas the optimized ones are pointer-based.

The neighborhood is built by pair-wise exchange operations (known as the 2exchange neighborhood) which is a standard way for permutation problems. The number of global iterations of the local search is set to 10000. The considered instances are the Taillard instances proposed in [START_REF] Taillard | Robust tabu search for the quadratic assignment problem[END_REF].

Experiments have been carried on top of an Intel Core i7 970 3.2Ghz with a GTX 480 card (15 multiprocessors with 32 cores). For measuring the acceleration factors, only a single core has been considered using the Intel i7 turbo mode (3.46Ghz). The average time has been measured in seconds for 30 runs. Results are reported in Table 1 for the different versions. From the instance tai40a, both GPU versions start to give positive accelerations (from ×1.6 to ×1.7). The poor performance for small instances can be explained by the fact that since the neighborhood is relatively small, the number of threads per block is not enough to fully cover the memory access latency. However, as long as the instance size increases, the acceleration factor grows accordingly (e.g. from ×5.4 to ×5.8 for tai80a). Finally, efficient speed-ups are obtained for the instance tai150b. They vary between ×13.4 and ×15.1.

A thorough examination of the acceleration factors highlights that they are quite similar for the two GPU versions. The performance difference which exists is certainly due to the CPU version provided by ParadisEO-MO. Indeed, regarding the two CPU versions, initially, there is already a slight performance gap regarding the execution time. It varies between 11% and 19% according to the instance. This gap can be explained by the overhead caused by the creation of generic objects in ParadisEO whereas the optimized version on CPU is a pure pointer-based implementation. This can also explain the performance difference between the two different GPU counterparts in which the same phenomenon occurs. However, for such transparent exploitation, the obtained results are really convincing. A conclusion of the experiments indicates that the performance results of the GPU version provided by ParadisEO are not so far from the GPU pointer-based one.

In this paper, we have presented a step towards a ParadisEO framework for the reusable design and implementation of the GPU-based parallel metaheuristics. The focus is set on S-metaheuristics and the iteration-level parallel exploration of the neighborhood of a solution. We have revisited the design and implementation of this last model in ParadisEO-MO to allow its efficient execution and its transparent use on GPU.

The implementation in ParadisEO-MO using CUDA has been experimented on the QAP and compared to the same implementation realized outside Par-adisEO. The preliminary results show that the performance gap which occurs between the two implementations is not really important. Indeed, for such flexibility and easiness of reuse at implementation, the obtained results with ParadisEO-MO-GPU are really promising (up to ×13). We are strongly convinced that the overall performance provided by ParadisEO-MO-GPU will be much better for larger neighborhoods or other problems requiring more computational calculations (see [START_REF] Luong | Parallel hybrid evolutionary algorithms on gpu[END_REF] and [START_REF] Luong | Local Search Algorithms on Graphics Processing Units. A Case Study: the Permutation Perceptron Problem[END_REF] for previous test cases).

The first release of ParadisEO-MO on GPU is available on the ParadisEO website 1 . Tutorials and documentation are provided to facilitate its reuse. This release is dedicated to parallel S-metaheuristics based on the iteration-level parallel model. In the short run, this release will be first extended to the algorithmic (multi-start) and solution-level parallel models. Second, it will be extended to other problem representations such as discrete representation and other solution methods. Third, it will be validated on a wider range of problems. In the long run, ParadisEO will be revisited following the same roadmap for evolutionary algorithms on GPU.

3. 2

 2 The ParadisEO-MO framework ParadisEO-MO is part of ParadisEO dedicated to S-Metaheuristics such as Hill Climbing, Simulated Annealing, Tabu Search, ILS, etc. ParadisEO [2] is a framework dedicated to the reusable design of parallel hybrid metaheuristics by providing a broad range of features including evolutionary algorithms (ParadisEO-EO), local search methods (ParadisEO-MO), parallel and distributed models (ParadisEO-PEO), different hybridization mechanisms, etc. ParadisEO is a C++

4. 1

 1 Architecture ofParadisEO-MO-GPU is a new framework which is a coupling between ParadisEO-MO and CUDA. It aims at deploying the S-Metaheuristics on GPU in a generic way. It is composed by a set of new C++ abstract and predefined classes that allows an easy and transparent development of S-metaheuristics on GPU accelerators. The architecture of ParadisEO-MO-GPU is layered as it is illustrated in Fig.1.

Fig. 1 .

 1 Fig. 1. A layered architecture of ParadisEO-MO-GPU.

Fig. 2 .

 2 Fig. 2. The parallel generation and evaluation of a neighborhood provided in ParadisEO-MO-GPU.

Table 1 .

 1 Measures in terms of efficiency of a ParadisEO-MO-GPU implementation with an optimized version made outside the platform.

	Instance	ParadisEO-MO Optimized version Perf. gap CPU GPU Acc. CPU GPU Acc. CPU GPU
	tai30a	0.8 0.9 ×0.9 0.7 0.7	×1.0	13% 23%
	tai40a	1.9 1.2 ×1.6 1.7 1.0	×1.7	11% 17%
	tai50a	3.6 1.6 ×2.2 3.0 1.3	×2.3	17% 19%
	tai60a	6.0 2.0 ×3.0 5.0 1.6	×3.1	17% 20%
	tai80a 15.1 2.8 ×5.4 12.3 2.1	×5.8	19% 25%
	tai100a 31.7 3.8 ×8.3 26.0 2.8	×9.2	18% 27%
	tai150b 119.7 8.9 ×13.4 98.1 6.5	×15.1 18% 27%

http://paradiseo.gforge.inria.fr


GPU-based Approaches for Multiobjective Local Search Algorithms. A Case Study: the Flowshop Scheduling Problem

Multiobjective local search algorithms are efficient methods to solve complex problems in science and industry. Even if these heuristics allow to significantly reduce the computational time of the solution search space exploration, this latter cost remains exorbitant when very large problem instances are to be solved. As a result, the use of GPU computing has been recently revealed as an efficient way to accelerate the search process. This paper presents a new methodology to design and implement efficiently GPU-based multiobjective local search algorithms. The experimental results show that the approach is promising especially for large problem instances.

Introduction

Real-world optimization problems are often complex and NP-hard, their modeling is continuously evolving in terms of constraints and objectives, and their resolution is time-consuming. Although near-optimal algorithms such as metaheuristics allow to reduce the temporal complexity of their resolution, they are unsatisfactory to tackle large problems.

Nowadays, GPU computing is recognized as a powerful way to achieve highperformance on long-running scientific applications [START_REF] Ryoo | Program optimization carving for gpu computing[END_REF]. Designing multiobjective local search (MLS) algorithms for solving real-world optimization problems are good challenges for GPU computing. However, only few research works related to evolutionary algorithms on GPU for monoobjective optimization exist [START_REF] Li | An efficient fine-grained parallel genetic algorithm based on gpu-accelerated[END_REF][START_REF] Chitty | A data parallel approach to genetic programming using programmable graphics hardware[END_REF][START_REF] Wong | Parallel evolutionary algorithms on consumer-level graphics processing unit[END_REF][START_REF] Fok | Evolutionary computing on consumer graphics hardware[END_REF]. Regarding existing works on LS algorithms on GPU, to the best of our knowledge, only coarse-grained implementations are provided in the context of the multi-start tabu search for monoobjective optimization [START_REF] Zhu | Simd tabu search with graphics hardware acceleration on the quadratic assignment problem[END_REF][START_REF] Janiak | Tabu search on gpu[END_REF].

Indeed, a finer-grained model such that the parallel exploration of the neighborhood on GPU is not immediate and several challenges persist and are particular related to the characteristics and underlined issues of the GPU architecture and the MLS algorithms. The major issues are the efficient distribution of data processing between CPU and GPU, the thread synchronization, the optimization of data transfer between the different memories, the capacity constraints of these memories, etc.

The main objective of this paper is to deal with such issues for the re-design of parallel MLS algorithms to allow solving of large scale optimization problems on GPU architectures. In this paper, we contribute with the first results of multiobjective local search algorithms on GPU. More exactly, we propose some new GPU-based approaches for building the parallel exploration of the neighborhood on GPU in a multiobjective context. These approaches are based on a decomposition of the GPU hierarchy allowing a clear separation between generic and problem-dependent LS features. Several challenges are dealt with: (1) the distribution of the search process among the CPU and the GPU minimizing the data transfer between them; (2) finding the efficient mapping of the hierarchical LS parallel models on the hierarchical GPU; (3) using efficiently the coalescing and texture memory in the context of MLS algorithms.

To validate the approaches presented in this paper, the flowshop scheduling problem (FSP) [START_REF] Taillard | Benchmarks for basic scheduling problems[END_REF] have been considered and implemented on GPU.

The remainder of the paper is organized as follows: Section 2 highlights the principles of MLS methods and their parallel models. In Section 3, for a better understanding of the difficulties of using the GPU architecture, its characteristics are described according to a decomposition of the GPU hierarchy. Section 4 presents generic concepts for designing parallel MLS methods on GPU. In Section 5, on the one hand, efficient mappings between state-of-the-art LS structures and GPU threads model are performed. On the other hand, a depth look on the GPU memory management adapted to MLS heuristics is depicted. Section 6 reports the performance results obtained for the implemented problem mentioned above. Finally, a discussion and some conclusions of this work are drawn in Section 7.

Parallel Multiobjective Local Search Algorithms

Multiobjective Local Search algorithms

The existing LS methods that intend to find an approximation of the Pareto optimal set of a multiobjective optimization problem fall into two categories: scalar approaches and Pareto approaches.

Approaches belonging to the first class contains the approaches that transform a multiobjective problem into a monoobjective one or a set of such problems. Many proposed algorithms in the literature are scalar approaches. Among these methods one can find the aggregation methods, the weighted metrics, the -constraint methods . . . A review of these methods is given in [START_REF] Talbi | Metaheuristics: From design to implementation[END_REF].

The second class consists in defining the acceptance of the LS according to a dominance relationship such as the Pareto dominance. The idea of Pareto approaches is to maintain an archive of non-dominated solutions, to explore the neighborhood of the solutions contained in the archive and to update the archive with the visited solutions. A complete description of the different algorithms can be found in [START_REF] Paquete | Stochastic Local Search Algorithms for Multiobjective Combinatorial Optimization: Methods And Analysis[END_REF].

Parallel Models of Local Search Algorithms

For non-trivial problems, executing the iterative process of a MLS on large neighborhoods requires a large amount of computational resources. In general, evaluating a fitness function for each solution is frequently the most costly operation of the MLS. Consequently, a variety of algorithmic issues are being studied to design efficient MLS heuristics. Parallelism arises naturally when dealing with a neighborhood, since each of the solutions belonging to it is an independent unit. Parallel design and implementation of metaheuristics have been studied as well on different architectures [START_REF] Alba | 4. Metaheuristics and Parallelism[END_REF][START_REF] Melab | Grid computing for parallel bioinspired algorithms[END_REF].

Basically, three major parallel models for LS heuristics can be distinguished: solution-level, iteration-level and algorithmic-level.

• Solution-level Parallel Model. A focus is made on the parallel evaluation of a single solution. That model is particularly interesting when the evaluation function can be itself parallelized as it is CPU time-consuming and/or IO intensive. In that case, the function can be viewed as an aggregation of a given number of partial functions. • Iteration-level Parallel Model. This model is a low-level Master-Worker model that does not alter the behavior of the heuristic. Exploration and evaluation of the neighborhood are made in parallel. At the beginning of each iteration, the master duplicates the current solution between parallel nodes. Each of them manages some candidates and the results are returned back to the master. • Algorithmic-level Parallel Model. Several LS algorithms are simultaneously launched for computing better and robust solutions. They may be heterogeneous or homogeneous, independent or cooperative, start from the same or different solution(s), configured with the same or different parameters.

The solution-level model is problem-dependent and does not present many generic concepts. In this paper, we focus only on the fine-grained problemindependent model: the iteration-level. Indeed, unlike the algorithmic-level, the iteration-level can be seen as an acceleration model which does not change the semantics of the algorithm.

GPU Computing for Metaheuristics

Driven by the demand for high-definition 3D graphics on personal computers, GPUs have evolved into a highly parallel, multithreaded and many-core environment. Indeed, this architecture provides tremendous computational horsepower and very high memory bandwidth compared to traditional CPUs. Since more transistors are devoted to data processing rather than data caching and flow control, GPU is specialized for compute-intensive and highly parallel computation. A complete review of GPU architectures can be found in [START_REF] Ryoo | Program optimization carving for gpu computing[END_REF].

The adaptation of MLS algorithms on GPU requires to take into account at the same time the characteristics and underlined issues of the GPU architecture and the LS parallel models. Mapping existing parallel models to the GPU in an efficient way involves a clear understanding of GPU characteristics. In this section, we propose a decomposition of the GPU adapted to the parallel iterationlevel model allowing to identify the different challenges that must be dealt with.

General GPU Model

In general-purpose computing on graphics processing units, the CPU is considered as a host and the GPU is used as a device coprocessor. This way, each GPU has its own memory and processing elements that are separate from the host computer. Data must be transferred between the memory space of the host and the memory of GPU during the execution of programs.

Each processor device on GPU supports the single program multiple data (SPMD) model, i.e. multiple autonomous processors simultaneously execute the same program on different data. For achieving this, the concept of kernel is defined. The kernel is a function callable from the host and executed on the specified device simultaneously by several processors in parallel.

Regarding the iteration-level parallel model, since the evaluation of neighboring candidates is often the most time-consuming part of MLSs, it must be done in parallel on GPU. Therefore, according the the Master-Worker paradigm, a kernel is associated with the evaluation of the neighborhood and the CPU controls the whole sequential part of the LS process.

However, memory transfer from the CPU to the device memory is a synchronous operation which is time consuming. Indeed, bus bandwidth and latency between the CPU and the GPU can significantly decrease the performance of the search. As a result, one of the challenges is to optimize the data transfer from CPU to GPU.

Parallelism Control: GPU threads model

The kernel handling is dependent of the general-purpose language. For instance, CUDA or OpenCL are parallel computing environments which provide an application programming interface for GPU architectures [13] [START_REF] Group | OpenCL 1.0 Quick Reference Card[END_REF]. Indeed, these toolkits introduce a model of threads which provides an easy abstraction for SIMD architecture. The concept of a GPU thread does not have exactly the same meaning as a CPU thread. Compared to CPU threads, GPU threads are lightweight. That means that changing the context between two threads is not a costly operation.

Regarding their spatial organization, threads are organized within so called thread blocks. A kernel is executed by multiple equally threaded blocks. Blocks can be organized into a one-dimensional or two-dimensional grid of thread blocks, and threads inside a block are grouped in a similar way. All the threads belonging to the same thread block will be assigned as a group to a single multiprocessor. Thereby, a unique id can be deduced for each thread to perform computation on different data.

Regarding MLS algorithms, a move which represents a particular neighbor candidate solution can also be associated with a unique id. However, according to the solution representation of the problem, finding a corresponding id for each move is not straightforward. As a consequence, another challenging issue is to find an efficient mapping between GPU threads and LS moves.

Memory Management: Kernel Management

From a hardware point of view, graphics cards consist of streaming multiprocessors, each with processing units, registers and on-chip memory. Since multiprocessors are used according to the SPMD model, threads share the same code and have access to different memory areas.

Communication between the CPU host and its device is done through the global memory. Since this memory is not cached and its access is slow, one needs to minimize accesses to global memory (read/write operations) and reuse data within the local multiprocessor memories. Graphics cards provide also read-only texture memory to accelerate operations such as 2D or 3D mapping. Constant memory is read only from kernels and is hardware optimized for the case where all threads read the same location. Shared memory is a fast memory located on the multiprocessors and shared by threads of each thread block. This memory area provides a way for threads to communicate within the same block. Registers among streaming processors are private to an individual thread, they constitute fast access memory.

The memory management on GPU is problem-specific. A clear understanding of the characteristics described above is required to provide an efficient implementation of the evaluation function. In other words, a last challenge consists in finding the association of the different LS structures with the different available memories to obtain the best performances.

Design of Local Search Algorithms on GPU

In this section, the focus is on the re-design of the iteration-level parallel model. Designing parallel LS model is a great challenge as nowadays there is no generic GPU-based MLS algorithms to the best of our knowledge.

A General Model for LS Algorithms

According to the SPMD model, multiple autonomous processors simultaneously execute the same program at independent points. Therefore, the mapping between the LS iteration-level parallel model and the GPU model becomes quiet natural. First, the CPU sends the number of expected running threads to the GPU, then candidate neighbors are evaluated on GPU, and finally newly evaluated solutions are returned back to the host. This model can be seen as a cooperative model where the GPU is used as a coprocessor in a synchronous manner. The resource-consuming part i.e. the generation and evaluation kernel, is calculated by the GPU and the rest is handled by the CPU.

Generation of the Neighborhood

As quoted above, CPU/GPU communication might be a major bottleneck in the performance of GPU applications. One of the major challenges is to optimize the data transfer between the CPU and the GPU. In other words, one of the issues which remains is to say where the neighborhood must be generated. For doing that, basically there are two different approaches:

• Generation of the neighborhood on CPU and its evaluation on GPU. At each iteration of the LS process, the neighborhood is generated on the CPU side and its associated structure storing the solutions is copied on GPU. This approach is the most straightforward since the mapping between a thread and a neighbor is direct. It has been widely used in the context of evolutionary algorithms on GPU [START_REF] Li | An efficient fine-grained parallel genetic algorithm based on gpu-accelerated[END_REF][START_REF] Chitty | A data parallel approach to genetic programming using programmable graphics hardware[END_REF][START_REF] Wong | Parallel evolutionary algorithms on consumer-level graphics processing unit[END_REF][START_REF] Fok | Evolutionary computing on consumer graphics hardware[END_REF]. • Generation of the neighborhood and its evaluation on GPU. In the second approach, the neighborhood is generated on GPU. This generation is performed in a dynamic manner which implies that no explicit structure needs to be allocated. For doing that, a neighbor is considered as a slight variation of the candidate solution which generates the neighborhood. Thereby, only the representation of this candidate solution must be copied from the CPU to the GPU. The advantage of such approach is to drastically reduce the data transfers since the whole neighborhood does not have to be copied. However, as previously said, finding a mapping between a thread and a neighbor might be challenging.

Even if the first approach is the most natural one, applying it on MLS on GPU will naturally result in a lot of data transfers for large neighborhoods, leading to a great loss of performance. That is the reason why, in the rest of the paper we will consider the second approach. An analysis of the data transfers for the two different approaches might be done but it is not the scope of this paper.

The Proposed GPU-based Algorithm

Adapting traditional MLS methods to GPU is not a straightforward task because the hierarchical memory management on GPU has to be handled. We propose (see algorithm 1) a methodology to adapt MLS methods on GPU in a generic way. The given template is applicable to most of scalar and Pareto approaches.

First of all, memory allocations on GPU are made: data inputs and candidate solution of the problem must be allocated (lines 4 and 5). Since GPUs require massive computations with predictable memory accesses, a structure has to be allocated for storing the results of the evaluation of each neighbor (neighborhood fitnesses structure) at different addresses (line 6). In the case of Pareto approaches, this structure can represent different objective vectors. Additional solution structures which are problem-dependent can also be allocated to facilitate the computation of the evaluation function (line 7). Second, problem data inputs, initial candidate solution and additional structures associated with this solution have to be copied on the GPU (lines 8 to 10). It is important to notice that problem data inputs (e.g. a matrix in TSP) are a read-only structure and never change during all the execution of LS algorithms. Therefore, their associated memory are copied only once during all the execution. Third, comes the for each neighbor in parallel on the GPU kernel do 13:

Complete or delta evaluation of the candidate solution 14:

Insert the resulting fitness into the neighborhood fitnesses structure 15:

end for 16:

Copy the neighborhood fitnesses structure on CPU host memory 17:

Specific MLS solution selection strategy on the neighborhood fitnesses structure 18:

Specific MLS post-treatment 19:

Copy the chosen solution on GPU device memory 20:

Copy additional solution structures on GPU device memory 21: until a stopping criterion satisfied parallel iteration-level, in which each neighboring solution is generated, evaluated and copied into the neighborhood fitnesses structure (from lines 12 to 15). Fourth, since the order in which candidate neighbors are evaluated is undefined, the neighborhood fitnesses structure has to be copied to the host CPU (line 16). Then a specific LS solution selection strategy is applied to this structure (lines 17 and 18). For instance, the archiving of non-dominated solutions might be done in the case of Pareto approaches. Finally, after a new candidate has been selected, this latter and its additional structures are copied to the GPU (lines 19 and 20). The process is repeated until a stopping criterion is satisfied.

Neighborhood Generation and Memory Management

Efficient Mappings of Neighborhood Structures on GPU

Since the generation of the neighborhood is done on GPU to reduce multiple data transfers, the issue is to say which solution must be handled by which thread. The answer is dependent of the solution representation. Indeed, the neighborhood structure strongly depends on the target optimization problem representation. In the following, we provide a methodology to deal with different structures of the literature.

Binary Representation: A solution is coded as a vector of bits. The neighborhood representation for binary problems is based on the Hamming distance where a given solution is obtained by flipping one bit of the solution.

A mapping between LS neighborhood encoding and GPU threads is quiet trivial. Indeed, on the one hand, for a binary vector of size n, the size of the neighborhood is exactly n. On the other hand, threads are provided with a unique id. As a result, a IN → IN mapping is straightforward.

Discrete Vector Representation: This is an extension of binary encoding using a given alphabet Σ. In this representation, each variable takes its value over the alphabet Σ. Assume that the cardinality of the alphabet Σ is k, the size of the neighborhood is (k -1) × n for a discrete vector of size n.

Let id be the identity of the thread corresponding to a given solution of the neighborhood. Compared to the initial candidate solution, id/(k -1) represents the position which differs from the initial solution and id%(k -1) is the available value from the ordered alphabet Σ (both using zero-index based numbering). As a consequence, a IN → IN mapping is possible.

Permutation Representation: Building a neighborhood by pairwise exchange operations is a standard way for permutation problems. Unlike the previous representations, a mapping between a neighbor and a GPU thread is not straightforward. Indeed, on the one hand, a neighbor is composed by two element indexes. On the other hand, threads are identified by a unique id. As a result, a IN → IN× IN mapping has to be considered to transform one index into two ones. In a similar way, a IN × IN → IN mapping is required to transform two indexes into one.

Proposition 1. Two-to-one index transformation

Given i and j the indexes of two elements to be exchanged in the permutation representation, the corresponding index f (i, j) in the neighborhood representation is equal to i × (n -1) + (j -1) -i×(i+1)

2

, where n is the permutation size.

Proposition 2. One-to-two index transformation

Given f (i, j) the index of the element in the neighborhood representation, the

corresponding index i is equal to n -2 - √ 8×(m-f (i,j)-1)+1-1 2
and j is equal to f (i, j) -i × (n -1) + i×(i+1) 2 + 1 in the permutation representation, where n is the permutation size and m the neighborhood size.

The proofs of these two index transformations can be found in [START_REF] Luong | Large neighborhood for local search algorithms[END_REF]. Notice that for different neighborhood operators such as 2-opt or the insertion operator, a slight variation of these mappings is easily applicable.

Memory Management of Local Search Algorithms on GPU

In this section, the focus is on the memory management. Understanding the GPU memory organization and issues is useful to provide an efficient implementation of parallel MLS heuristics.

Memory Coalescing Issues:

In GPGPU paradigm, each block of threads is split into SIMD groups of threads called warps. At any clock cycle, each processor of the multiprocessor selects a half-warp (16 threads) that is ready to execute the same instruction on different data. Global memory is conceptually organized into a sequence of 128-byte segments. The number of memory transactions performed for a half-warp will be the number of segments having the same addresses than those used by that half-warp.

For more efficiency, global memory accesses must be coalesced, which means that a memory request performed by consecutive threads in a half-warp is associated with precisely one segment. The requirement is that threads of the same warp must read global memory in an ordered pattern. If per-thread memory accesses for a single half-warp constitute a contiguous range of addresses, accesses will be coalesced into a single memory transaction. Otherwise, accessing scattered locations results in memory divergence and requires the processor to perform one memory transaction per thread. The performance penalty for noncoalesced memory accesses varies according to the size of the data structure.

Regarding LS evaluation kernels, coalescing is problem-dependent and difficult when global memory accesses have a data-dependent unstructured pattern (especially for permutation representation). As a result, non-coalesced memory accesses imply many memory transactions and it can lead to a significant performance decrease for LS methods.

Texture Memory: Optimizing the performance of GPU applications often involves optimizing data accesses which includes the appropriate use of the various GPU memory spaces. The use of texture memory is a solution for reducing memory transactions due to non-coalesced accesses. Texture memory provides a surprising aggregation of capabilities including the ability to cache global memory. Therefore, texture memory can be seen as a relaxed mechanism for the thread processors to access global memory because the coalescing requirements do not apply to texture memory accesses.

The use of texture memory is well adapted for LS algorithms since cached texture data is laid out to give best performance for 1D/2D access patterns. The best performance will be achieved when the threads of a warp read locations that are close together from a spatial locality perspective. Since optimization problem inputs are generally 2D matrices or 1D solution vectors, LS structures can be bound to texture memory. The use of textures in place of global memory accesses is a completely mechanical transformation. Details of texture coordinate clamping and filtering is given in [13].

Kernel management: Table 1 summarizes the kernel memory management in accordance with the different LS structures. The inputs of the problem (e.g. matrix in TSP) and the solution which generates the neighborhood are associated with the texture memory. The fitnesses structure which stores the obtained results for each neighbor is declared as global memory. Indeed, since only one writing operation per thread is performed at each iteration, this structure is 6 Application to the Flowshop Scheduling Problem

Configuration

To validate the approaches presented in this paper, the FSP have been implemented on GPU. This problem is one of the most well-known scheduling problems. The problem can be presented as a set of n jobs J 1 , J 2 , . . . , Jn to be scheduled on m machines. Each job J i is composed of m consecutive tasks t i1 , . . . , t im , where t ij represents the jth task of the job J i requiring the machine M j . To each task t ij is associated a processing time p ij , and to each job J i a release time r i and a due date d i (deadline of the job) are given.

For the following experiments, three objectives are used in scheduling tasks on different machines:

-Makespan (total completion time); max{C i |i ∈ [1...n]} -Total tardiness; n i=1 max(0, C i -d i )
-Number of jobs delayed with regard to their due date d i where s ij represents the time at which the task t ij is scheduled and C i = s im + p im represents the completion time of job J i .

The problem has been implemented using a permutation representation and the neighborhood is based on a standard insertion operator. The incremental evaluation function has a time complexity of O(n) and a space complexity of O(m × n). The considered instances are the Taillard instances extended by Liefooghe in a multiobjective context [START_REF] Liefooghe | Combinatorial optimization of stochastic multi-objective problems: An application to the flow-shop scheduling problem[END_REF].

As the iteration-level parallel model does not change the semantics of the sequential algorithm, the effectiveness in terms of quality of solutions is not addressed here. Only average execution times and acceleration factors are reported in comparison with a single-core CPU. The objective is to evaluate the impact of a GPU-based implementation in terms of efficiency.

Experiments have been implemented on top of two different configurations. The GPU desktop cards have a different number of cores (cadenced between 600 and 700 Mhz) which determines the number of active threads being executed. The number of global iterations of the each LS is set to 10000 which corresponds to a realistic scenario in accordance with the algorithm convergence. For each algorithm, a single-core CPU implementation, a CPU-GPU, and a CPU-GPU version using texture memory (GP U tex ) are considered for each configuration. To build the CPU test code, the g++ compiler has been used with the -O2 optimization flag and SSE instructions. The average time has been measured in seconds for 50 runs. For an ease of reading, the standard deviation is not represented in the following tables. Regarding the results, there is no difference of the quality of the solutions provided by both CPU and GPU.

Aggregated Tabu search

For the first experiment, a tabu search based on an aggregation (or weighted) method is used for the generation of Pareto solutions. Thereby, the FSP is transformed into a monoobjective problem by combining the various objective functions into a single one in a linear way. The results are shown in Table 2. The associated standard deviation values are close to zero. From the instance 50-10, GPU versions start to give positive accelerations for both configurations (from ×2.1 to ×3.8). Indeed, the slow speed-ups for small instances can be explained by the fact that since the neighborhood is relatively small, the number of threads per block is not enough to fully cover the memory access latency. As long as the instance size increases, the acceleration factor grows accordingly. For example, from a larger instance such as 100-10, the provided speed-ups are much better (from ×7.2 to ×12.6).

For each instance, in a general manner, the use of texture memory allows to provide additional acceleration. However, constraints of memory alignment in latest G200 and G400 series are relaxed in comparison with the previous cards (e.g. G80 and G90 series). As a consequence, programs running on the used cards get a better global memory performance and the benefits of using the texture memory are less evident.

Finally, efficient speed-ups are obtained for the instance 200-20. They vary between ×9.3 and ×16. As a consequence, parallelization on top of GPU provides an efficient way for handling large neighborhoods.

Pareto Local Search algorithms

For the next experiments, the PLS-1 (Pareto local search proposed by Paquete et al [START_REF] Paquete | Stochastic Local Search Algorithms for Multiobjective Combinatorial Optimization: Methods And Analysis[END_REF]) has been considered. At each iteration, PLS-1 selects a non-dominated solution from the unbounded archive and explores its neighborhood in an exhaustive way. The termination condition is done when all the solutions contained in the archive have been visited. For the experiments, a restart mechanism is performed while the number of iterations has not reached 10000. The solutions archiving on CPU is done in a general way where the solution to insert is compared with each of those contained in the archive.

A first algorithm P LS have been considered where only the dominating neighbors are added to the archive. The results of the experiments are provided in Table 3. Kolmogorov-Smirnov statistical tests can be made to check the normal distribution of the different results. In comparison with Table 2, similar observations can be made regarding the performance results where the maximal speed-up reaches the value ×15.7 for the biggest instance. A look at the average number of non-dominated solutions obtained by the algorithm shows that this number is rather low whatever the instance size. Therefore, this may explained why the performance results of the Pareto algorithm are similar to the aggregated tabu search. To emphasize this point, a second algorithm P LS ≺ has been considered. In this version, all the non-dominated neighbors are added to the archive. Table 4 reports the different measurements. For the instance 100-10, in comparison with the previous table, one can clearly start to see the impact of the number of non-dominated solutions in terms of performance. Indeed, the acceleration factors vary between ×3.1 to ×5.1. The performance results are significantly reduced in comparison with the analog instance for P LS . Finally, the speed-up still grows with the size increase until reaching the value ×7.1 for the last instance.

This global performance loss can be explained by the fact that the number of non-dominated solutions is more important. Indeed, one can see that P LS ≺ is more time-consuming than its counterpart P LS . Therefore, the time spent on the archiving of solutions may be significant. Furthermore, this step is performed exclusively on CPU.

Table 5 confirms this point by showing a profiling of the different steps of the different algorithms. As one can see, for P LS ≺ , the time dedicated to the solutions archiving on CPU is more important than its counterpart P LS . As a consequence, less time is dedicated to the parallel generation and evaluation on GPU. The conclusion of these experiments indicates that since the solutions archiving is significant in the PLS process, it may be significant to proceed this archiving on GPU.

Discussion and Conclusion

High-performance computing based on the use of computational GPUs is recently revealed to be a good way to get at hand such computational power. However, the exploitation of parallel models is not trivial and many issues related to the GPU memory hierarchical management of this architecture have to be considered. To the best of our knowledge, GPU-based parallel LS approaches have never been widely investigated.

In this paper, efficient mapping of the iteration-level parallel model on the hierarchical GPU has been proposed. In the cooperation layer, the CPU manages the whole MLS process and let the GPU be used as a coprocessor dedicated to intensive calculations. Then, efficient mappings between neighborhood candidate solutions and GPU threads are necessary to generate the neighborhood on GPU and thus reduce CPU/GPU transfers. Finally, memory management is applied to the evaluation function kernel.

Apart from being generic, we proved the effectiveness of our methodology by making extensive experiments. Applying such mechanism with an efficient memory management allows to provide promising speed-ups (up to ×16). We strongly believe that the overall performance could be better for other multiobjective optimization problems requiring 1) more computational calculations and 2) less resources in terms of memory (linear time complexity and quadratic space complexity for each neighbor of the FSP).

The re-design of the parallel MLS iteration-level model on GPU fits well for deterministic scalar and Pareto approaches. However, few other MLS algorithms only partially explore neighborhoods and take the first improving local neighbor that is detected. Applied to LS algorithm such as multiobjective simulated annealing, this model needs to be re-thought.

Furthermore, in the case of Pareto approaches, the experimental results have shown a global performance loss with the increase of non-dominated solutions. For other multiobjective optimization problems whose objectives are uncorrelated, this number could be huge, leading to a serious performance decrease. Even if some archiving techniques could be applied to bound the archive size, this would not completely solve the issue at all. As a consequence, for being complete, the next step of our method is to provide a SIMD parallel archiving on GPU. This way, it will allow to significantly enhance the performance of the provided algorithms. However, performing such step is challenging since it requires to ensure additional synchronizations, non-concurrent writings and to manage some dynamic allocations on the GPU.

  Algorithm 1 Multiobjective Local Search Template on GPU

	1: Choose an initial solution
	2: Evaluate the solution
	3: Specific MLS initializations
	4: Allocate problem data inputs on GPU device memory
	5: Allocate a solution on GPU device memory
	6: Allocate a neighborhood fitnesses structure on GPU device memory
	7: Allocate additional solution structures on GPU device memory
	8: Copy problem data inputs on GPU device memory
	9: Copy the solution on GPU device memory
	10: Copy additional solution structures on GPU device memory
	11: repeat
	12:

Table 1 .

 1 Summary of the different memories used in the evaluation function. Declared variables for the computation of the evaluation function of each neighbor are automatically associated with registers by the compiler. Additional complex structures which are private to a neighbor will reside in local memory. Regarding the shared memory, since it is local to each threads block, its use is not addressed in the generic parallel iteration-level model. However, according to the specific problem, the shared memory might be used in the solution-level model where the evaluation function is divided into partial functions.

	Type of memory	LS structure
	Texture memory data inputs, solution representation
	Global memory	fitnesses structure
	Registers	additional variables
	Local memory	additional structures
	not part of intensive calculation.

Table 2 .

 2 Time measurements for the tabu search.

			Xeon 3Ghz	Core i7 3.2Ghz
	Instance		GTX 285 240 cores		GTX 480 480 cores
		CPU	GPU	GP UT ex CPU	GPU	GP UT ex
	20-10	1.1	3.9×0.3	3.7×0.4	0.9	1.9×0.5	1.7×0.6
	20-20	2.3	7.1×0.3	6.6×0.4	1.9	3.3×0.6	6×0.7
	50-10	19.8 9.4×2.1	8.9×2.2	16.5	5.0×3.3	4.3×3.8
	50-20	38.0 17.4×2.2 16.3×2.3 31.9	9.1×3.5	7.8×4.1
	100-10 170.8 23.6×7.2 20.9×8.2 144.5 12.6×11.5 11.5×12.6
	100-20 321.1 44.1×7.3 38.1×8.4 270.7 23.3×11.6 20.9×12.9
	200-10 1417.4 159.4×8.9 144.3×9.8 1189.7 81.88×14.5 77.2×15.4
	200-20 2644.1 284.4×9.3 263.9×10.0 2220.7 147.8×15.0 139.0×16.0

Table 3 .

 3 Time measurements for P LS .

			Xeon 3Ghz	Core i7 3.2Ghz
	Instance		GTX 285 240 cores		GTX 480 480 cores	# ≺ solutions
		CPU	GPU GP UT ex CPU	GPU	GP UT ex
	20-10	1.1	5.3×0.2	3.7×0.3	1.0	1.9×0.5	1.7×0.6	11
	20-20	2.5	8.4×0.3	6.6×0.4	1.9	3.4×0.6	3.0×0.6	13
	50-10	19.6 10.9×1.8 8.9×2.2 16.6	4.9×3.4	4.3×3.8	19
	50-20	39.7 18.9×2.1 16.5×2.4 32.7	9.1×3.6	7.8×4.2	20
	100-10 167.0 25.3×6.6 21.2×7.9 139.5 12.8×10.9 11.7×12.0	31
	100-20 329.8 45.8×7.2 40.4×8.1 274.9 23.5×11.7 21.1×13.1	34
	200-10 1391.5 161.8×8.6 145.1×9.6 1171.5 82.5×14.2 77.9×15.1	61
	200-20 2707.8 307.7×8.8 285.9×9.4 2196.3 148.4×14.8 139.7×15.7	71

Table 4 .

 4 Time measurements for P LS ≺ .

			Xeon 3Ghz	Core i7 3.2Ghz
	Instance		GTX 285 240 cores		GTX 480 480 cores	# ≺ solutions
		CPU	GPU GP UT ex CPU	GPU GP UT ex
	20-10	1.1	5.4×0.2	3.9×0.3	1.2	2.0×0.6	1.8×0.7	77
	20-20	2.6	8.5×0.3	6.7×0.4	2.0	3.4×0.6	3.0×0.7	83
	50-10	25.7 15.1×1.7 13.1×1.9 20.5 8.9×2.3 8.2×2.5	396
	50-20	44.5 24.7×1.8 22.3×2.2 39.5 14.1×2.8 13.2×3.0	596
	100-10 220.1 71.0×3.1 66.9×3.3 253.3 51.7×4.9 49.2×5.1	1350
	100-20 386.8 96.7×4.0 91.3×4.2 312.5 62.5×5.0 59.0×5.3	1530
	200-10 1631.2 362.5×4.5 346.3×4.7 1361.5 223.2×6.1 217.7×6.4	1597
	200-20 2998.8 588.0×5.1 566.1×5.3 2672.6 398.9×6.7 378.0×7.1	2061

Table 5 .

 5 Profiling of the two Pareto algorithms.

	Instance	PLS Evaluation Archiving LS process Evaluation Archiving LS process PLS ≺
	20-10	99.2%	0.7%	0.1%	95.2%	4.7%	0.1%
	20-20	99.4%	0.5%	0.1%	98.2%	1.7%	0.1%
	50-10	99.5%	0.4%	0.1%	89.7%	10.2%	0.1%
	50-20	99.7%	0.2%	0.1%	94.2%	5.7%	0.1%
	100-10	99.8%	0.1%	0.1%	88.2%	11.8%	0.1%
	100-20	99.85%	0.1%	0.05%	93.6%	6.35%	0.05%
	200-10	99.9%	0.05%	0.05%	92.3%	7.65%	0.05%
	200-20	99.9%	0.05%	0.05%	94.2%	5.75%	0.05%


Thé Van Luong

Nouredine Melab

El-Ghazali Talbi
email: el-ghazali.talbi]@lifl.fr

GPU-based Multi-start Local Search Algorithms

Keywords:

In practice, combinatorial optimization problems are complex and computationally time-intensive. Local search algorithms are powerful heuristics which allow to significantly reduce the computation time cost of the solution exploration space. In these algorithms, the multistart model may improve the quality and the robustness of the obtained solutions. However, solving large size and time-intensive optimization problems with this model requires a large amount of computational resources. GPU computing is recently revealed as a powerful way to harness these resources. In this paper, the focus is on the multi-start model for local search algorithms on GPU. We address its re-design, implementation and associated issues related to the GPU execution context. The preliminary results demonstrate the effectiveness of the proposed approaches and their capabilities to exploit the GPU architecture.

Introduction

Over the last years, interest in metaheuristics (generic heuristics) has risen considerably in the field of optimization. Indeed, plenty of hard problems in a wide range of areas including logistics, telecommunications, biology, etc., have been modeled and tackled successfully with metaheuristics. Local search (LS) algorithms are a class of metaheuristics which handle with a single solution iteratively improved by exploring its neighborhood in the solution space. Different parallel models have been proposed in the literature for the design and implementation of LSs [START_REF] Talbi | Metaheuristics: From design to implementation[END_REF]. The multi-start model consists in executing in parallel many LSs in an independent/cooperative manner. This mechanism may provide more effective, diversified and robust solutions.

Nevertheless, although LS methods have provided very powerful search algorithms, problems in practice are becoming more and more complex and CPU time-intensive and their resolution requires to harness more and more computational resources. In parallel, the recent advances in hardware architecture allow to provide such required tremendous computational power through GPU infrastructures. This new emerging technology is indeed believed to be extremely useful to speed up many complex algorithms. However, the exploitation of such computational infrastructures in metaheuristics is not straightforward. Indeed, several scientific challenges mainly related to the hierarchical memory management or to the execution context have to be faced. The major issues are the efficient distribution of data processing between the CPU and the GPU, the thread synchronization, the optimization of data transfer between the different memories, the capacity constraints of these memories, etc. The main objective of our research work is to deal with such issues for the re-design of parallel metaheuristics models to allow solving of large scale optimization problems on GPU architectures. In [START_REF] Luong | Local search algorithms on graphics processing units. a case study: the permutation perceptron problem[END_REF][START_REF] Luong | Large neighborhood for local search algorithms[END_REF], we have proposed to re-design the parallel evaluation of the neighborhood model for LSs on GPU. To go on this way, the main objective of this paper is to deal with the well-known multi-start model on GPU architectures where many LSs are executed in parallel.

We deal with the entire re-design of the multi-start model on GPU by taking into account the particular features related to both the LS process and the GPU computing. More exactly, we provide two different general schemes for building efficient multi-start LSs on GPU. The first scheme combines the multistart model with the parallel evaluation of the neighborhood on GPU previously mentioned above. In the second scheme, the search process of each LS algorithm is fully distributed on GPU. The advantage of the full distribution of the search process on GPU is to reduce CPU/GPU memory copy latency. We will essentially focus on this approach throughout this paper.

Despite the fact that the second scheme for the multi-start model has already been applied in some previous works in the context of the tabu search on GPU [START_REF] Zhu | Simd tabu search with graphics hardware acceleration on the quadratic assignment problem[END_REF][START_REF] Janiak | Tabu search on gpu[END_REF], to the best of our knowledge, it has never been widely investigated in terms of 1) reproducibility for any other LS algorithm and 2) memory management. Indeed, the contribution of this paper is to provide a general methodology for the design of multi-start LSs on GPU applicable to any class of LS algorithms such as hill climbing, tabu search or simulated annealing. Furthermore, a particular focus is made on finding efficient associations between the different available memories and the data commonly used in the multi-start LS algorithms.

The remainder of the paper is organized as follows: on the hand, Section 2 highlights the principles of LS parallel models. On the other hand, a brief review of the GPU architecture is also depicted. Section 3 presents a methodology for the design and the implementation of parallel multi-start LS methods on GPU. The performance results obtained for the associated implementations are reported in Section 4. Finally, a discussion and some conclusions of this work are drawn in Section 5.

Parallel Local Search Algorithms and GPU Computing

Parallel Models of LS Algorithms

For non-trivial problems, executing the iterative process of a simple LS on large neighborhoods requires a large amount of computational resources. Consequently, a variety of algorithmic issues are being studied to design efficient LS heuristics. Parallelism arises naturally when dealing with a neighborhood, since each of the solutions belonging to it is an independent unit. Due to this, the performance of LS algorithms is particularly improved when running in parallel. Parallel design and implementation of metaheuristics have been studied as well on different architectures [START_REF] Alba | 4. Metaheuristics and Parallelism[END_REF][START_REF] Zomaya | Sequential and parallel meta-heuristics for solving the single row routing problem[END_REF][START_REF] Melab | Grid computing for parallel bioinspired algorithms[END_REF].

Basically, three major parallel models for LS heuristics can be distinguished: solution-level, iteration-level and algorithmic-level.

• Solution-level Parallel Model. A focus is made on the parallel evaluation of a single solution. Problem-dependent operations performed on solutions are parallelized. In that case, the function can be viewed as an aggregation of a given number of partial functions. The solution-level model is problem-dependent and does not present many generic concepts. In this paper, we will focus on the multi-start model which is an instantiation of the algorithmic-level model where LS algorithms are all homogeneous.

GPU Computing

GPUs have evolved into a highly parallel, multithreaded and many-core environment. Indeed, since more transistors are devoted to data processing rather than data caching and flow control, GPU is specialized for compute-intensive and highly parallel computation. A complete review of GPU architecture can be found in [START_REF] Ryoo | Program optimization carving for gpu computing[END_REF].

In general-purpose computing on graphics processing units, the CPU is considered as a host and the GPU is used as a device coprocessor. This way, each GPU has its own memory and processing elements that are separate from the host computer. Memory transfer from the CPU to the GPU device memory is a (a)synchronous operation which is time consuming. Bus bandwidth and latency between the CPU and the GPU can significantly decrease the performance of the search, so data transfers must be minimized.

Each processor device on GPU supports the single program multiple data (SPMD) model, i.e. multiple processors simultaneously execute the same program on different data. For achieving this, the concept of kernel is defined. The kernel is a function callable from the host and executed on the specified device by several processors in parallel. This kernel handling is dependent of the general-purpose language. For instance, CUDA [10] or OpenCL [START_REF] Group | OpenCL 1.0 Quick Reference Card[END_REF] are parallel computing environments which provide an application programming interface. These toolkits introduce a model of threads which provides an easy abstraction for single-instruction and multipledata (SIMD) architecture.

Regarding their spatial organization, threads are organized within so called thread blocks. A kernel is executed by multiple equally threaded blocks. Blocks can be organized into a one-dimensional or two-dimensional grid of thread blocks, and threads inside a block are grouped in a similar way. All the threads belonging to the same thread block will be assigned as a group to a single multiprocessor, while different thread blocks can be assigned to different multiprocessors.

From a hardware point of view, graphics cards consist of streaming multiprocessors, each with processing units, registers and on-chip memory. Since multiprocessors are used according to the SPMD model, threads share the same code and have access to different memory areas. Basically, the communication between the CPU host and its device is done through the global memory.

Design and implementation of multi-start local search algorithms on GPU

With the recent advances in parallel computing particularly based on GPU computing, the multi-start model has to be re-visited from the design and implementation points of view. In this section, we propose multiple deployment schemes of the multi-start model for LS algorithms on GPU.

Multi-start Local Search Algorithms Based on the Iteration-level

In [START_REF] Luong | Local search algorithms on graphics processing units. a case study: the permutation perceptron problem[END_REF][START_REF] Luong | Large neighborhood for local search algorithms[END_REF], we have proposed the design and the implementation of the parallel evaluation of the neighborhood (iteration-level) model for a single LS on GPU. That is the reason why, a natural way for designing multi-start LSs on GPU based on the iteration-level is to iterate the whole process (i.e. the execution of a single LS on GPU) to deal with as many LSs as needed (see Fig. 1). Indeed, in general, evaluating a fitness function for each neighbor is frequently the most costly operation of the LS. Therefore, in this scheme, task distribution is clearly defined: the CPU manages the whole sequential LS process for each LS algorithm and the GPU is dedicated only to the parallel evaluation of solutions. Algorithm 2 gives the template of this model. The reader is referred to [START_REF] Luong | Local search algorithms on graphics processing units. a case study: the permutation perceptron problem[END_REF][START_REF] Luong | Large neighborhood for local search algorithms[END_REF] for more details about the original algorithm. Basically, for each LS, the CPU first sends the number of expected neighbors to be evaluated to the GPU and then these solutions are processed on GPU. Regarding the kernel thread organization, as quoted above, a GPU is organized following the SPMD model, meaning that each GPU thread associated with one neighbor executes the same evaluation function kernel. Finally, results of the evaluation function are returned back to the host via the global memory. Choose an initial solution 8:

Evaluate the solution 9:

Specific LS initializations 10: end for 11: repeat 12:

for m = 1 to #local searches do 13:

Copy the solution on GPU memory 14:

Copy additional solution structures on GPU memory 15:

for each neighbor in parallel on GPU do 16:

Incremental evaluation of the candidate solution 17:

Insert the resulting fitness into the neighborhood fitnesses structure 18:

end for 19:

Copy back the neighborhood fitnesses structure on CPU memory 20:

Specific LS solution selection strategy on the neighborhood fitnesses structure 21:

Specific LS post-treatment 22:

end for 23:

Possible cooperation between the different solutions 24: until a stopping criterion satisfied This way, the GPU is used as a coprocessor in a synchronous manner. The time-consuming part i.e. the incremental evaluation kernel is calculated by the GPU and the rest is handled by the CPU. The advantage of this scheme resides in its highly parallel structure (i.e. an important number of generated neighbors to handle), leading to a significant multiprocessors occupancy of the GPU. However, depending on the number of LS algorithms, the main drawback of this scheme is that copying operations from the CPU to the GPU can become frequent and thus can lead to a significant performance decrease.

Design of Multi-start Local Search Algorithms Based on the Algorithmic-level

A natural way for designing multi-start LSs on GPU is to parallelize the whole LS process on GPU by associating one GPU thread with one LS. This way, the main advantage of this approach is to minimize the data transfers between the host CPU memory and the GPU. Figure 2 illustrates this idea of this full distribution (algorithmic-level). In the rest of this paper, we will focus on this approach. The details of the algorithm are given in Algorithm 1. First of all, at initialization stage, memory allocations on GPU are made: data inputs of the problem must be allocated and copied on GPU (lines 1 and 2). It is important to notice that problem data inputs (e.g. a matrix in the traveling salesman problem [START_REF] Burkard | The travelling salesman problem on permuted monge matrices[END_REF]) are a read-only structure and never change during all the execution of LS algorithms. Therefore, their associated memory is copied only once during all the execution. Second, a certain number of solutions corresponding to each LS must be allocated on GPU (line 3). Additional solution structures which are problem-dependent can also be allocated to facilitate the computation of incremental evaluation (line 4). Third, during the initialization of the different LS algorithms on GPU, each solution is generated and evaluated (from lines 5 to 9). Fourth, comes the algorithmic-level, in which the iteration process of each LS is performed in parallel on GPU (from lines 11 to 17). Since each neighbor is evaluated in a sequential manner on GPU, unlike the iteration-level scheme, there is no need to allocate and manipulate any neighborhood fitness structure. Fifth, an exchange of the best-so-far solutions could be made to accelerate the search process (line 18). In that case, operations on the global memory may be considered. Finally, the process is repeated until a stopping criterion is satisfied. Choose an initial solution 7:

Evaluate the solution 8:

Specific LS initializations 9: end for 10: repeat 11:

for each LS in parallel on GPU do 12:

for each neighbor do 13:

Incremental evaluation of the candidate solution 14:

Specific LS solution selection strategy 15:

end for 16:

Specific LS post-treatment 17:

end for 18:

Possible cooperation between the different solutions 19: until a stopping criterion satisfied

Memory Management of Multi-start Local Search Algorithms

on the Algorithmic-level Memory Coalescing Issues. When an application is executed on GPU, each block of threads is split into SIMD groups of threads called warps. At any clock cycle, each processor of the multiprocessor selects a half-warp (16 threads) that is ready to execute the same instruction on different data. Global memory is conceptually organized into a sequence of 128-byte segments. The number of memory transactions performed for a half-warp will be the number of segments having the same addresses than those used by that half-warp. Fig. 3 illustrates an example of the memory management layer for a simple vector addition.

Fig. 3. An example of kernel execution for vector addition.

For more efficiency, global memory accesses must be coalesced, which means that a memory request performed by consecutive threads in a half-warp is associated with precisely one segment. The requirement is that threads of the same warp must read global memory in an ordered pattern. If per-thread memory accesses for a single half-warp constitute a contiguous range of addresses, accesses will be coalesced into a single memory transaction. In the example of vector addition, memory accesses to the vectors a and b are fully coalesced, since threads with consecutive thread indices access contiguous words.

Otherwise, accessing scattered locations results in memory divergence and requires the processor to perform one memory transaction per thread. The performance penalty for non-coalesced memory accesses varies according to the size of the data structure. Regarding LS structures, coalescing is difficult when global memory accesses have a data-dependent unstructured pattern (especially for a permutation representation). As a result, non-coalesced memory accesses imply many memory transactions and it can lead to a significant performance decrease for LS methods.

Memory Organization. Optimizing the performance of GPU applications often involves optimizing data accesses which includes the appropriate use of the various GPU memory spaces. For instance, the use of texture memory is a solution for reducing memory transactions due to non-coalesced accesses. Texture memory provides a surprising aggregation of capabilities including the ability to cache global memory (separate from register, global, and shared memory). Regarding the data management on the different GPU memories, the following observations can be made whatever the used multi-start LS algorithm:

• Global memory: For each running LS on GPU (one thread), its associated solution is stored on the global memory. The same goes on for additional solution structures. This way, it ensures a global visibility among the different threads (LSs) during the entire search process for a possible cooperation. In a general way, all the data in combinatorial problems could be also associated with the global memory. However, as previously said, non-coalesced memory accesses may lead to a performance decrease. Therefore, the texture memory might be preferred since it can be seen as a relaxed mechanism for the threads to access the global memory. Indeed, the coalescing requirements do not apply to texture memory accesses. • Texture memory: This read-only memory is adapted to LS algorithms since the problem inputs do not change during the execution of the algorithm.

In most of optimization problems, problem inputs do not often require a large amount of allocated space memory. As a consequence, these structures can take advantage of the 8KB cache per multiprocessor of texture units. Indeed, minimizing the number of times that data goes through cache can increase significantly the efficiency of algorithms [START_REF] Bader | A cache-aware parallel implementation of the pushrelabel network flow algorithm and experimental evaluation of the gap relabeling heuristic[END_REF]. Moreover, cached texture data is laid out to give best performance for structures with 1D/2D access patterns such as matrices. The use of textures in place of global memory accesses is a completely mechanical transformation. Details of texture coordinate clamping and filtering is given in [10,[START_REF] Nickolls | Scalable parallel programming with cuda[END_REF]. • Constant memory: This memory is read only from kernels and is hardware optimized for the case where all threads read the same location. It might be used when the calculation of the evaluation function requires a common lookup table for all solutions (e.g. a decoder table for an indirect encoding on the job shop scheduling problem [START_REF] Dell'amico | Applying tabu search to the job-shop scheduling problem[END_REF]). • Shared memory: The shared memory is a fast memory located on the multiprocessors and shared by threads of each thread block. Since this memory area provides a way for threads to communicate within the same block, it might be used with the global memory in the context of a possible cooperation between different LS algorithms. In the case of the multi-start LS model, the type of shared information is the best-so-far solution found at each iteration of the search process. • Registers: Among streaming processors, they are partitioned among the threads running on it and they constitute fast access memory. In the kernel code, each declared variable is automatically put into registers. • Local memory: In a general way, additional structures such as declared array will reside in local memory. In fact, local memory resides in the global memory allocated by the compiler and its visibility is local to a thread (a LS).

Table 1 summarizes the kernel memory management in accordance with the different LS components. For the management of random numbers in SA, effi- cient techniques are provided in many books such as [START_REF]GPU Gems 3. Chapter 37: Efficient Random Number Generation and Application Using CUDA[END_REF] to implement random generators on GPU. For deterministic multi-start LSs based on HC or TS, the random initialization of solutions might be done on CPU and then they can be copied on the GPU via the global memory to perform the LS process. This way, it ensures that the obtained results are the same as a multi-start LS performed on a traditional CPU. Regarding the management of the tabu list on GPU, since the list is particular to a TS execution, a natural mapping is to associate a tabu list to the local memory. However, since this memory has a limited size, large tabu lists should be associated with the global memory instead.

Experiments

To validate our approach, the multi-start model has been implemented on the quadratic assignment problem (QAP) on GPU using CUDA. The QAP arises in many applications such as facility location or data analysis. Let A = (a ij ) and B = (b ij ) be n × n matrices of positive integers. Finding a solution of the QAP is equivalent to finding a permutation π = (1, 2, . . . , n) that minimizes the objective function:

z(π) = n i=1 n j=1 a ij b π(i)π(j)
The problem has been implemented using a permutation representation. The chosen neighborhood for all the experiments is based on a 2-exchange operator ( n×(n-1) 2 neighbors). The incremental evaluation function has a time complexity of O(n). The considered instances are the Taillard instances proposed in [START_REF] Taillard | Robust taboo search for the quadratic assignment problem[END_REF]. They are uniformly generated and are well-known for their difficulty.

The used configuration is an Intel Xeon 3GHz 2 cores with a GTX 280 (30 multiprocessors). From an implementation point of view, to build the CPU test code, the g++ compiler has been used with the -O2 optimization flag and SSE instructions. The specific parameters for each single LS algorithm are given in Table 2.

Table 2. Used parameters for each particular LS.

Tabu search

Simulated annealing geometric cooling schedule tabu list size: tl = n×(n-1) 16 initial temperature: T0 = 10000 threshold: thr = 1 iterations: iters = 10000 ratio: r = 0.9 iterations: iters = n×(n-1) 2 equilibrium state: T < thr

Measures of the Efficiency of Multi-start Algorithms Based on the Algorithmic-level

In the next experiments, the effectiveness in terms of quality of solutions is not addressed here. Only execution times and acceleration factors are reported in comparison with a mono-core CPU. The objective is to evaluate the impact of a GPU implementation of multi-start algorithms based on the algorithmic-level (i.e. the full distribution of the search process on GPU) in terms of efficiency.

For each multi-start algorithm, a standalone mono-core CPU implementation, a pure GPU one, and a GPU version using texture memory (GP U tex ) are considered. The number of LS algorithms of the multi-start model is set to 4096 which corresponds to a realistic scenario in accordance with the algorithm convergence.

The average time has been measured in seconds for 30 runs. The standard deviation is not represented since its value is very low for each measured instance.

The obtained results are reported in Table 3 for the different LS multi-start algorithms on GPU. Regarding the acceleration for a pure implementation on GPU based on HC (HC GPU), it varies between ×1.7 for the instance tai30a to ×3.6 for the last instance. In comparison with a pure CPU implementation, the obtained accel-eration factors are positive but not impressive. Indeed, due to high misaligned accesses to global memories (flows and distances in QAP), non-coalescing memory reduces the performance of the implementation. Binding texture on global memory allows to overcome the problem (HC GPUTex). Indeed, from the instance tai30a, using texture memory starts providing significant acceleration factors (×5.4). GPU keeps accelerating the LS process as long as the size grows and the best results are obtained for the instance tai100a (×7.8).

Regarding the performance for the other multi-start algorithms (TS and SA based), similar observations can be made. Indeed, on the hand, the obtained speed-ups for the texture version of multi-start algorithms based on TS vary between ×6.1 to ×10.8. And on the other hand, they vary from ×5.7 to ×12.0 for the multi-start algorithms based on SA. In a general manner, the performance variation obtained with the different algorithms on GPU is in accordance with the algorithm complexity (Complexity(SA) >= Complexity(T S) >> Complexity(HC)).

The point to highlight in these experiments is that organizing data into cache such as texture memory clearly allows to improve the speed-ups in comparison with a standard GPU version where inputs are stored in the global memory.

Measures of the Efficiency of Large GPU-based Implementations

Another experiment consists in measuring the impact in terms of efficiency by varying the number of LSs in the multi-start based on the algorithmic-level. In addition, we propose to compare this approach with the multi-start based on the iteration-level model (parallel evaluation of the neighborhood on GPU) presented in Section 3.1. For doing this, we propose to deal with the instance tai50a with the same parameters used before in the context of multi-start methods based on TS. The obtained results are depicted in Fig. 4 for the texture optimization.

For the algorithmic-level, one can notice that it starts providing a positive acceleration of ×1.7 from a number of 512 LSs (one thread per LS). From 1024 LSs, the acceleration factors are drastically improved until reaching ×8.7 for 4096 LSs. After that, the speed-up keeps improving slowly with the size increase. However, as one can see in Fig. 5, no significant difference can be made in terms of the quality of the solutions obtained for more than 168384 LSs. Therefore, since the execution is already time-consuming, it might not be relevant to perform more LSs.

Regarding a small number of running LSs, from 1 to 256 LSs, the multi-start for the algorithmic-level is clearly inefficient. This can be explained by the fact that since the number of threads is relatively small, the number of threads per block is not enough to fully cover the memory access latency.

Unlike the previous model, for the multi-start based on the iteration-level, the obtained speed-ups are quiet regular (from ×4.4 to ×5.1) whatever the number of running LSs. Indeed, since one thread is associated with one neighbor ( n×(n-1) 2 neighbors), during the kernel execution, there is enough threads to keep the GPU multiprocessors busy. However, as one can notice, the maximal performance of Fig. 4. Measures of the efficiency of the two multi-start approaches using the texture memory algorithmic-level approach in comparison with the iteration-level by varying the number of tabu searches (instance tai50a). this scheme is quiet limited because of the multiple data copies between the CPU and the GPU (see [START_REF] Luong | Local search algorithms on graphics processing units. a case study: the permutation perceptron problem[END_REF] for an analysis of data transfers).

Discussion and Conclusion

Parallel metaheuristics such as the multi-start model allow to improve the effectiveness and robustness in optimization problems. Their exploitation for solving real-world problems is possible only by using a great computational power. High-performance computing based on GPU accelerators is recently revealed as an efficient way to use the huge amount of resources at disposal. However, the exploitation of the multi-start model is not trivial and many issues related to the context execution and to the memory hierarchy of this architecture have to be considered.

In this paper, we have proposed a guideline to design and implement general GPU-based multi-start LS algorithms. The different concepts addressed throughout this paper takes into account popular LS algorithms such as HC, SA or TS. The designed and implemented approaches have been experimentally validated on a combinatorial optimization problem. To the best of our best of knowledge, multi-start parallel LS approaches have never been widely investigated so far.

The idea of our methodology is based on two natural schemes which exploit the GPU in a different manner. In the first scheme, the multi-start model is combined with the parallel evaluation of the neighborhood. The advantage of this is to maximize the GPU in terms of multiprocessor occupancy. However, the performance of this scheme is limited due to the data transfers between the CPU and the GPU. To deal with this issue, we have particularly focused on the full distribution of the search process on GPU with the appropriate use of memory. Applying such mechanism with an efficient memory management allows to provide significant speed-ups (up to ×12). However, this second scheme could also present some performance limitations when dealing with a small number of LS executions.

In a general manner, the two proposed schemes are complementary and their use strongly depends of the number of LSs to deal with. It would be interesting to test the performance of our approaches with some combinatorial optimization problems involving the use of different memories such as the constant and the shared memory.

Another perspective of this work is to combine the multi-start on GPU with a pure multi-core approach. Indeed, since this model has a high degree of parallelism, the CPU cores can also work in parallel in an independent manner. Moreover, since nowadays the actual configurations have 4 and 8 cores, instead of waiting the results back from the GPU, this computational power should be well-exploited in parallel to provide additional accelerations.

Fig. 1 .

 1 Fig.1. Multi-start LS algorithms based on the parallel evaluation of the neighborhood on GPU (iteration-level). In this scheme, one thread is associated with one neighbor.

Fig. 2 .

 2 Fig. 2. Multi-start LS algorithms based on the full distribution of LSs on GPU (algorithmic-level). One thread is associated with one local search.

Algorithm 2

 2 Multi-start local search algorithms template on GPU based on the algorithmic-level model 1: Allocate problem data inputs on GPU memory 2: Copy problem data inputs on GPU memory 3: Allocate #local searches solutions on GPU memory 4: Allocate #local searches additional solution structures on GPU memory 5: for each LS in parallel on GPU do 6:

Fig. 5 .

 5 Fig. 5. Measures of the quality of the solutions for the multi-start model based on the algorithmic-level (tai50a). The average fitness is reported for 30 executions where each point represents a certain number of LSs.

  • Iteration-level Parallel Model. This model is a low-level Master-Worker model that does not alter the behavior of the heuristic. Exploration and evaluation of the neighborhood are made in parallel. At the beginning of each iteration, each parallel node manages some candidates and the results are returned back to the master. An efficient execution is often obtained particularly when the evaluation of each solution is costly.• Algorithmic-level Parallel Model. Several LS algorithms are simultaneously launched for computing better and robust solutions. They may be heterogeneous or homogeneous, independent or cooperative, start from the same or different solution(s), configured with the same or different parameters.

Table 1 .

 1 Kernel memory management. Summary of the different memories used in the multi-start LS algorithms on GPU.

	Type of memory	LS structure
	Texture memory	problem data inputs
	Global memory candidate solutions, additional candidate solution structures
	Shared memory	possible solutions to exchange
	Registers	additional LS variables
	Local memory	additional LS structures
	Constant memory	additional problem lookup tables

Table 3 .

 3 Measures of the efficiency of the algorithmic-level on the QAP. The average time is reported in seconds for 30 executions, the number of LSs is fixed to 4096.

		tai30a	tai40a	tai50a	tai60a	tai80a	tai100a
	HC CPU	5.48	17.18	44.56	88.32	302.43	810.39
	HC GPU	3.19×1.7 7.44×2.3 15.79×2.8 30.06×2.9 90.45×3.3 224.51×3.6
	HC GPUTex 1.02×5.4 2.96×5.8 6.69×6.7 12.52×7.1 41.65×7.3 103.59×7.8
	TS CPU	335.57	725.39	1539.60	2439.86	6097.61	13004.76
	TS GPU 105.12×3.2 207.12×3.5 414.50×3.7 655.32×3.7 1544.32×3.9 3222.01×4.0
	TS GPUTex 55.12×6.1 105.65×6.9 176.29×8.7 262.31×9.3 588.33×10.4 1207.77×10.8
	SA CPU	412.64	874.44	1672.63	2699.89	6807.88	13960.69
	SA GPU 115.32×3.6 223.65×3.9 422.32×4.0 677.28×4.0 1578.21×4.3 3121.28×4.5
	SA GPUTex 72.25×5.7 135.21×6.5 205.74×8.1 278.88×9.7 609.78×11.2 1161.52×12.0