Optimization of the Multi-Facility Location Problem Using Widely Available Ofﬁce Software

: Multi-facility location problem is a type of task often solved (not only) in logistics. It consists in ﬁnding the optimal location of the required number of centers for a given number of points. One of the possible solutions is to use the principle of the genetic algorithm. The Solver add-in, which uses the evolutionary method, is available in the Excel ofﬁce software. It was used to solve the benchmark in 4 levels of difﬁculty (from 5 centers for 25 points to 20 centers for 100 points), and one task from practice. The obtained results were compared with the results obtained by the metaheuristic simulated annealing method. It was found that the results obtained by the evolutionary method are sufﬁciently accurate. Their accuracy depends on the complexity of the task and the performance of the HW used. The advantage of the proposed solution is easy availability and minimal requirements for user knowledge.


Introduction
In this work we deal with the problem of finding the coordinates of c centers for b points.A typical example is finding the optimal location for c central warehouses that will serve b branches.In the literature, this problem is called the multi-facility location problem (MFLP) and there are a number of ways to solve this type of task.This problem is often solved in the field of logistics, but there are also other applications.
The article is organized into five sections.In this introduction, our motivation and the objectives of this article are presented.The first section also covers the problem formulation and the literature review.Section 2 presents the details concerning the approach and methods used.Section 3 shows the results of the experiments.Section 4 follows with the discussion about the results achieved.Finally, Section 5 concludes the article.

Literature Review
Various methods for solving facility or multi-facility location problem are described in the literature.

Facility Location Problems Solved in Excel Environment
Ref. [1] describes an Excel-based Decision Support System for Facility Location Problems.The aim is to select a subset of locations from a set of candidate locations, and to determine which customer locations will be served by which facility, to optimize an objective function that is based on the distances (or the costs) between the facilities and the demands of customer locations they serve.
Use of Solver Add-in in logistics is presented in [2].The aim of the paper is to present a solution of facility location problem using the Solver add-in in Excel.In the case study discussed in the paper, the company's central warehouse location was selected based on the classic location theory, which addresses the need to minimize the cost of transport.The mathematical model of the exercise is based on the Euclidean coordinates metrics.
Article [3] concentrates on facility location problem, focusing on the linear programming model and a genetic algorithm in the location problem analysis and analytical method.In the linear programming model, because the given complex table calculation method is too complicated and the workload is very large, the Excel software is proposed to solve the location problem, which can greatly improve the efficiency of enterprise facility location problem.In addition, a genetic algorithm based on MATLAB toolbox is applied to another type of facility location problem, which provides a referential method for location decision under different conditions and different facilities.

In the Area of the Rectangular Distance Multi-Facility Location Problem
A nonlinear approximation method was developed by in [4], where any number of linear and (or) nonlinear constraints defining a convex feasible region can be included.
Ref. [5] presents a new method that, as it states, efficiently handles the rectilinear distance problem MFLP having large clusters, i.e., where several new facilities are located together at one point.This paper states and proves a new necessary and sufficient optimality condition.This condition transforms the problem of computing a descent direction into a constrained linear least-squares problem.The latter problem is solved by a relaxation method that takes advantage of its special structure.The new technique is incorporated into the direct search method.
Ref. [6] develops a dual problem for the constrained multi-facility minimum location problems involving mixed norms.General optimality conditions were obtained providing new algorithms which are decomposition methods based on the concept of partial inverse of a multifunction.
As an alternative to linear programming, coordinate descent was presented in [7] as a simple approach which sometimes finds optimal locations.By deleting the term that shows the relationship between new facilities in objective function, the problem is converted to some single facility location problems to which we can apply median conditions.The first coordinate we choose is the first variable and the second coordinate is the second variable, and so on.It is continued until we obtain the same vector by coordinate descent that we have obtained previously by coordinate descent, at which point we stop.

In the Area of the Euclidean Multi-Facility Location Problem
In [8], the authors proposed a pseudo-gradient technique that classifies the new facilities into distinct categories based on their coincidence with other facilities in order to derive a descent method for solving Euclidean MFLP.
Author in [9] also approached Euclidean MFLP by using a general e-subgradient method in which search directions are generated based on the subdifferential of the objective function over a neighborhood of the current iterate.
In [10] developed an algorithm which, from any initial point, generates a sequence of points that converges to the closed convex set of optimal solutions to the problem Euclidean MFLP.
For the multi-facility location problem with no constraints on the location of the new facilities, authors in [11] derived some sufficient conditions for the coincidence of facilities that are valid in a general symmetric metric.
Xue et al. [12] suggested the use of polynomial-time interior point algorithm to solve the MFLP problem.They presented a procedure in which an approximate optimum to Euclidean MFLP can be recovered by solving a sequence of linear equations, each associated with an iterate of the interior point algorithm used to solve the dual problem.

In the Area of the Rectangular Distance Minimax Location Problem
Morris [13] introduced this problem with linear constraints which (a) limit the new facilities location and (b) enforce upper bounds on the distances between new and existing facilities and between new facilities.He used dual variables that provide information about the complete range of new facility locations which satisfies the MiniMax criterion.
Authors in [14] presented a method involving the numerical integration of ordinary differential equations and was computationally superior to methods using nonlinear programming.

Solutions for Other Models
A specialized simplex-based algorithm was derived in [15] for solving rectangular multiproduct MFLP.
Authors in [16] presented an iterative solution for MFLP on sphere.The procedure involved the approximation of the domain of objective function which in the limit approaches to that of the original objective function.In [17], authors considered Euclidean, squared Euclidean and the great circle distances.They formulated an algorithm and investigated its convergence properties.
1.1.6.Some Heuristic and Metaheuristic Methods Authors [18] introduced a simple heuristic for solving MFLP with Euclidean distance.This procedure locates each of the new facilities in a temporary location at each step and locates the next new facility according to the facilities located so far.After all n new facilities are located in this manner the process is repeated and the readjustment process is continued until no further movements occur during a complete round of adjustment evaluations.
In [19] solved the MFLP with ant-colony optimization metaheuristic when the distances are rectangular and Euclidean.This algorithm produces optimal solutions for problem instances of up to 20 new facilities.

Solutions Using Genetic Algorithm
In the study "A Novel Nondominated Sorting Simplified Swarm Optimization for Multi-Stage Capacitated Facility Location Problems (CFLP) with Multiple Quantitative and Qualitative Objectives," [20], a novel solution based on simplified swarm optimization (SSO) and a nondominated sorting technique is proposed to provide Pareto-optimal solutions for enhancing search efficiency and solution quality.To yield feasible solutions, three repairer mechanisms, namely, random repair, cost-based, and utility-based mechanisms, are proposed to enhance the search efficiency and diversity of each population.A fuzzy analytic hierarchy process is used to calculate the weight of qualitative objectives.To evaluate the efficiency and effectiveness of the proposed algorithm, extensive experiments are conducted on benchmark and newly generated instances of the four stages of CFLPs.Then, results are compared with those of the nondominated sorting genetic algorithm-II, multi-objective SSO, and multi-objective particle swarm optimization reported from the literature.The computational results demonstrate that the proposed algorithm is highly competitive and performs well in terms of solution quality and computational time.The Pareto set in the investigated type of facility location problems leads to solutions that may better support decision-making.
In article [21], a multi-objective multi-facility model of green closed-loop supply chain under uncertain environment is presented.In this model, the proposed green closed-loop supply chain considers three classes in the case of the leading chain and three classes in terms of the recursive chain.
In the study "A Hybrid Genetic Algorithm for Multi-emergency Medical Service Center Location-allocation Problem in Disaster Response," [22] the final patient mortality risk value (injury severity) caused by both initial mortality risk value and travel distance (travel time) is considered to determine the location-allocation of temporary emergency medical service centers.
According to the article "Multi-facilities Location and Allocation Problem of Three-Echelon Supply Chain Based on an Improved Genetic Algorithm," [23] a new mathematics model on multi-facilities location and allocation problem (MLAP) of three-echelon supply chain was built, which mainly took some logistics cost such as inventory cost, carrying cost, transportation cost into consideration.Then, an improved Genetic Algorithm was proposed to solve the MLAP.In this algorithm, the real encoding method was used to encode the solution directly, and meanwhile the crossover and mutation operator were improved.Finally, a numerical experiment showed clearly the feasibility of the IGA method and the effectiveness of the model.
In [24] the multi-objective facility layout problem is defined in the literature as an extension of the famous quadratic assignment problem (QAP).Most previous mathematical models tried to combine both the quantitative and the qualitative objectives into a single objective by using weighting factors.This paper introduces a multi-objective mathematical model and solves it using the revised Strength Pareto Evolutionary Algorithm.The purpose of work is to find an efficient set of solutions "Pareto optimal set" which could be introduced to the decision maker to select the best alternative, while considering conflicting and non-commensurate objectives.A computer program is developed to define the mathematical model, code candidate solutions into genetic form, and use Evolutionary Multi-Objective Optimization algorithms to find the efficient set of solutions.

Applications
The MFLP and related problems have a broad range of applications in various domains such as logistics, transportation, planning and scheduling, supply chains and many others.Foltin et al. [25] show future logistics applications in the military; they discuss the planning process and its subsequent implementation to ensure logistical support to deployed units and propose discrete event simulation for professional training.
Another example where the MFLP is used is the Tactical Decision Support System (TDSS).This system was developed to support commanders in their decision-making processes [26].This system is composed of a number of models of military tactics for optimization of the operation task at hand, see, e.g., [27][28][29][30][31][32].
From the above, it is clear that the authors use complex proprietary methods and software to solve the problem, which are complicated and inaccessible for the ordinary manager.

Main Aim of this Work
The article discusses the possibility of using the commercial software Excel for solving MFLP and describes the problems associated with it.Plane coordinates and calculation of distances in Euclidean space are used in the calculations.The calculated data must therefore be considered approximate, serving mainly as a basis for further decisions.
In order to verify its applicability, we solve five MFLP benchmark tasks with different levels of difficulty.The benchmark tasks are designed so that it is possible to calculate the optimal solution (with the exception of the last task, which is created based on a realworld problem).We then compare the optimal solution with the result obtained by calculation in the Excel application environment and calculation by the metaheuristic method (Simulated Annealing).The tasks are solved on two computers with different performance parameters in order to verify the influence of the hardware configuration on the accuracy of the calculation.Both methods are also used to solve a practical example based on real data.
Excel software and its Solver add-in are commercial products.This work does not examine the algorithms they use.From the point of view of this research, it is a "black box" for which the inputs, settings of basic parameters and outputs are known.
The proposed procedure uses available commercial software (Excel) for solving Multi-facility location problem, as it is described in Section 1.3.This solution has not yet been described in the scientific literature.The use of the proposed procedure will allow the required calculations to be performed on a personal computer with commonly available office software.No special competencies are required; user-level knowledge is sufficient.The aim of this work is to verify the proposed procedure on test examples and to find the problems that this simplified approach brings.

Problem Formulation
The aim is to find a solution to the following problem.There are b points where the coordinates of their positions are known.There is a requirement to find the coordinates of the positions of the c centers.Each point will be assigned to exactly one center which is closest to it.Variable c is an integer greater than or equal to 1 and b is an integer greater than or equal to c.
The objective is to find the positions of the centers so that the sum of the distances between all the centers and their assigned points will be minimal.
A practical example is a company with a number of existing branches which needs to set up a number of centers that will provide services to these branches.The locations of the branches are known; the centers will be built as new.The location of these centers and the branches they will serve are not specified.
The mathematical formulation of this problem is as follows.Let   ,  , … ,  be a set of points ( 1); the coordinates of these points are known.Let   ,  , … ,  be a set of centers (  ; the coordinates of these points are not known as they are the subject of the optimization. The Euclidean distance between a center and any point in two-dimensional space is calculated according to Formula (1).In general, any space, number of dimensions and distance function can be used.
where   ,  is the distance between center  and point  ;  ,  are coordinates of the center  , and  ,  are coordinates of point  .
The objective function of the MFLP problem is shown in Formula (2).It corresponds to the sum of distances between all points and their assigned centers which are closest to them.
The optimization goal is in Formula (3); the aim is to find the locations of c centers so that the objective function D is minimal.In the case of two-dimensional space, there are 2c continuous optimization variables.

𝑀𝑖𝑛𝑖𝑚𝑖𝑧𝑒 𝐷
(3) Software configuration.System type 64-bit operating system, Windows 10 Enterprise LSTC, MS Excel 2016, 64 bit, part of the Microsoft Office Suite, with the Solver addin installed, which performs calculations.In the basic version it allows only one method of calculation using genetic algorithms.The Excel environment is used to enter data, calculate intermediate results, and display the result of the entire process.

Benchmark Instances
To verify the applicability of the proposed solution, both from the accuracy and optimization time points of view, five benchmark instances differing in difficulty are proposed.
The first four instances (labeled A, B, C, D) are designed in such a way that the optimal solution can be easily determined.The points in these instances are arranged in groups of five, four at the vertices of the square, the fifth at its center.For each group, a center exists.The optimal position of the center for such a group is identical with the position of the fifth point in the center of the square.Instances A and D are shown in Figure 1 for illustration of the principle; blue circles represent points; orange circles represent centers in their optimal position.The fifth blue point in the center is covered by an orange circle.The last instance (labeled E) is created based on a real-world problem using real geographic data.There are 27 points and 4 centers to deploy.The coordinates of points are recorded in Table 1 and their visualization is shown in Figure 2. In contrast to benchmark tasks, for which it is possible to manually calculate the optimal position of the centers thanks to artificially designed points, for this type of task the optimal solution cannot be calculated in advance.The basic parameters and the optimal solution (if available) of all the benchmark instances are shown in Table 2 and in Figure 3.For each benchmark task, 100 optimizations were performed both in the Excel-Solver software as well as using the metaheuristic algorithm (see Section 2.4) and thus, 100 different solutions are obtained.
Each task was solved in 4 ways:  by Excel-Solver evolutionary method on PC1;  by Excel-Solver evolutionary method on PC2;  by metaheuristic algorithm on PC1;  by metaheuristic algorithm on PC2.
For each task, method and hardware configuration, the following parameters and outputs calculated from all individual solutions are recorded:  minimal solution;  average solution;  standard deviation;  difference between the optimal and the minimal solutions in percent;  difference between the optimal and the average solutions in percent;  histogram, classes are from 1 percent to 10 percent deviation;

Evolutionary Method
This method is used for computing in Excel software with the Solver add-in.In the available help and information from the manufacturer, the method is noted as evolutionary method, with the explanation that it is used to solve non-smooth problems.Convergence values, mutation frequency, base file size, random number, and maximum time without enhancement can also be set for this method.What the effect is of changes in these parameters on the calculation is not stated; therefore, default values were left for all calculations.Default values of parameters are listed in Table 3.

Metaheuristic Algorithm
The metaheuristic algorithm used to complement the Excel Solver is based on the Simulated Annealing (SA) principle adapted for the MFLP problem.This choice was based on the previous experiences of the authors as this algorithm proved to be very successful in other position optimization problems (see [23] for more details).The principle is inspired by annealing in metallurgy, where this process is used to reduce the defects of material by way of heating and controlled cooling.
Algorithm At the beginning (point 1), the initial deployment of centers is generated using the random number generator with uniform distribution.The value of every variable  is generated in range limited by the minimum and maximum value of all the points in the corresponding axes x or y.In point 2, this solution is evaluated using the objective function in Formula (2).
The algorithm works in iterations.During an iteration (points 4 to 18), the value of temperature T is constant.In each iteration, a number of transformations of the current solution  into solution  is performed (see point 7).This transformed solution replaces the original (points 10 to 13) with some probability given by the Metropolis criterion (4).A better solution is always replaced.However, the key is accepting even worse solutions, thus expanding the search space explored for the global optimum.In such cases, the probability depends on the difference between both solutions and the current temperature; the smaller the difference and the bigger the temperature, the bigger the probability.An iteration is terminated (point 6) when the maximum number of transformations or replacements is performed.
The transformation of a current solution is a simple process.One variable (i.e., x or y coordinate of one of the centers) is randomly selected.Then this variable is changed using the Formula (5).The size of the change depends on the current temperature; the bigger the temperature, the bigger the change.This causes the wide exploration of the state space at the beginning of the optimization and the refinement of the solution towards the end.

𝐶
RandN , where  is the variable selected to be changed in solution ;  is the changed variable in transformed solution  ; RandN is a random number generator with normal distribution with mean  and standard deviation ;  is the maximum range in which the variable can change given by the minimum and maximum coordinates in the corresponding axis x or y.
The algorithm is terminated (point 4) when the current temperature decreases below the lower threshold  .The best solution found during the process is stored (see points 14 to 16) and returned at the end of the algorithm (point 19).

Results
In this section, results obtained by the benchmark calculations are analyzed.For a better understanding, the terms and abbreviations used are explained in Tables 4 and 5. Tables 6-9 are four tables of the results of the calculations of benchmark instances A, B, C, D obtained by evolutionary and metaheuristic methods on PC1 and PC2.The results are also presented in the form of histograms and graphs (Section 3.1).The best solution for benchmark E is recorded in Section 3.2.

Optimal
Optimal solution for an instance.As the optimal solution for instance E is not known, the minimal distance computed by metaheuristic algorithm is used.

Minimal
Minimal value from set of 100 results

Average
Average value from set of 100 results Dif.optimal-min Difference between optimal and minimal value (%) Dif.optimal-avg Difference between optimal and average value (%) Standard dev.
Standard deviation for set of 100 results

Avg. time
Average time of computing one result in set of 100 results     4) and graphs with the best solution (Figure 5) are presented in this section only for benchmark instance D (20 centers for 100 points), because for this most difficult benchmark the illustration of the differences in the methods and hardware used is most obvious.
In each method 100 different results were obtained for this benchmark; they can be statistically evaluated using a histogram.The intervals are designed in the range of 101, 102 … 110 percent of the optimal value.The histogram shows the frequency of occurrence of the obtained values in a given interval.From the graph it is possible to assess how accurate the method is.For the metaheuristic method on PC1 and PC2, the histograms are the same.

Best Solution for Benchmark E
Coordinates of centers are in Table 10 and are calculated using MPC2 method and configuration, as this method found the best solution.The illustration of this best solution is shown in Figure 6.

Discussion
This work is focused on verifying the possibility of using Excel to solve MFLP-type tasks.The quality of the result depends on the HW used and the complexity of the task.Because the software works with the generated random number during the calculations, 100 calculations were performed for each of the variants, in order to partially eliminate this fact.
The results obtained using the metaheuristic algorithm are balanced in all benchmarks on both PCs; deviations from the optimal value are in the order of hundredths of percent.Thus, the method is verified and considered suitable for obtaining the reference value of the result in the case of benchmark E.
When using the evolutionary method, the influence of the complexity of the benchmark and the HW used is obvious.The EPC1 graph (Figure 7) compares two coefficients (Dif.optimal-min, Dif.optimal-avg-see Table 4), on which it is possible to present an influence of the complexity.Both values increase with the benchmark complexity.From the higher growth of the value of Dif.optimal-avg it can be concluded that for this type of task it is important to conduct more calculations and select the best; the value of Dif.optimalmin increases less.Thus, more calculations may partially eliminate inaccuracy.The values Standard deviation, Dif.optimal-avg and Dif.optimal-min can be used for illustration of the influence of HW used and method too.For all values, a smaller value is better.Figure 8     The advantage of performing a larger number of calculations on more powerful HW when using Solver is also evident from the histograms of benchmark D on Figure 4.The metaheuristic method has the same result on both PCs; for the evolutionary method there is an obvious difference in the quality of the results obtained on PC1 and PC2.

Conclusions
The obtained results show that the use of the Excel-Solver software for solving MFLP tasks is possible.The accuracy of the results depends on the complexity of the task and the HW used.Improving the results can be achieved by repeating the calculations and finding the lowest value in the result set.Even for the most demanding task from the used benchmarks, a deviation of the best solution from the optimal one of less than 1 percent can be achieved.In a real environment, it will not be possible to place the device exactly in the calculated locations, so these results can be considered accurate enough.Very good results have been achieved on the example of practice.
The use of the proposed method is especially suitable in situations where finding a solution is required and access to sophisticated methods for solving MFLP-type tasks is not possible.In education, it offers to acquaint students with the simple possibility of solving this type of problem and thus enable a wider use of optimization procedures in their future practice.

Figure 2 .
Figure 2. Visualization of benchmark instance E.


graph with coordinates for result for the minimal solution.

Figure 6 .
Figure 6.Representation of the best solution on map, left-only centers, right-centers and connected points.

Figure 7 .
Figure 7. Influence of the complexity of the benchmark.
presents decreasing value of Standard deviation, Figure 9 presents decreasing value of Dif.optimal-avg and Figure 10 decreasing value of Dif.optimal-min.The best visibility of influence is mainly for benchmark D.

Figure 8 .
Figure 8.Standard deviation for different configurations.

Table 1 .
Coordinates of points for instance E.

Instance Number of Points Number of Centers Optimal Solution
Figure 3. Graphical comparison of benchmark complexity.

Table 3 .
Parameters used for optimization in the Excel Solver.

Table 4 .
Definitions of terms used.

Table 5 .
Definitions of abbreviations used.

Table 10 .
Best solution found for benchmark E.