Paper deep dive
Linear Proposal Operators and Stochastic Search Geometry in SOMA and Differential Evolution
Vojtěch Novák, Ivan Zelinka
Intelligence
Status: succeeded | Model: Gemma-4-26B-A4B | Prompt: intel-v1 | Confidence: 94%
Last extracted: 8/3/2026, 2:47:43 AM
Summary
This paper introduces an operator-selection factorization to analyze the stochastic search geometry of the Self-Organizing Migrating Algorithm (SOMA) and Differential Evolution (DE). It derives closed-form expressions for proposal moments, covariance, and active dimensionality under Bernoulli perturbation masks. The study demonstrates that canonical SOMA proposals are affine/linear in an augmented state and proposes geometry-controlled and rotation-aware SOMA variants that outperform canonical SOMA on the BBOB benchmark.
Entities (11)
Relation Signals (8)
Vojtěch Novák → authored → Linear Proposal Operators and Stochastic Search Geometry in SOMA and Differential Evolution
confidence 99% · Linear Proposal Operators and Stochastic Search Geometry in SOMA and Differential Evolution Vojtěch Novák
Ivan Zelinka → affiliatedwith → VSB - Technical University of Ostrava
confidence 95% · Ivan Zelinka Department of Computer Science, Faculty of Electrical Engineering and Computer Science, VSB - Technical University of Ostrava
Vojtěch Novák → affiliatedwith → VSB - Technical University of Ostrava
confidence 95% · Vojtěch Novák vojtech.novak.st1@vsb.cz ... VSB - Technical University of Ostrava
Self-Organizing Migrating Algorithm → analyzedwith → operator-selection factorization
confidence 95% · This paper introduces an operator--selection factorization that separates objective-independent variation from boundary repair and fitness-dependent selection, and uses it to study the proposal geometry of the Self-Organizing Migrating Algorithm (SOMA)
Differential Evolution → analyzedwith → operator-selection factorization
confidence 95% · uses it to study the proposal geometry of the Self-Organizing Migrating Algorithm (SOMA) and Differential Evolution (DE)
Self-Organizing Migrating Algorithm → hasproperty → affine proposal
confidence 92% · The canonical SOMA proposal is shown to be affine in the search space and exactly linear in an augmented migrant--leader state.
Geometry-Controlled SOMA → →
Cypher Suggestions (0)
No Cypher suggestions yet.
Abstract
Abstract:Swarm and evolutionary algorithms are usually analyzed as complete procedural systems in which nonlinear selection, replacement, and adaptation obscure simpler structure within candidate generation. This paper introduces an operator--selection factorization that separates objective-independent variation from boundary repair and fitness-dependent selection, and uses it to study the proposal geometry of the Self-Organizing Migrating Algorithm (SOMA) and Differential Evolution (DE). The canonical SOMA proposal is shown to be affine in the search space and exactly linear in an augmented migrant--leader state. In leader-relative coordinates, the resulting operator provides a direct interpretation of interpolation, projection, overshooting, and coordinate masking. Under Bernoulli perturbation masks, we derive closed-form expressions for the proposal mean, covariance, expected squared step length, expected squared distance from the leader, active dimensionality, and coordinate coverage. For canonical DE/rand/1/bin, we derive the finite-population moments of differential mutation and characterize the additional covariance and coordinate dependence induced by forced-coordinate binomial crossover. Exact enumeration and Monte Carlo experiments verify the analytical identities and quantify the effects of mask conditioning, boundary repair, and fitness-based selection. The analysis further motivates geometry-controlled and rotation-aware SOMA variants, together with an adaptive population-reducing extension of iSOMA. Experiments on the complete noiseless BBOB benchmark show that these operator-guided variants substantially improve upon canonical SOMA and are competitive with established DE methods in several dimension--budget regimes. The results demonstrate how proposal-level operator analysis can support both the interpretation and design of population-based optimizers.
Tags
Links
- Source: https://arxiv.org/abs/2607.29228v1
- Canonical: https://arxiv.org/abs/2607.29228v1
Trouble viewing inline? Open PDF directly →
Full Text
100,512 characters extracted from source content.
Expand or collapse full text
Linear Proposal Operators and Stochastic Search Geometry in SOMA and Differential Evolution Vojtěch Novák vojtech.novak.st1@vsb.cz Ivan Zelinka Department of Computer Science, Faculty of Electrical Engineering and Computer Science, VSB - Technical University of Ostrava, Ostrava, Czech Republic IT4Innovations National Supercomputing Center, VSB - Technical University of Ostrava, 708 00 Ostrava, Czech Republic Department of Informatics and Statistics, Marine Research Institute, Klaipeda University, Lithuania Abstract Swarm and evolutionary algorithms are usually analyzed as complete procedural systems in which nonlinear selection, replacement, and adaptation obscure simpler structure within candidate generation. This paper introduces an operator–selection factorization that separates objective-independent variation from boundary repair and fitness-dependent selection, and uses it to study the proposal geometry of the Self-Organizing Migrating Algorithm (SOMA) and Differential Evolution (DE). The canonical SOMA proposal is shown to be affine in the search space and exactly linear in an augmented migrant–leader state. In leader-relative coordinates, the resulting operator provides a direct interpretation of interpolation, projection, overshooting, and coordinate masking. Under Bernoulli perturbation masks, we derive closed-form expressions for the proposal mean, covariance, expected squared step length, expected squared distance from the leader, active dimensionality, and coordinate coverage. For canonical DE/rand/1/bin, we derive the finite-population moments of differential mutation and characterize the additional covariance and coordinate dependence induced by forced-coordinate binomial crossover. Exact enumeration and Monte Carlo experiments verify the analytical identities and quantify the effects of mask conditioning, boundary repair, and fitness-based selection. The analysis further motivates geometry-controlled and rotation-aware SOMA variants, together with an adaptive population-reducing extension of iSOMA. Experiments on the complete noiseless BBOB benchmark show that these operator-guided variants substantially improve upon canonical SOMA and are competitive with established DE methods in several dimension–budget regimes. The results demonstrate how proposal-level operator analysis can support both the interpretation and design of population-based optimizers. keywords: Evolutionary computation , Swarm intelligence , Linear operators , Self-Organizing Migrating Algorithm , Differential Evolution 1 Introduction Population-based metaheuristics constitute an important class of methods for black-box, non-convex, and derivative-free optimization. They maintain a set of candidate solutions and repeatedly generate, evaluate, and retain new search points without requiring analytical derivatives of the objective function. Differential Evolution (DE) [1] and the Self-Organizing Migrating Algorithm (SOMA) [2, 3] are representative continuous-domain methods. DE constructs trial vectors from scaled differences between population members, whereas SOMA samples candidate points along masked paths directed toward selected targets. Despite their different terminology and population dynamics, both methods combine an objective-independent geometric proposal mechanism with objective-dependent decisions concerning guidance, survival, and replacement. The distinction between variation and selection is fundamental in evolutionary computation [4, 5]. In practical algorithms, however, these roles are usually embedded in complete procedural iterations that additionally contain random index sampling, mask generation, objective evaluation, ranking, memory updates, boundary handling, and conditional replacement. When the complete iteration is treated as a single mathematical map, discontinuous operations such as ranking, argmin , and greedy accept–reject decisions can dominate its apparent structure and obscure simpler algebraic regularities within candidate generation. Existing theory has examined evolutionary algorithms through runtime and drift analysis, Markov chains, stochastic processes, and population-dynamics models [6]. DE has been surveyed extensively from both algorithmic and theoretical perspectives [7, 8, 9], while longitudinal benchmark evidence also highlights the increasing importance of rotated nonseparable functions and adaptive, coordinate-robust search mechanisms [10]. More focused probabilistic studies have derived expectation vectors and covariance matrices for differential mutation strategies [11], examined the contour-fitting behavior induced by population differences [12], and analyzed the influence of crossover on DE dynamics [13]. These results establish that meaningful distributional structure can be isolated within individual DE operators. The present work complements them by deriving an explicit augmented-state linear representation, retaining exact finite-population sampling without replacement, and propagating the resulting moments through forced-coordinate binomial crossover while keeping repair and selection as separate transformations. For SOMA, empirical studies have clarified how the perturbation parameter PRT affects movement, coordinate activation, parameter-space coverage, diversity, and convergence [14], while broader accounts describe the historical development and strategy variants of the algorithm [2, 3]. Nevertheless, an explicit proposal-level comparison of the exact algebraic and stochastic geometry of canonical SOMA and DE remains limited. In particular, complete population dynamics are difficult to analyze without assumptions on the objective function because leader choice and survivor replacement are fitness-dependent. This paper addresses that gap through an operator–selection factorization. A complete update is decomposed into conceptually distinct maps, current state and random choices→raw proposals→ℛΩfeasible candidates→fnext population state,current state and random choices\; \;V\;\;raw proposals\; \;R_ \;\;feasible candidates\; \;S_f\;\;next population state, where V denotes variation, ℛΩR_ denotes an optional boundary-repair transformation, and fS_f denotes fitness-dependent selection. Conditional on sampled indices, masks, parameters, and auxiliary population states, V may be exactly linear or affine even though fS_f is nonlinear and potentially discontinuous. Boundary repair is retained as a separate transformation because clipping, reflection, resampling, and related mechanisms can alter the mean, covariance, support, and dependence structure of the raw proposal distribution [15]. For SOMA and DE/rand/1/bin, the factorization exposes conditional linear or affine proposal rules whose stochastic moments can be analyzed before boundary repair and fitness-dependent selection. This provides a common basis for comparing their masking, covariance, and parameter-controlled search geometry. The novelty claimed here is neither the general distinction between variation and selection nor the first use of matrices in evolutionary optimization. Linear algebra has previously been used to design new DE transformations, such as Antisymmetric and Reversible Differential Evolution [16]. In contrast, the present work uses affine and linear operators to analyze unchanged canonical proposal rules and then extends the deterministic representation to their stochastic proposal geometry. A more detailed positioning relative to the closest DE and SOMA studies is given in Section 2.3 and Table 2. The main contributions are as follows: 1. We formulate canonical SOMA and DE/rand/1/bin candidate generation as conditional affine or augmented linear proposal operators, explicitly separated from boundary repair and fitness-dependent selection. 2. For SOMA, we derive the target-relative operator and closed-form proposal moments under independent Bernoulli masks and under masks conditioned to contain at least one active coordinate. We additionally quantify expected step lengths, expected squared target distances, active dimensionality, and coordinate coverage. 3. For DE, we derive the finite-population mean and covariance of the DE/rand/1/bin mutant under ordered donor sampling without replacement and combine these results with the exact moments and cross-coordinate dependence induced by forced-coordinate binomial crossover. 4. We characterize how PRT, the SOMA path parameter t, the differential weight F, and the crossover rate CRCR control proposal scale, dimensionality, covariance, contraction, and orientation. 5. We use the resulting geometric quantities to construct Geometry-Controlled SOMA, Rotation-Aware SOMA, and an experimental iSOMA extension combining rotation-aware masking, success-history adaptation, and linear population-size reduction. Their performance is evaluated on the complete noiseless BBOB benchmark across three dimensions and two evaluation budgets. 6. Exact enumeration and Monte Carlo experiments verify the analytical identities and quantify the distortions introduced by mask conditioning, boundary repair, and fitness-based selection. Section 2 introduces SOMA and DE and defines the scope of the canonical variants considered in the analysis. The following sections develop the SOMA operator factorization and stochastic proposal geometry and then derive the corresponding finite-population results for DE. The analytical quantities are subsequently used to construct and evaluate operator-guided SOMA variants on the noiseless BBOB benchmark. The final sections discuss the empirical implications, limitations, and directions for further development. The appendices provide detailed numerical verification of the analytical identities and complete function-wise convergence profiles. 2 SOMA and Differential Evolution: Algorithmic Background and Positioning 2.1 Self-Organizing Migrating Algorithm SOMA was introduced as a population-based optimizer inspired by the competitive–cooperative behavior of a group of agents. Its defining mechanism is migration rather than reproduction in the conventional genetic-algorithm sense: individuals evaluate a finite set of candidate positions along paths directed toward selected targets. In the canonical continuous formulation, the principal controls are the population size, the migration Step, the PathLength, and the perturbation probability PRT [2, 3]. The associated binary perturbation vector determines which coordinates participate in a displacement and hence restricts an individual proposal to an active coordinate subspace. Different SOMA strategies primarily modify target construction and update semantics. In canonical All-to-One, every non-leader individual migrates toward the current best individual. All-to-All evaluates migrations toward multiple population members; All-to-Random replaces the deterministic best target by a sampled member; and cluster-based strategies select targets within locally identified subpopulations. Discrete, multi-objective, constrained, hybrid, and parallel variants modify the representation, selection, or implementation layers while retaining the migration principle [2, 3]. A further distinction concerns the temporal use of the perturbation mask. Some formulations reuse one perturbation vector throughout a complete migration path, whereas others regenerate the mask at individual path points. A shared mask couples all path proposals through one random coordinate-selection operator; regenerated masks produce path-point-specific random operators. The proposal moments derived later therefore state explicitly whether masks are independent, reused, or conditioned to contain at least one active coordinate. The present analysis concerns the canonical continuous, leader-directed, static-origin proposal rule. Every path proposal is generated from the same migration-origin point, and the best evaluated path point is selected only after the path has been sampled. Immediate-update variants instead replace the migration origin as soon as an improvement is found and therefore generate a sequential, state-dependent composition of operators. Modern variants can add group-based target selection, reverse-ordered path evaluation, progressive search-space narrowing, and stagnation-triggered replacement, as illustrated by iSOMA [17]. Such mechanisms remain compatible with the operator–selection viewpoint, but their exact representation requires an augmented state containing the variables that determine the current target, bounds, history, and update condition. Table 1 summarizes representative established SOMA strategies in the terminology used by the operator framework. The three operator-guided variants proposed in this paper are better introduced later, together with their definitions and experimental motivation, rather than in the background taxonomy. Table 1: Representative established SOMA strategies interpreted through the proposal-operator framework. The table is a structural taxonomy rather than an exhaustive historical catalogue. Strategy or variant Target construction Perturbation or search geometry Update rule Operator interpretation Canonical All-to-One Global best individual Binary diagonal PRT mask and discretized path Best point selected after evaluating the path Family of affine proposals sharing one migrant–leader state. All-to-All Every other population member Binary diagonal PRT mask Pairwise path evaluation Collection of pair-conditioned affine operators with a changing target. All-to-Random Randomly sampled population member Binary diagonal PRT mask Pathwise selection Random-anchor affine operator; target randomness contributes an additional covariance component. Cluster-based SOMA Leader selected within an inferred cluster Usually coordinate-masked, locally anchored movement Cluster-conditioned migration and selection Mixture of locally conditioned proposal operators indexed by cluster membership. iSOMA [17] Leaders and migrants selected from sampled groups Reverse-ordered jumps and progressively narrowed search space Immediate improvement update and stagnation-triggered replacement State-, bound-, and history-dependent composition requiring an augmented algorithmic state. 2.2 Differential Evolution DE is a population-based evolutionary optimizer for continuous parameter spaces. In the canonical DE/rand/1/bin strategy, each population member acts as a target vector. Three distinct donor indices are sampled, a mutant is constructed by adding a scaled population difference to a randomly selected base vector, and binomial crossover combines the mutant with the target. The standard crossover implementation forces at least one coordinate to be inherited from the mutant. Finally, greedy one-to-one selection retains the trial vector when it is no worse than its target [1, 7]. The strategy notation DE/x/y/z records the rule used to choose the base vector, the number of difference vectors, and the crossover type. Thus, DE/rand/1/bin uses a random base, one scaled difference, and binomial crossover. The behavior of canonical DE is governed mainly by the population size, the differential weight F, the crossover rate CRCR, and the mutation strategy. The differential weight scales population-derived directions, while crossover determines how much of the mutant is transferred to the trial. The state-of-the-art and updated surveys of DE document the subsequent development of alternative base-vector rules, multiple-difference strategies, adaptive and self-adaptive parameter control [18, 19, 20], ensembles of mutation operators, population size adaptation, hybridization, and extensions to constrained, multi-objective, large-scale, and uncertain optimization [7, 8]. Although these variants differ substantially at the algorithmic level, many retain the same modular structure: donor construction, recombination, objective evaluation, and target–trial selection. An influential adaptive lineage begins with SHADE, which stores successful values of the differential weight and crossover rate in historical memories [21]. L-SHADE augments this mechanism with linear population-size reduction [22], while iL-SHADE introduces further modifications to improve its single-objective real-parameter performance [23]. These methods provide the adaptation background for the DE reference algorithm and for the success-history and population-reduction components used later in the experimental iL-SHOMA-RA method. Theoretical analyses of DE have addressed convergence, invariance, differential mutation, crossover, population diversity, and population dynamics [9]. Proposal generation is comparatively amenable to algebraic and probabilistic analysis because, after donor indices and crossover choices are fixed, the trial is a linear combination of population members. Greedy replacement is harder to characterize in a problem-independent manner because the acceptance event depends on the objective function. This separation motivates the present focus on exact finite-population mutant moments and the additional covariance and dependence introduced by binomial crossover, while treating selection and boundary repair as subsequent transformations. 2.3 Common proposal–selection structure and relation to prior work SOMA and DE differ in how their proposal endpoints are constructed. Canonical SOMA uses a selected target and samples a discretized, coordinate-masked path from the migrant toward that target. Canonical DE/rand/1/bin constructs a random mutant from three sampled population members and then uses a coordinate mask to combine it with the target vector. Conditional on the required state and random choices, both proposal mechanisms are affine or linear; their complete updates remain nonlinear because target selection, boundary repair, and survivor selection depend on fitness or feasibility. The scope of the exact derivations is therefore deliberately narrower than a complete population-dynamics theory. For SOMA, the formulas apply to the raw continuous, leader-directed, static-origin proposal under the stated mask model. For DE, they apply to canonical DE/rand/1/bin under ordered donor sampling without replacement and standard forced-coordinate binomial crossover. Pairwise, random-target, clustered, adaptive, or immediate-update variants can be represented conditionally by augmenting the state with the variables that determine their target, bounds, history, or update event, but those variants do not share identical unconditional dynamics. Table 2 positions the present analysis relative to the most closely related studies. Table 2: Positioning of the present analysis relative to closely related work. Study Primary objective and methodology Relation to the present work DE surveys [7, 8, 9] Review DE strategies and theoretical results concerning mutation, crossover, parameter control, convergence, invariance, diversity, and population dynamics. Establish the broader DE context. The present work develops an explicit augmented-state representation and exact finite-population proposal moments for canonical DE/rand/1/bin, while treating repair and selection as separate transformations. Guo and Yang [24]; Caraffini and Neri [25] Develop covariance-eigenvector and rotation-invariant crossover mechanisms for DE and study their behavior under rotated problems. Establish the closest representation-oriented precedent. The present work applies a population-derived basis to the perturbation mask of the unchanged leader-directed SOMA proposal and derives its operator interpretation. Opara and Arabas [11, 12] Derive distributional properties of differential mutation strategies and formalize the contour-fitting behavior of DE. Provide the closest probabilistic DE analysis. The present work additionally uses an augmented linear representation, exact finite-population donor sampling, and forced-coordinate crossover moments, with repair and selection treated separately. Tomczak et al. [16] Design ADE and RevDE candidate-generation rules by applying invertible linear transformations to triplets of population members and analyze reversibility and eigenvalue behavior. Use linear algebra for algorithm synthesis. The present work instead uses linear algebra to analyze unchanged SOMA and DE rules and extends the deterministic representation to stochastic proposal geometry. Pluháček et al. [14] Empirically study the influence of PRT on SOMA trajectories, parameter-space coverage, population diversity, and convergence. Provide the main empirical SOMA baseline. The present work derives proposal-level moments and geometric quantities that quantify mask activation under explicitly stated probabilistic models. Davendra and Zelinka; Skanderová [2, 3] Describe the historical development, migration strategies, variants, applications, and implementation choices of SOMA. Establish the SOMA-family context. The present analysis is restricted to the canonical continuous proposal and develops exact operator-level stochastic geometry rather than a general algorithmic survey. Diep et al. [17] Introduce iSOMA with group-based organization, reverse-ordered jumps, immediate updating, narrowed search bounds, stagnation handling, and a drone path-planning application. Represent a modern state-dependent SOMA variant. The present framework identifies the additional state variables required to express such a variant as a conditional composition of proposal and selection operators. Present work Factorizes canonical SOMA and DE/rand/1/bin into raw proposal, repair, and fitness-dependent selection stages and derives exact operator forms and stochastic moments. Combines conditional linear representation, finite-population sampling, random-mask analysis, and explicit examination of distortions introduced by repair and selection. 3 Operator–Selection Factorization of SOMA Consider the minimization of an objective function f:Ω⊆ℝD→ℝf: ^D . At migration cycle m, SOMA maintains a population (m)=1(m),…,N(m)X^(m)=\ x^(m)_1,…, x^(m)_N\ and selects a leader according to (m)∈argmin∈(m)f(). L^(m)∈ *arg\,min x∈ X^(m)f( x). (1) The leader-selection mechanism may vary across SOMA strategies, but it is an objective-dependent operation. To analyze the variation stage, fix one migration cycle and one migrant and write the migration origin as x and the selected target as L. Let t1,…,tJt_1,…,t_J denote the path parameters determined by Step and PathLength, and let Πj=diag(πj,1,…,πj,D) _j=diag( _j,1,…, _j,D) be the perturbation mask at path point j. The raw proposal is j=+tjj(−). z_j= x+t_j _j( L- x). (2) The algorithmic meaning of the path, mask, and static-origin assumption was specified in Section 2.1. The purpose of the present section is to separate this proposal equation from repair and pathwise survivor selection and to derive its exact affine and augmented linear representations. Figure 1 visualizes the deterministic fully-active case =I =I, in which each SOMA proposal path reduces to a discretized line segment from the migrant toward the current leader. This figure is intended to illustrate the geometric structure of the proposal rule. Stochastic coordinate masking and its induced proposal distribution are examined separately in Figure 3. Figure 1: Discrete SOMA proposal paths and selected replacements on the Three-Hump Camel function for three successive migration cycles. For each non-leader individual, candidate points are generated along a leader-directed linear path using a fully active perturbation mask =I =I, so that the path points satisfy (tj)=+tj(−) z(t_j)= x+t_j( L- x). The colored circular markers show the discretized proposal points associated with the path parameters tjt_j, and the dashed segments indicate the corresponding leader-directed proposal lines. The perturbation mask determines which coordinates participate in the movement and serves the purpose of exploration. If πj,d=1 _j,d=1, coordinate d is displaced along the migrant–leader direction. If πj,d=0 _j,d=0, that coordinate remains equal to its value in the migration-origin point. The path parameter controls the geometry of the active coordinates. Values 0<tj<10<t_j<1 interpolate between the migrant and the leader, tj=1t_j=1 copies the active leader coordinates exactly, and tj>1t_j>1 generates an overshooting proposal beyond the leader. This overshooting mechanism allows SOMA to investigate regions that are not contained in the line segment joining the two individuals. Although the proposals are indexed by j, they should not generally be interpreted as a recursively generated trajectory. In the standard SOMA path construction, every proposal j z_j is obtained from the same migration-origin point x. In particular, j+1 z_j+1 is not normally generated by applying the next path transformation to j z_j. The term linear proposal path therefore refers to a discrete collection of geometrically ordered proposals, each produced by an operator acting on the same migrant–leader state. This static-origin interpretation applies to the canonical pathwise-update formulation analyzed in this paper. Adaptive or immediate-update variants may replace the migration origin as soon as an improving point is found, in which case the migration becomes a genuinely sequential, state-dependent composition of proposal operators rather than a set of proposals sharing one fixed origin. After the proposals have been evaluated, the next position of the migrant is selected from =,1,…,J,+∈argmin∈f().C= x, z_1,…, z_J, x^+∈ *arg\,min y f( y). (3) Including the original migrant in C ensures that the selected position is not worse than the migration origin under exact objective evaluation. Equations (2) and (3) describe two mathematically different operations. The first is a parameterized geometric proposal rule. The second is an objective-dependent decision rule involving fitness evaluation and an argmin operation. Consequently, the proposal rule may possess an exact linear representation even though the complete migration remains nonlinear. 3.1 Exact linear representation of the proposal rule For a fixed path parameter t and mask , define the augmented migrant–leader state =[]∈ℝ2D. s= bmatrix x\\ L bmatrix ^2D. The proposal rule can then be written as []=[I−tt0I]⏟ℳ(t,)[]. bmatrix z\\ L bmatrix= bmatrixI-t &t \\ 0&I bmatrix_ (t, ) bmatrix x\\ L bmatrix. (4) In the original migrant coordinates, Eq. (2) is affine because the leader acts as a translation anchor. After augmenting the state with the leader position, the same proposal becomes an exact linear matrix–vector product. The coefficients of ℳ(t,)M(t, ) depend only on the path parameter and perturbation mask. They do not depend on the numerical values of x or L, and they contain no objective value, gradient, rank, or fitness-conditioned switch. Once the migrant, leader, path parameter, and mask have been specified, the proposal is therefore produced by an objective-independent linear operator. This statement does not imply that the proposal itself is independent of previous objective information. The leader appearing in the input state was selected using fitness values, and the migrant is the result of previous selection decisions. The exact decoupling concerns the coefficients of the proposal operator: objective information enters through the current state, but not through the matrix that maps this state to a path proposal. At the population level, one SOMA migration follows the standard variation–selection logic of evolutionary computation. Fitness-based selection first determines the leader that guides the search. The resulting migrant–leader pair then generates a path of candidate solutions through an exact linear variation operator. Finally, the generated candidates are evaluated and selection determines which point survives into the next population state. Figure 2 summarizes this decomposition. Current migrant-leader state (,)( x, L) from population at migration X(m)X^(m) Linear proposal path ℳ(tj,j)M(t_j, _j) Evaluate candidates and retain best Updated migrant / next population state Variation exact linear generation of path candidates Selection survival / replacement under fitness (m)=,1,…,JC^(m)=\x,z_1,…,z_J\ Figure 2: Variation–selection decomposition of one SOMA migration. The current population first undergoes selection to identify the leader that defines the search direction. Conditional on the migrant, leader, path parameter, and perturbation mask, the candidate points along the path are then generated by an exact objective-independent linear variation operator. The resulting candidate set is finally subjected to selection, where objective-function evaluations determine which candidate survives into the next population state. This separates the exactly representable proposal mechanism from the objective-dependent selection pressure exerted before and after variation. The decomposition is also consistent with the variation–selection organization underlying evolutionary computation [4, 6]. In the biological interpretation associated with the Modern Evolutionary Synthesis, mechanisms that generate heritable variation are conceptually distinguished from differential survival and reproduction under selection [26]. Evolutionary algorithms retain this abstract organization without reproducing the details of biological evolution: variation constructs candidate solutions, whereas fitness-based selection changes their representation in the subsequent population. Within this interpretation, the SOMA proposal path is a structured leader-guided variation mechanism, and the comparison of its path points exerts selection pressure. The contribution of the present analysis is to show that the internal variation mechanism can possess an exact linear operator form even though the surrounding selection process is nonlinear. The biological analogy concerns this separation of roles; it does not imply that biological variation itself is linear or that SOMA constitutes a literal model of biological evolution. 3.2 Extension to random, pairwise, and state-dependent targets The canonical expression can be embedded in a broader class of state-conditioned anchored affine proposals, =+θ(,,ω)(θ(,,ω)−), z= x+ A_θ( X, H,ω) ( L_θ( X, H,ω)- x ), (5) where X is the current population, H denotes optional algorithmic history, ω collects random auxiliary choices, θ L_θ is a target-construction rule, and θ A_θ is the proposal operator. Canonical SOMA is recovered by θ=t A_θ=t and by taking θ L_θ as the current leader. All-to-All uses the same conditional operator with a pair-specific target; All-to-Random makes the target itself random; cluster-based methods condition it on cluster membership; and adaptive variants allow both terms to depend on H. For a random target L, the conditional moments derived below remain valid given L. The unconditional proposal covariance then satisfies Cov()=[Cov(∣)]+Cov([∣]),Cov( z)=E_ L [Cov( z L) ]+Cov_ L (E[ z L] ), (6) which separates mask-induced uncertainty from uncertainty caused by target selection. Equation (5) is not used to claim that all SOMA variants share identical dynamics. Rather, it identifies the minimum additional conditioning needed to transfer the present proposal-level analysis to pairwise, clustered, adaptive, or bounded-search variants. 4 Geometry of the SOMA Operator The SOMA proposal becomes especially transparent when it is written relative to the current leader. Define the leader-relative displacement =−. e= x- L. Since the proposal rule is =+t(−), z= x+t ( L- x), the corresponding proposal displacement relative to the leader is +=−=(I−t). e^+= z- L=(I-t ) e. (7) Equation (7) is the key geometric form of the SOMA operator. It shows that the proposal mechanism does not require a more complicated nonlinear interpretation: relative to the leader, SOMA simply rescales the current displacement coordinate-wise. Inactive coordinates remain unchanged, while active coordinates are multiplied by the factor 1−t1-t. For an active coordinate, Eq. (7) reduces to ed+=(1−t)ed.e_d^+=(1-t)e_d. The path parameter t therefore has an immediate geometric meaning. It determines whether the proposal interpolates toward the leader, projects onto the leader coordinate, or overshoots beyond it. The set ℰ0=(,)∈ℝ2D:=E_0=\( x, L) ^2D: x= L\ forms a consensus manifold of the proposal rule: if the migrant already coincides with the leader, then the proposal leaves this state unchanged. This is a geometric property of the proposal operator itself and should not be interpreted as a convergence theorem for the complete SOMA optimizer. Table 3: Coordinate-wise geometry of an active SOMA proposal. The quantities are expressed relative to the current leader. Path parameter Multiplier 1−t1-t Geometric effect t=0t=0 11 no movement 0<t<10<t<1 (0,1)(0,1) interpolation toward the leader t=1t=1 0 projection onto the leader coordinate 1<t<21<t<2 (−1,0)(-1,0) overshooting with reduced distance t=2t=2 −1-1 reflection with unchanged distance t>2t>2 less than −1-1 overshooting with increased distance Table 3 shows that overshooting is not necessarily expansive. For 1<t<21<t<2, the proposal crosses the leader but remains closer to it than the migration origin. Only values t>2t>2 increase the magnitude of an active leader-relative coordinate. If K=tr()K=tr( ) denotes the number of active coordinates, then the proposal rescales K leader-relative directions by the common factor 1−t1-t. Consequently, the volume change associated with one proposal operator is detℳ(t,)=(1−t)K. (t, )=(1-t)^K. (8) For 0<t<20<t<2, every active relative coordinate contracts in magnitude. At t=1t=1, the operator becomes singular whenever at least one coordinate is active, because those active directions are mapped exactly onto the leader. For t>2t>2, the active relative directions expand. The determinant in Eq. (8) describes the deformation induced by a single proposal operator. It should not be interpreted as repeated contraction of the standard SOMA trajectory, because the path points are generated from a common migration origin rather than by iterative composition of the same state update. 4.1 Bernoulli perturbation masks The perturbation mask can be interpreted as a random coordinate-selection operator. Assume that its diagonal entries are independent Bernoulli variables, πd∼Bernoulli(p), _d (p), where p∈[0,1]p∈[0,1] corresponds to the SOMA perturbation parameter PRT under the independent-mask model [2]. For the migrant–leader displacement =− r= L- x, the random proposal is =+t. z= x+t r. Its conditional first and second moments are [∣,,t] [ z x, L,t] =+pt, = x+pt r, (9) Cov(∣,,t) ( z x, L,t) =p(1−p)t2diag(r12,…,rD2). =p(1-p)t^2diag (r_1^2,…,r_D^2 ). (10) The expected proposal therefore lies on the migrant–leader line, with ptpt acting as an effective path parameter. The covariance is diagonal because the mask coordinates are independent, and its variance in coordinate d is proportional to rd2r_d^2. Hence the proposal uncertainty is largest in coordinates where the migrant and leader differ most, while coordinates in which they agree have zero variance. This induced covariance should not be confused with the covariance maintained by CMA-ES [27]. CMA-ES explicitly learns and updates a generally dense sampling covariance from selected search steps. SOMA stores no covariance matrix; Eq. (10) is only a statistical consequence of its existing masked proposal rule. A compact measure of leader-directed movement is the normalized expected squared distance [‖−‖2∣,,t]‖−‖2=1−pt(2−t). E[\| z- L\|^2 x, L,t]\| x- L\|^2=1-pt(2-t). (11) For p>0p>0, this quantity is smaller than one when 0<t<20<t<2, equals one at t=2t=2, and exceeds one for t>2t>2. Thus, overshooting with 1<t<21<t<2 still reduces the expected distance from the leader. The expected squared step length from the migrant is, similarly, pt2‖2pt^2\| r\|^2. Figure 3: Bernoulli-mask SOMA proposal geometry for fixed migrant x, leader L, and path parameter t. The gray markers show the exact mask-generated support, while the colored markers and ellipsoids show the proposal means and covariance summaries for three values of p=PRTp= PRT. Changing p redistributes probability over the same support, translates the mean along the migrant–leader direction, and scales the covariance by p(1−p) p(1-p). The inset compares the corner probability masses for p=0.2p=0.2 and p=0.8p=0.8. The ellipsoids are second-moment summaries rather than support boundaries or Gaussian confidence regions. The covariance scale is governed by p(1−p)p(1-p), which is maximal at p=1/2p=1/2. Increasing PRT therefore does not monotonically increase mask-induced variability: small p activates few coordinates, intermediate p produces the largest variance, and p=1p=1 eliminates mask randomness. The number of active coordinates satisfies K=tr()∼Binomial(D,p),[K]=Dp,Pr(K=0)=(1−p)D.K=tr( ) (D,p), [K]=Dp, (K=0)=(1-p)^D. (12) Thus, PRT also controls the expected dimensionality of an individual proposal. Implementations that reject all-zero masks instead use the corresponding distribution conditioned on K>0K>0. If independent masks are regenerated at J path points, coordinate coverage is summarized by Pr(one coordinate is covered) (one coordinate is covered) =1−(1−p)J, =1-(1-p)^J, (13) [uncovered coordinates] [uncovered coordinates] =D(1−p)J, =D(1-p)^J, Pr(all coordinates are covered) (all coordinates are covered) =[1−(1−p)J]D. = [1-(1-p)^J ]^D. These expressions do not apply when a single mask is reused throughout the path or when coordinate activations are correlated. Figure 3 shows that the support of the three-dimensional proposal distribution consists of the eight corners of the axis-aligned box between the migrant and the fully active proposal. For 0<p<10<p<1, changing PRT leaves these locations unchanged but alters their probability masses, the mean, and the covariance. The cases p and 1−p1-p produce ellipsoids of equal size because they share the same factor p(1−p)p(1-p), although their centers differ. The ellipsoids remain aligned with the Cartesian axes because the independent Bernoulli mask introduces no cross-coordinate covariance. 5 Linear Operator Structure of Differential Evolution Using the DE/rand/1/bin notation introduced in Section 2.2, let i x_i be the target vector and let r1,r2,r3r_1,r_2,r_3 be distinct sampled donor indices, also distinct from i. The mutant and trial vectors are =r1+F(r2−r3),i=(I−C)i+C, v= x_r_1+F( x_r_2- x_r_3), u_i=(I-C) x_i+C v, (14) where C is the diagonal binomial-crossover mask. The background section introduced the procedural roles of mutation, crossover, and greedy replacement; here the same operations are represented as conditional linear maps and their finite-population proposal moments are derived. Conditional on the sampled parent indices, crossover mask, and value of F, Eq. (14) is an exact objective-independent linear transformation of the stacked population state. For example, [⋅i⋅r1⋅r2⋅r3]=[I−F−F0I0000I0000I]⏟⋅ℳ⋅DE[⋅i⋅r1⋅r2⋅r3]. bmatrix u· i\ x·r_1\\ x·r_2\\ x·r_3 bmatrix= bmatrixI- C& C&F C&-F C\\ 0&I&0&0\\ 0&0&I&0\\ 0&0&0&I bmatrix· ·DE bmatrix x· i\\ x·r_1\\ x·r_2\\ x·r_3 bmatrix. (15) The lower block rows merely preserve the sampled donor vectors, in the same way that the SOMA operator preserves the leader. The similarity between DE and SOMA becomes even more explicit in leader-directed DE strategies. For example, the DE/current-to-best/1/bin proposal can be written as ⋅i=⋅i+[F(⋅best−⋅i)+F(⋅r1−⋅r2)]. u· i= x· i+ C [F( x·best- x· i)+F( x·r_1- x·r_2) ]. (16) The first difference term in Eq. (16) has exactly the same leader-attraction structure as the SOMA proposal t(−)t ( L- x), with the correspondence t=Ft=F, = = C, and =best L= x_best. DE adds a second population-difference term that introduces a direction derived from the current population distribution. SOMA and leader-directed DE can therefore be viewed as members of a broader class of masked affine difference operators, with DE combining directed attraction and differential exploration. The complete DE update is nevertheless nonlinear because the trial vector is accepted through greedy fitness-based replacement: i+=i,f(i)≤f(i),i,otherwise. x_i^+= cases u_i,&f( u_i)≤ f( x_i),\\ x_i,&otherwise. cases (17) Thus, DE possesses an exact linear variation operator conditional on its sampled indices and crossover mask, followed by a nonlinear selection operator. 5.1 Stochastic proposal geometry of Differential Evolution The similarity between SOMA and Differential Evolution extends beyond their deterministic operator forms. Consider the DE/rand/1/bin trial vector i=i+(−i),=r1+F(r2−r3), u_i= x_i+ C( v- x_i), v= x_r_1+F( x_r_2- x_r_3), (18) where C is the diagonal binomial-crossover mask. If its coordinates are initially modeled as independent Bernoulli variables with activation probability q=CRq=CR, then, conditional on the target and mutant vectors, [i∣] [ u_i v] =i+q(−i), = x_i+q( v- x_i), (19) Cov(i∣) ( u_i v) =q(1−q)diag((−i)2). =q(1-q)diag (( v- x_i)^2 ). (20) Conditional on v, the DE crossover geometry is therefore equivalent to the Bernoulli-mask SOMA geometry with the mutant acting as the fully active endpoint. The exact support contains at most 2D2^D corners of the axis-aligned box joining i x_i and v. Changing CRCR redistributes probability over this support and translates the conditional mean, but does not alter the support itself for 0<CR<10<CR<1. Figure 4: Influence of the differential weight F and crossover rate CRCR on the donor-marginalized proposal geometry of DE/rand/1/bin. Both panels use the same population and target vector i x_i. Gray markers denote the current population, blue points represent possible mutant vectors =r1+F(r2−r3) v= x_r_1+F( x_r_2- x_r_3), and red points represent the corresponding trial vectors after standard forced-coordinate binomial crossover. The blue and red ellipsoids summarize the mutant and trial covariance, with centers v=[] μ_v=E[ v] and u=[i] μ_u=E[ u_i], respectively. (a) Smaller F reduces differential displacement, while smaller CRCR retains more target coordinates and concentrates the trial distribution near i x_i. (b) Larger F amplifies population differences, while larger CRCR transfers more mutant coordinates into the trial, producing a broader distribution that more closely follows the mutant geometry. The covariance ellipsoids are second-moment summaries of the finite proposal distributions rather than Gaussian confidence regions. The principal difference is that the DE endpoint v is itself random because it is constructed from sampled population members. Let v=[],v=Cov() m_v=E[ v], _v=Cov( v) denote its mean and covariance under donor sampling. Applying the law of total covariance gives Cov(i)=q2v+q(1−q)diag[diag(v)+(v−i)2].Cov( u_i)=q^2 _v+q(1-q)diag [diag( _v)+( m_v- x_i)^2 ]. (21) The first term represents uncertainty induced by mutation and donor selection, whereas the second is the coordinate-aligned covariance introduced by binomial crossover. Unlike the conditional crossover covariance, the mutation-induced component may be dense and may rotate the principal axes of the total trial distribution. Standard binomial crossover additionally forces one randomly selected coordinate to be inherited from the mutant. The numerical visualization therefore enumerates this forced coordinate explicitly rather than treating all crossover-mask entries as independent. 6 Operator-Guided SOMA Design and Evaluation 6.1 Algorithms and experimental protocol The preceding analysis provides quantities that can be used directly for algorithm design. In particular, the perturbation probability controls the expected active dimensionality, its interaction with the path parameter determines the expected contraction toward the leader, and the diagonal proposal covariance reveals an explicit dependence on the coordinate system. We use these observations to construct Geometry-Controlled SOMA and Rotation-Aware SOMA. We additionally combine rotation-aware masking with success-history adaptation and population-size reduction in an experimental extension of iSOMA. For a migrant x, leader L, path parameter t, and Bernoulli mask m∈0,1Dm∈\0,1\^D, the standard proposal is z=x+tm⋅(L−x).z=x+t\,m·(L-x). (22) Geometry-Controlled SOMA. Geometry-Controlled SOMA (GC-SOMA) replaces direct specification of the perturbation probability p and path parameter t by an expected active dimension k⋆k and a desired normalized squared leader distance ρ⋆ρ . From the proposal moments derived in Section 4, p=k⋆D,ρ⋆=1−pt(2−t),t=1±1−1−ρ⋆p,0<p≤1,1−p≤ρ∗≤1.p= k D, ρ =1-p\,t(2-t), t=1± 1- 1-ρ p, 0<p≤ 1, 1-p≤ρ^*≤ 1. (23) The plus and minus branches generate overshooting and interpolating proposals, respectively. The implementation varies k⋆k and ρ⋆ρ with the consumed evaluation budget and evaluates several nearby contraction targets during each migration. GC-SOMA therefore retains a discretized SOMA-like path while expressing its parameters through intended geometric effects. Rotation-Aware SOMA. Rotation-Aware SOMA (RA-SOMA) addresses the axis alignment of the standard Bernoulli mask. Let Q contain the eigenvectors of the current population covariance matrix. The displacement is transformed into this basis, masked, and transformed back: z=x+tQ[m⋅Q(L−x)].z=x+tQ [m· Q T(L-x) ]. (24) Equivalently, the diagonal mask is replaced by the correlated operator Qdiag(m)Qdiag(m)Q T. The basis is recomputed during the optimization, so the active proposal directions follow the orientation represented by the current population rather than the fixed coordinate axes. Covariance-derived coordinate systems have previously been used to reduce coordinate dependence in evolutionary search. In DE, eigenvector-based crossover applies recombination in a basis estimated from the population covariance [24], and rotationally invariant DE operators have been studied systematically [25]. More broadly, invariance can materially affect performance on ill-conditioned and nonseparable problems [28]. RA-SOMA transfers this principle specifically to the SOMA perturbation mask: the leader-directed proposal remains unchanged in form, but coordinate activation is performed in a population-derived basis. Experimental iL-SHOMA-RA. The third method is an experimental composite optimizer constructed on top of iSOMA. It retains the rank-based migrant selection, tournament-based leader selection, decreasing iSOMA path, first-improvement acceptance, and stagnation restart. The success-history and population-reduction components are motivated by the SHADE and L-SHADE design principles [21, 22]. Their transfer to a rotation-aware SOMA mechanism nevertheless produces a new composite optimizer rather than a variant or implementation of iL-SHADE. Three mechanisms are added: 1. the perturbation mask is applied in the current population covariance basis as in RA-SOMA; 2. the perturbation probability and a multiplicative path scale are sampled from success-history memories; and 3. the population size is reduced linearly during the run. For memory index r, the sampled parameters have the form p∼(αMp,r+(1−α)piSOMA,σp2),s∼Cauchy(Ms,r,γs),t=stiSOMA.p \! (α M_p,r+(1-α)p_iSOMA, _p^2 ), s (M_s,r, _s), t=s\,t_iSOMA. (25) Only parameter pairs associated with accepted improvements update the memories. The perturbation probability uses an improvement-weighted arithmetic mean, whereas the path scale uses an improvement-weighted Lehmer mean. The population size follows N(q)=maxNmin,round[N0+(Nmin−N0)qB],N(q)= \N_ ,round\! [N_0+(N_ -N_0) qB ] \, (26) where q is the number of consumed evaluations and B is the evaluation budget. We refer to this composite method as iL-SHOMA-RA. Because several mechanisms are introduced simultaneously, it is treated as an exploratory optimizer rather than as evidence for the effect of any single component. Reference methods and benchmark design. The comparison includes canonical All-to-One SOMA, the supplied iSOMA implementation, SciPy differential evolution using the best1bin strategy [29], and the pyade-python implementation of iL-SHADE [23, 30]. GC-SOMA, RA-SOMA, iSOMA, iL-SHOMA-RA, SciPy-DE, and iL-SHADE use an initial population of 4D4D. Canonical SOMA uses its conventional population setting of 10D10D. The matched SOMA variants and SciPy-DE receive a common initial population within each experimental block. For iL-SHADE, the initial population size is overridden to 4D4D, while population initialization remains internal to the PyADE implementation. Following established benchmarking principles [31, 32], the algorithms were evaluated on the complete noiseless BBOB suite f1f_1–f24f_24, comprising separable functions, functions with low or moderate conditioning, highly conditioned unimodal functions, multimodal functions with adequate global structure, and multimodal functions with weak global structure [33, 34, 35, 36]. Dimensions D∈5,10,20D∈\5,10,20\ and budgets B∈500D,2500DB∈\500D,2500D\ were used. Each function was evaluated on five BBOB instances with ten independent optimizer repetitions per instance. Each algorithm therefore completed 24×5×10=120024× 5× 10=1200 runs for every dimension–budget configuration. The primary performance value is the final absolute error Δf=f(xbest)−fopt. f=f(x_best)-f_opt. (27) Repeated runs were first aggregated by their median within each function–instance block. Average ranks were then calculated over the resulting 120 blocks, with lower ranks indicating better performance. The Friedman test assessed overall differences. Pairwise comparisons used the Wilcoxon signed-rank test with Holm correction [37, 38]. Anytime performance was measured by the normalized absolute-error area under the best-so-far curve, Aabs=1B∑q=1BminΔfqmax(Δf1,ε),1.A_abs= 1B _q=1^B \ f_q ( f_1, ),1 \. (28) Lower values indicate that small errors were obtained earlier and maintained over a larger proportion of the evaluation budget. Final error remains the primary endpoint; AabsA_abs, target attainment, and convergence curves describe the behavior during the run. 6.2 Numerical results Overall ranking Table 4 reports the average final-error ranks. Canonical SOMA was ranked last or nearly last in every configuration, whereas all operator-guided and improved SOMA variants were substantially more competitive. The Friedman test rejected the null hypothesis of equal rank distributions in every dimension–budget configuration, with all overall p-values below 10−3310^-33. Table 4: Average ranks over the 120 BBOB function–instance blocks. Lower values are better. The best result in each column is shown in bold. D=5D=5 D=10D=10 D=20D=20 Algorithm 500D500D 2500D2500D 500D500D 2500D2500D 500D500D 2500D2500D Canonical-SOMA 6.833 5.888 6.933 6.588 6.842 6.596 GC-SOMA 4.067 4.112 4.467 4.362 4.537 4.521 RA-SOMA 3.175 2.596 3.892 3.233 3.929 3.996 iSOMA 4.408 4.204 4.088 3.917 3.325 3.533 iL-SHOMA-RA 3.663 3.725 3.246 3.983 2.675 3.667 SciPy-DE 2.683 3.650 2.900 3.117 4.771 3.379 iL-SHADE 3.171 3.825 2.475 2.800 1.921 2.308 RA-SOMA produced the best overall rank at D=5D=5 and 2500D2500D evaluations. At D=10D=10, RA-SOMA remained competitive with the DE algorithms under the larger budget, but iL-SHADE obtained the best rank at both budgets. At D=20D=20, iL-SHADE was the strongest general method, although iL-SHOMA-RA ranked second under the restricted 500D500D budget. iL-SHOMA-RA versus iSOMA Table 5 isolates whether the composite method improves its iSOMA base. Under 500D500D evaluations, iL-SHOMA-RA significantly outperformed iSOMA at all three dimensions. The advantage remained significant at D=5D=5 under 2500D2500D, but not at D=10D=10 or D=20D=20. The combined method therefore improves the limited-budget behavior of iSOMA more consistently than its longer-budget final precision. Table 5: Paired comparison of iL-SHOMA-RA against iSOMA over 120 function–instance blocks. W/T/L gives wins, ties, and losses of iL-SHOMA-RA. The ratio is the median blockwise final-error ratio; values below one favor iL-SHOMA-RA. D Budget W/T/L Error ratio pHolmp_Holm 5 500D500D 75/13/32 0.810 1.35×10−31.35× 10^-3 5 2500D2500D 66/18/36 0.889 3.80×10−43.80× 10^-4 10 500D500D 79/11/30 0.866 4.28×10−44.28× 10^-4 10 2500D2500D 57/16/47 1.000 1.001.00 20 500D500D 77/11/32 0.824 5.55×10−45.55× 10^-4 20 2500D2500D 50/11/59 1.000 0.8580.858 The improvement was concentrated on landscapes for which oriented search directions are expected to be useful. On the highly conditioned unimodal group under 500D500D, iL-SHOMA-RA obtained W/T/L counts of 20/0/520/0/5, 23/0/223/0/2, and 25/0/025/0/0 for D=5,10,20D=5,10,20, respectively. The corresponding corrected p-values were 0.02030.0203, 0.002630.00263, and 1.25×10−61.25× 10^-6. At D=20D=20, its median blockwise error was approximately 9.3%9.3\% of that of iSOMA on this group. Significant advantages were also observed on the structured multimodal group at D=10D=10 and D=20D=20 under 500D500D. No similarly consistent advantage occurred on separable functions or on multimodal functions with weak global structure. Low-dimensional performance of RA-SOMA The best RA-SOMA configuration was not only highly ranked but also statistically better than both DE references. At D=5D=5 and 2500D2500D, RA-SOMA won 75 of 120 blocks against iL-SHADE and 72 against SciPy-DE. Its median blockwise error was approximately 32%32\% of the iL-SHADE error and 61%61\% of the SciPy-DE error. Table 6: Pairwise final-error comparisons for RA-SOMA at D=5D=5 and 2500D2500D evaluations. W/T/L is reported from the perspective of RA-SOMA. Comparator W/T/L Error ratio pHolmp_Holm SciPy-DE 72/21/27 0.605 7.96×10−67.96× 10^-6 iL-SHADE 75/15/30 0.324 6.07×10−36.07× 10^-3 At D=10D=10, no statistically significant difference was detected between RA-SOMA and the DE references under the larger budget, whereas RA-SOMA was inferior to iL-SHADE at D=20D=20. The benefit of the population covariance basis is therefore substantial but dimension-dependent. Anytime behavior Table 7 counts the functions on which an algorithm obtained the lowest median normalized absolute-error AUC. SciPy-DE dominated this measure at D=5D=5, whereas iL-SHOMA-RA became the leading anytime method as dimension increased. It had the best AUC on 17 of 24 functions at D=20D=20 and 500D500D, and on 11 functions at D=20D=20 and 2500D2500D. Table 7: Algorithm with the lowest median normalized absolute-error AUC on the largest number of BBOB functions in each configuration. D Budget Algorithm Function wins 5 500D500D SciPy-DE 17/24 5 2500D2500D SciPy-DE 14/24 10 500D500D iL-SHOMA-RA 9/24 10 2500D2500D SciPy-DE 12/24 20 500D500D iL-SHOMA-RA 17/24 20 2500D2500D iL-SHOMA-RA 11/24 The low AUC of iL-SHOMA-RA was not explained solely by rapid progress to a poor plateau. Table 8 shows selected success rates for the strict target Δf≤10−8 f≤ 10^-8. RA-SOMA had the highest attainment rate in the low-dimensional, high-budget configuration. iL-SHOMA-RA had the highest rate at D=10D=10 and D=20D=20 under the restricted budget, exceeding both iSOMA and iL-SHADE. Table 8: Percentage of runs attaining Δf≤10−8 f≤ 10^-8 in selected configurations. The best result in each column is shown in bold. Algorithm D=5, 2500D=5,\;2500D D=10, 500D=10,\;500D D=20, 500D=20,\;500D Canonical-SOMA 10.0 2.0 0.8 GC-SOMA 20.6 5.4 4.9 RA-SOMA 52.4 5.0 3.4 iSOMA 27.8 9.8 9.2 iL-SHOMA-RA 22.7 13.8 10.1 SciPy-DE 42.1 6.7 3.2 iL-SHADE 21.6 12.2 7.0 Geometry control GC-SOMA consistently improved upon canonical SOMA but did not obtain the best aggregate rank in any configuration. This result separates interpretability from optimization dominance: the analytical formulas provide direct control over expected active dimensionality and contraction, but the particular budget-dependent schedule used here is not necessarily optimal. GC-SOMA should therefore be interpreted as a proof of concept for geometry-based parameterization rather than as a fully tuned state-of-the-art optimizer. Complete function-wise convergence profiles are provided in Appendix B for every BBOB function, dimension, and evaluation budget. The profiles retain the individual function scales rather than pooling absolute errors across heterogeneous objectives and therefore provide a trajectory-level complement to the average ranks, paired tests, normalized AUC values, and target-attainment summaries. They also illustrate that the relative advantages of the algorithms are landscape- and budget-dependent rather than uniform across the benchmark. The target-attainment and convergence summaries follow the general anytime-performance perspective used in black-box optimization benchmarking [39]. The normalized absolute-error AUC used here is an additional within-run progress measure and is not intended as a replacement for the standard target-based COCO indicators. 7 Discussion The operator–selection factorization proved useful not only as an analytical description but also as a design principle. The SOMA analysis exposed two actionable properties: the indirect relation between algorithmic parameters and proposal geometry, and the coordinate dependence introduced by axis-aligned masking. GC-SOMA and RA-SOMA address these properties in different ways, showing that proposal-level quantities can be converted into explicit algorithmic controls. Their nonuniform performance also confirms that proposal geometry constrains the candidates supplied to selection but does not alone determine the resulting population dynamics. The performance of the RA-SOMA configuration is consistent with a practical benefit from population-aligned masking, although the present comparison does not isolate the basis transformation from the accompanying population, path, and repair settings. Its population-derived basis aligns masked movement with directions represented by the current population, which is particularly plausible on conditioned and nonseparable landscapes. The weaker behavior observed as dimension increased suggests, however, that a full empirical covariance basis is not uniformly reliable. With a finite and progressively selected population, its eigendirections may become noisy or overly restrictive. More robust alternatives include shrinkage estimates, low-rank subspaces, mixtures of coordinate and eigenvector masking, and less frequent or incremental basis updates. Related work on efficient eigenvector-based crossover indicates that rank-one updates may reduce the cost of repeatedly reconstructing such bases [40]. GC-SOMA illustrates a separate distinction between interpretability and adaptation quality. Expected active dimensionality and leader contraction are meaningful control variables, but specifying them geometrically does not by itself determine an effective schedule. Their values may need to respond to selection success, population diversity, stagnation, and the remaining evaluation budget. GC-SOMA should therefore be regarded as evidence that the derived quantities are operationally usable, rather than as evidence that the particular schedule examined here is optimal. The composite iL-SHOMA-RA method appears primarily to improve evaluation efficiency rather than uniformly improve eventual precision. Its combination of rotation-aware masking, success-history adaptation, and population-size reduction can concentrate evaluations on productive proposals, but the present design does not identify which component is responsible for the observed behavior. Component-level conclusions therefore require an ablation comparing the complete method with variants that introduce each mechanism separately and in selected combinations. Such experiments would also clarify whether the mechanisms are complementary or whether one accounts for most of the gain. The DE methods remain an important reference because they show that analytically motivated SOMA modifications do not imply general dominance over established adaptive optimizers. More broadly, the results emphasize that raw proposal geometry interacts with target construction, population sampling, boundary repair, acceptance rules, adaptation, and population management. The current evidence is also limited to noiseless, unconstrained BBOB functions and dimensions up to 2020. Blockwise ranks, paired tests, target attainment, and convergence profiles should remain the primary evidence because absolute errors are not directly comparable across heterogeneous functions, while the normalized absolute-error AUC is only a within-run progress measure. Computational overhead may also matter when objective evaluations are inexpensive. Future work should therefore examine larger-dimensional, constrained, and noisy problems, including variational quantum landscapes in which sampling noise can alter apparent minima and optimizer rankings [41, 42, 43]. These extensions should determine when proposal-level geometry is sufficient and when it must be coupled to stronger adaptation, selection, and population-management mechanisms. 8 Conclusion This paper introduced an operator–selection factorization for the analysis of population-based optimization algorithms. The central observation is that the nonlinearity of a complete optimization iteration does not preclude exact linear or affine structure within its variation stage. Candidate generation can be isolated as an objective-independent proposal operator once the required population state and random choices have been fixed, while boundary repair, leader selection, survivor replacement, ranking, and adaptation remain separate nonlinear transformations. For SOMA, the canonical path proposal was represented exactly as a linear transformation of an augmented migrant–leader state. In leader-relative coordinates, the operator reduces to coordinate-wise rescaling and exposes the geometric roles of interpolation, projection, overshooting, and masking. Treating the perturbation mask as a Bernoulli random operator yielded closed-form expressions for the proposal mean, covariance, expected step length, expected distance from the leader, active dimensionality, and coordinate coverage. The analysis also clarified how conditioning the mask on nonzero activation introduces cross-coordinate dependence. For DE/rand/1/bin, the same viewpoint produced finite-population moments for differential mutation and separated mutation-induced covariance from the additional uncertainty and dependence introduced by forced-coordinate binomial crossover. These analytical quantities were then used as design variables rather than only descriptive measures. Geometry-Controlled SOMA parameterized proposals through intended active dimensionality and leader contraction, while Rotation-Aware SOMA applied perturbation masks in a population-derived covariance basis. An experimental extension of iSOMA additionally combined rotation-aware masking with success-history adaptation and linear population-size reduction. On the complete noiseless BBOB benchmark, all operator-guided variants substantially improved upon canonical SOMA. Rotation-Aware SOMA was particularly effective in the low-dimensional, high-budget regime, whereas the composite iSOMA extension showed its clearest advantages in restricted-budget and anytime performance. The improvements were nevertheless dependent on dimension, function group, and evaluation budget, and the proposed methods did not uniformly dominate the DE reference algorithms. The results therefore support proposal-level geometry as a useful connection between algorithm analysis and algorithm design, while also showing that proposal structure alone does not determine complete optimizer behavior. Selection, target construction, boundary handling, adaptation, and population management determine how the available geometry is exploited. Future work should examine component-level ablations of the composite method, more robust low-rank or regularized rotation models, adaptive geometry-control schedules, and extensions to larger-dimensional, noisy, and constrained optimization problems. Acknowledgements This project has received funding from the Research Council of Lithuania (LMTLT), agreement No. P-ITP-24-9. This research was also supported by research grants SGS No. SP2026/063 of VSB-Technical University of Ostrava, Czech Republic. Data availability The source code, numerical-verification experiments, BBOB benchmark workflow, processed results, and figure-generation scripts are available at https://github.com/VojtechNovak/SOMA-operator. Further details are available through the corresponding author. Appendix A Numerical Verification and Practical Distortions The analytical results were evaluated using exact enumeration whenever the finite support was computationally manageable and Monte Carlo simulation otherwise. All pseudorandom experiments used the fixed seed 2026031920260319. Each Monte Carlo moment estimate used 150 000150\,000 samples, the coordinate-coverage study used 60 00060\,000 independent trials, and the boundary and selection studies used 150 000150\,000 and 50 00050\,000 samples, respectively. Covariance matrices were calculated using population normalization 1/n1/n, matching the probabilistic definition used in the analytical derivations. For SOMA, the independent-mask experiments covered dimensions D∈3,10,30D∈\3,10,30\, activation probabilities p∈0.1,0.3,0.5,0.8p∈\0.1,0.3,0.5,0.8\, and path parameters t∈0.5,1,1.5,2.5t∈\0.5,1,1.5,2.5\, with five independently generated migrant–leader pairs for each parameter combination. The conditioned mask study used D∈3,10,30D∈\3,10,30\ and p∈0.05,0.2,0.5p∈\0.05,0.2,0.5\. For DE, ordered and distinct donor indices were sampled uniformly without replacement from the target-excluded population. Standard binomial crossover was implemented with one uniformly selected forced coordinate jrandj_rand. The experiments were designed to verify distributional identities rather than to compare optimizer rankings. Boundary repair and fitness-based selection were then examined separately to quantify how these operations distort the raw proposal geometry. A.1 Verification of SOMA proposal moments Under independent Bernoulli masks, the empirical proposal moments closely matched the analytical expressions []=+pt(−)E[ z]= x+pt( L- x) and Cov()=p(1−p)t2diag((−)2).Cov( z)=p(1-p)t^2diag (( L- x)^2 ). Across 240240 primary parameter cases, the mean relative error of the empirical proposal mean was 1.20×10−31.20× 10^-3, and the mean relative Frobenius error of the covariance was 6.50×10−36.50× 10^-3. The corresponding maximum errors were 5.82×10−35.82× 10^-3 and 1.85×10−21.85× 10^-2. The expected squared step length and expected squared leader distance agreed with their analytical values with mean relative errors of 1.81×10−31.81× 10^-3 and 7.09×10−47.09× 10^-4, respectively. A separate convergence experiment increased the sample count from 10310^3 to 1.5×1051.5× 10^5. Over this range, the mean error decreased from 2.14×10−22.14× 10^-2 to 1.01×10−31.01× 10^-3, while the covariance error decreased from 4.63×10−24.63× 10^-2 to 5.76×10−35.76× 10^-3. These results are consistent with ordinary Monte Carlo convergence and provide no indication of a systematic discrepancy between the procedural and analytical proposal distributions. A.2 Conditioning SOMA masks on nonzero activation When all-zero masks are rejected, the mask distribution is no longer a product of independent Bernoulli variables. Writing A=1−(1−p)D,α=pA,A=1-(1-p)^D, α= pA, the conditional mask mean is [πd∣K>0]=αE[ _d K>0]=α, and distinct coordinates satisfy Cov(πd,πe∣K>0)=−p2(1−p)DA2,d≠e.Cov( _d, _e K>0)=- p^2(1-p)^DA^2, d≠ e. Consequently, for Γ=Cov(∣K>0) =Cov( π K>0), [∣K>0]=+tα(−),E[ z K>0]= x+tα( L- x), and Cov(∣K>0)=t2Diag(−)ΓDiag(−).Cov( z K>0)=t^2Diag( L- x) ( L- x). Across nine conditioned-mask cases, the empirical activation probabilities, means, and covariances agreed with these formulas. The mean relative error of the proposal mean was 6.60×10−46.60× 10^-4, and the covariance error was 6.23×10−36.23× 10^-3, with a maximum of 1.17×10−21.17× 10^-2. The induced dependence is substantial in low dimension and at small p. For example, at D=3D=3 and p=0.05p=0.05, conditioning increases the marginal activation probability from 0.050.05 to 0.350570.35057, while the theoretical off-diagonal mask covariance is −0.10537-0.10537; the empirical value was −0.10526-0.10526. By contrast, at D=30D=30 and p=0.5p=0.5, rejection of the all-zero mask has a negligible effect because (1−p)D(1-p)^D is extremely small. The coordinate-coverage simulations also confirmed the path-level expressions. With independently regenerated masks, the probability of covering a coordinate at least once was 1−(1−p)J1-(1-p)^J, whereas reusing one mask left the coverage probability approximately equal to p, independently of J. For example, with D=30D=30 and p=0.2p=0.2, the probability that all coordinates were covered increased from essentially zero at J=5J=5 to 0.7060.706 at J=20J=20 and 0.99960.9996 at J=50J=50. The corresponding theoretical values were 6.72×10−66.72× 10^-6, 0.70620.7062, and 0.999570.99957. A.3 Finite-population DE mutation For a target-excluded donor pool of size M=N−1M=N-1, let −i μ_-i and −i S_-i denote the donor-pool mean and population covariance, with normalization 1/M1/M. Uniform sampling of three ordered, distinct donors gives []=−i,E[ v]= μ_-i, and Cov()=(1+2F2M−1)−i.Cov( v)= (1+2F^2 MM-1 ) S_-i. The formula was tested for N∈6,10,20,50N∈\6,10,20,50\, D∈3,10D∈\3,10\, and F∈0.2,0.5,0.8,1.0F∈\0.2,0.5,0.8,1.0\, using two target indices whenever applicable. Exact enumeration over all ordered donor triples matched the analytical moments to floating-point precision: the mean scaled error was 2.48×10−162.48× 10^-16, and the covariance relative error was 1.71×10−151.71× 10^-15. Monte Carlo estimates produced mean scaled and covariance errors of 2.30×10−32.30× 10^-3 and 4.37×10−34.37× 10^-3, respectively. The finite-population factor is visible even when the donor-pool covariance is fixed. For N=6N=6, corresponding to M=5M=5, the covariance multiplier increases from 1.101.10 at F=0.2F=0.2 to 3.503.50 at F=1F=1. As N grows, the factor approaches 1+2F21+2F^2, recovering the large-population limit while retaining the finite-sampling correction. A.4 Exact forced-coordinate DE crossover The forced coordinate in standard DE/rand/1/bin changes both the marginal activation probability and the dependence structure of the crossover mask. For dimension D and crossover rate q=CRq=CR, define q⋆=q+1−qD.q_ =q+ 1-qD. After marginalizing over jrandj_rand, the mask covariance matrix Γ has entries Γdd=q⋆(1−q⋆),Γde=−(1−q)2D2,d≠e. _d=q_ (1-q_ ), _de=- (1-q)^2D^2, d≠ e. For a fixed target and mutant with =−i δ= v- x_i, the exact conditional moments are [i∣]=i+q⋆E[ u_i v]= x_i+q_ δ and Cov(i∣)=Diag()ΓDiag().Cov( u_i v)=Diag( δ) ( δ). Exact enumeration was performed for D∈2,3,5,10D∈\2,3,5,10\ and CR∈0,0.1,0.5,0.9,1CR∈\0,0.1,0.5,0.9,1\. The exact mean and covariance errors were at floating-point precision, while Monte Carlo mean and covariance errors averaged 7.91×10−47.91× 10^-4 and 3.64×10−33.64× 10^-3. The discrepancy from the independent-mask approximation is largest at small dimension and low crossover rate. For example, at D=3D=3 and CR=0.1CR=0.1, the effective activation probability is q⋆=0.4q_ =0.4, and the off-diagonal mask covariance is −0.09-0.09. Thus, treating the mask entries as independent Bernoulli variables with activation probability 0.10.1 substantially mischaracterizes the standard crossover distribution. A.5 Combined DE trial moments and parameter effects Let =[]−i. δ=E[ v]- x_i. Combining donor randomness with forced crossover yields [i]=i+q⋆E[ u_i]= x_i+q_ δ and Cov(i)=q⋆2Cov()+Γ⊙[Cov()+⊤,]Cov( u_i)=q_ ^2Cov( v)+ [Cov( v)+ δ δ , ] where ⊙ denotes the Hadamard product. The combined expression was tested in 2727 cases covering three population–dimension pairs, F∈0.3,0.6,0.9F∈\0.3,0.6,0.9\, and CR∈0.2,0.5,0.9CR∈\0.2,0.5,0.9\. Exact enumeration again matched the closed form to numerical precision, with covariance error below 9.01×10−159.01× 10^-15. The Monte Carlo trial-mean scaled error averaged 2.20×10−32.20× 10^-3, and the covariance error averaged 6.45×10−36.45× 10^-3, with a maximum of 1.13×10−21.13× 10^-2. The parameter study confirms that F primarily controls mutation spread, whereas CRCR controls how strongly that geometry is transferred to the trial distribution. In the N=30,D=10N=30,D=10 population, at CR=0.5CR=0.5, increasing F from 0.30.3 to 0.90.9 increased the trial covariance trace from 9.229.22 to 19.5519.55. At F=0.9F=0.9, increasing CRCR from 0.20.2 to 0.90.9 increased the trace from 10.2510.25 to 31.1431.14. The principal-axis angle between mutant and trial covariance decreased from 23.82∘23.82 to 1.91∘1.91 over the same crossover-rate change, showing that high crossover rates make the trial geometry closely inherit the mutant orientation. The closed-form proposal moments describe raw proposals before feasibility repair. To quantify the resulting limitation, raw SOMA and DE proposals were generated near the boundary of Ω=[−1,1]D =[-1,1]^D and repaired by clipping or repeated reflection. Both operators changed the proposal means and covariances, and the effect increased rapidly with the proportion of infeasible samples. For SOMA with D=10D=10, p=0.8p=0.8, and t=3t=3, all sampled raw proposals required repair. Clipping retained only 27.1%27.1\% of the raw covariance trace, while reflection retained 8.7%8.7\%. For DE with D=10D=10, F=1F=1, and CR=0.9CR=0.9, 95.6%95.6\% of trials required repair; clipping and reflection retained 63.0%63.0\% and 44.8%44.8\% of the raw covariance trace, respectively. At less aggressive settings, the distortion was considerably smaller. These results confirm that repair cannot be treated as a negligible implementation detail in high-spread regimes. Finally, raw proposal distributions were compared with the distributions remaining after standard fitness-based selection on the sphere function and on a rotated quadratic objective. For DE, greedy target–trial selection produced acceptance rates between 0.4370.437 and 0.9620.962. The selected covariance trace ranged from 42.3%42.3\% to 96.0%96.0\% of the raw trial trace, depending on dimension, objective, F, and CRCR. The higher-spread setting (F,CR)=(0.9,0.9)(F,CR)=(0.9,0.9) generally produced a larger discrepancy between pre-selection and post-selection covariance. The SOMA experiment selected the best point from a discretized leader-directed proposal family. Because the controlled leader was placed closer to the minimizer, at least one proposal was selected in almost every trial. Even in this favorable setting, selection substantially reshaped the distribution: the selected-to-raw covariance-trace ratio ranged from 0.0810.081 to 1.4101.410. Selection therefore need not act as a uniform contraction of the proposal cloud; it can suppress, preserve, or amplify different directions depending on the objective and candidate family. A.6 Summary The exact-enumeration results verify the finite-population DE donor and forced-crossover formulas to floating-point precision, while the Monte Carlo studies validate the SOMA and DE moment expressions with covariance errors generally below 2%2\%. More importantly, the conditioned-mask, boundary, and selection experiments identify where the raw linear-operator description must be qualified. Rejecting all-zero masks and forcing a DE crossover coordinate introduce cross-coordinate dependence, while repair and selection can substantially alter both proposal means and covariance. The operator formulas therefore provide an exact description of the raw variation stage under their stated assumptions, but not of the complete feasible or selected update. Appendix B Function-wise convergence profiles This appendix reports function-wise convergence profiles for all 24 noiseless BBOB functions, organized by the five standard function groups and the two evaluation budgets. Each figure combines the results for D∈5,10,20D∈\5,10,20\ for one function group and budget. The curves show the median best-so-far absolute error over all instances and independent repetitions, while the shaded regions indicate the corresponding interquartile ranges. These profiles complement the aggregate ranks, statistical comparisons, normalized AUC values, and target-attainment results reported in Section 6.2 by showing when performance differences emerge and whether they persist on individual functions. Figure 5: Convergence profiles for BBOB Group 1 (Separable functions) under a total budget of 500D500D function evaluations (500500 FE/D). The three vertically stacked blocks correspond to D=5D=5, D=10D=10, and D=20D=20. Curves show the median best-so-far absolute error; shaded bands show the interquartile range over instances and repeats. Figure 6: Convergence profiles for BBOB Group 2 (Low or moderate conditioning) under a total budget of 500D500D function evaluations (500500 FE/D). The three vertically stacked blocks correspond to D=5D=5, D=10D=10, and D=20D=20. Curves show the median best-so-far absolute error; shaded bands show the interquartile range over instances and repeats. Figure 7: Convergence profiles for BBOB Group 3 (High conditioning and unimodal) under a total budget of 500D500D function evaluations (500500 FE/D). The three vertically stacked blocks correspond to D=5D=5, D=10D=10, and D=20D=20. Curves show the median best-so-far absolute error; shaded bands show the interquartile range over instances and repeats. Figure 8: Convergence profiles for BBOB Group 4 (Multimodal with adequate global structure) under a total budget of 500D500D function evaluations (500500 FE/D). The three vertically stacked blocks correspond to D=5D=5, D=10D=10, and D=20D=20. Curves show the median best-so-far absolute error; shaded bands show the interquartile range over instances and repeats. Figure 9: Convergence profiles for BBOB Group 5 (Multimodal with weak global structure) under a total budget of 500D500D function evaluations (500500 FE/D). The three vertically stacked blocks correspond to D=5D=5, D=10D=10, and D=20D=20. Curves show the median best-so-far absolute error; shaded bands show the interquartile range over instances and repeats. Figure 10: Convergence profiles for BBOB Group 1 (Separable functions) under a total budget of 2500D2500D function evaluations (25002500 FE/D). The three vertically stacked blocks correspond to D=5D=5, D=10D=10, and D=20D=20. Curves show the median best-so-far absolute error; shaded bands show the interquartile range over instances and repeats. Figure 11: Convergence profiles for BBOB Group 2 (Low or moderate conditioning) under a total budget of 2500D2500D function evaluations (25002500 FE/D). The three vertically stacked blocks correspond to D=5D=5, D=10D=10, and D=20D=20. Curves show the median best-so-far absolute error; shaded bands show the interquartile range over instances and repeats. Figure 12: Convergence profiles for BBOB Group 3 (High conditioning and unimodal) under a total budget of 2500D2500D function evaluations (25002500 FE/D). The three vertically stacked blocks correspond to D=5D=5, D=10D=10, and D=20D=20. Curves show the median best-so-far absolute error; shaded bands show the interquartile range over instances and repeats. Figure 13: Convergence profiles for BBOB Group 4 (Multimodal with adequate global structure) under a total budget of 2500D2500D function evaluations (25002500 FE/D). The three vertically stacked blocks correspond to D=5D=5, D=10D=10, and D=20D=20. Curves show the median best-so-far absolute error; shaded bands show the interquartile range over instances and repeats. Figure 14: Convergence profiles for BBOB Group 5 (Multimodal with weak global structure) under a total budget of 2500D2500D function evaluations (25002500 FE/D). The three vertically stacked blocks correspond to D=5D=5, D=10D=10, and D=20D=20. Curves show the median best-so-far absolute error; shaded bands show the interquartile range over instances and repeats. Appendix C Implementation and parameter settings The experiments used the noiseless COCO/BBOB functions f1f_1–f24f_24, dimensions D∈5,10,20D∈\5,10,20\, instances 11–55, and budgets B∈500D,2500DB∈\500D,2500D\. Ten independent repetitions were performed per function–instance block. All pseudorandom seeds were derived deterministically from the base seed 2026072320260723, the dimension, budget, function, instance, and repetition using NumPy SeedSequence. Initial populations were sampled uniformly within the COCO bounds. GC-SOMA, RA-SOMA, iSOMA, iL-SHOMA-RA, and SciPy-DE used population size 4D4D; the first four and SciPy-DE received the same initial population within each block, while iL-SHADE used its internally generated seeded population. Canonical SOMA used a separately seeded population of size 10D10D. Initial-population evaluations were included in the budget, and a common objective wrapper enforced an exact maximum of B true function evaluations. Table 9: Optimizer configurations used in the BBOB experiments. Masks were regenerated independently at every proposal point and were permitted to be all zero. Method Configuration Canonical-SOMA All-to-One with a fixed leader during each migration and a common origin for all path points. Population 10D10D, PRT=0.10PRT=0.10, Step=0.11Step=0.11, and PathLength=3.0PathLength=3.0, giving t∈0.11,0.22,…,2.97t∈\0.11,0.22,…,2.97\. Each migrant retained its best path point only when it improved upon the origin. Boundary violations were clipped. GC-SOMA Population 4D4D, identity masking basis, six proposals per migrant, and reflection repair. For consumed-budget fraction s=q/Bs=q/B, p=0.35+0.55s,ρc=max0.85−0.70s, 1−p+10−6.p=0.35+0.55s, 17.00024pt _c= \0.85-0.70s,\,1-p+10^-6\. For δj∈linspace(0.12,−0.12,6) _j (0.12,-0.12,6), ρj=clip(ρc+δj,max1−p+10−6,0.02, 0.98), _j=clip ( _c+ _j,\, \1-p+10^-6,0.02\,\,0.98 ), and tj=1±1−1−ρjp.t_j=1± 1- 1- _jp. The overshooting branch was used for s<0.55s<0.55, and the interpolating branch thereafter. The best improving proposal was retained. RA-SOMA Population 4D4D, PRT=0.30PRT=0.30, and six equally spaced path values t∈0.4,0.8,1.2,1.6,2.0,2.4t∈\0.4,0.8,1.2,1.6,2.0,2.4\. At the start of each migration, the eigenvectors of the current population covariance formed the masking basis. The covariance used normalization 1/N1/N and regularization 10−12(trΣ/D)I10^-12(tr /D)I, or 10−12I10^-12I when the trace vanished. Reflection repair and best-path-point acceptance were used. iSOMA The supplied ISOMAOptimizer implementation was loaded without modification. Population 4D4D, Njump=10N_jump=10, Step=0.3Step=0.3, MaxMigration=10000MaxMigration=10000, m=10m=10, n=5n=5, and k=15k=15. The maximum function evaluations and COCO bounds were passed directly to the implementation. Its seeded numpy.random.rand initialization was matched to the common 4D4D population. iL-SHOMA-RA Initial population N0=4DN_0=4D, Nmin=max4,D=DN_ = \4,D\=D, with Njump=10N_jump=10, Step=0.3Step=0.3, and initial selection parameters m=10m=10, n=5n=5, k=15k=15; their proportions were preserved as the population decreased. The covariance basis was recomputed each migration. Six-element memories were initialized to Mp=0.5M_p=0.5 and Ms=1M_s=1. For a randomly selected memory entry, p=clip((0.7Mp+0.3[0.1+0.9q/B], 0.12), 1/D, 1),p=clip\! (N\! (0.7M_p+0.3[0.1+0.9q/B],\,0.1^2 ),\,1/D,\,1 ), while the path scale was sampled from Cauchy(Ms,0.1)Cauchy(M_s,0.1), restricted to 0.05<s≤20.05<s≤ 2. For jump j=1,…,10j=1,…,10, t=s(10−j+1)0.3t=s(10-j+1)0.3. First improvement was accepted. Successful p values updated the memory by an improvement-weighted arithmetic mean and successful scales by an improvement-weighted Lehmer mean. The population was reduced linearly while retaining the best individuals. After more than 50N50N unsuccessful migrant attempts without improvement, 10%10\% of the non-best population was resampled uniformly. Boundary violations were clipped. SciPy-DE SciPy differential_evolution with population 4D4D supplied through init, strategy best1bin, differential weight F=0.8F=0.8, crossover rate CR=0.9CR=0.9, immediate updating, one worker, no polishing, and tol=atol=0tol=atol=0. The SciPy population multiplier was 44, and the iteration limit was ⌈(B−4D)/(4D)⌉ (B-4D)/(4D) ; the common objective wrapper enforced the exact budget. iL-SHADE PyADE ilshade.get_default_params(D) was used, with only the objective, COCO bounds, random seed, budget B, and initial population size 4D4D overridden. Population initialization remained internal to PyADE. Callbacks and optional arguments were disabled, and no fallback implementation was used. The reported objective value was the absolute BBOB error Δf=maxf(x)−fopt,0 f= \f(x)-f_opt,0\. Target attainment was evaluated at 102,101,100,10−1,…,10−8.10^2,10^1,10^0,10^-1,…,10^-8. Best-so-far histories were retained at every true evaluation; the published convergence files used the union of 201 linearly spaced and 201 geometrically spaced evaluation points. The numerical floor used for logarithms and the normalized AUC was 10−30010^-300. References [1] R. Storn, K. Price, Differential evolution—a simple and efficient heuristic for global optimization over continuous spaces, Journal of Global Optimization 11 (4) (1997) 341–359. doi:10.1023/A:1008202821328. [2] D. Davendra, I. Zelinka, Self-Organizing Migrating Algorithm: Methodology and Implementation, Springer International Publishing, Cham, 2016. doi:10.1007/978-3-319-28161-2. [3] L. Skanderova, Self-organizing migrating algorithm: Review, improvements and comparison, Artificial Intelligence Review 56 (1) (2023) 101–172. doi:10.1007/s10462-022-10167-8. [4] A. E. Eiben, J. E. Smith, Introduction to Evolutionary Computing, 2nd Edition, Springer, Berlin, Heidelberg, 2015. doi:10.1007/978-3-662-44874-8. [5] I. Zelinka, A survey on evolutionary algorithms dynamics and its complexity—mutual relations, past, present and future, Swarm and Evolutionary Computation 25 (2015) 2–14. doi:10.1016/j.swevo.2015.06.002. [6] B. Doerr, F. Neumann, A survey on recent progress in the theory of evolutionary algorithms for discrete optimization, ACM Transactions on Evolutionary Learning and Optimization 1 (4) (2021) 1–43. doi:10.1145/3472304. [7] S. Das, P. N. Suganthan, Differential evolution: A survey of the state-of-the-art, IEEE Transactions on Evolutionary Computation 15 (1) (2011) 4–31. doi:10.1109/TEVC.2010.2059031. [8] S. Das, S. S. Mullick, P. N. Suganthan, Recent advances in differential evolution—an updated survey, Swarm and Evolutionary Computation 27 (2016) 1–30. doi:10.1016/j.swevo.2016.01.004. [9] K. R. Opara, J. Arabas, Differential evolution: A survey of theoretical analyses, Swarm and Evolutionary Computation 44 (2019) 546–558. doi:10.1016/j.swevo.2018.06.010. [10] V. Novák, T. Bezděk, I. Zelinka, S. Das, M. Beseda, A longitudinal analysis of the CEC single-objective competitions (2010–2024) and implications for variational quantum optimization, Swarm and Evolutionary Computation (2026) 102469doi:10.1016/j.swevo.2026.102469. [11] K. R. Opara, J. Arabas, Comparison of mutation strategies in differential evolution—a probabilistic perspective, Swarm and Evolutionary Computation 39 (2018) 53–69. doi:10.1016/j.swevo.2017.12.007. [12] K. R. Opara, J. Arabas, The contour fitting property of differential mutation, Swarm and Evolutionary Computation 50 (2019) 100441. doi:10.1016/j.swevo.2018.09.001. [13] D. Zaharie, Influence of crossover on the behavior of differential evolution algorithms, Applied Soft Computing 9 (3) (2009) 1126–1138. doi:10.1016/j.asoc.2009.02.012. [14] M. Pluháček, A. Kazikova, T. Kadavy, A. Viktorin, R. Senkerik, Explaining SOMA: The relation of stochastic perturbation to population diversity and parameter space coverage, in: Proceedings of the Genetic and Evolutionary Computation Conference Companion, Association for Computing Machinery, New York, NY, USA, 2021, p. 1944–1952. doi:10.1145/3449726.3463211. [15] R. Biedrzycki, J. Arabas, D. Jagodziński, Bound constraints handling in differential evolution: An experimental study, Swarm and Evolutionary Computation 50 (2019) 100453. doi:10.1016/j.swevo.2018.10.004. [16] J. M. Tomczak, E. Węglarz-Tomczak, A. E. Eiben, Differential evolution with reversible linear transformations, in: Proceedings of the 2020 Genetic and Evolutionary Computation Conference Companion, Association for Computing Machinery, New York, NY, USA, 2020, p. 205–206, extended version available as arXiv:2002.02869. doi:10.1145/3377929.3389972. [17] Q. B. Diep, T. C. Truong, S. Das, I. Zelinka, Self-organizing migrating algorithm with narrowing search space strategy for robot path planning, Applied Soft Computing 116 (2022) 108270. doi:10.1016/j.asoc.2021.108270. [18] J. Brest, S. Greiner, B. Bošković, M. Mernik, V. Žumer, Self-adapting control parameters in differential evolution: A comparative study on numerical benchmark problems, IEEE Transactions on Evolutionary Computation 10 (6) (2006) 646–657. doi:10.1109/TEVC.2006.872133. [19] A. K. Qin, V. L. Huang, P. N. Suganthan, Differential evolution algorithm with strategy adaptation for global numerical optimization, IEEE Transactions on Evolutionary Computation 13 (2) (2009) 398–417. doi:10.1109/TEVC.2008.927706. [20] J. Zhang, A. C. Sanderson, JADE: Adaptive differential evolution with optional external archive, IEEE Transactions on Evolutionary Computation 13 (5) (2009) 945–958. doi:10.1109/TEVC.2009.2014613. [21] R. Tanabe, A. S. Fukunaga, Success-history based parameter adaptation for differential evolution, in: 2013 IEEE Congress on Evolutionary Computation, IEEE, 2013, p. 71–78. doi:10.1109/CEC.2013.6557555. [22] R. Tanabe, A. S. Fukunaga, Improving the search performance of SHADE using linear population size reduction, in: 2014 IEEE Congress on Evolutionary Computation, IEEE, 2014, p. 1658–1665. doi:10.1109/CEC.2014.6900380. [23] J. Brest, M. Sepesy Maučec, B. Bošković, iL-SHADE: Improved L-SHADE algorithm for single-objective real-parameter optimization, in: 2016 IEEE Congress on Evolutionary Computation, IEEE, 2016, p. 1188–1195. doi:10.1109/CEC.2016.7743922. [24] S.-M. Guo, C.-C. Yang, Enhancing differential evolution utilizing eigenvector-based crossover operator, IEEE Transactions on Evolutionary Computation 19 (1) (2015) 31–49. doi:10.1109/TEVC.2013.2297160. [25] F. Caraffini, F. Neri, A study on rotation invariance in differential evolution, Swarm and Evolutionary Computation 50 (2019) 100436. doi:10.1016/j.swevo.2018.08.013. [26] J. Huxley, Evolution: The Modern Synthesis, George Allen & Unwin, London, 1942. [27] N. Hansen, A. Ostermeier, Completely derandomized self-adaptation in evolution strategies, Evolutionary Computation 9 (2) (2001) 159–195. doi:10.1162/106365601750190398. [28] N. Hansen, R. Ros, N. Mauny, M. Schoenauer, A. Auger, Impacts of invariance in search: When CMA-ES and PSO face ill-conditioned and non-separable problems, Applied Soft Computing 11 (8) (2011) 5755–5769. doi:10.1016/j.asoc.2011.03.001. [29] P. Virtanen, R. Gommers, T. E. Oliphant, M. Haberland, T. Reddy, D. Cournapeau, E. Burovski, P. Peterson, W. Weckesser, J. Bright, S. J. van der Walt, M. Brett, J. Wilson, K. J. Millman, P. van Mulbregt, SciPy 1.0 Contributors, SciPy 1.0: Fundamental algorithms for scientific computing in python, Nature Methods 17 (2020) 261–272. doi:10.1038/s41592-019-0686-2. [30] xKuZz, PyADE: Python advanced differential evolution algorithms library, PyPI package, version 1.1, accessed 31 July 2026 (2019). URL https://pypi.org/project/pyade-python/ [31] J. J. Moré, S. M. Wild, Benchmarking derivative-free optimization algorithms, SIAM Journal on Optimization 20 (1) (2009) 172–191. doi:10.1137/080724083. [32] V. Beiranvand, W. Hare, Y. Lucet, Best practices for comparing optimization algorithms, Optimization and Engineering 18 (4) (2017) 815–848. doi:10.1007/s11081-017-9366-1. [33] N. Hansen, A. Auger, R. Ros, O. Mersmann, T. Tušar, D. Brockhoff, COCO: A platform for comparing continuous optimizers in a black-box setting, Optimization Methods and Software 36 (1) (2021) 114–144. doi:10.1080/10556788.2020.1808977. [34] S. Finck, N. Hansen, R. Ros, A. Auger, Real-parameter black-box optimization benchmarking 2009: Noiseless functions definitions, Tech. Rep. R-6829, INRIA, updated version as of February 2019 (2009). URL https://inria.hal.science/inria-00362633v2/document [35] N. Hansen, A. Auger, S. Finck, R. Ros, Real-parameter black-box optimization benchmarking 2009: Experimental setup, Tech. Rep. R-6828, INRIA (2009). URL https://inria.hal.science/inria-00362649v2/document [36] M. Jamil, X.-S. Yang, A literature survey of benchmark functions for global optimisation problems, International Journal of Mathematical Modelling and Numerical Optimisation 4 (2) (2013) 150–194. doi:10.1504/IJMMNO.2013.055204. [37] J. Derrac, S. García, D. Molina, F. Herrera, A practical tutorial on the use of nonparametric statistical tests as a methodology for comparing evolutionary and swarm intelligence algorithms, Swarm and Evolutionary Computation 1 (1) (2011) 3–18. doi:10.1016/j.swevo.2011.02.002. [38] S. García, A. Fernández, J. Luengo, F. Herrera, Advanced nonparametric tests for multiple comparisons in the design of experiments in computational intelligence and data mining: Experimental analysis of power, Information Sciences 180 (10) (2010) 2044–2064. doi:10.1016/j.ins.2009.12.010. [39] N. Hansen, A. Auger, D. Brockhoff, T. Tušar, Anytime performance assessment in blackbox optimization benchmarking, IEEE Transactions on Evolutionary Computation 26 (6) (2022) 1293–1305. doi:10.1109/TEVC.2022.3210897. [40] T. J. Choi, An efficient eigenvector-based crossover for differential evolution: Simplifying with rank-one updates, AIMS Mathematics 10 (2) (2025) 3500–3522. doi:10.3934/math.2025162. [41] S. Illésová, V. Novák, T. Bezděk, C. Possel, M. Beseda, Numerical optimization strategies for the variational hamiltonian ansatz in noisy quantum environments (2025). arXiv:2505.22398, doi:10.48550/arXiv.2505.22398. [42] V. Novák, S. Illésová, T. Bezděk, I. Zelinka, M. Beseda, Reliable optimization under noise in quantum variational algorithms (2025). arXiv:2511.08289, doi:10.48550/arXiv.2511.08289. [43] V. Novák, I. Zelinka, V. Snášel, Optimization strategies for variational quantum algorithms in noisy landscapes, Evolutionary IntelligenceAccepted for publication; preprint available as arXiv:2506.01715 (2026). arXiv:2506.01715.