## 1 Introduction  
In the nanophotonics field, wavefront design is one of the timely ways of enabling expected functionalities for a flat and thin metasurface. These flat functional metasurfaces provide many applications such as beam deflection or focusing, etc., in which high efficiencies devices are required. The inverse design method based on forward and adjoint computations [1; 2; 3; 4; 5; 6; 7; 8], has proven to be one of the most efficient in practice for achieving high performance optimized structures. Like any gradient-based method, the adjoint-based-inverse-design method relies on the computation of the gradient of a objective function with respect to the design parameters. The gradient of the figure of merit (FOM) is used to modify the design parameters in order to increase the figure of merit (FOM). The key point, hence the main advantage of the adjoint-based strategy, is that it allows an evaluation of the gradient of the FOM with a cost independent of the number of the design parameters. In addition, this gradient vector is used, at each iteration, to simultaneously update all the design parameters.The simultaneous-pushing process of design parameters, in each iteration, lacks randomness and this may trap the gradient-based algorithm in basins of local minimum. In linear problems with a simply convex shape merit function, using the gradient information, it is relatively easy for some gradient-based algorithms to reach an optimal solution. Therefore these gradient-based algorithms are eminently suitable for finding an optimal solution when the starting point is close to a local minimum. However,in photonics, one often faces constrained and non linear problems stated in an infinite and indeterminate design space. Therefore, the conventional gradient-based optimization method requires multiple restarts initial points and multiple runs, hoping that a best single solution may be discovered at the end of these runs. This process may require higher computing memory capacity. On the other hand, global optimizers, continue to emerge, providing researchers a large choice of tools, in several disciplines. Some of them, enable an efficient random sampling of the design space, in each iteration, in order to detect and always maintain several better solutions during the optimization. This is the case for some metaheuristic algorithms [9; 10; 11; 12; 13; 14; 15]. Referring to the research-strategy, optimization algorithms are mainly based on four approaches: complete, random, heuristic, and metaheuristic. In the complete strategy, all possible combinations are evaluated in order to target the best solution. This strategy allows finding optimal solutions but it is computationally costly, which makes it impractical for high dimension problems. The random search is a way to attempt to overcome the drawback raised by the complete research-strategy. This strategy considers many initial candidates as possible solutions. In the best case, an optimal solution can be found earlier among these candidates, in the optimization process. However, this random strategy may behave as a complete method in its worst case. In heuristic algorithms, heuristic information which highly depends on the problem at hand is enforced in order to indicate certain ways to find optimal solutions. This strategy cannot be generalized since each problem involves (requires) a specific heuristic information. In the metaheuristic algorithm, few assumptions about the optimization problem being solved are made. This makes it possible to apply them for a variety of problems contrary to the heuristic algorithms. These methods range from Genetic algorithm [9], particle swarm optimization [10], ant colony optimization [11], differential evolution (DE) [12], grey wolves, whales and Harris Hawks [13; 14; 15], and the slime mould algorithm [16]. In general, they can be performed based on two approaches: a single-solution-based approach and a population-based approach.