## 5 Issues for Numerical Solutions  
In general, numerical solving techniques are independent of the formulation of an optimization problem. However, a proper solver to address an optimization problem can be closely related to its characteristics in such a way that the structure of an optimization problem can suggest convenient solving algorithms. Deterministic approaches that are based on gradient calculations impose rigorous conditions over a mathematical formulation structure, so the original problem might be reformulated for a convenient representation to apply a specific numerical solving technique as described in [27]. The main advantage of a deterministic solver is that once all algorithm requirements are fulfilled, it is possible to guarantee a locally optimal solution. Another option is the use of a heuristic technique that does not require any reformulation and can be applied directly to the formulated problem [28]. There are multiple options to select a solver algorithm that can fit an optimization problem structure. The development of new algorithms is still a topic of interest, as pointed out in [28]. This work uses the Differential Evolution (DE) algorithm proposed in [29] along with the feasibility rules by Deb as a constraint manager [30], and the use of Pareto-dominance approach for solving the presented optimization problem. This research selected DE as a first option because it is a popular heuristic algorithm extensively treated in the literature and is an excellent alternative to solve optimization problems formulated from real-world design problems described in [31]. However, other heuristic algorithms or any combination with gradient-based algorithms can be used to solve the optimization problem. Among other advantages, DE has a a few parameters compared to other heuristic algorithms and its implementation is fairly simple. DE heuristic algorithm is a continuous stochastic search method that is based on evolution strategies [29]. The implementation of this heuristic method is represented with a scheme in Figure 6 and a detailed computation routine is shown in Algorithm 1. An initial population is generated using a uniform random distribution within the bounds that are described in Equation (7) or within the parameter space A as in Figure 5, where the elements of the generated candidates correspond to the design variables. The external memory data base (EMdb) corresponds to a variable that stores all the feasible non-dominated solutions, which this data base is initialized with an empty value. The next step is to evaluate all the initial population candidates in the objective functions, inequality, and equality constraints of Equations (1)-(7). The crowding distance of an arbitrary point in the search space's domain corresponds to an estimation of the distance from the nearest neighbors along each of the objectives to measure the population's density. The higher the crowding distance, the higher is its isolation. Therefore, when the crowding criterion in Figure 6 is active, the algorithm prioritize the search around isolated solutions to raise the number of candidates, as pointed out in [28]. A crowding criterion selects the percentage of the total number of iterations computed using individuals from the current population or EMdb. Therefore, a crowding parameter (\(DE_{CRW}\)) tunes the algorithm's behavior to explore the search space or to exploit the search near to a found candidate that remains distant from other found solutions. DE algorithm is based on the use of two operators, which are the mutation and crossover. These two operators are used to search in the parameter space, and they are applied for every design variable. The Crossover operator in Figure 6 determines the probability of selecting a value that can be computed using the mutation operator or the value from a previous iteration as in line 15 of the pseudo-code in Algorithm 1.