First Mathematical Runtime Analyses of Multi-Objective Evolutionary Algorithms for Multi-Valued Decision Variables

1 Introduction

Evolutionary algorithms are widely used for solving complex optimization problems, in particular, for those with multiple conflicting objectives. Multi-objective evolutionary algorithms (MOEAs) address such problems by exploring the different trade-offs between the objectives in a single run. Well-known MOEAs such as NSGA-II, NSGA-III, SPEA2, and SMS-EMOA are widely used in both research and applications. Their empirical success has motivated a considerable body of theoretical works aimed at understanding how MOEAs work and why they perform well in practice. These include, among others, guarantees on the time required to cover the entire Pareto front [Laumanns et al., 2002; Zheng and Doerr, 2023; Wietheger and Doerr, 2023; Zheng and Doerr, 2024; Ren et al., 2024], analyses on the benefits of crossover operators [Doerr and Qu, 2023b; Dang et al., 2024; Opris, 2026a], the role of diversity preservation and archiving mechanisms [Covantes Osuna et al., 2020; Bian et al., 2024; Ren et al., 2025], the impact of non-elitist strategies [Zheng et al., 2024; Li et al., 2025b; Bian et al., 2025; Bendahi et al., 2025], and approximation guarantees for situations in which the full Pareto front is too costly to be computed [Horoba and Neumann, 2008; Zheng and Doerr, 2025; Deng et al., 2025; Alghouass et al., 2025; Li et al., 2025a].

∗Corresponding author.

The vast majority of the theoretical works (in discrete optimization) have focused on problems defined on binary decision spaces. To the best of our knowledge, only two runtime analyses of MOEAs move beyond binary decision spaces, and they do so by considering the other extreme of unbounded integer-valued search domains [Rudolph, 2023; Doerr et al., 2025b]. Among them, Rudolph [2023] analyzed subproblems of a multi-objective problem in integer search spaces, but the analysis remains a single-objective analysis. Following this work, Doerr et al. [2025b] studied evolutionary algorithms on unbounded integer domains, focusing on how different mutation strength distributions affect the runtime. Till now, no mathematical runtime analyses exist for MOEAs dealing with decision variables that take a finite number r of values that is larger than two.

Our Contributions: In this work, we undertake a first attempt to analyze MOEAs for problems with multi-valued decision variables with mathematical means. We first extend the classical bi-objective ONEMINMAX benchmark to an rvalued domain, r > 2. This benchmark preserves the simplicity and structural clarity of its binary counterpart, but (naturally) has a larger Pareto front. We then consider the classic SEMO algorithm with unit-strength local mutation, which can be regarded as the natural multi-valued analogue of onebit mutation. For this algorithm, we prove an upper bound of O(n2 r2 log n) and a lower bound of Ω(n2 r(r+log n)) on the expected number of iterations until the population covers the entire Pareto front. The dependence on r in these bounds is moderate. Although the search space has size rn , the Pareto front of our benchmark grows only linearly with r. Since the runtime is measured by the time to cover this front, its dependence on r is tied more to the linearly growing front than to the exponentially growing search space. The (small) gap between these bounds stems from the population dynamics, which (as this work shows) are much more complex in the multi-valued setting.

We further consider a very mild modification of the tiebreaking rule in the selection of the next population, namely that we accept only strictly better solutions (that is, in case of equal function values, we keep the old solution and discard the new one). This modification allows for a more tractable analysis of the population dynamics, and allows us to prove an improved upper bound of O(n2 r(r + logn)), which matches our lower bound (that is valid also for this modification of the algorithm).

We do not see any reason why the modified algorithm should be faster than the original one. Our (simple) experimental results also could not detect a significant difference (rather the modified algorithm appeared slightly slower). To check that this observation is not specific to r-valued ONEMINMAX, where all solutions are Pareto optimal, we additionally introduce and experimentally study an r-valued analogue of LOTZ. This benchmark has the same Paretofront size, but, among many other differences, most solutions are not Pareto-optimal. Again, we observe no noticeable runtime difference between the original and the modified algorithm. For these reasons, we conjecture that the O(n2 r(r + logn)) bound, proven for the modified SEMO on r-valued ONEMINMAX and matching our lower bound, is also the true runtime for the original SEMO. Proving this conjecture, unfortunately, seems out of reach with the current methods of this field, and might be a worthwhile challenge for the future, that could give deep insights into the more complex population dynamics in multi-valued search spaces.

2 RelatedWork

The theory of randomized search heuristics is a relatively young subarea of artificial intelligence. While first mathematical runtime analyses for single-objective algorithms appeared already in the 1990s, say [M¨uhlenbein, 1992; B¨ack, 1993], the first such works for multi-objective algorithms took around another ten years [Laumanns et al., 2002] and the truly practically relevant algorithms could only be analyzed very recently [Li et al., 2016; Zheng and Doerr, 2023; Wietheger and Doerr, 2023; Bian et al., 2025; Ren et al., 2024].

As stated already in the introduction, there are no theoretical works analyzing how an MOEA solves a problem defined over decision variables taking a finite number r > 2 of values. Closest to our work are the two recent papers [Rudolph, 2023; Doerr et al., 2025b] that analyze discrete multi-objective problems defined on Zn . Only the second work considers a true MOEA, namely the SEMO regarded also in this work and the related GSEMO. The main difficulty in this work, however, is dealing with the unboundedness of the set of values taken by the decision variables. Consequently, the methods employed there are of little help for our problem setting.

In contrast, in the theory of single-objective evolutionary algorithms a few sporadic works on multi-valued optimization exist. Doerr et al. [2011] analyzed the (1+1) EA with uniform mutation strength, that is, resetting a component to a random other value, on linear functions over {0, 1, 2}n . They proved an expected runtime of O(n logn), matching the well-known bound for the binary case. The result was then extended by Doerr and Pohl [2012], where they studied r-valued linear functions and derived an upper bound of O(rn logn) + O(r3n log logn) and a lower bound of Ω(rn log n).

Subsequent work considered different and more complex mutation operators and in return only regarded the simpler ONEMAX problem. K¨otzing et al. [2015] analyzed the (1+1) EA with unit-strength mutation on dynamically changing ONEMAX where each dimension has r different possible values. Doerr et al. [2018] considered the (1+1) EA with three different mutation strengths and proved asymptotically tight runtimes for a variety of ONEMAX-type test functions over an alphabet of size r. They further analyzed the performance of randomized local search, a variant of the (1 + 1) EA using local mutation, with a selfadjusting mutation strength and proved an expected optimization time of Θ(n(logn + log r)), which is optimal among all comparison-based search heuristics. Harder et al.[2024] studied the (1 + 1) EA with different mutation operators on the unbounded integer search space Zn . All these works are hardly comparable to ours since they regard single-objective optimization via an algorithm with parent and offspring population equal to one. Hence the population dynamics are less interesting here.

In [Qian et al., 2018a; Qian et al., 2018b], submodular functions with multi-valued discrete domains are studied. Both works formulate a constrained single-objective problem as a bi-objective problem. They then use problem-specific algorithms and algorithms similar to classic MOEAs to solve these. However, since their objective is to find a single solution with a particular constraint satisfaction (equivalent to a particular value of the second objective), they do not analyze how the full Pareto front is computed or approximated, but only how long it takes until such a particular solution is contained in the population. For that reason, we could not derive deeper insights on multi-objective optimization from these works, and feel that they should be viewed as singleobjective optimization works.

More recently, estimation-of-distribution algorithms (EDAs) have been analyzed in multi-valued domains, see, e.g., [Ben Jedidia et al., 2024; Adak and Witt, 2024]. Due to the very different algorithmic concept of an EDA, these works are hardly comparable to ours.

Overall, these results demonstrate that runtime analysis in multi-valued decision spaces is feasible in the singleobjective setting, but typically requires more involved techniques than in the binary case.

3 Preliminaries

Throughout the paper, for all a ≤ b, we write [a..b] := {i ∈ Z | a ≤ i ≤ b} to denote the interval of integers between a and b.

3.1 Multi-Objective Optimization

This paper considers multi-objective optimization problems defined on the search space [0..r - 1]n , where r > 2. Such a problem is given by f : [0..r - 1]n → Rm. Our goal is to maximize the objectives. For two solutions x, y ∈ [0..r-1]n , we say that x weakly dominates y if fi(x) ≥ fi(y) holds for all i ∈ [1..m], denoted as x ⪰ y. If, in addition, there exists at least one objective j such that fj(x) > fj(y), then x is said to dominate y, denoted as x ≻ y. A solution is called Pareto optimal if it is not dominated by any other solution in the search space. The set of all Pareto-optimal solutions is referred to as the Pareto set. The Pareto front consists of the objective vectors corresponding to these solutions. As common in the evolutionary computation theory community [Zhou et al., 2019; Doerr and Neumann, 2020], we measure the runtime by the number of function evaluations for the population to cover the entire Pareto front. We often use asymptotic notation (big-Oh notation, Landau symbols). If so, then the underlying variable is the problem size n. We allow other quantities to be functionally dependent on n, for example, the number r of different values taken by the decision variables. This is natural, for example, in rooted tree problems in graphs on n vertices, we have n - 1 decision variables (describing the predecessor of a non-root node on the path to the source), and each of them can take any other node as value (so there are n - 1 different values).

3.2 The Multi-Valued Algorithm

The simple evolutionary multi-objective optimizer (SEMO) was the algorithm first used in a mathematical runtime analysis of an MOEA [Laumanns et al., 2002]. It is still the algorithm best understood from a theoretical perspective. Starting from a single random solution, the algorithm repeatedly selects a random parent from the population, generates an offspring via mutation, and keeps as next population a largest pair-wise non-dominated subset of the old population and the offspring, giving preference to the offspring in case of ties.

In the classic binary setting, the SEMO employs one-bit mutation, which flips a single bit chosen uniformly at random from the n bit positions. In the multi-valued setting considered in this work, we employ a unit-strength local mutation operator, which has, to the best of our knowledge, first been proposed and argued for in the single-objective work [Doerr et al., 2018]. This operator, as one of several alternatives, was also studied in the recent work on multi-objective optimization for unbounded integer-valued problems [Doerr et al., 2025b]. When applying this mutation operator to an individual x ∈ [0..r - 1]n , one position i ∈ [1..n] is chosen uniformly at random, and the value xi of this position is increased or decreased by 1 with equal probability. A mutation that would lead to a value outside the allowed interval [0..r - 1] is considered infeasible and discarded (that is, the parent of the mutation operation is returned as offspring). The pseudocode of the algorithm is given in Algorithm 1.

Algorithm 1: The r-valued SEMO with unit-strength local mutation to maximize f : [0..r - 1]n → Rm.

• 1 Initialize x ∈ [0..r -1]n uniformly at random and set P0 ← {x}; 2 for t = 0, 1, 2, . . . do 3 Select one individual x uniformly at random from Pt; 4 Generate y via choosing one entry of x uniformly at random and modifying it by ±1 with equal probability; 5 if there is no z ∈ Pt such that z ≻ y then 6 Pt+1 ← {z ∈ Pt | y ⪰̸ z} ∪ {y}; 7 else 8 Pt+1 ← Pt;

3.3 The Multi-Valued Benchmark

The ONEMINMAX problem is the best studied benchmark in the theory of MOEAs. Its two objectives are the number of zeros in the bit string and the number of ones. This is equivalent to saying that the objectives are the sum of the bit values and the sum of the inverses of the bit values. We therefore extend the binary ONEMINMAX benchmark to the multi-valued search space [0..r-1]n by letting the first objective be the sum of the values of the entries of x and the second objective be the sum of the inverses. This choice of the objective order makes the definition closer to the COCZ benchmark, but has no effect on the theoretical results. Moreover, the same objective order was used in [Zheng, 2026], where the corresponding binary problem was called ONEMAXMIN. Since in the single-objective multi-valued setting, Adak and Witt [2024] called this first objective G-ONEMAX, we refer to our multi-valued extension of the ONEMINMAX benchmark as the generalized ONEMAXMIN problem, denoted by G-ONEMAXMIN.

Definition 1. Let n ∈ N and r ∈ N≥2. The G-ONEMAXMIN function f(x) = (f1(x), f2(x)) : [0..r-1]n -→R2 is defined by

The following lemma characterizes the Pareto set and the Pareto front of G-ONEMAXMIN. As in the binary case, every solution in the search space is Pareto optimal, but the size of the Pareto front also increases with r. Due to space limitations, all proofs are omitted. A complete version is available on arXiv [Li et al., 2026].

Lemma 2. The Pareto set of G-ONEMAXMIN is

The Pareto front is

whose size is nr - n + 1.

As a second benchmark in the experimental section, we use a multi-valued analogue of the classic LOTZ benchmark. With very few solutions Pareto-optimal, it can be seen as a problem complementary to G-ONEMAXMIN. The first objective follows the idea of the r-valued LEADINGONES benchmark [Ben Jedidia et al., 2024], but refines it by also giving value-progress credit to the first position after the certified prefix. Hence a complete prefix of entries equal to r -1 contributes fully, the first subsequent position contributes its current value, and all later positions contribute nothing. Symmetrically, the second objective measures the same progress from the right towards value 0. Values lying strictly between the active leading position and the active trailing position are deliberately ignored: they are neither part of the prefix of r - 1 values nor part of the suffix of 0 values. This construction leaves only a narrow set of Pareto-optimal solutions, while preserving a Pareto front of size n(r - 1) + 1, as in GONEMAXMIN.

Definition 3. Let n ∈ N and r ∈ N≥2. The G-LOTZ function f(x) = (G-LO(x), G-TZ(x)) : [0..r - 1]n -→ R2 is defined by

For r = 2, this definition coincides with the usual binary LOTZ benchmark. The following lemma characterizes the Pareto-optimal solutions and Pareto front of this benchmark.

Lemma 4. For each k ∈ [0..n(r-1)], write q = ⌊k/(r-1)⌋ and a = k - q(r - 1). Let

The Pareto set of G-LOTZ is the set of solutions

and its Pareto front is

F ∗ G-LOTZ = {(k, n(r - 1) - k) | k ∈ [0..n(r - 1)]}. Hence both have size n(r - 1) + 1.

4 Runtime Analysis

4.1 Mathematical Tools

Before presenting the mathematical runtime analyses, we first introduce several tools that will be used later. To bound the objective value of the initial solution, we employ the following Chernoff bound.

Theorem 5 ([Chernoff, 1952]). Let X1, . . . ,Xn be independent random variables taking values in [0, r - 1]. Let X = ∑ni=1 Xi. Let δ ∈ [0, 1]. Then

The following tail bound for sums of independent geometric random variables provides strong guarantees for analyzing population dynamics in our lower-bound analysis.

Theorem 6 ([Witt, 2014]). Let k ∈ N≥1, and let {Di}i∈[k] be independent geometric random variables with respective positive success probabilities (pi)i∈[k]. Let T ∗ := ∑i∈[k] Di, s := ∑i∈[k] 1 p2i , and pmin := min{pi | i ∈ [k]}. Then for all λ ∈ R≥0, we have

Drift analysis is a fundamental tool in the runtime analysis of randomized search heuristics [He and Yao, 2001]. The multiplicative drift theorem [Doerr et al., 2012] applies when the expected progress is proportional to the current distance to the target.

Theorem 7 ([Doerr et al., 2012]). Let X(0) ,X(1) , . . . be a random process taking values in S := {0} ∪ [smin,∞) ⊆ R. Assume that X(0) = s0 with probability one. Assume that there is a δ > 0 such that for all t ≥ 0 and all s ∈ S with Pr[X(t) = s] > 0 we have

Then T := min{t ≥ 0 | X(t) = 0} satisfies

4.2 Upper Bound for the Runtime of the SEMO on G-ONEMAXMIN

We first bound the size of a pair-wise non-dominated set w.r.t. G-ONEMAXMIN.

Lemma 8. Let P ⊆ [0..r - 1]n be any set of solutions that none weakly dominates the other w.r.t. G-ONEMAXMIN. Then |P | ≤ nr - n + 1.

The following theorem provides an upper bound on the expected number of iterations until Algorithm 1 covers the whole Pareto front of G-ONEMAXMIN. In contrast to the binary ONEMINMAX problem, where the objective value directly determines the number of flippable bits, solutions with the same objective value in our setting have a much more diverse structure. Hence a key part of the proof is to estimate the number of entries of a solution that can be mutated to derive an improvement in one objective. More precisely, once the population contains an objective value, a neighboring value can be created by selecting a corresponding individual and mutating one suitable non-boundary coordinate in the correct direction. Summing the resulting waiting times over all missing front values yields the upper bound.

Theorem 9. Consider Algorithm 1 optimizing GONEMAXMIN. Then the expected number of iterations until the population covers the whole Pareto front is O(n2 r2 log n).

4.3 Lower Bound for the SEMOAlgorithm

To derive a lower bound on the expected runtime of the SEMO, a technical difficulty is controlling the population size before the Pareto front is fully covered. Following the proof strategy introduced in [Doerr et al., 2025a], we consider a delayed SEMO that simplifies the analysis in this phase. The delayed SEMO is identical to the original SEMO (Algorithm 1) except for the parent selection in line 3, where it proceeds as follows. In each iteration, a value i ∈ [0..nr - n] is chosen uniformly at random. If there is no individual z ∈ P with f1(z) = i, the iteration is skipped; otherwise, the unique individual x with f1(x) = i is selected and the algorithm proceeds exactly as in line 4 of Algorithm 1. This modification may introduce iterations that do not change the state of the algorithm, but it has no influence otherwise. Hence it can only slow down the process. In particular, structural statements like the shape of the population at a particular stopping time are not affected. The advantage of working with the delayed SEMO is the more regular way parents are selected: In each iteration, the probability to choose a parent with some existing objective value is exactly the reciprocal of the size of the Pareto front. Since we will introduce another variant of the SEMO algorithm later on, we present the pseudocode of the delayed SEMO algorithm here to clearly distinguish and refer to the different algorithms. See Algorithm 2.

Algorithm 2: The delayed r-valued SEMO with unitstrength local mutation to maximize f : [0..r-1]n → Rm.

1 Initialize x ∈ [0..r -1]n uniformly at random and set P0 ← {x}; 2 for t = 0, 1, 2, . . . do 3 Choose i ∈ [0..nr -n] uniformly at random; 4 if there is no x ∈ Pt such that f1(x) = i then 5 Continue; 6 else 7 Select one individual x such that f1(x) = i; 8 Generate y via choosing one entry of x uniformly at random and modifying it by ±1 with equal probability; 9 if there is no z ∈ Pt such that z ≻ y then 10 Pt+1 ← {z ∈ Pt | y ⪰̸ z} ∪ {y}; 11 else 12 Pt+1 ← Pt;

We use Algorithm 2 for analysis until the size of the population has grown to linear order Θ(n(r - 1)) and then switch back to the original SEMO to then analyze the time required to cover the Pareto front. In Lemma 10, we show that the size of the population grows to Θ(n(r-1)) within O(n2 (r-1)2 ) iterations with high probability. The key idea of the proof is that as long as the population is still small, there always exists a missing neighboring objective value, which can be generated with a constant probability conditioned on selecting a suitable individual. Repeatedly filling such gaps, the size of the population increases steadily until it reaches the desired order. This result provides a lower bound on the population size that will be crucial in the final lower-bound analysis.

Lemma 10. Consider Algorithm 2 optimizing GONEMAXMIN. Then with probability 1 - exp(-Θ(n)), the population contains at least n(r - 1)/130 individuals after at most 8n2 (r - 1)2 /129 iterations.

For a population P , we call the two extreme objective vectors (0, n(r - 1)) and (n(r - 1), 0) the Pareto borders. We measure the distance of P to these borders by

Hence dPF (P) = 0 is necessary for covering the full Pareto front. The following lemma shows that, with high probability, the distance of the population to the Pareto borders is still considerable within n2 (r - 1)2 /16 iterations. The proof idea is as follows. By a Chernoff bound, the random initial individual lies in a central region of the objective space with high probability, and is therefore far from both Pareto borders. Afterwards, unit-strength local mutation can move the population towards a fixed Pareto border only gradually, since each mutation can change the relevant objective value by at most one. Moreover, such progress requires selecting an individual with the minimum value in the corresponding objective. Combining these with a concentration argument for the resulting waiting times shows that the population cannot approach either Pareto border too quickly.

Lemma 11. Consider Algorithm 2 optimizing GONEMAXMIN. Then with probability at least 1-exp(-Θ(n)), in all iterations t ∈ [0..n2 (r-1)2 /16], the distance of the population to the Pareto borders is at least n(r - 1)/8.

As mentioned above, we switch back to the analysis of the original SEMO algorithm. Conditioning on a sufficiently large population but a considerable distance to the Pareto borders, we analyze the expected time required to reach the Pareto borders. Since reaching one of the borders is necessary for covering the full Pareto front, it suffices to analyze the time needed to reach one fixed border, say the border with first objective value 0; the other case is symmetric. The proof splits this progress into a far-from-border phase and a closeto-border phase. In both phases, we upper-bound the probability of decreasing the remaining distance by one, which yields lower bounds on the required waiting times.

Lemma 12. Consider Algorithm 1 optimizing GONEMAXMIN, starting with a population size of at least n(r-1) 130 and the distance to the Pareto borders of at least n(r - 1)/8. Then, the expected number of iterations until the population covers the whole Pareto front is Ω(n2 r(r + log n)).

Combining the above three lemmas, we derive the following theorem.

Theorem 13. Consider Algorithm 1 optimizing GONEMAXMIN. Then the expected number of iterations until the population covers the whole Pareto front is Ω(n2 r(r + log n)).

4.4 A Tighter Bound for a Variant of the SEMO Algorithm

We note that the two bounds proven so far are not absolutely tight, but differ by a factor of O(logn). Due to the complicated population dynamics in a run of an MOEA, tight runtime bounds are generally hard to obtain, and in fact, the only ones in this field we are aware of are those in [Laumanns et al., 2002; Doerr and Zheng, 2021; Bossek and Sudholt, 2024; Doerr and Qu, 2023a; Doerr et al., 2025a; Opris, 2026b]. To try to shed more light on the true runtime of the multi-valued SEMO, we now consider a mild modification of the SEMO. This variant is not introduced as an artificial algorithmic improvement. It differs from the original algorithm only in the tie-breaking rule in line 5 of Algorithm 1, where now we accept only offspring that are not even weakly dominated. In other words, if a solution with some objective value is already present in the population, then it is not replaced by a later solution with equal objective value. Since equal objective vectors represent the same point on the Pareto front, this modification affects only the representative solution stored for a front point, but not which objective vectors are already covered. We note that this is exactly the algorithm proposed in the original work introducing the SEMO [Laumanns et al., 2002], and only later the research community gradually adopted the version presented in Algorithm 1. To avoid any ambiguity, the pseudocode of the algorithm regarded in this subsection can be found in Algorithm 3.

Algorithm 3: The variant of the r-valued SEMO that does not accept solutions which are weakly dominated. As before, we use unit-strength local mutation and maximize f : [0..r -1]n →Rm.

1 Initialize x ∈ [0..r -1]n uniformly at random and set P0 ← {x}; 2 for t = 0, 1, 2, . . . do 3 Select one individual x uniformly at random from Pt; 4 Generate y via choosing one entry of x uniformly at random and modifying it by ±1 with equal probability; 5 if there is no z ∈ Pt such that z ⪰ y then 6 Pt+1 ← {z ∈ Pt | y ⪰̸ z} ∪ {y}; 7 else 8 Pt+1 ← Pt;

For Algorithm 3, we can prove an upper bound that agrees with the previously shown lower bound (which is valid here as well). The proof builds on an idea from [Doerr et al., 2018], where a potential function was used to analyze how randomized local search optimizes a variant of the r-valued single-objective ONEMAX problem. In our setting, the modified tie-breaking rule removes neutral replacements. Once a front point is covered by the population, its representative can no longer be changed by another solution with the same objective value. The analysis does not have to track neutral replacements among the many solutions with the same objective value. This is in sharp contrast to Algorithm 1, where such replacements may change the internal structure of a front point without changing its objective value. These neutral changes appear harmless algorithmically, but they obscure the drift of natural potential functions. For Algorithm 3, we can therefore define a potential function that decreases only when real progress towards one of the two extreme Pareto-optimal solutions, 0n and (r -1)n , is made. A multiplicative-drift argument then gives the desired upper bound. Since we do not believe that the small difference of accepting or not solutions with equal objective value changes the asymptotic runtime (which will be supported by our experiments), we see this result as a strong indication for the true runtime of the r-valued SEMO.

Theorem 14. Consider optimizing G-ONEMAXMIN via Algorithm 1 with the modification that it accepts only strictly better solutions. Then the expected number of iterations until the population covers the whole Pareto front is Θ(n2 r(r + log n)).

5 Experiments

Above, we used Algorithm 3 as a model for the SEMO which we could analyze well enough to obtain asymptotically tight bounds for the runtime. Our strong belief was that the small difference in the selection procedure has little influence on the performance. We now support this belief with experiments on both G-ONEMAXMIN and G-LOTZ. The reason for including G-LOTZ is that G-ONEMAXMIN has the special property that every solution is Pareto optimal. Consequently, neutral replacements are particularly frequent and may appear more harmless than in problems with dominated regions. GLOTZ removes this special feature while keeping the same Pareto-front size, and therefore provides a check for our comparison of the two tie-breaking rules.

For G-ONEMAXMIN, we conduct two sets of experiments. In the first set, we fix r = 4 and consider the problem sizes n ∈ {20, 40, 60, 80, 100}; in the second, we fix n = 100 and consider the r-values r ∈ {2, 3, 4, 5}. For each setting, we conduct 100 independent runs of both algorithms and record the first iteration after which the population covers the full Pareto front. We repeat the same experimental setup for G-LOTZ, which has the same front size n(r - 1) + 1 but only one Pareto-optimal solution for each objective vector. This allows us to test whether our empirical conclusions are robust beyond the special G-ONEMAXMIN case, where all solutions are Pareto optimal.

Figures 1 and 2 display the observed runtimes (averages and standard deviations). The results indicate not greater differences between the original SEMO algorithm (Algorithm 1) and its variant Algorithm 3; possibly the original SEMO performs slightly better. Figures 3 and 4 show the corresponding G-LOTZ results. Again Algorithm 3 behaves similarly to the original SEMO algorithm. Across both benchmarks and both parameter regimes, we observe no systematic effect of accepting or rejecting neutral replacements. Hence our experiments support our belief that the small variation in the selection does not drastically change the algorithm behavior, and more specifically, that our improved upper bound on the runtime shown for Algorithm 3 is likely to hold also for the SEMO.

n

r

n

r

Figure 4: The mean (with standard deviations) number of iterations of Algorithms 1 and 3 for solving G-LOTZ with n = 100 and r ∈ [2..5] in 100 independent runs.

6 Conclusion

This paper conducted the first mathematical runtime analyses of MOEAs for multi-valued decision variables. We extended the classic ONEMINMAX problem and the SEMO multiobjective algorithm from the binary setting to finite multivalued domains. We proved an upper bound of O(n2 r2 log n) on the expected number of iterations required to compute the entire Pareto front. We also proved a near-tight lower bound of Ω(n2 r(r+log n)). Both agree with the known Θ(n2 log n) runtime for the binary case.

These results suggest that existing MOEAs, when equipped with natural extensions of the mutation operator to multi-valued domains, perform well on r-valued problems. The moderate increase in runtime as r grows is to be expected, among others, because the size of the Pareto front of our optimization problem grows linearly with r.

In contrast, the current mathematical analysis techniques appear insufficient to completely understand such multivalued optimization processes. At present, we can determine the true asymptotic order of magnitude of the runtime only for the variant of the SEMO that accepts only strictly dominating solutions. Here we prove a runtime of Θ(n2 r(r + log n)), showing a linear increase of the runtime with r when r is at most logarithmic in n. Both from our understanding of the two algorithms and our experiments, we conjecture that both algorithms have the same asymptotic runtime on the G-ONEMAXMIN benchmarks, which then would be Θ(n2 r(r + log n)). Proving such a result formally, however, appears to be challenging and might need considerably new analysis methods.

Acknowledgments

This work was supported by Guangdong Basic and Applied Basic Research Foundation (Grant No. 2025A1515011936), National Natural Science Foundation of China (Grant No. 62306086), and Xinjiang Tianshan Innovative Research Team (2025D14009). This research benefited from the support of the FMJH Program PGMO.