## 2 Methods
### CUDA Implementation of Greedy Agglomerative Procedures
We examined the following algorithms: (a) Lloyd: the ALA algorithm in the multi-start mode; (b) j-means: j-means algorithm (regular search in SWAP\({}_{1}\) neighborhood in combination with the ALA algorithm) in the multi-start mode; (c) AGGL\({}_{r}\): randomized search in the AGGL\({}_{r}\) neighborhood, \(r=\overline{1},p\); (d) SWAP\({}_{r}\): local search in SWAP, neighborhoods, \(r=\overline{1},p\) (only the best result for \(r=\overline{1},p\) is given); (e-g) GH-VNS1, GH-VNS2, GH-VNS3: Variable Neighborhood Search algorithms with neighborhoods formed by application of AGGL\({}_{r}\)(0) procedure, see [52]; (h) GA-1 POINT: genetic algorithm with a standard 1-point crossover; (i) GA-UNIFORM: genetic algorithm with a standard 1-point crossover and uniform random mutation [56]; (i) Aggl-EA: our new algorithm. For all datasets, 30 attempts were made to run each of the algorithms (see Tables 1 and 1, Tables 2-15 in Appendix A). In Table 1, we present the results of our new algorithm Aggl-EA(0) and the best of listed known algorithms, i.e., known algorithms which provided the best average or median values of the objective function (1) after 30 runs. For the smallest and simplest test problems (see Tables 1 and 4 in Appendix A), several algorithms including Aggl-EA() and search in the _AGGL_ neighborhoods resulted in the same objective function value (standard deviation of the result is 0.0) which is probably the global minimum. Nevertheless, the Aggl-EA() algorithm outperforms both genetic algorithms and simplest local search algorithms such as Lloyd() or j-means which do not reach this minimum value in all attempts. For more complex test problems, the comparative efficiency varies. However, there were no test problems for which our new algorithm demonstrated a statistically significant disadvantage in comparison with the best of known algorithms. The best (minimum), worst (maximum), average, and median values of the achieved objective function were averaged after 30 runs of each algorithm with fixed time limitation. The best average and median values of the objective Function (1) are underlined. We compared our new Aggl-EA (0) algorithm with known examined algorithm having the best median or average results. The significance of advantage/disadvantage of Aggl-EA() algorithm was estimated with the \(t\)-test [106; 107] and non-parametric Wilcoxon rank sum test [108; 109]. In the analysis of algorithm efficiency, the researchers often estimate the computational expenses as the number of the objective function calculations. The ALA procedure is the most computationally expensive part of all examined algorithms. In the ALA algorithm, the objective Function (1) is never calculated directly. The program code profiling results show that Step 2 of Algorithm 1 and Step 3 of Algorithm 2 occupy more than 99% of computational resources (processor time). These steps were implemented on GPUs. Moreover, we use the same implementation of the ALA procedure and BasicAggl() algorithm for all examined algorithms. Thus, the consumed time could be proportional with the number of iterations of the ALA algorithm (its Step 2) and the most time-consuming iterations of Algorithm 2 (its Step 3). However, the time consumption of these iterations depends on the number of centers which successively decreases during the work of the BasicAggl (0) algorithm. Thus, for the p-median problems, we used the astronomical time as the most informative unit for the computational expenses. \begin{table}
\begin{tabular}{c c c c c c} \hline \hline \multirow{2}{*}{**Algorithm or Neighborhood**} & \multicolumn{4}{c}{**Achieved Objective Function Values (1) Summarized after 30 Runs**} \\ \cline{2-5}  & **Min (Record)** & **Max (Worst)** & **Average** & **Median** & **Std. Dev** \\ \hline  & BIRCH3 dataset.