**1. Introduction**

Combinatorial optimization problems appear naturally in the areas of engineering and science. From the research point of view, these problems present interesting challenges in the areas of operations research, computational complexity, and algorithm theory. Examples of combinatorial problems are found in, scheduling problems [1,2], transport [2], machine learning [3], facility layout design [4], logistics [5], allocation resources [6,7], routing problems [8,9], robotics applications [10], civil engineering problem [11–13], engineering design problem [14], fault diagnosis of machinery [15], and social sustainability of infrastructure projects [16], among others. Combinatorial optimization algorithms should explore the solutions space to find optimal solutions. When the space is large, this cannot be fully explored. Combinatorial optimization algorithms address this difficulty, reducing the effective space size, and exploring the search space efficiently. Among these algorithms, metaheuristics have been a good approximation to obtain adequate solutions. However, addressing large instances and requiring almost real-time solutions for some cases, motivates lines of research that strengthen the methods that address these problems.

One way to classify metaheuristics is according to the search space in which they work. In that sense, we have metaheuristics that work in continuous spaces, discrete and mixed spaces [17]. An important line of inspiration for metaheuristic algorithms are natural phenomena, many of which develop in a continuous space. Examples of metaheuristics inspired by natural phenomena in continuous spaces are: Optimization of particle swarm [18], black hole [19], cuckoo search (CS) [20], Bat algorithm [21], Algorithm of fire fly [22], fruit fly [23], artificial fish swarm [24], gravitational search algorithm [25], among others. The design of discrete versions of these algorithms carries important challenges because they must try to preserve the intensification and diversification properties [17].

With the aim of modifying and improving the results of a metaheuristic algorithm, the main tool available to them is tuning their control parameters. These types of modifications correspond to small modifications of the algorithm since the operating mechanism of the algorithm is not altered. However, many optimization problems cannot be efficiently addressed by an algorithm through modification of its parameters and require deep modifications that alter its mechanism of operation [17,26]. A strategy that has strengthened the results of metaheuristic algorithms has been the hybridization of these with techniques that come from the same and the other areas. Hybridization strategies involve deep changes to the original algorithm. The main hybridization proposals found in the literature are the following: (i) matheuristics, where mathematical programming is combined with metaheuristic techniques [27], (ii) hybrid heuristics, hybridization between different metaheuristic methods [28], (iii) simheuristics, that combine simulation and metaheuristics [29], and (iv) the hybridization between machine learning and metaheuristics. The latest line of research exploring the integration between the areas of machine-learning and metaheuristic algorithms is an emerging line of research in the areas of computing, mathematics, and operational research. In developing the state-of-the-art, we find that hybridization occurs primarily for two purposes. The first, with the goal that metaheuristics help machine-learning algorithms improve their results (for example, [30,31]). In the second intention, machine-learning techniques are used to strengthen metaheuristic algorithms (for example, [32,33]). The details of the hybridization forms are specified in Section 2.

In this article, inspired by the research lines mentioned above. A hybrid algorithm is designed, which explores the application of a machine-learning algorithm in a discrete operator to allow continuous metaheuristics to address combinatorial optimization problems. We will apply our hybrid algorithm to the discrete problem of the design of counterfort retaining walls. The contributions of this work are detailed below:


The structure of the remaining of the paper is as follows: a state-of-the-art of hybridizing metaheuristics with machine learning is developed in Section 2. In Section 3 the optimization problem, the variables involved, and the restrictions are defined. The discrete k-means algorithm is detailed in Section 4. The experiments and results obtained are shown in Section 5. Finally, in Section 6 the conclusions and new lines of research are summarized.

#### **2. Hybridizing Metaheuristics with Machine Learning**

Metaheuristics form a wide family of algorithms. These algorithms are classified as incomplete optimization techniques and are usually inspired by natural or social phenomena [17,35]. The main objective of these is to solve problems of high computational complexity and they own the property of not having to deeply alter their optimization mechanism when the problem to be solved is modified. On the other hand, machine-learning techniques correspond to algorithms capable of learning from a dataset [36]. If we make a classification according to the way of learning, there are three main categories: Supervised learning, unsupervised learning, and learning by reinforcement. Usually, these algorithms are used in problems of regression, classification, transformation, dimensionality reduction, time series, anomaly detection, and computational vision, among others.

In the state-of-the-art of algorithms that integrate machine-learning techniques with metaheuristic algorithms, we have found two main approaches [26]. In the first approach, machine-learning techniques are used in order to improve the quality of the solutions and convergence rates obtained by the metaheuristic algorithms. The second approach uses metaheuristic algorithms to improve the performance of machine-learning techniques [26]. Usually, the metaheuristic is responsible for solving more efficiently an optimization problem related to the machine-learning technique. Adapted and extended from [26], in Figure 1, we have proposed a general scheme of techniques where machine learning and metaheuristics are combined.

**Figure 1.** General scheme: Combining Machine Learning and Metaheuristics.

When we analyze the integration mechanisms of the first approach, we identify two lines of research. In the first line, machine-learning techniques are used as metamodels in order to select different metaheuristics, therein choosing the most appropriate for each instance. The second line aims to use specific operators that make use of machine-learning algorithms and subsequently this specific operator is integrated into a metaheuristic [26].

In the general integration of machine-learning algorithms on metaheuristic techniques we find three main groups: algorithm selection, hyper-heuristics, and cooperative strategies. In algorithm selection, the goal is to choose from a set of algorithms and considering the associated characteristics for each instance of the problem, an algorithm that works best for similar instances. Another way to approach the problem is to automate the design of heuristic or metaheuristic methods to tackle a wide range of problems, we call this hyperheuristic approximation. Finally, the objective of cooperative strategies is to mix algorithms in a parallel or sequential way to obtain more robust methods. The cooperation mechanism can be complete, i.e., sharing the complete solution, or partial when only part of the solution is shared. In [32], the problem of docking scheduling in massive terminals

was addressed by algorithm selection techniques. In [37], the nurse training problem through a tensor-based hyperheuristic algorithm was addressed. Finally, a distributed framework that uses agents was proposed in [38]. In this framework, each agen<sup>t</sup> corresponds to a metaheuristic and the framework can adapt through direct cooperation. This framework was applied to the problem of permutation flow stores.

Additionally, a metaheuristic incorporates operators that allow strengthening its performance. Examples of these operators are initialization operators, solution perturbation, population management, binarization, parameter configuration, and local search operators, among others [26]. Specific integrations explore the machine-learning application in some of these operators [26]. In the design of binary versions of algorithms that work naturally in continuous spaces, we find binarization operators in [2]. These binary operators use unsupervised learning techniques to perform the binarization process. In [39], the concept of percentile was explored in the process of generating binary algorithms. In addition, in [1], the Apache spark big data framework was applied to manage the population size of solutions to improve convergence times and the quality of results. Another interesting research line was found in the adjustment of metaheuristic parameters. In [40], the parameter setting of a chess classification system was implemented. Based on decision trees, and using fuzzy logic, a semi-automatic parameter setting algorithm was designed in [41]. Usually, in the initiation of solutions, a random mechanism is used. However, using machine-learning techniques, algorithms have been developed that improve the initiation stage of a solution. Applied to the problem of designing a weighted circle in [42], a case-based reasoning technique was used to start a genetic algorithm. In [43] Hopfield neural networks were applied to initiate solutions of a genetic algorithm. This initiation was applied to an economic dispatch problem.

To improve the quality of machine-learning algorithms, we see that metaheuristics have a broad contribution. Contributions are found in problems of feature selection, classification, clustering, extraction of characteristics, among others. In the identification of breast cancer, in [44] a genetic algorithm was designed that improves the performance of support vector machine (SVM) when applied in image analysis. The genetic algorithm was used in the characteristic extraction stage from the images. In [45] a multiverse algorithm was used to adjust the parameters of an SVM classifier. An improved monarch butterfly algorithm was used in [46] with the goal to optimize the weights of a feed-forward neural network. In [47], a geotechnical problem was addressed. In this article, a firefly algorithm was integrated with the least-squares support vector machine technique. When considering regression problems, there are again several cases where metaheuristics have contributed. For example, in predicting the strength of high-performance concrete, in [48] a metaheuristic was used to improve the least-squares technique. In [49], share price prediction was improved by integrating metaheuristics into a neural network. Again, in the prediction of the share price of construction companies in Taiwan, in [30], a sliding-window metaheuristic-optimized machine-learning algorithm was designed, which integrates a metaheuristic in its learning process. Optimization of the parameters of a least-squares regression was improved, in [50] using a firefly algorithm. We also found metaheuristic applications for unsupervised learning techniques. In particular, there are significant contributions from metaheuristics to clustering techniques. For example, in [51], we found an algorithm application that integrates a metaheuristic with a kernel intuitionistic fuzzy c-means. This algorithm was applied to different data sets. Another integration found for clustering techniques is the search for centroids. This problem is to find the set of centroids that best group the data under a certain metric. This is an NP-hard type problem and therefore it is natural to approach it with a metaheuristic. In [52] an algorithm was proposed that uses a bee colony to solve a problem of grouping energy efficiency data obtained from a wireless sensor network. In [53], a metaheuristic was used to find centroids in planning the transportation of employees from an oil platform via helicopters. In the case of neural networks, metaheuristics have contributed to different types of integration. In [54] a metaheuristic algorithm is used to automatically find the appropriate number of neurons in a feed-forward neural network. In [55] a metaheuristic algorithm allows determining the dropout probability in convolutional neural

networks. The cuckoo search and bat algorithm algorithms were used in [56] to adjust the weights of neural networks of feed-forward type. In [57] an algorithm was designed using metaheuristics to optimize convolutional neural networks. An application to long term short memory (LSTM) neural network training was studied in [58]. LSTM networks were applied in healthcare analysis.

Inspired by the lines of research detailed above, this work proposes a hybrid algorithm in which the unsupervised k-means learning technique is used to obtain binary versions of the cuckoo search optimization algorithm. This hybrid algorithm was used to solve the problem of the counterfort retaining walls. The k-means technique has been widely used and in recent studies, it has been applied in [59] to bioinformatics for detecting gene expression profile, image segmentation for pest detection [60] in agriculture, and brain tumor identification [61], among others. Particularly the k-means technique has been previously applied in obtaining binary versions of continuous metaheuristics and used to solve the multidimensional knapsack problem [33] and the set covering problem [5] which are NP-hard problems. In the case of the cuckoo search, hybrid versions of this algorithm have been applied in [62] to numerical optimization problems in engineering. In [63] cuckoo search was applied to benchmark functions and engineering design problems. An algorithm with reinforced learning was used in [64] to solve a logistic distribution problem. In [65] an improved version of the cuckoo search was applied to drone location.

#### **3. Problem Definition**

In this Section, we will detail the buttressed earth-retaining walls optimization problem. In Section 3.1 we will give the definition of the optimization problem. Later in Section 3.2, we will detail the design variables. Then the design parameters will be described in Section 3.3, and finally, in Section 3.4, we define the constraints.

#### *3.1. Optimization Problem*

The optimization problem considers two objective functions, which will be addressed independently. The first function corresponds to the cost (*pi*) of the wall construction, expressed in euros. The second one considers the CO2 equivalent emissions units (*ei*). These construction units correspond to formwork, materials, excavation and earth-fill. The cost and emission functions are based ona1m wide strip [34]. The emission and cost values were obtained from [34,66] and are shown in Table 1. Then, in a general way, our optimization problem is defined according to Equation (1).

$$O(\mathbf{x}) = \sum\_{i=1}^{r} a\_i \mathbf{x}\_i \tag{1}$$

where *x* is the set of decision variables and *ai* ∈ {*<sup>c</sup>*,*<sup>e</sup>*}, corresponds to cost or emissions. Additionally, the optimization problem is subject to a set of restrictions determined by the ultimate (ULS) and service (SLS) state limits.

#### *3.2. Problem Design Variables*

In the design of the buttressed retaining wall to study, three groups of variables are defined. The geometric variables, the concrete and steel grade, and the passive reinforcement of the wall. In total there are 32 variables where each of these has a discrete amount of possibilities. In the group of concrete and steel grades. Concrete HA-25 to HA-50 is considered in discrete intervals of 5 MPa. In the case of steel, the B500S and B400S types are considered. In the group of geometric variables, there is the thickness of the footing (*c*), thickness of the stem (*b*), the length of the toe (*p*), the thickness of the buttresses (*ec*), the length of the heel (*t*) and the distance between buttresses (*d*). The last group of variables that correspond to the passive reinforcement of the wall. There are 24 variables shown in Figures 2 and 3. The diameter and the number of bars define the reinforcement. Three reinforcement flexural bars defined as A1, A2 and A3 contribute in the main bending of the stem. The vertical reinforcement of foundation in the rear side of the stem is given by A4, up to a height L1. The secondary longitudinal reinforcement is given by A5 for shrinkage and thermal effects in the stem. The longitudinal reinforcement of the buttress is given by A6. The area of reinforcement bracket from the bottom of the buttress is given by A7 and A8. The upper and bottom heel reinforcement are defined by A9 and A11 and the shear reinforcement in the footing by A12. The longitudinal effects in the toe are defined by A10. The set of combinations of the values for the 32 variables shown in Table 2, constitutes the solution space.


**Table 1.** Unit breakdown of emissions and cost.

**Figure 2.** Reinforcement variables for the design of earth-retaining walls. *Source:* [67].

**Figure 3.** Earth-retaining buttressed wall. Floor cross-section. *Source:* [67].



(1) Number of values depends on the height.

#### *3.3. Problem Design Parameters*

The data that will remain fixed in the optimization process will be called problem design parameters. The main design parameters are shown in Figure 4. The height of the wall and the depth of the soil in front of the wall are denoted by H and H2, respectively. The maximum bearing pressure is the soil foundation ultimate bearing capacity divided by the bearing capacity factor of safety. The maximum support pressure for the operating conditions is represented by *σ*, the base-friction coefficient is *μ* and the backfill slope at the top of the rod is *β*. The density, the internal friction angle and the friction angle that determine the earth's pressure angle are given respectively by P (*γ*, *φ*, *δ*). The roughness between the wall and the fill is determined by a fraction of *φ*. The values of the problem design parameters are shown in Table 3.


**Table 3.** Problem design parameters values.

**Figure 4.** Problem design parameters. *Source:* [67].

#### *3.4. Problem Constraints*

The feasibility of the structure is verified in accordance with the Spanish technical standards defined in [68] and the recommendations detailed in [69]. The flexural and shear limit states are checked. The structure is verified in accordance with the approach specified in [70]. To check the structure limit states, a uniform surface load at the top of the fill is considered [71]. For the active earth pressure calculation the surface loads and fill are considered. The major forces in wall analysis consider wall weight, heel backfill load, surface load, earth pressure, weight at the front toe, and passive resistance against the toe. Buttresses receive a load equivalent the product of the distance between the buttresses by the pressure distribution in the stem.

The structural model considers that the upper part of the stem works as a cantilever, while the lower part of the stem is strongly coerced by two elements: the footing and the lower part of the buttress, located at the rear of the stem. Calculation bending moments are taken in the midsection between the buttresses and are given by *B*1 and *B*2 described in Equations (2) and (3) respectively.

$$B\_1 = -0.03p\_1d(H-\mathcal{c})\tag{2}$$

$$B\_2 = -0.0075p\_1d(H-c)\tag{3}$$

*B*1 is the bending moment at the connection of the stem to the footing, *B*2 is the maximum bending moment on the stem and *p*1 represents the pressure over the slab on the upper side of the footing. When the spacing of the buttresses is less than 50% of the height, Equation (4) defines the shear resistance (*s*) at the connection of the plate to the footing. For an accurate estimate of the moments in each section of the stem, as a result of the vertical bending stress in the stem, the trapezoidal pressure distribution is considered [70]. 50% of the maximum pressure at the top of the foundation is taken as the maximum value. Taking into account the vertical bending moment in the upper quarter part of the stem may be insignificant due to the involvement of [70].

$$s = 0.4p\_1d\tag{4}$$

Verification of the bending stress in the T-shaped horizontal cross-section is made considering the effective width, according to [72]. Mechanical bending and shear capacity are evaluated using the equations expressed in [71]. In this manual, the construction limit states of the EHE Structural Concrete Code are considered. The checks against sliding, overturning, and soil stresses, are carried out taking into account the effect of the buttresses, and are given in Equations (5)–(7). In the overturning check, Equation (5) ensures that the favorable overturning moments are sufficiently greater than the unfavorable overturning moments. In Equation (6), *Bo f* is the total favorable overturning moment. In Equation (7), *Bou* is the unfavorable total overturn moment, and the overturning safety factor is *γto* and is taken as 1.8 for frequent events. Equation (8) represents the reaction of soil against sliding. As *μ* is the base-friction coefficient, *N* corresponds to the total sum of the ground and wall weights located at the heel and toe, and *Ep* determines the passive resistance against the toe defined by Equation (9).

$$B\_{of} - \gamma\_{to} B\_{on} \stackrel{>}{\geq} 0\tag{5}$$

$$B\_{of} = N'(\frac{B}{2} - e\_p) - E\_p(H\_l - h') \tag{6}$$

$$B\_{\rm vol} = E\_h \* h\_\varepsilon - E\_\nu(\frac{B}{2} - f) \tag{7}$$

$$R = N'\mu + E\_p \tag{8}$$

$$E\_p = \frac{1}{2\gamma (H\_t^2 - (H\_t - c))^2} \frac{(1 + \sin(\phi))}{(1 - \sin(\phi))}\tag{9}$$

#### **4. The K-Means Discrete Algorithm**

The first stage of the algorithm corresponds to generating valid solutions randomly. In the generation procedure, it is first validated if all the variables that make up the solution have been started. If not, the variables are started randomly. Once all the variables are generated, the next step is to validate if the solution obtained is feasible. In the event that it is not feasible, all the variables are cleaned, and they have generated again. The detail of the initiation procedure is shown in Figure 5. Once we have a valid set of solutions, continuous metaheuristics are used in this case CS, to produce a new solution in continuous space. The CS algorithm will be described in Section 4.1. Subsequently, we apply the k-means operator to transform continuous movements into transition probabilities. The k-means operator will be detailed in Section 4.2. The result of the transition probabilities generated by k-means are used by the discretization operator to generate a new solution. The discretization operator is detailed in Section 4.3. Finally, the new solution is validated and, in the case that it complies with the restrictions, it is compared with the best solution obtained. In the case that the new one is superior, it is replaced. The detail flow chart of the solution is shown in Figure 6. , p r,

**Figure 5.** Solution initiation procedure.

**Figure 6.** The discrete k-means algorithm flow chart.

#### *4.1. Cuckoo Search Algorithm*

The phenomenon of cuckoo species, which lay their eggs in the nests of other bird species, has inspired the CS algorithm. Such is the level of sophistication of cuckoo birds that in some cases even the colors and patterns of the eggs of the chosen host species are mimicked. In the analogy, an egg corresponds to a solution. The concept behind the analogy is to use the best solutions (cuckoos) with the aim of replacing those that do not perform well. The CS algorithm uses three basic rules:


The algorithm pseudo-code is shown in Algorithm 1.

**Algorithm 1** Cuckoo search algorithm


#### *4.2. k-Means Operator*

The goal of the k-means operator is to group the different solutions obtained by executing continuous metaheuristics in this case CS. However, this operator can be applied to any continuous swarm intelligence metaheuristics. When we consider solutions as particles, we understand the position of the particle as the location of the solution in the search space, while the velocity represents the vector of the transition of the solution from the *t* iteration to the *t* + 1 iteration.

In the k-means operator (kmeanOp), the k-means clustering technique is used to make transitions in a discrete space. The proposal is based on using the movements generated by the CS metaheuristic in each dimension for all the particles. Let *x*(*t*) be a solution in iteration t, then <sup>Δ</sup>*<sup>i</sup>*(*x*(*t*)) represents the magnitude of the displacement <sup>Δ</sup>(*x*(*t*)) in the i-th position, considering the iterations *t* and *t* + 1. Once all the displacements have been calculated, they are grouped using the magnitude of the displacement |Δ*<sup>i</sup>*(*x*(*t*))|, we do not consider the sign. This grouping is done using the k-means technique where k represents the number of groups used. Finally, we propose a generic function *Ptr* shown in Equation (10) with the objective of assigning a transition probability to each displacement.

$$P\_{\rm tr}: \mathbb{Z}/\mathbb{k}\mathbb{Z} \to [0, 1] \tag{10}$$

Then using the function *Ptr*, a probability is assigned to each group obtained from the clustering process. In this article, we use the linear function given in Equation (11). Where Clust (*x<sup>i</sup>*) indicates the location of the group to which <sup>Δ</sup>*<sup>i</sup>*(*x*) belongs. The coefficient *α* represents the transition probability and the coefficient *β* models the transition separation for the different groups. Both parameters must be estimated. The pseudo-code of the discretization procedure is shown in Algorithm 2.

$$P\_{tr}(\mathbf{x}^i) = P\_{tr}(\text{Clust}(\mathbf{x}^i)) = \boldsymbol{\mathfrak{a}} + \beta \text{Clust}(\mathbf{x}^i)\mathbf{a} \tag{11}$$

#### **Algorithm 2** k-means operator

1: **Function** kmeanOp(*lx*(*t*), *lx*(*t* + 1))2:**Input***<sup>x</sup>*(*t*),*x*(*<sup>t</sup>* + 1)

 3: **Output** lTranProb(t+1)

4: *<sup>l</sup>*Δ*<sup>i</sup>*(*x*(*<sup>t</sup>* + 1)) ← getDelta(*lx*(*t*), *lx*(*t* + 1))


#### *4.3. Discretization Operator*

As input the operator receives the transition probabilities (*lTranProb*(*t* + 1)) obtained in the previous stage and the list of solutions *lx*(*t*). For each solution *x*(*t*) ∈ *lx*(*t*), we consider each component *x<sup>i</sup>*(*t*) and evaluate if it should be modified according to its transition probability. The transition probability is compared with a random number *r*1. In case the change is selected, the movement can increase the value (+1) or decrease it (−1). Finally, the selected value is compared with the best value obtained by the algorithm and remains with the minimum of both. The pseudo-code of the discrete procedure is shown in Algorithm 3.

**Algorithm 3** Discretization operator.

```
1: Function DiscOp(lTranProb(t + 1), lx(t))
2: Input lTranProb(t + 1)
3: Output x(t + 1)–Where x(t+1) is discrete.
4: movement = 0
5: for xi ∈ x(t) ∈ lx(t) do
6: if r1 > 0.5 then
 7: movement = 1
8: else
9: movement = −1
10: end if
11: xi = max(1,min(xibest,xi + movement))
12: end for
```
#### **5. Results and Discussion**

This Section aims to describe the results obtained by the hybrid algorithm when applying it to the counterfort retaining walls problem. In the first stage, the contribution of the k-means operator in the discretization process will be determined. This is described in Section 5.2. Later, in a second stage, we will compare our proposal with another algorithm that has solved the problem, Section 5.3. For each problem, we make 30 independent runs. The value 30 is widely accepted for statistical conclusions [73]. The Wilcoxon signed-rank [74] method was the test selected to determine if the difference is statistically significant. The *p*-value used was 0.05. The selection of the test is based on the methodology proposed in [73]. In this methodology, the Shapiro–Wilk or Kolmogorov—Smirnov-Lilliefors normality test is applied first. In the event that one of the populations is not normal, and both populations have the same number of points, the Wilcoxon signed-rank is suggested to verify the difference. In our case, the Wilcoxon test was applied to the entire population of instances considering all executions and comparing the results of both algorithms. For the execution of the instances, a laptop with a Core i7-4770 as a processor and 16GB in RAM has been used. The algorithm was programmed in Python 3.6.

#### *5.1. Parameter Settings*

To calibrate the parameters, heights 8 and 12 were used. These heights were selected with the intention of considering walls of different complexity. 8 represents small size walls and 12 represents larger size walls where the satisfaction of stability restrictions is much more difficult. Subsequently, each configuration was executed 5 times for each height considering all the configurations proposed in Table 4. The value 5 was chosen with the intention of being able to execute all the combinations in a limited time. In Table 4, the range column shows the scanned values to perform the CS adjustment. These ranges were inspired by previous studies in which the k-means technique was applied to solve the multidimensional knapsack problem [33] and the problem of the set covering [5]. The flow chart of the parameter settings is shown in Figure 7. For more detail on the method, reference [33] can be consulted.

To determine the best configuration, the method proposed in [33] was used. In this method, 4 measurements are defined that I know are shown in Equations (12)–(15). The BestGlobalValue represents the best value obtained considering all the configurations. The BestLocalValue, corresponds to the best value obtained in the evaluated configuration. The minGlobalTime is the minimum convergence time considering all settings. Finally, the convergenceLocalTime represents the average convergence time for a configuration. Each of the measures defined, the closer to 1, the better performance that indicator is having. On the other hand, the closer to 0 the worse performance. Since there are 4 measurements to be able to make the evaluation, we incorporate them into a radar chart and calculate their area. As a consequence of the measurement definition, the largest area corresponds to the configuration that has the best performance considering the 4 defined measures.

1. The deviation of the best local value obtained in five executions compared with the best global value:

$$bSolution = 1 - \frac{BestGlobalValue - BestlocalValue}{BestGlobalValue} \tag{12}$$

2. The deviation of the worst value obtained in five executions compared with the best global value:

$$NewSolution = 1 - \frac{BestGlobalValue - WordLocalValue}{BestGlobalValue} \tag{13}$$

3. The deviation of the average value obtained in five executions compared with the best global value:

$$aSolution = 1 - \frac{BestGlobalValue - - -AverageLocalValue}{BestGlobalValue} \tag{14}$$

4. The convergence time for the average value in each experiment is normalized according to Equation (15).

$$nTime = 1 - \frac{convergenceLocalTime - minGlobalTime}{maxGlobalTime - minGlobalTime} \tag{15}$$

**Figure 7.** Flow chart of parameter setting.


**Table 4.** Parameter setting for the Cuckoo search algorithm.

#### *5.2. Random Operator*

This experiment aims to quantify the contribution of the k-means operator in the process of optimizing the cost of the wall. To achieve this objective, a random operator was designed which replaces the k-means operator described in Section 4.2. This random operator behaves very similarly to that of k-means but instead of assigning transition probabilities per cluster, it assigns a fixed transition probability with a value of 0.5. In the experiment, wall height values between 6 and 14 m are considered. The k-means operator is compared with the random operator through the best value and the average obtained over 30 runs for each height. Later we will use the violin charts to compare both distributions of results and finally the Wilcoxon signed-rank test to validate that the difference is significant.

The results are shown in Table 5 and Figure 8 for cost optimization. In the case of emission optimization, the results are detailed in Figure 9 and in Table 6. In the case of cost optimization when analyzing the best value indicator, we observe that for all heights the k-means operator obtains better results compared to the random operator. In the case of a height of 6 and 7 m, this difference does not reach 2%. However, as the wall grows, the difference increases, reaching in the case of a height of 14 (m), 26%. The average indicator has similar behavior to that previously analyzed. The k-means operator has better performance in all cases. When applying the Wilcoxon test, it indicates that the difference is significant. When analyzing the distributions shown in Figure 8, we observe that the interquartile ranges obtained by the k-means operator are completely displaced to values close to zero with respect to the results obtained by the random operator. On the other hand, the dispersion of the distribution is greater in the solutions obtained by the random operator. In the case of optimization of CO2 emissions, when analyzing the best value indicator, we observe that the k-means operator obtains better results at all heights. However, the maximum difference from the random operator does not exceed 10%. The average indicator is again higher in the case of k-means, the Wilcoxon test indicating that the difference is significant. When analyzing the violin plots of the emissions shown in Figure 9, we see that the interquartile range obtains better quality results in the case of the k-means discretization. However, the result is not as remarkable as in the case of cost optimization.


**Table 5.** Comparison between random and k-means operators in cost optimization.

**Figure 8.** Cost violin plots, comparison between k-means and random operators.



**Figure 9.** Emission Violin plots, comparison between k-means and random operators.

Figure 10a,b show the convergence diagrams for cost optimization for the algorithm using k-means and random, respectively. The heights chosen were 6, 9, 12 and 14. From Figure 10a it is

observed that heights 6 and 9 have a better convergence than in the case of 12 and 14. The same effect is observed for the case of the random operator, Figure 10b, although in the latter case it is not so notorious. When we make a comparison between the graphs, we observe for the random case the slope stabilizes much earlier than in the case of k-means, as a consequence of this the latter obtains better results.

**Figure 10.** K-means and Random convergence chart for cost optimization.
