Ringmaster LMO: Asynchronous Linear Minimization Oracle Momentum Method
arXiv:2605.18174v1 [cs.LG] 18 May 2026
Abdurakhmon Sadiev KAUST [email protected] Ivan Ilin KAUST [email protected]
Artavazd Maranjyan KAUST [email protected] Peter Richtárik KAUST [email protected]
Abstract Muon has recently emerged as a strong alternative to AdamW for training neural networks, with encouraging large-scale pretraining results and growing evidence that matrix-structured updates can be faster in practice. Yet Muon, and more generally Linear Minimization Oracle (LMO) based methods, are typically used synchronously. This is problematic in heterogeneous distributed systems, where workers complete gradient computations at different speeds and synchronous training must repeatedly wait for slower workers. In this work, we introduce Ringmaster LMO, an asynchronous LMO-based momentum method for unconstrained stochastic nonconvex optimization. Our method builds on the delay-thresholding idea of Ringmaster ASGD. For SGD-type methods, Ringmaster ASGD achieves optimal time complexity by discarding overly stale gradients. Ringmaster LMO extends this mechanism to general LMO-based updates. We establish convergence guarantees under generalized (L0 , L1 )-smoothness and further develop a parameter-agnostic variant with decreasing stepsizes and adaptive delay thresholds. Finally, we translate our iteration guarantees into time complexity bounds under heterogeneous worker computation times. In the classical Euclidean smooth setting, these bounds recover the optimal time complexity of Ringmaster ASGD. Experiments on stochastic quadratic problems and NanoChat language-model pretraining show that the advantages of Ringmaster LMO grow with system heterogeneity and that the method outperforms strong synchronous and asynchronous baselines.
1
Introduction
For more than a decade, Adam [Kingma and Ba, 2014] has been the default choice for training large neural networks, with its decoupled weight-decay variant AdamW [Loshchilov and Hutter, 2019] more recently becoming the default variant. For matrix-valued layers, however, these methods treat each weight matrix as a vector, rather than explicitly using its matrix structure. Recently, Jordan et al. [2024b] proposed the Muon optimizer, showing that using the matrix structure of hidden-layer weights can be beneficial. Muon forms the update direction from the matrix gradient via approximate orthogonalization. This matrix-structured update has led to strong empirical results, including the nanoGPT speedrun where Muon improves on the AdamW baseline [Jordan et al., 2024a]. Recent pretraining benchmarks also report that Muon can be substantially more compute-efficient than AdamW [Shah et al., 2025, Liu et al., 2025, Wen et al., 2026]. More generally, Muon fits into Code: github.com/vectozavr/ringmaster-lmo.
Preprint.
Table 1: High-level comparison with the closest theoretical works. All listed methods can be viewed through the LMO-based update perspective; the norm column reports the geometry covered by the available analysis, and ‘Euclidean’ indicates results specialized to normalized SGD. The asynchronous column indicates whether the method is analyzed in the delayed-gradient asynchronous setting, the parameter-agnostic column indicates whether the method avoids tuning the main schedule from problem parameters, and the generalized smoothness column indicates whether the theory covers our generalized (L0 , L1 )-smoothness model (Assumption 2.3). Generalized Parameter(L0 , L1 )Method Norm Asynchronous agnostic smoothness Existing LMO-based methods [Pethick et al., 2025a,b] NSGD-M
[Hübler et al., 2024] Ringmaster Asynchronous NSGD-M [Wu and Luo, 2026] Ringmaster LMO [new]
(Algorithm 1)
Any norm
✘
✘
✔
Euclidean
✘
✔
✔
Euclidean
✔
✘
✘
Any norm
✔
✔
✔
the framework of linear minimization oracle (LMO)-based optimizers introduced by Pethick et al. [2025a]. This viewpoint makes the geometry of the update explicit: the choice of norm determines the update direction, and for matrix parameters the spectral-norm ball recovers the orthogonalized direction used by Muon. Despite this progress, LMO-based optimizers are still typically used within synchronous training pipelines. When several workers compute gradients in parallel, the usual synchronous implementation must wait for the slowest worker before forming the next update. Asynchronous methods, going back at least to HOGWILD! [Recht et al., 2011], avoid this synchronization barrier by updating the model whenever a worker returns a gradient. The price is delay: a gradient arriving at iteration k may have been computed at an older iterate with delay δk ≥ 0. Recent analyses have clarified the effect of such delays for asynchronous SGD, and delay-thresholding variants can achieve optimal time complexity under heterogeneous worker speeds [Tyurin and Richtárik, 2023, Maranjyan et al., 2025d, Maranjyan and Richtárik, 2026]. These results, however, are built for SGD-type updates in Euclidean setting, and do not directly apply to LMO-based optimizers. This leaves the central question of this paper: Can asynchronous LMO-based methods recover the optimal time complexity known for Euclidean setting. We show that this can be done. Building on the delay-thresholding mechanism of Ringmaster ASGD [Maranjyan et al., 2025d] and the LMO-based momentum framework of Pethick et al. [2025a], we introduce an asynchronous LMO-based method for unconstrained stochastic nonconvex optimization. To the best of our knowledge, this is the first asynchronous LMO-based method with time complexity guarantees under heterogeneous worker computation times. In the classical Euclidean smooth setting, the bound recovers the optimal time complexity of Ringmaster ASGD. 1.1
Contributions
Our contributions are the following: • First asynchronous LMO-based method. We propose Ringmaster LMO (Algorithm 1), an asynchronous LMO-based method that combines momentum with delay-thresholding. To the best of our knowledge, this is the first asynchronous LMO-based method with time complexity guarantees under heterogeneous worker computation times. The method extends the Ringmaster ASGD delay-control principle from SGD-type updates to LMObased updates; see Section 3. 2
• Convergence under generalized smoothness. The analysis in Sections 4.1 and 4.2 establishes convergence under generalized (L0 , L1 )-smoothness [Zhang et al., 2020b,a], which includes the classical smooth case as L1 = 0. The resulting rates recover the standard stochastic LMO-type dependence while allowing delayed gradients generated by asynchronous workers; see Theorems 4.1 and 4.4. • Parameter-agnostic scheduling. Section 4.2 gives a variant with decreasing stepsizes and adaptive delay thresholds. This removes the need to tune the main schedule using problem parameters and, in particular, provides a parameter-agnostic extension that was not present in the original Ringmaster ASGD analysis; see Theorem 4.4. • time complexity guarantees. Theorems 4.3, 4.6 and E.6 translate the iteration complexity into time complexity under both the fixed computation model and the universal time-varying computation model. In the standard smooth Euclidean setting, these bounds recover the same optimal time complexity as Ringmaster ASGD. • Experiments. Section 5 compares the proposed method against synchronous and asynchronous baselines. On both stochastic quadratic problems and NanoChat pretraining, the proposed method is competitive in mild heterogeneity and its advantage grows as workerspeed heterogeneity increases. Table 1 places our method in the context of the closest theoretical works. 1.2
Related work
We review the prior work most closely related to our setting. LMO-based and matrix-structured optimizers. The empirical success of Muon [Jordan et al., 2024b,a] has motivated a line of work seeking to explain and generalize its matrix-structured update. Bernstein and Newhouse [2024] relate Muon to older adaptive and preconditioned methods, including Shampoo [Gupta et al., 2018], by emphasizing the role of the underlying norm. Pethick et al. [2025a] give a broader optimization framework in which the update direction is chosen by a linear minimization oracle over a norm ball. LMOs are classical in Frank-Wolfe and conditional-gradient methods [Frank and Wolfe, 1956, Hazan, 2008, Clarkson, 2010, Jaggi, 2013], but the Scion framework uses them for unconstrained stochastic optimization and recovers several normalized or structuredgradient methods through the choice of norm. In parallel, Kovalev [2025] interprets gradient orthogonalization through non-Euclidean trust-region optimization. More broadly, the concept of spectral descent and orthogonalized updates has historical roots in deep learning [Carlson et al., 2015a,b, 2016, Tuddenham et al., 2022]. Subsequent work has refined this theory under generalized smoothness [Riabinin et al., 2025, Pethick et al., 2025b], incorporated second-order momentum variants [Khirirat et al., 2025], analyzed inexact Newton-Schulz orthogonalization [Shulgin et al., 2026], and adapted the LMO viewpoint to layerwise neural-network training [Riabinin et al., 2025]. These works are synchronous: they do not analyze delayed gradients produced by asynchronous workers, nor do they provide time complexity guarantees under heterogeneous computation times. Asynchronous stochastic optimization. Asynchronous stochastic optimization removes the synchronization barrier by updating the model whenever a worker returns a gradient [Tsitsiklis et al., 1986, Agarwal and Duchi, 2011, Recht et al., 2011, Maranjyan, 2025]. The main analytical difficulty is delay: the gradient used at iteration k may have been computed at an earlier iterate xk−δk . Classical analyses control this error through bounded-delay assumptions [Feyzmahdavian et al., 2016, Lian et al., 2015, Arjevani et al., 2020], while more recent work uses delay-adaptive stepsizes [Mishchenko et al., 2022, Koloskova et al., 2022, Wu et al., 2022]. For heterogeneous worker speeds, the relevant objective is time complexity rather than only counting iterations. Tyurin and Richtárik [2023] identified optimal time complexities under a fixed computation model. This lower-bound perspective was later broadened to cover time-varying worker behavior [Tyurin, 2024] as well as variance-reduced parallel stochastic optimization [Tovmasyan et al., 2026]. Ringmaster ASGD [Maranjyan et al., 2025d] achieves these guarantees for asynchronous SGD by accepting only gradients whose delay is below a threshold. Maranjyan and Richtárik [2026] and Mahran et al. [2026] address the corresponding gap for data- and system-heterogeneous settings. These results are Euclidean SGD-type methods; they do not cover LMO-based directions or the non-Euclidean geometries needed to recover Muon-like updates. 3
Closest baselines and parameter agnosticism. The closest asynchronous baseline to our method is Wu and Luo [2026], which combines the delay-thresholding mechanism of Maranjyan et al. [2025d] with normalized SGD and momentum under heavy-tailed noise. Their method corresponds to the Euclidean-norm case of our update and therefore does not cover general LMOs, matrix-structured directions, or the (L0 , L1 )-smoothness model used here. Generalized (L0 , L1 )-smoothness was introduced as a way to capture losses whose local smoothness can grow with the gradient norm [Zhang et al., 2020b,a], and it is also the setting of recent synchronous LMO analyses [Riabinin et al., 2025, Pethick et al., 2025b, Khirirat et al., 2025]. Separately, Hübler et al. [2024] develop parameter-agnostic normalized SGD with momentum under relaxed smoothness in the synchronous Euclidean setting. Our parameter-agnostic variant brings this idea to asynchronous LMO-based optimization: the same schedule controls both the momentum averaging and the delay threshold, removing the need to tune the threshold from problem parameters. Overall, our work connects these three lines by extending delay-thresholding from Euclidean SGD to general LMO-based momentum methods and by translating the resulting iteration bounds into time complexity guarantees.
2
Preliminaries and setup
This section fixes the optimization model, the LMO geometry, the distributed computation model, and the assumptions used in the analysis. 2.1
Optimization model and LMO geometry
We study the nonconvex stochastic optimization problem minimize {f (x) := Eξ∼D [f (x; ξ)]} , x∈Rd
where x ∈ Rd is the model parameter, D is the data distribution, ξ ∼ D is a data sample, and f (x; ξ) is the corresponding loss. Let ∥·∥ be a norm on Rd . Following Pethick et al. [2025a], the update direction is defined through a linear minimization oracle (LMO) over the unit ball of this norm: lmo(y) ∈ arg min ⟨y, u⟩ u : ∥u∥≤1
for y ̸= 0,
lmo(0) := 0.
When the minimizer is not unique, any minimizer may be chosen. This definition implies that the magnitude of the update is always bounded because ∥lmo(y)∥ ≤ 1. 2.2
Function and oracle assumptions
We use the following assumptions throughout the convergence analysis. Assumption 2.1 (Lower boundedness). There exists f ∗ > −∞ such that f (x) ≥ f ∗ for all x ∈ Rd . We define ∆0 := f (x0 ) − f ∗ , where x0 is the starting point of the optimization methods. Assumption 2.2 (Unbiased stochastic gradients with bounded variance). For every ξ, the function f (x; ξ) is differentiable with respect to its first argument x. Moreover, the stochastic gradients are unbiased and have bounded variance σ 2 ≥ 0, that is, Eξ∼D [∇f (x; ξ)] = ∇f (x), ∀x ∈ Rd , h i 2 Eξ∼D ∥∇f (x; ξ) − ∇f (x)∥2 ≤ σ 2 , ∀x ∈ Rd .
Assumption 2.3 ((L0 , L1 )–smoothness). There exist constants L0 , L1 ≥ 0 such that, for all x, y ∈ Rd , ∥∇f (x) − ∇f (y)∥∗ ≤ sup (L0 + L1 ∥∇f (θx + (1 − θ)y)∥∗ ) ∥x − y∥ . θ∈[0,1]
Here ∥·∥∗ denotes the dual norm associated with ∥·∥, i.e. ∥y∥∗ :=
sup u : ∥u∥≤1
⟨y, u⟩ .
Assumption 2.3 reduces to the standard L0 -smoothness condition with respect to ∥·∥ when L1 = 0. When L1 > 0, the effective local smoothness constant is allowed to grow with the gradient norm, 4
which motivates generalized smoothness models for objectives arising in modern machine learning [Chen et al., 2023, Pethick et al., 2025b]. Assumption 2.2 controls stochastic-gradient noise in the Euclidean norm, whereas Assumption 2.3 is stated in the norm induced by ∥·∥ and its dual. We denote by ∥z∥∗ z̸=0 ∥z∥2
ρ := sup
the corresponding norm-equivalence constant, so that ∥z∥∗ ≤ ρ ∥z∥2 for all z ∈ Rd . 2.3
Distributed computation model
We consider a distributed setting with n workers. We write [n] := {1, . . . , n}. Each worker i has access to the same data distribution D. To model heterogeneous worker speeds, we use the fixed computation model of Mishchenko et al. [2022]. We extend the analysis to arbitrarily changing computation times in Section E, while keeping the fixed model in the main text for a clearer presentation. Assumption 2.4 (Fixed computation model). For each worker i ∈ [n], computing a stochastic gradient ∇f (x; ξ) takes τi > 0 seconds. Without loss of generality, we assume τ1 ≤ τ2 ≤ · · · ≤ τn . For m ∈ [n], define !−1 m 1 X 1 Hm := , m i=1 τi the harmonic mean of the first m computation times τ1 , . . . , τm . Since τ1 ≤ · · · ≤ τn , the sequence Hm is non-decreasing in m. Throughout this paper, we assume zero-latency communication between the workers and the server, adopting the standard modeling assumption used in prior asynchronous optimization literature [Mishchenko et al., 2022, Koloskova et al., 2022, Tyurin and Richtárik, 2023, Maranjyan et al., 2025d]. While communication costs are undeniably important in practice, modeling them properly requires the adoption of specific mechanisms like compression or local steps. Incorporating those techniques introduces a separate set of theoretical challenges that are largely orthogonal to our primary objective: understanding asynchronous dynamics in the presence of heterogeneous computation speeds.
3
Ringmaster LMO method
We now describe the asynchronous LMO-based method. A central server maintains the iterate and momentum, while n workers repeatedly compute stochastic gradients. Whenever worker i finishes, it returns a gradient computed at the iterate it last received from the server. If the server accepts that gradient, it first updates the momentum and then applies the LMO direction: (1 − α0 )m0 + αinit g0 , k = 0, mk+1 = (1 − αk )mk + αk gk , k ≥ 1, (1) xk+1 = xk + ηk lmo(mk+1 ). ik Here gk = ∇f xk−δk ; ξk−δ is the gradient that arrives from worker ik at server iteration k, δk is its k delay, αk is the momentum parameter, and ηk is the stepsize. Because gk may have been computed at the stale iterate xk−δk rather than at xk , large delays can degrade the update. To control this effect, we use the delay-thresholding principle of Ringmaster ASGD [Maranjyan et al., 2025d]: the server accepts an arriving gradient only if δk < Rk . The threshold sequence Rk controls the freshnessthroughput tradeoff: smaller thresholds are more conservative, while larger thresholds allow more asynchronous parallelism. Allowing Rk to vary with k will be important for the parameter-agnostic variant in Section 4.2. Algorithm 1 gives the full description of the method. Algorithm 1 discards stale gradients only after they arrive, which may appear to waste computation. In practice, the server can often determine in advance that an in-flight computation will miss the current threshold and interrupt that worker early instead of waiting for the stale gradient to finish; 5
Algorithm 1 Ringmaster LMO 1: Input: initial iterate x0 ∈ Rd ; step-size sequence ηk > 0; delay-threshold sequence Rk > 0;
momentum parameter sequence αk ∈ (0, 1] with α0 = 1 and αinit ∈ (0, 1]; number of model updates K ∈ N 2: Each worker i ∈ [n] starts computing ∇f (x0 ; ξ0i ) in parallel 3: k ← 0 4: while k < K do ik 5: Receive gk = ∇f xk−δk ; ξk−δ from worker ik ▷ δk is the delay of the arriving gradient k 6: if δk < Rk then (1 − α0 )m0 + αinit g0 , k = 0, 7: mk+1 ← ▷ m0 need not be set since α0 = 1 (1 − αk )mk + αk gk , k ≥ 1. 8: xk+1 ← xk + ηk lmo(mk+1 ) 9: The server sends xk+1 to worker ik ik 10: Worker ik starts computing ∇f xk+1 ; ξk+1 11: k ←k+1 12: else 13: Discard the stale gradient gk 14: The server sends xk to worker ik 15: Worker ik starts computing ∇f xk ; ξkik 16: end if 17: end while
this computation-stopping implementation was also proposed in Ringmaster ASGD [Maranjyan et al., 2025d]. One can go further when worker runtimes are random rather than fixed [Maranjyan et al., 2025a]: the runtime distributions can be learned online and used to make better assignment or interruption decisions on the fly [Maranjyan et al., 2025b]. These systems-level refinements are complementary to our main objective, so we keep the simpler fixed-time model with delay-based acceptance and rejection in the main text. Special cases. If every accepted gradient is fresh, then (1) reduces to the synchronous LMO-based momentum method [Pethick et al., 2025a, Kovalev, 2025]. In Algorithm 1, this is enforced by setting Rk = 1 for all k, which means that only gradients with δk = 0 are accepted. If the LMO is defined with respect to the Euclidean norm, then the update reduces to normalized SGD with momentum [Hazan et al., 2015, Cutkosky and Mehta, 2020]: xk+1 = xk − ηk
mk+1 ∥mk+1 ∥2
instead of (1). In this case, Algorithm 1 reduces to the method studied by Wu and Luo [2026]. Their analysis, however, is restricted to the Euclidean setting and classical smoothness (L1 = 0), and it does not cover the parameter-agnostic threshold schedule studied here. Other induced methods. Because the method is defined through an LMO, changing the underlying norm immediately yields asynchronous variants of other normalized optimizers, as discussed by Pethick et al. [2025a]. For example, choosing the spectral norm recovers Muon [Jordan et al., 2024b], while choosing the max norm ∥·∥∞ recovers SignSGD or Signum [Bernstein and Newhouse, 2024]. See Pethick et al. [2025a] for a more comprehensive discussion of which norm choices recover which methods.
4
Theoretical analysis
We begin with the fixed-parameter regime before turning to the parameter-agnostic case. 4.1
Fixed-parameter version
We first state the corresponding iteration-complexity guarantee of Ringmaster LMO (Algorithm 1). 6
Theorem 4.1 (Fixed-parameter iteration complexity; proof in Section B). Suppose Assumptions 2.1– 2.3 hold. Run Algorithm 1 with α0 = 1. Set the momentum parameter to the constant α for all k ≥ 1, and use the constant delay threshold R for all k ≥ 0, where ) ( 1/2 1/2 1 ∆0 L0 , ∀k ≥ 1, Rk = R = , ∀k ≥ 0, αk = αinit = α = min 1, α ρσK 1/2 with the convention that, when σ = 0, we take α = 1, and choose the constant stepsize as follows: • If L1 = 0, take (r ηk = η = min
3/4
∆0 ∆0 , 1/4 L0 K L (ρσ)1/2 K 3/4
) ,
∀k ≥ 0.
3/4
)
0
• If L1 > 0, take (r ηk = η = min
1/2
1/2
∆0 ∆0 1 ∆0 L 0 , , , L0 K 8L1 8L1 ρσK 1/2 L1/4 (ρσ)1/2 K 3/4 0
,
∀k ≥ 0.
Then, after K ≥ 1 iterations, the iterates of Algorithm 1 satisfy ! r 1/2 1/4 1/4 L 1 ∆0 L1 ∆0 ρσ L0 ∆0 (ρσ)1/2 L0 ∆0 ρσ + + 1/2 min E [∥∇f (xk )∥∗ ] ≤ O + 1/2 . + K K k∈{0,...,K−1} K 1/4 K L K 1/2 0
In the case L1 = 0, the terms containing L1 in the preceding bound are zero. Consequently, the number of iterations needed to reach an ε-stationary point satisfies L 1 ∆0 L0 ∆0 (ρσ)2 L21 ∆0 (ρσ)2 (ρσ)2 L0 ∆0 Kε = O + + + + . ε2 ε ε4 L0 ε2 ε2 Again, when L1 = 0, the terms containing L1 are omitted. Discussion In the standard smooth Euclidean case (ρ = 1 and L1 = 0), our rate matches the lower bounds proved by Arjevani et al. [2023] for first-order stochastic smooth non-convex optimization. When L1 > 0, the deterministic part of our complexity matches the best known rates under generalized smoothness and agrees with lower bounds for gradient descent [Zhang et al., 2020b, Crawshaw et al., 2022]1 . GGNC [Pethick et al., 2025b] has the same complexity in the non-stochastic case. In the stochastic case, however, our bounds are tighter because of the optimal parameter selection. Gluon [Riabinin et al., 2025] has the same deterministic rate when each layer shares the same generalized smoothness constants. However, their analysis assumes asymmetric generalized smoothness, which is less general than symmetric generalized smoothness [Chen et al., 2023]. Moreover, our stochastic complexity bounds are tighter than those of GGNC and Gluon. In the ℓ1 -norm case, our bounds match those of SignSGD [Crawshaw et al., 2022], although that work assumes coordinate-wise generalized smoothness. Khirirat et al. [2025] proved better rates O(ε−3 ) for LMO-MVR and LMO-SOM, but their algorithms assume access either to multiple stochastic gradient queries with shared randomness or to Hessian-vector products. These assumptions are restrictive and, especially in the former case, difficult to adapt to the asynchronous setting. Finally, in the bounded-variance case, our method recovers the rate of RANSGD-M from Wu and Luo [2026], since in the Euclidean setting our method reduces to RANSGD-M. Time complexity The time complexity under Assumption 2.4 follows directly from Maranjyan et al. [2025d]. Their Ringmaster ASGD method and our Ringmaster LMO differ only in how they form updates; this changes the iteration complexity, but the same time complexity analysis applies. More concretely, we need the following lemma. 1 These comparisons are limited to algorithm-specific lower bounds, since lower bounds for the full class of generalized
(L0 , L1 )-smooth functions remain open.
7
Lemma 4.2 (Duration of R updates [Maranjyan et al., 2025d, Lemma 4.1]). Let the delay threshold of Algorithm 1 be fixed as Rk ≡ R. Under Assumption 2.4, the time needed to complete any R consecutive updates is at most R t(R) := 2 min Hm 1 + . m m∈[n] With this lemma we can get the time complexity of our method with the following theorem. Theorem 4.3 (Fixed-parameter total time complexity; proof in Section D.1). Under the assumptions and parameter choices of Theorem 4.1, under Assumption 2.4, the time complexity of reaching an ε-stationary point satisfies √ L1 ∆0 L0 ∆0 L21 ∆0 (ρσ)2 (ρσ)2 L0 ∆0 L0 ∆0 (ρσ)2 + + + . Tε = O min Hm + + ε2 ε ε mε4 mL0 ε2 mε2 m∈[n] Proof sketch. Let Kε denote the iteration bound from Theorem 4.1, and recall that in the fixedparameter regime the delay threshold is constant with Rk ≡ R = 1/α. We partition the first Kε updates into ⌈Kε/R⌉ consecutive blocks of length at most R and apply Lemma 4.2 to each block. Therefore, Tε ≤ ⌈Kε/R⌉ t(R). Substituting the choices of α and the bound on Kε from Theorem 4.1 yields the claimed time complexity bound. Discussion In the standard smooth Euclidean case (ρ = 1 and L1 = 0), our bound reduces to the optimal time complexity scaling identified by Tyurin and Richtárik [2023]2 . In particular, it recovers the same optimal time complexity as Ringmaster ASGD [Maranjyan et al., 2025d] and RANSGD-M [Wu and Luo, 2026]. Thus, extending delay thresholding from normalized SGD to general LMO-based updates does not incur any loss in the classical Euclidean regime. 4.2
Parameter-agnostic case
We now turn to the parameter-agnostic case and begin with its iteration complexity. Theorem 4.4 (Parameter-agnostic iteration complexity; proof in Section C). Suppose Assumptions 2.1–2.3 hold. Run Algorithm 1 with Rk = 1/αk , αinit = α0 = 1, and αk = k −1/2 for k ≥ 1. In the first two cases, let η > 0 be any constant and define 2 2
2 2
Ψ(L0 , L1 ) := eL1 η ∆0/η + ρσ + eL1 η L0 η. Choose the stepsizes as follows, and let Ψ denote the corresponding value in the final bound: η ∆0/η + ρσ + L η. • If L1 = 0, take ηk = (k+1) 3/4 and Ψ = Ψ(L0 , 0) = 0 η • If L1 > 0, take ηk = 17(k+1) 3/4 and Ψ = Ψ(L0 , L1 ). 1 L0/L1 , which is the • If L1 > 0 is known, take ηk = 17L (k+1) 3/4 and Ψ = L1 ∆0 + ρσ + 1 preceding envelope with η = 1/L1 up to universal constants.
Then, after K ≥ 1 iterations, min k∈{0,...,K−1}
E [∥∇f (xk )∥∗ ] ≤ O
Ψ log K K 1/4
.
Consequently, the number of iterations needed to reach an ε-stationary point satisfies 4 ! Ψ4 1 Kε = O log . ε4 ε 2 The proved lower bounds hold in the case when 0 < ε ≤ c′
bound from Theorem 4.3.
8
√
L0 ∆0 for some constant c′ > 1, and matches with our
Discussion In the Euclidean case, our parameter-agnostic bound matches that of NSGD-M from Hübler et al. [2024], while working under symmetric generalized smoothness rather than their asymmetric variant; for twice-differentiable objectives, the two notions are equivalent. When L1 is 2 2 unknown, our bound contains the factor eL1 η , which is unavoidable and is predicted by the lower bounds of Hübler et al. [2024]. The parameter-agnostic analysis of Gluon [Riabinin et al., 2025] does not contain this factor, but it relies on more restrictive assumptions. Compared with existing LMO-based results, our guarantee does not require prior knowledge of either K or L1 . In particular, Khirirat et al. [2025] analyze several LMO-based momentum variants, but their parameter choices depend on the time horizon and generalized smoothness constants. Finally, relative to Wu and Luo [2026], we provide, to the best of our knowledge, the first parameter-agnostic guarantee for an asynchronous method under generalized smoothness, and more broadly the first convergence guarantee for asynchronous optimization in this setting. Time complexity We now analyze the time complexity of Algorithm 1 in the parameter-agnostic case, where the threshold Rk changes over time according to the square-root schedule in the theorem above. This case is slightly more difficult than the fixed-threshold case below. Lemma 4.5 (Time complexity for a square-root delay threshold; proof in Section D.2). Consider √ Algorithm 1 with delay-threshold sequence given in Theorem 4.4, i.e. R0 = 1 and Rk = k for k ≥ 1. Under Assumption 2.4, for any K ≥ 1, the time needed to complete the first K iterations satisfies √ K T (K) = O min Hm K+ . m m∈[n] Proof sketch. We group iterations by the integer threshold level Kr := {k ∈ {0, . . . , K − 1} | ⌊Rk ⌋ = r}. Since Rk is non-decreasing, each Kr is a contiguous block. The fixed-threshold bound from Lemma 4.2 can therefore be applied blockwise, after splitting Kr into groups of at most r updates. √ For the square-root rule R0 = 1 and Rk = k for k ≥ 1, we have |Kr | ≤ 2r + 2, so each√threshold level contributes only a constant number of such groups. Summing over r ≤ ⌊RK−1 ⌋ ≤ K gives the stated bound. Combining Theorem 4.4 with Lemma 4.5 yields the total time complexity for reaching an ε-stationary point. Theorem 4.6 (Total time complexity). Under the assumptions and parameter choices of Theorem 4.4, let Kε be the number of iterations needed to reach an ε-stationary point. Then, under Assumption 2.4, the corresponding time complexity Tε := T (Kε ) satisfies 2 Ψ4 Ψ e Tε = O min Hm + , ε2 mε4 m∈[n] where Hm is the harmonic mean of the first m computation times τ1 , . . . , τm , Ψ is the corresponding e hides logarithmic factors. quantity from Theorem 4.4, and O(·) Proof. By Theorem 4.4,
√
e Kε = O
Tε = O
min Hm
m∈[n]
Ψ2 ε2
. Applying Lemma 4.5 with K = Kε gives
p 2 4 Kε e min Hm Ψ + Ψ Kε + =O . m ε2 mε4 m∈[n]
Discussion. Up to logarithmic factors, this bound has the same worker-heterogeneity dependence as Theorem 4.3, namely minm∈[n] Hm · + m· , so it preserves the delay-thresholding benefit of adapting to the faster workers. Its advantage is parameter agnosticism: the square-root schedule attains this guarantee without knowing K or L1 . To the best of our knowledge, no existing asynchronous method provides a comparable parameter-agnostic time complexity guarantee. 9
5
Experiments
Reproducibility. The implementation and code for reproducing our experiments are available at github.com/vectozavr/ringmaster-lmo. Overview. We test the main empirical prediction of the theory: delay thresholding should be most useful when asynchronous workers have heterogeneous runtimes. Choosing the spectral-norm LMO in Algorithm 1 recovers Muon, so our two tested methods are Ringmaster Muon and its parameteragnostic variant Parameter-agnostic Ringmaster Muon. We compare them with two asynchronous Muon baselines, Rennala Muon and Delay-adaptive Muon, on a stochastic quadratic benchmark and on NanoChat language-model pretraining. Methods and tuning. For all fixed-parameter methods, we use the same internal Muon configuration: momentum coefficient β = 0.95, five Newton–Schulz iterations, and Nesterov lookahead. Thus the tuning only concerns asynchronous control parameters such as the stepsize scale η, delay threshold, or baseline-specific refresh parameter. For Parameter-agnostic Ringmaster Muon, we use the schedule motivated by Theorem 4.4, η ηk = , α0 = 1, αk = k −1/2 (k ≥ 1), Rk = max{1, ⌊1/αk ⌋}, (k + 1)3/4 so the only tuned parameter is the scale η. Throughout the experimental grids, η denotes the tuned stepsize scale; for delay-adaptive baselines, it is the nominal scale before applying the method-specific delay correction. All methods are tuned independently for each benchmark and delay regime under the same simulated runtime budget. Runtime model and compute resources. All experiments use a server–worker simulator instead of a physical multi-GPU deployment. The simulator follows the asynchronous protocol in Algorithm 1: a worker receives the current model, computes one stochastic gradient, returns after its assigned runtime, and is then either accepted or rejected according to the method’s asynchronous rule. The horizontal axis in every plot is this simulated runtime, not the host machine’s elapsed wall-clock time. Indexing workers by increasing deterministic runtime as i = 0, . . . , n − 1, we use the following three delay profiles: √ gilin = 1 + i. (2) gisim = 1, gisub = 1 + i, The first profile is nearly homogeneous, the second introduces moderate heterogeneity, and the third creates a strong straggler tail. We add a small independent perturbation to worker runtimes so that the similar-delay regime is not exactly tied across workers. The synthetic quadratic experiments were run on a MacBook Pro with an Apple M1 Max chip and 64GB RAM using Python 3.12. The NanoChat gradients were computed on one NVIDIA A100-PCIE-40GB GPU using bfloat16 autocast and torch.compile; the multi-worker dynamics were simulated on top of these measured gradient times. Synthetic quadratic benchmark. Our synthetic benchmark is the standard worst-case Nesterov tridiagonal quadratic in dimension d = 1729: 1 ⊤ 1 1 x Ax − b⊤ x, A = tridiag(−1, 2, −1), b = − e1 . 2 4 4 √ We initialize at x0 = d e1 , compute the exact minimizer by solving Ax = b, and report the objective t ⋆ gap f (x ) − f as a function of simulated runtime. The stochastic oracle is f (x) =
g(x, ξ) = ∇f (x) + ξ1, where ξ ∼ N (0, 0.012 ) is a scalar Gaussian perturbation shared across all coordinates. For this benchmark, worker i has deterministic base runtime gi , where gi is chosen from the three profiles in (2). We simulate n = 6174 workers and tune each method over a runtime horizon of 2000 simulated seconds, with one trial per hyperparameter setting. The search grids are • Ringmaster Muon: η ∈ {5−6 , 5−5 , . . . , 51 } and R ∈ {1, 2, 4, 6, 8, 16, 32}; • Parameter-agnostic Ringmaster Muon: η ∈ {5−4 , 5−3 , . . . , 53 }; 10
Nesterov’s Quadratic Muon Async Training (τi ≈ const)
101
Nesterov’s Quadratic Muon Async Training (τi ∝
√
Nesterov’s Quadratic Muon Async Training (τi ∝ i)
i)
Ringmaster Muon: γ = 0.2, R = 2 PA Ringmaster Muon: η = 5 Rennala Muon: γ = 0.2, B = 4 Delay-Adaptive Muon: γ = 0.0016
102
101
f (xt) − f inf
f (xt) − f inf
101
100
100
10−1
100
10−1
10−2
Ringmaster Muon: γ = 0.2, R = 1 PA Ringmaster Muon: η = 5 Rennala Muon: γ = 0.2, B = 1 Delay-Adaptive Muon: γ = 0.04
102
f (xt) − f inf
10
Ringmaster Muon: γ = 0.04, R = 2 PA Ringmaster Muon: η = 1 Rennala Muon: γ = 0.2, B = 32 Delay-Adaptive Muon: γ = 1.28e − 05
2
10−1
10−2 10−2
10−3
10−3 0
250
500
750
1000
1250
1500
1750
2000
0
250
Runtime (seconds)
(a) Nearly homogeneous delays
500
750
1000
1250
1500
1750
2000
0
250
500
750
Runtime (seconds)
1000
1250
1500
1750
2000
Runtime (seconds)
(b) Sublinear delays
(c) Linear delays
Figure 1: Comparison on the stochastic tridiagonal quadratic objective under similar, sublinear, and linear worker-delay regimes. The advantage of Ringmaster Muon appears as the delay heterogeneity becomes stronger. Table 2: Measured stochastic-gradient times for the NanoChat setup used in our experiments: a 6-layer, 192-hidden model with sequence length 2048 and device batch size 4, run on a single GPU with torch.compile enabled and mixed precision where supported. Each entry reports the mean and standard deviation over 100 synchronized gradient steps after 3 warm-up iterations. The table shows substantial cross-device variability, with mean step times ranging from 14.18 ms on an H100 to 215.20 ms on a GTX 1080 Ti, which supports the heterogeneous delay ranges used in our simulated runtime model. GPU Precision Mean time (ms) Std. (ms) NVIDIA H100-PCIE-80GB NVIDIA A100-PCIE-40GB NVIDIA RTX 4090 24GB NVIDIA RTX 3090 24GB NVIDIA V100-PCIE-16GB NVIDIA T4 16GB NVIDIA GTX 1080 Ti 11GB
bf16 bf16 bf16 bf16 fp16 fp16 fp32
14.18 22.74 26.53 38.42 64.15 128.60 215.20
7.32 11.64 13.11 15.10 20.30 35.50 40.15
• Rennala Muon: η ∈ {5−6 , 5−5 , . . . , 51 } and B ∈ {1, 2, 4, 6, 8, 16, 32};
• Delay-adaptive Muon: η ∈ {5−8 , 5−7 , . . . , 5−1 }.
Figure 1 shows the expected tradeoff. In the similar-delay regime, where workers are nearly homogeneous, Ringmaster Muon slightly underperforms the strongest baselines because discarding updates provides little benefit and can waste gradient evaluations. As heterogeneity increases, the value of thresholding becomes clearer: Ringmaster Muon performs best in the sublinear-delay regime and remains highly competitive in the linear-delay regime, where it matches Rennala Muon. NanoChat benchmark. Our language-modeling experiments use an architecture-faithful implementation of NanoChat [Karpathy, 2025]. The model has context length 2048, window pattern SSSL, vocabulary size 8192, 6 transformer blocks, hidden size 192, and grouped-query attention with 3 query heads and 3 key-value heads. Each block uses RMS normalization, rotary positional embeddings, value embeddings with gating, learned scalar coefficients for residual connections, and the MLP nonlinearity x 7→ ReLU(x)2 . The tokenizer is trained once using an 8192-wordpiece BPE vocabulary on the Karpathy climbmix-400b-shuffle shards while preserving the GPT-style regex split pattern of the reference implementation. The data pipeline prepends a BOS token to each document and packs documents into length-2049 streams using a best-fit heuristic, so each stochastic gradient processes fully packed sequences of length 2048. Optimization uses device batch size 4. The optimization variable concatenates all trainable tensors into a single flattened representation. Within this representation, Muon orthogonalization is applied blockwise only to matrix-valued tensors inside transformer blocks; embeddings, output weights, and scalar residual parameters use the standard momentum recursion without orthogonalization. To avoid the overhead of a separate evaluation pass, we report the training loss of the most recently computed minibatch from the worker that last returned a gradient. 11
NanoChat Muon Async Training (¿i ¼ const)
NanoChat Muon Async Training (¿i /
Ringmaster Muon: ° = 1; R = 8
Ringmaster Muon: ° = 1; R = 12
Ringmaster Muon: ° = 1; R = 12
Rennala Muon: ° = 1; B = 2 Delay-Adaptive Muon: ° = 1
6 × 100
4 × 10
0
101
PA Ringmaster Muon: ´ = 5 Rennala Muon: ° = 1; B = 2 Delay-Adaptive Muon: ° = 1
6 × 100
4 × 100 0
100
200
300
400
500
Ringmaster Muon: ° = 1; R = 2
9 × 100
Ringmaster Muon: ° = 1; R = 8
Latest Minibatch Loss
10
NanoChat Muon Async Training (¿i / i)
i)
Ringmaster Muon: ° = 1; R = 4
PA Ringmaster Muon: ´ = 5 1
p
Ringmaster Muon: ° = 1; R = 2
Ringmaster Muon: ° = 1; R = 4
Latest Minibatch Loss
Latest Minibatch Loss
Ringmaster Muon: ° = 1; R = 2
Ringmaster Muon: ° = 1; R = 4 Ringmaster Muon: ° = 1; R = 8
8 × 100
Ringmaster Muon: ° = 1; R = 12 PA Ringmaster Muon: ´ = 5 Rennala Muon: ° = 1; B = 2
7 × 100
Delay-Adaptive Muon: ° = 1
6 × 100
5 × 100
4 × 100 0
100
Runtime (seconds)
200
300
400
500
0
100
Runtime (seconds)
(a) n = 32, nearly homogeneous delays
200
300
400
500
Runtime (seconds)
(b) n = 32, sublinear delays
(c) n = 32, linear delays
Figure 2: NanoChat training loss versus simulated runtime for a 6-layer, 192-hidden, contextlength-2048 model with 32 simulated workers. The three panels correspond to nearly homogeneous, sublinear, and linear delay profiles, with delays anchored to the measured A100 mean gradient time. Ringmaster Muon is competitive in the nearly homogeneous regime and achieves the lowest loss as heterogeneity increases, with the clearest margin in the linear-delay setting. NanoChat Muon Async Training (¿i ¼ const)
NanoChat Muon Async Training (¿i /
i)
Rennala Muon: ° = 1; B = 2 Delay-Adaptive Muon: ° = 1
101
Latest Minibatch Loss
PA Ringmaster Muon: ´ = 5
Delay-Adaptive Muon: ° = 1
4 × 100
6 × 100
4 × 100 100
200
Ringmaster Muon: ° = 1; R = 2
Rennala Muon: ° = 1; B = 2
6 × 100
0
NanoChat Muon Async Training (¿i / i)
101
PA Ringmaster Muon: ´ = 5
Rennala Muon: ° = 1; B = 2 Delay-Adaptive Muon: ° = 1
101
p
Ringmaster Muon: ° = 1; R = 2
PA Ringmaster Muon: ´ = 5
Latest Minibatch Loss
Latest Minibatch Loss
Ringmaster Muon: ° = 1; R = 2
300
400
500
Runtime (seconds)
(a) n = 8, nearly homogeneous delays
6 × 100
4 × 100 0
100
200
300
400
500
0
100
Runtime (seconds)
(b) n = 8, sublinear delays
200
300
400
500
Runtime (seconds)
(c) n = 8, linear delays
Figure 3: NanoChat training loss versus simulated runtime for the same model with 8 simulated workers. The ordering is preserved at the smaller worker count, again with a larger gap in the more heterogeneous delay regimes. To anchor the simulated runtime to measured hardware performance, we first profile one stochasticgradient computation for the NanoChat model on the A100. After 3 warm-up steps, we average 100 synchronized gradient computations and obtain a mean step time of τb = 22.74 ms, as reported in Table 2. The NanoChat simulator then converts the dimensionless profiles in (2) into milliseconds by assigning worker i the deterministic base runtime τi = τbgi , where gi is chosen as gisim , gisub , or gilin depending on the delay regime. Finally, we inject half-normal runtime noise with standard deviation equal to 5% of this deterministic base runtime. For n = 32, the similar-delay regime is centered at 22.74 ms per worker; before adding noise, the sublinear regime ranges from 22.74 ms to 128.62 ms and the linear regime ranges from 22.74 ms to 727.57 ms. For n = 8, the corresponding upper bounds are 64.31 ms and 181.89 ms. These simulated ranges are consistent with the cross-device variability in Table 2, while remaining tied to a measured gradient time for the exact NanoChat setup used in the experiments. We evaluate NanoChat with n = 32 and n = 8 simulated workers. Hyperparameters are tuned over a runtime horizon of 500 simulated seconds with one trial per configuration. The search grids are • Ringmaster Muon: η ∈ {5−2 , 5−1 , 1, 5} and R ∈ {1, 2, 4, 8, 16}; • Parameter-agnostic Ringmaster Muon: η ∈ {5−2 , 5−1 , 1, 5, 25}; • Rennala Muon: η = 5 and B ∈ {1, 2, 4, 8, 16}; • Delay-adaptive Muon: η ∈ {5−2 , 5−1 , 1, 5}. Figure 2 shows that Ringmaster Muon achieves the lowest training loss across all three 32-worker regimes. In the nearly homogeneous setting, the margin over the strongest asynchronous baselines is modest, consistent with the limited straggler effect. As heterogeneity increases, the gain widens, with the clearest separation in the linear-delay regime. Figure 3 shows the same qualitative ordering for 12
8 workers. Together with the quadratic benchmark, these results support the main empirical claim: delay thresholding is most useful in the heterogeneous, straggler-dominated settings targeted by the theory, and it does not require changing the underlying Muon update.
6
Conclusion
We introduced Ringmaster LMO, an asynchronous LMO-based momentum method for stochastic nonconvex optimization, together with a parameter-agnostic variant. Our analysis establishes iteration and time complexity guarantees under generalized (L0 , L1 )-smoothness and shows that, in the classical Euclidean smooth setting, the time bound recovers the optimal time complexity known for Ringmaster ASGD. Experiments on stochastic quadratics and NanoChat pretraining further show that the benefits of asynchrony become more pronounced as worker-speed heterogeneity increases. An important limitation of our analysis is that it assumes a homogeneous data distribution across workers; extending both the method and the theory to data-heterogeneous distributed settings is a natural next step. Another limitation is that lower bounds for the asynchronous time complexity of general LMO-based methods in arbitrary norms remain open. This question is challenging because matching lower bounds are not known even for the corresponding synchronous case. Developing such lower bounds for the parameter-agnostic case appears even harder.
Acknowledgments and Disclosure of Funding The research reported in this publication was supported by funding from King Abdullah University of Science and Technology (KAUST): i) KAUST Baseline Research Scheme, ii) CRG Grant ORFS-CRG12-2024-6460, and iii) Center of Excellence for Generative AI, under award number 5940.
References Alekh Agarwal and John C Duchi. Distributed delayed stochastic optimization. Advances in Neural Information Processing Systems, 24, 2011. (Cited on page 3.) Yossi Arjevani, Ohad Shamir, and Nathan Srebro. A tight convergence analysis for stochastic gradient descent with delayed updates. In Algorithmic Learning Theory, pages 111–132. PMLR, 2020. (Cited on page 3.) Yossi Arjevani, Yair Carmon, John C Duchi, Dylan J Foster, Nathan Srebro, and Blake Woodworth. Lower bounds for non-convex stochastic optimization. Mathematical Programming, 199(1):165–214, 2023. (Cited on page 7.) Jeremy Bernstein and Laker Newhouse. Old optimizer, new norm: An anthology. arXiv:2409.20325, 2024. (Cited on pages 3 and 6.)
arXiv preprint
David Carlson, Volkan Cevher, and Lawrence Carin. Stochastic Spectral Descent for Restricted Boltzmann Machines. In Guy Lebanon and S. V. N. Vishwanathan, editors, Proceedings of the Eighteenth International Conference on Artificial Intelligence and Statistics, volume 38 of Proceedings of Machine Learning Research, pages 111–119, San Diego, California, USA, 09–12 May 2015a. PMLR. URL https://proceedings.mlr. press/v38/carlson15.html. (Cited on page 3.) David Carlson, Edo Collins, Ya-Ping Hsieh, Lawrence Carin, and Volkan Cevher. Preconditioned spectral descent for deep learning. In C. Cortes, N. Lawrence, D. Lee, M. Sugiyama, and R. Garnett, editors, Advances in Neural Information Processing Systems, volume 28. Curran Associates, Inc., 2015b. URL https://proceedings.neurips.cc/paper_files/paper/2015/file/ f50a6c02a3fc5a3a5d4d9391f05f3efc-Paper.pdf. (Cited on page 3.) David Carlson, Ya-Ping Hsieh, Edo Collins, Lawrence Carin, and Volkan Cevher. Stochastic spectral descent for discrete graphical models. IEEE Journal of Selected Topics in Signal Processing, 10(2):296–311, 2016. doi: 10.1109/JSTSP.2015.2505684. (Cited on page 3.) Ziyi Chen, Yi Zhou, Yingbin Liang, and Zhaosong Lu. Generalized-smooth nonconvex optimization is as efficient as smooth nonconvex optimization. In International Conference on Machine Learning, pages 5396–5427. PMLR, 2023. (Cited on pages 5, 7, and 17.)
13
Kenneth L. Clarkson. Coresets, sparse greedy approximation, and the Frank-Wolfe algorithm. ACM Trans. Algorithms, 6(4), September 2010. ISSN 1549-6325. doi: 10.1145/1824777.1824783. URL https: //doi.org/10.1145/1824777.1824783. (Cited on page 3.) Michael Crawshaw, Mingrui Liu, Francesco Orabona, Wei Zhang, and Zhenxun Zhuang. Robustness to unbounded smoothness of generalized signsgd. Advances in Neural Information Processing Systems, 35: 9955–9968, 2022. (Cited on page 7.) Ashok Cutkosky and Harsh Mehta. Momentum improves normalized SGD. In International Conference on Machine Learning, pages 2260–2268. PMLR, 2020. (Cited on page 6.) Sanghamitra Dutta, Gauri Joshi, Soumyadip Ghosh, Parijat Dube, and Priya Nagpurkar. Slow and stale gradients can win the race: Error-runtime trade-offs in distributed SGD. In International Conference on Artificial Intelligence and Statistics, pages 803–812. PMLR, 2018. (Cited on page 36.) Hamid Reza Feyzmahdavian, Arda Aytekin, and Mikael Johansson. An asynchronous mini-batch algorithm for regularized stochastic optimization. IEEE Transactions on Automatic Control, 61(12):3740–3754, 2016. (Cited on page 3.) Marguerite Frank and Philip Wolfe. An algorithm for quadratic programming. Naval Research Logistics Quarterly, 3(1-2):95–110, 1956. doi: https://doi.org/10.1002/nav.3800030109. URL https://onlinelibrary. wiley.com/doi/abs/10.1002/nav.3800030109. (Cited on page 3.) Eduard Gorbunov, Nazarii Tupitsa, Sayantan Choudhury, Alen Aliev, Peter Richtárik, Samuel Horváth, and Martin Takáč. Methods for convex (L0 , L1 )-smooth optimization: Clipping, acceleration, and adaptivity. In The Thirteenth International Conference on Learning Representations, 2025. (Cited on page 17.) Vineet Gupta, Tomer Koren, and Yoram Singer. Shampoo: Preconditioned stochastic tensor optimization. In International Conference on Machine Learning, pages 1842–1850. PMLR, 2018. (Cited on page 3.) Elad Hazan. Sparse approximate solutions to semidefinite programs. In LATIN 2008: Theoretical Informatics, pages 306–316, Berlin, Heidelberg, 2008. Springer Berlin Heidelberg. ISBN 978-3-540-78773-0. (Cited on page 3.) Elad Hazan, Kfir Levy, and Shai Shalev-Shwartz. Beyond convexity: Stochastic quasi-convex optimization. In C. Cortes, N. Lawrence, D. Lee, M. Sugiyama, and R. Garnett, editors, Advances in Neural Information Processing Systems, volume 28. Curran Associates, Inc., 2015. URL https://proceedings.neurips. cc/paper_files/paper/2015/file/934815ad542a4a7c5e8a2dfa04fea9f5-Paper.pdf. (Cited on page 6.) Florian Hübler, Junchi Yang, Xiang Li, and Niao He. Parameter-agnostic optimization under relaxed smoothness. In Sanjoy Dasgupta, Stephan Mandt, and Yingzhen Li, editors, Proceedings of The 27th International Conference on Artificial Intelligence and Statistics, volume 238 of Proceedings of Machine Learning Research, pages 4861–4869. PMLR, 02–04 May 2024. URL https://proceedings.mlr.press/v238/hubler24a. html. (Cited on pages 2, 4, 9, and 17.) Martin Jaggi. Revisiting Frank-Wolfe: Projection-free sparse convex optimization. In Sanjoy Dasgupta and David McAllester, editors, Proceedings of the 30th International Conference on Machine Learning, volume 28 of Proceedings of Machine Learning Research, pages 427–435, Atlanta, Georgia, USA, 17–19 Jun 2013. PMLR. URL https://proceedings.mlr.press/v28/jaggi13.html. (Cited on page 3.) Keller Jordan, Jeremy Bernstein, Brendan Rappazzo, @fernbear.bsky.social, Boza Vlado, You Jiacheng, Franz Cesista, Braden Koszarsky, and @Grad62304977. modded-nanogpt: Speedrunning the NanoGPT baseline, 2024a. URL https://github.com/KellerJordan/modded-nanogpt. (Cited on pages 1 and 3.) Keller Jordan, Yuchen Jin, Vlado Boza, Jiacheng You, Franz Cesista, Laker Newhouse, and Jeremy Bernstein. Muon: An optimizer for hidden layers in neural networks, 2024b. URL https://kellerjordan.github. io/posts/muon/. (Cited on pages 1, 3, and 6.) Andrej Karpathy. nanochat: The best ChatGPT that $100 can buy, 2025. URL https://github.com/ karpathy/nanochat. (Cited on page 11.) Sarit Khirirat, Abdurakhmon Sadiev, Yury Demidovich, and Peter Richtárik. Better LMO-based momentum methods with second-order information. arXiv preprint arXiv:2512.13227, 2025. (Cited on pages 3, 4, 7, 9, 17, and 18.) Diederik P Kingma and Jimmy Ba. Adam: A method for stochastic optimization. arXiv preprint arXiv:1412.6980, 2014. (Cited on page 1.)
14
Anastasiia Koloskova, Sebastian U Stich, and Martin Jaggi. Sharper convergence guarantees for asynchronous SGD for distributed and federated learning. Advances in Neural Information Processing Systems, 35: 17202–17215, 2022. (Cited on pages 3 and 5.) Dmitry Kovalev. Understanding gradient orthogonalization for deep learning via non-euclidean trust-region optimization. arXiv preprint arXiv:2503.12645, 2025. (Cited on pages 3 and 6.) Tian Li, Anit Kumar Sahu, Manzil Zaheer, Maziar Sanjabi, Ameet Talwalkar, and Virginia Smith. Federated optimization in heterogeneous networks. Proceedings of Machine Learning and Systems, 2:429–450, 2020. (Cited on page 36.) Xiangru Lian, Yijun Huang, Yuncheng Li, and Ji Liu. Asynchronous parallel stochastic gradient for nonconvex optimization. Advances in Neural Information Processing Systems, 28, 2015. (Cited on page 3.) Jingyuan Liu, Jianlin Su, Xingcheng Yao, Zhejun Jiang, Guokun Lai, Yulun Du, Yidao Qin, Weixin Xu, Enzhe Lu, Junjie Yan, et al. Muon is scalable for LLM training. arXiv preprint arXiv:2502.16982, 2025. (Cited on page 1.) Ilya Loshchilov and Frank Hutter. Decoupled weight decay regularization. In International Conference on Learning Representations, 2019. (Cited on page 1.) Ammar Mahran, Artavazd Maranjyan, and Peter Richtárik. Rescaled asynchronous SGD: Optimal distributed optimization under data and system heterogeneity. arXiv preprint arXiv:2605.13434, 2026. (Cited on page 3.) Artavazd Maranjyan. First Provably Optimal Asynchronous SGD for Homogeneous and Heterogeneous Data. PhD thesis, King Abdullah University of Science and Technology, 2025. (Cited on page 3.) Artavazd Maranjyan and Peter Richtárik. Ringleader ASGD: The first asynchronous SGD with optimal time complexity under data heterogeneity. In The Fourteenth International Conference on Learning Representations, 2026. URL https://openreview.net/forum?id=5wqTal0EuC. (Cited on pages 2 and 3.) Artavazd Maranjyan, Omar Shaikh Omar, and Peter Richtárik. Mindflayer SGD: Efficient parallel SGD in the presence of heterogeneous and random worker compute times. In The 41st Conference on Uncertainty in Artificial Intelligence, 2025a. (Cited on page 6.) Artavazd Maranjyan, El Mehdi Saad, Peter Richtárik, and Francesco Orabona. ATA: Adaptive task allocation for efficient resource management in distributed machine learning. In International Conference on Machine Learning, 2025b. (Cited on page 6.) Artavazd Maranjyan, Mher Safaryan, and Peter Richtárik. GradSkip: Communication-accelerated local gradient methods with better computational complexity. Transactions on Machine Learning Research, 2025c. ISSN 2835-8856. URL https://openreview.net/forum?id=6R3fRqFfhn. (Cited on page 36.) Artavazd Maranjyan, Alexander Tyurin, and Peter Richtárik. Ringmaster ASGD: The first asynchronous SGD with optimal time complexity. In International Conference on Machine Learning, 2025d. (Cited on pages 2, 3, 4, 5, 6, 7, 8, 33, 36, and 39.) Konstantin Mishchenko, Francis Bach, Mathieu Even, and Blake E Woodworth. Asynchronous SGD beats minibatch SGD under arbitrary delays. Advances in Neural Information Processing Systems, 35:420–433, 2022. (Cited on pages 3 and 5.) Thomas Pethick, Wanyun Xie, Kimon Antonakopoulos, Zhenyu Zhu, Antonio Silveti-Falls, and Volkan Cevher. Training deep learning models with norm-constrained LMOs. In Aarti Singh, Maryam Fazel, Daniel Hsu, Simon Lacoste-Julien, Felix Berkenkamp, Tegan Maharaj, Kiri Wagstaff, and Jerry Zhu, editors, Proceedings of the 42nd International Conference on Machine Learning, volume 267 of Proceedings of Machine Learning Research, pages 49069–49104. PMLR, 13–19 Jul 2025a. URL https://proceedings.mlr.press/v267/ pethick25a.html. (Cited on pages 2, 3, 4, and 6.) Thomas Pethick, Wanyun Xie, Mete Erdogan, Kimon Antonakopoulos, Tony Silveti-Falls, and Volkan Cevher. Generalized gradient norm clipping & non-euclidean (L0 , L1 )-smoothness. In The Thirty-ninth Annual Conference on Neural Information Processing Systems, 2025b. URL https://openreview.net/forum? id=rMdf8jhLR7. (Cited on pages 2, 3, 4, 5, and 7.) Benjamin Recht, Christopher Re, Stephen Wright, and Feng Niu. HOGWILD!: A lock-free approach to parallelizing stochastic gradient descent. Advances in Neural Information Processing Systems, 24, 2011. (Cited on pages 2 and 3.)
15
Artem Riabinin, Egor Shulgin, Kaja Gruntkowska, and Peter Richtárik. Gluon: Making muon & scion great again! (bridging theory and practice of LMO-based optimizers for LLMs). arXiv preprint arXiv:2505.13416, 2025. (Cited on pages 3, 4, 7, and 9.) Ishaan Shah, Anthony M Polloreno, Karl Stratos, Philip Monk, Adarsh Chaluvaraju, Andrew Hojel, Andrew Ma, Anil Thomas, Ashish Tanwer, Darsh J Shah, et al. Practical efficiency of Muon for pretraining. arXiv preprint arXiv:2505.02222, 2025. (Cited on page 1.) Egor Shulgin, Sultan AlRashed, Peter Richtárik, and Francesco Orabona. Beyond the ideal: Analyzing the inexact Muon update. In The 29th International Conference on Artificial Intelligence and Statistics, 2026. URL https://openreview.net/forum?id=IBRMWPBouf. (Cited on page 3.) Zhirayr Tovmasyan, Artavazd Maranjyan, and Peter Richtárik. Rennala MVR: Improved time complexity for parallel stochastic optimization via momentum-based variance reduction. arXiv preprint arXiv:2605.08871, 2026. (Cited on page 3.) John Tsitsiklis, Dimitri Bertsekas, and Michael Athans. Distributed asynchronous deterministic and stochastic gradient optimization algorithms. IEEE Transactions on Automatic Control, 31(9):803–812, 1986. (Cited on page 3.) Mark Tuddenham, Adam Prügel-Bennett, and Jonathan Hare. Orthogonalising gradients to speed up neural network optimisation. arXiv preprint arXiv:2202.07052, 2022. (Cited on page 3.) Alexander Tyurin. Tight time complexities in parallel stochastic optimization with arbitrary computation dynamics. arXiv preprint arXiv:2408.04929, 2024. (Cited on pages 3 and 36.) Alexander Tyurin and Peter Richtárik. Optimal time complexities of parallel stochastic optimization methods under a fixed computation model. In Advances in Neural Information Processing Systems, volume 36, pages 16515–16577, 2023. (Cited on pages 2, 3, 5, and 8.) Kaiyue Wen, David Leo Wright Hall, Tengyu Ma, and Percy Liang. Fantastic pretraining optimizers and where to find them. In The Fourteenth International Conference on Learning Representations, 2026. URL https://openreview.net/forum?id=2J51qUZ0iG. (Cited on page 1.) Xuyang Wu, Sindri Magnusson, Hamid Reza Feyzmahdavian, and Mikael Johansson. Delay-adaptive step-sizes for asynchronous learning. In Kamalika Chaudhuri, Stefanie Jegelka, Le Song, Csaba Szepesvari, Gang Niu, and Sivan Sabato, editors, Proceedings of the 39th International Conference on Machine Learning, volume 162 of Proceedings of Machine Learning Research, pages 24093–24113. PMLR, 17–23 Jul 2022. URL https://proceedings.mlr.press/v162/wu22g.html. (Cited on page 3.) Yidong Wu and Luo Luo. Optimal asynchronous stochastic nonconvex optimization under heavy-tailed noise. arXiv preprint arXiv:2601.19379, 2026. (Cited on pages 2, 4, 6, 7, 8, and 9.) Bohang Zhang, Jikai Jin, Cong Fang, and Liwei Wang. Improved analysis of clipping algorithms for non-convex optimization. Advances in Neural Information Processing Systems, 33:15511–15521, 2020a. (Cited on pages 3 and 4.) Jingzhao Zhang, Tianxing He, Suvrit Sra, and Ali Jadbabaie. Why gradient clipping accelerates training: A theoretical justification for adaptivity. In International Conference on Learning Representations, 2020b. (Cited on pages 3, 4, and 7.)
16
A
Auxiliary lemmas
We collect several auxiliary estimates used in the convergence proofs, including a summation/product bound, consequences of generalized (L0 , L1 )-smoothness, and a descent inequality for the LMO update. Lemma A.1 (Hübler et al. [2024]). Consider parameters q ∈ (0, 1) and p ≥ 0. Let a, b ∈ N such that b ≥ a ≥ 2. 1. The following inequality is true: b Y 1 −q 1−q 1−q (1 − t ) ≤ exp a −b 1−q t=a 2. If p ≤ q, then the following inequality holds: 1−q b t X Y a − (a − 1)1−q −p −q q−p . t (1 − τ ) ≤ (a − 1) exp 1−q t=a τ =a 1
(3)
(4)
1
1−q }, then the following inequality holds: 3. If the initial index satisfies a ≥ max{p 1−q , ( p−q 2 ) b b X Y 1 t−p (1 − τ −q ) ≤ 2 exp (b + 1)q−p . (5) 1 − q t=a τ =t+1
In particular, for p ≤ 1, these requirements are satisfied for all valid a, b.
Throughout this section, we assume the objective function is symmetric (L0 , L1 )-smooth. To facilitate the proof of our main convergence theorem, we first establish an auxiliary lemma detailing several of its useful properties. Lemma A.2 (Chen et al. [2023]). Let f satisfy Assumption 2.1 (lower boundedness) and Assumption 2.3 (generalized smoothness). Then, for all x, y ∈ Rd , the following holds: 1 2 f (y) ≤ f (x) + ⟨∇f (x), y − x⟩ + (L0 + L1 ∥∇f (x)∥∗ ) exp (L1 ∥x − y∥) ∥x − y∥ , (6) 2 ∥∇f (x) − ∇f (y)∥∗ ≤ (L0 + L1 ∥∇f (x)∥∗ ) exp (L1 ∥x − y∥) ∥x − y∥ , (7) 2
∥∇f (x)∥∗ ≤ 4 (L0 + L1 ∥∇f (x)∥∗ ) (f (x) − f∗ ).
(8)
Proof. Inequalities (6) and (7) were originally proven by Chen et al. [2023] for the Euclidean case, and subsequently extended to the non-Euclidean setting by Khirirat et al. [2025]. Regarding inequality (8), the Euclidean case was established by Gorbunov et al. [2025]. Here, we α∥∇f (x)∥∗ lmo(∇f (x)) extend this result to the non-Euclidean case. To do so, we define y = x+ L0 +L1 ∥∇f (x)∥ ∗ and apply (6): 1 2 f∗ ≤ f (y) ≤ f (x) + ⟨∇f (x), y − x⟩ + (L0 + L1 ∥∇f (x)∥∗ ) exp (L1 ∥x − y∥) ∥x − y∥ 2 α ∥∇f (x)∥∗ ≤ f (x) + ⟨∇f (x), lmo(∇f (x))⟩ L0 + L1 ∥∇f (x)∥∗ 2 α2 ∥∇f (x)∥∗ αL1 ∥∇f (x)∥∗ 1 + (L0 + L1 ∥∇f (x)∥∗ ) exp 2 (L0 + L1 ∥∇f (x)∥∗ )2 L0 + L1 ∥∇f (x)∥∗ 2
2
α ∥∇f (x)∥∗ ∥∇f (x)∥∗ α2 exp(α) + . L0 + L1 ∥∇f (x)∥∗ 2 L0 + L1 ∥∇f (x)∥∗ Setting α = 1/2, we obtain: ≤ f (x) −
2
f∗ ≤ f (x) −
2
∥∇f (x)∥∗ ∥∇f (x)∥∗ + 2(L0 + L1 ∥∇f (x)∥∗ ) 4(L0 + L1 ∥∇f (x)∥∗ ) 2
= f (x) −
∥∇f (x)∥∗ . 4(L0 + L1 ∥∇f (x)∥∗ )
17
The following lemma establishes a standard descent inequality for the update step xk+1 = xk + ηk lmo(mk+1 ). Lemma A.3 (Khirirat et al. [2025]). Under Assumption 2.1 (lower boundedness) Assumption 2.2 (bounded variance) and Assumption 2.3 (generalized smoothness), the iterates generated by Algorithm 1 satisfy: K−1 X k=0
ηk ϕk ∥∇f (xk )∥∗ ≤ ∆0 + 2
K−1 X k=0
K−1
ηk ∥êk+1 ∥∗ +
L0 X exp (L1 ηk ) ηk2 , 2
(9)
k=0
where we define the error term as êk+1 := mk+1 − ∇f (xk ) and the constant ϕk = 1 − L1 ηk 2 exp(L1 ηk ).
Proof. By Lemma A.2, we have: f (xk+1 ) ≤ f (xk ) + ⟨∇f (xk ), xk+1 − xk ⟩ 1 2 + (L0 + L1 ∥∇f (xk )∥∗ ) exp (L1 ∥xk+1 − xk ∥) ∥xk+1 − xk ∥ 2 1 ≤ f (xk ) + ηk ⟨∇f (xk ), lmo(mk+1 )⟩ + (L0 + L1 ∥∇f (xk )∥∗ ) exp (L1 ηk ) ηk2 2 = f (xk ) + ηk ⟨∇f (xk ) − mk+1 , lmo(mk+1 )⟩ + ηk ⟨mk+1 , lmo(mk+1 )⟩ 1 + (L0 + L1 ∥∇f (xk )∥∗ ) exp (L1 ηk ) ηk2 2 ≤ f (xk ) + ηk ∥∇f (xk ) − mk+1 ∥∗ − ηk ∥mk+1 ∥∗ 1 + (L0 + L1 ∥∇f (xk )∥∗ ) exp (L1 ηk ) ηk2 2 ≤ f (xk ) + 2ηk ∥∇f (xk ) − mk+1 ∥∗ − ηk ∥∇f (xk )∥∗ 1 + (L0 + L1 ∥∇f (xk )∥∗ ) exp (L1 ηk ) ηk2 . 2 In the steps above, the second inequality follows from the property ∥lmo(·)∥ ≤ 1. The third inequality applies Hölder’s inequality alongside the definition of the linear minimization oracle (which implies ⟨mk+1 , lmo(mk+1 )⟩ = − ∥mk+1 ∥∗ ). Finally, the last step utilizes the triangle inequality. Rearranging the terms to isolate ∥∇f (xk )∥∗ and summing over k = 0, . . . , K − 1 yields the desired result.
We remind that momentum updates defined as follows m1 = (1 − α0 )m0 + αinit g0 , for k = 0; mk+1 = (1 − αk )mk + αk gk , for k ≥ 1 . We define the error term as follows êk+1 := mk+1 − ∇f (xk ) .
Lemma A.4. Suppose Assumption 2.2 (bounded variance) and Assumption 2.3 (generalized smoothness). If stepsize ηk is non-increasing and momentum parameter αk defined as α0 = 1 and 18
αinit , αk ∈ (0, 1] for k ≥ 1, then the iterates of Algorithm 1 satisfy
E [∥êk+1 ∥∗ ] ≤
k Y k X
(1 − αj ) (L0 + L1 E [∥∇f (xt )∥∗ ]) exp(L1 ηt−1 )ηt−1
t=1 j=t
+
k k X Y
(1 − αj )αt (L0 + L1 E [∥∇f (xt )∥∗ ]) exp (L1 Rt ηt−δt ) Rt ηt−δt
t=1 j=t+1
v u k k k Y uX Y (1 − αj )2 α2 + (1 − αj )2 α2 + ρσ t t
t=1 j=t+1
+
k Y
init
j=1
(1 − αt )(1 − αinit )E [∥∇f (x0 )∥∗ ] .
t=1
Throughout the paper, an empty product is interpreted as 1 and an empty sum as 0.
Proof. Using momentum update, we derive for any k ≥ 1: êk+1 = mk+1 − ∇f (xk ) = (1 − αk )mk + αk gk − ∇f (xk ) = (1 − αk )(mk − ∇f (xk−1 )) + (1 − αk )(∇f (xk−1 ) − ∇f (xk )) + αk gk − αk ∇f (xk ) = (1 − αk )(mk − ∇f (xk−1 )) + (1 − αk )(∇f (xk−1 ) − ∇f (xk )) +αk (∇f (xk−δk ) − ∇f (xk )) + αk (gk − ∇f (xk−δk )) . Denoting Dk := ∇f (xk−1 ) − ∇f (xk ), Sk := ∇f (xk−δk ) − ∇f (xk ), and ek := gk − ∇f (xk−δk ), we obtain êk+1 = (1 − αk )êk + (1 − αk )Dk + αk Sk + αk ek =
k Y
(1 − αt )ê1 +
t=1
+
k k X Y
k k X Y
(1 − αj )(1 − αt )Dt
t=1 j=t+1
(1 − αj )αt St +
t=1 j=t+1
k k X Y
(1 − αj )αt et .
(10)
t=1 j=t+1
Now we consider ê1 at iteration k = 1. Using α0 = 1 and δ0 = 0 (which implies x0−δ0 = x0 ), we have ê1 = m1 − ∇f (x0 ) = (1 − α0 )m0 + αinit g0 − ∇f (x0 ) = αinit (g0 − ∇f (x0 )) − (1 − αinit )∇f (x0 ) = αinit e0 − (1 − αinit )∇f (x0 ) . 19
(11)
Plugging (11) into (10), we have êk+1 =
k Y
(1 − αt )αinit e0 −
t=1
+
k k X Y
k Y
(1 − αt )(1 − αinit )∇f (x0 )
t=1
(1 − αj )(1 − αt )Dt +
t=1 j=t+1
=
k k X Y
(1 − αj )(1 − αt )Dt +
t=1 j=t+1
{z
| +
}
:=T1 k k X Y
(1 − αj )αt et +
t=1 j=t+1
|
k Y
k k X Y
(1 − αj )αt St +
t=1 j=t+1
k k X Y
k k X Y
(1 − αj )αt et
t=1 j=t+1
(1 − αj )αt St
t=1 j=t+1
{z
|
}
:=T2 k Y
(1 − αt )αinit e0 −
t=1
{z
|
}
:=T3
(1 − αt )(1 − αinit )∇f (x0 ) .
t=1
{z
}
:=T4
To continue the proof, we need to bound the norm of terms T1 , T2 , T3 and T4 . By triangle inequality, we get ∥T1 ∥∗ ≤ =
k k X Y
(1 − αj )(1 − αt ) ∥Dt ∥∗
t=1 j=t+1 k Y k X
(1 − αj ) ∥∇f (xt−1 ) − ∇f (xt )∥∗
t=1 j=t
≤
k Y k X
(1 − αj ) (L0 + L1 ∥∇f (xt )∥∗ ) exp(L1 ηt−1 )ηt−1 .
t=1 j=t
We used Lemma A.2 in the last inequality. Bounding ∥T2 ∥∗ , we obtain ∥T2 ∥∗ ≤ =
k k X Y
(1 − αj )αt ∥St ∥∗
t=1 j=t+1 k k X Y
(1 − αj )αt ∥∇f (xt−δt ) − ∇f (xt )∥∗
t=1 j=t+1
≤ ≤
k k X Y
(1 − αj )αt (L0 + L1 ∥∇f (xt )∥∗ ) exp (L1 ∥xt − xt−δt ∥) ∥xt − xt−δt ∥
t=1 j=t+1 k k X Y
(1 − αj )αt (L0 + L1 ∥∇f (xt )∥∗ )
t=1 j=t+1
t−1 X
× exp L1
!
t−1 X
ητ lmo(mτ +1 )
ητ lmo(mτ +1 ) .
τ =t−δt
τ =t−δt
Assuming that stepsize ηk is non-increasing, we have ∥T2 ∥∗ ≤
k k Y X
(1 − αj )αt (L0 + L1 ∥∇f (xt )∥∗ ) exp (L1 Rt ηt−δt ) Rt ηt−δt .
t=1 j=t+1
Bounding ∥T3 ∥∗ , we obtain ∥T3 ∥∗ ≤ ρ
k k X Y
(1 − αj )αt et +
t=1 j=t+1
20
k Y
(1 − αt )αinit e0
t=1
. 2
Taking expectation, we have E [∥T3 ∥∗ ] ≤ ρE
k k X Y
(1 − αj )αt et +
t=1 j=t+1
k Y
(1 − αt )αinit e0
t=1
2 v u 2 u k k k Y Y u X ≤ ρu (1 − αj )αt et + (1 − αt )αinit e0 tE t=1 j=t+1
t=1
2
v u k k k h i Y h i uX Y 2 2 E ∥e ∥2 ≤ ρt (1 − αj )2 αt2 E ∥et ∥2 + (1 − αj )2 αinit 0 2 t=1 j=t+1
j=1
v u k k k Y uX Y (1 − αj )2 α2 + ≤ ρσ t (1 − αj )2 α2
init ,
t
t=1 j=t+1
j=1
where the second inequality uses Jensen’s inequality and the third uses iid ness. Bounding ∥T4 ∥∗ , we obtain ∥T4 ∥∗ ≤
k Y
(1 − αt )(1 − αinit ) ∥∇f (x0 )∥∗ .
t=1
Thus, we have E [∥êk+1 ∥∗ ] ≤
k Y k X
(1 − αj ) (L0 + L1 E [∥∇f (xt )∥∗ ]) exp(L1 ηt−1 )ηt−1
t=1 j=t
+
k k X Y
(1 − αj )αt (L0 + L1 E [∥∇f (xt )∥∗ ]) exp (L1 Rt ηt−δt ) Rt ηt−δt
t=1 j=t+1
v u k k k Y uX Y (1 − αj )2 α2 + (1 − αj )2 α2 +ρσ t t
t=1 j=t+1
+
k Y
init
j=1
(1 − αt )(1 − αinit )E [∥∇f (x0 )∥∗ ] .
t=1
21
B
Proof of Theorem 4.1
This section proves the fixed-parameter result and derives the parameter choices by balancing the terms in the convergence bound. We first restate the theorem and then optimize the constant momentum, threshold, and stepsize parameters. Let us first restate the theorem. Theorem 4.1 (Fixed-Parameter Iteration Complexity; Proof in Section B). Suppose Assumptions 2.1– 2.3 hold. Run Algorithm 1 with α0 = 1 and with the constant momentum parameter and delay threshold ( ) 1/2 1/2 ∆0 L0 1 αk = αinit = α = min 1, , ∀k ≥ 1, Rk = R = , ∀k ≥ 0, α ρσK 1/2 with the convention that, when σ = 0, we take α = 1, and choose the constant stepsize as follows: • If L1 = 0, take (r ηk = η = min
3/4
∆0 ∆0 , L0 K L1/4 (ρσ)1/2 K 3/4
) ,
∀k ≥ 0.
3/4
)
0
• If L1 > 0, take (r ηk = η = min
1/2
1/2
∆0 1 ∆0 L 0 ∆0 , , , L0 K 8L1 8L1 ρσK 1/2 L1/4 (ρσ)1/2 K 3/4
,
0
∀k ≥ 0.
Then, after K ≥ 1 iterations, the iterates of Algorithm 1 satisfy ! r 1/2 1/4 1/4 L0 ∆0 L 1 ∆0 L1 ∆0 ρσ L0 ∆0 (ρσ)1/2 ρσ + + 1/2 + 1/2 . min E [∥∇f (xk )∥∗ ] ≤ O + K K k∈{0,...,K−1} K 1/4 K L K 1/2 0
In the case L1 = 0, the terms containing L1 in the preceding bound are zero. Consequently, the number of iterations needed to reach an ε-stationary point satisfies L 1 ∆0 L0 ∆0 (ρσ)2 L21 ∆0 (ρσ)2 (ρσ)2 L0 ∆0 + + + + . Kε = O ε2 ε ε4 L0 ε2 ε2 Again, when L1 = 0, the terms containing L1 are omitted. Proof. We start with Lemma A.3, taking ηk = η and αk = 1 for k = 0 and αk = α for k ≥ 1, K−1 K−1 X X 1 1 0 ≤ ∆0 + 2η ∥êk+1 ∥∗ + η 2 L0 eL1 η K − η 1 − L1 ηeL1 η ∥∇f (xk )∥∗ . (12) 2 2 k=0
k=0
Taking expectation, we bound the second term from the last inequality. By Lemma A.4 (take αinit = α ), we have 2η
K−1 X k=0
E [∥êk+1 ∥∗ ] ≤ 2η
K−1 k XX k=0 t=1
(1 − α)k−t (L0 + L1 E [∥∇f (xt )∥∗ ]) eL1 η η
|
{z
}
:=T1
+ 2η
K−1 k XX k=0 t=1
(1 − α)k−t α(L0 + L1 E [∥∇f (xt )∥∗ ])eL1 Rη Rη
| + 2η
{z
K−1 X
v u k K−1 uX X ρσ t (1 − α)2(k−t) α2 + 2η (1 − α)k+1 E [∥∇f (x )∥ ] . 0
k=0
|
}
:=T2
t=0
∗
k=0
{z
}
:=T3
22
|
{z
:=T4
}
Next, we bound each T1 , T2 , T3 , T4 and sum them up. T1 = 2L0 η 2 eL1 η ≤ 2L0 η 2 eL1 η ≤ ≤
K−1 k XX k=0 t=1
(1 − α)k−t + 2L1 η 2 eL1 η
K−1 ∞ XX k=0 t=0
(1 − α)t + 2L1 η 2 eL1 η K−1 ∞ XX
2 L0 η 2 eL1 η K + 2L1 η 2 eL1 η α t=1 2 2 L0 η 2 eL1 η K + L1 η 2 eL1 η α α
k=0
K−1 X
K−1 k XX k=0 t=1
K−1 X K−1 X t=1 k=t
(1 − α)k−t E [∥∇f (xt )∥∗ ]
(1 − α)k−t E [∥∇f (xt )∥∗ ]
(1 − α)k E [∥∇f (xt )∥∗ ]
E [∥∇f (xk )∥∗ ] .
(13)
k=1
Now we bound T2 , By the same argument for bounding T1 , we obtain K−1
X 2 2 T2 ≤ L0 Rη 2 αeL1 Rη K + L1 Rη 2 αeL1 Rη E [∥∇f (xk )∥∗ ] . α α
(14)
k=1
Next, we bound T3 : v
T3 ≤ 2ηρσ ≤ 2ηρσ
u k K−1 X uX t
t=0
k=0
v u t
u∞ K−1 X X t=0
k=0
√
(1 − α)2(k−t) α2 (1 − α)t α2
≤ 2ηρσ αK.
(15)
Finally, we bound T4 : T4 ≤ 2η
K−1 X k=0
(1 − α)k+1 E [∥∇f (x0 )∥∗ ]
2η ≤ E [∥∇f (x0 )∥∗ ] . α Plugging (13), (14), (15), and (16) into (12), we have K−1 X 1 2 1 L1 η L1 η 0 ≤ ∆0 + η L0 e K − η 1 − L1 ηe E [∥∇f (xk )∥∗ ] 2 2 k=0
K−1
X 2 2 + L0 η 2 eL1 η K + L1 η 2 eL1 η E [∥∇f (xk )∥∗ ] α α k=1
K−1
X 2 2 + L0 Rη 2 αeL1 Rη K + L1 Rη 2 αeL1 Rη E [∥∇f (xk )∥∗ ] α α k=1
√
2η +2ηρσ αK + E [∥∇f (x0 )∥∗ ] . α 2 Multiplying by ηK , we have
0≤
√ 2∆0 4 4 4 + 4ρσ α + ηL0 eL1 η + eL1 η + RαeL1 Rη + E [∥∇f (x0 )∥∗ ] ηK α α αK K−1 2 X 1 2 2 − 1 − L1 ηeL1 η − L1 ηeL1 η − L1 RηαeL1 Rη E [∥∇f (xk )∥∗ ] . K 2 α α k=0
23
(16)
α and α ≤ 1, we have Assuming R = α1 and η ≤ 8L 1 K−1
√ 1 X 2∆0 4 L0 E [∥∇f (xk )∥∗ ] ≤ + 4ρσ α + E [∥∇f (x0 )∥∗ ] + 10η + 2ηL0 . K ηK αK α k=0
To select α, we minimize the upper bound w.r.t α with restriction that α < 1. Thus, we ∆
1/2
L
1/2
0 0 have α = min 1, ρσK 1/2
. To do so for η, we first plug in α in the upper bound and
α minimize w.r.t. η under constraints that η ≤ 8L . This allows us to set a stepsize as η = 1 q 3/4 1/2 1/2 ∆ ∆ L ∆0 1 0 0 0 . Also, by Lemma A.2 (inequality(8)) applied at min 1/2 L0 K , 8L1 , 8L ρσK 1/2 , 3/4 1/4 K
1
L0 (ρσ)
x0 , we have ( 2 ∥∇f (x)∥∗ ≤
8L0 ∆0 (8L1 ∆0 )2
n o 0 if ∥∇f (x)∥ ≤ L L1 ; ≤ max 8L ∆ , (8L ∆ )2 . 0 0 1 0 otherwise;
Therefore, we derive min k∈{0,...,K−1}
E [∥∇f (xk )∥∗ ] r
√
1/2
1/4
1/4
L0 ∆0 + L1 ∆0 L1 ∆ L ∆ (ρσ)1/2 ρσ ≤O + 1/2 0 ρσ + 0 0 1/4 + 1/2 1/2 K K K L0 K ! r 1/2 1/4 1/4 L0 ∆0 L1 ∆0 L1 ∆0 L0 ∆0 (ρσ)1/2 ρσ =O + + 1/2 + 1/2 . ρσ + K K K 1/4 K L0 K 1/2 L0 ∆0 + K
Therefore, iteration complexity is defined as L 1 ∆0 L21 ∆0 (ρσ)2 L0 ∆0 (ρσ)2 (ρσ)2 L0 ∆0 + + + + Kε = O . ε2 ε L0 ε2 ε4 ε2
24
!
C
Proof of Theorem 4.4
We next turn to the parameter-agnostic regime, where the momentum and threshold schedules vary with the iteration index. The proof follows the same general template as in the fixed-parameter case, but the time-varying parameters require more careful control of the accumulated error terms. Let us first restate the theorem before proving it. Theorem 4.4. Suppose Assumptions 2.1–2.3 hold. Run Algorithm 1 with Rk = 1/αk , αinit = α0 = 1, and αk = k −1/2 for k ≥ 1. In the first two cases, let η > 0 be any constant and define 2 2
2 2
Ψ(L0 , L1 ) := eL1 η ∆0/η + ρσ + eL1 η L0 η. Choose the stepsizes as follows, and let Ψ denote the corresponding value in the final bound: η ∆0/η + ρσ + L η. • If L1 = 0, take ηk = (k+1) 3/4 and Ψ = Ψ(L0 , 0) = 0 η • If L1 > 0, take ηk = 17(k+1) 3/4 and Ψ = Ψ(L0 , L1 ).
1 L0/L1 , which is the • If L1 > 0 is known, take ηk = 17L (k+1) 3/4 and Ψ = L1 ∆0 + ρσ + 1 preceding envelope with η = 1/L1 up to universal constants.
Then, after K ≥ 1 iterations,
Ψ log K min E [∥∇f (xk )∥∗ ] ≤ O . k∈{0,...,K−1} K 1/4 Consequently, the number of iterations needed to reach an ε-stationary point satisfies 4 ! 1 Ψ4 . log Kε = O ε4 ε Proof. We start with Lemma A.3 K−1 X k=0
ηk E [∥∇f (xk )∥∗ ] ≤ ∆0 +
K−1 K−1 X L0 X 2 L1 ηk ηk e +2 ηk E [∥êk+1 ∥∗ ] 2 k=0 k=0 | {z } | {z } :=T1
+
:=T2
K−1 X
L1 ηk2 eL1 ηk E [∥∇f (xk )∥∗ ] . 2 k=0 | {z } :=T3
We need to bound three terms T1 , T2 , T3 . Bounding T1 . The term T1 captures the accumulated error dependent purely on the stepsize η̂ sequence. To bound this, we substitute our chosen stepsize ηk = (k+1) 3/4 and upper-bound the resulting summation using an integral approximation. Specifically, we have: K−1
T1 =
L0 X η̂ 2 exp (L1 ηk ) 2 (k + 1)3/2 k=0
L0 η̂ 2 exp (L1 η̂) ≤ 2
K−1 X
1 1+ (k + 1)3/2 k=1 ! Z K L0 η̂ 2 exp (L1 η̂) 1 ≤ 1+ dt 2 t3/2 1 K ! L0 η̂ 2 exp (L1 η̂) 2 1 + −√ = 2 t 1 L0 η̂ 2 exp (L1 η̂) 2 = 3− √ . 2 K 25
!
Thus, by dropping the negative term, we derive the following conservative upper bound: T1 ≤
3 L0 η̂ 2 exp (L1 η̂) . 2
(17)
Bounding T2 . Next, we analyze T2 , which encapsulates the accumulated tracking error. By applying Lemma A.4 (set αinit = α0 ), we can expand this error. To facilitate a tight convergence analysis, we systematically decompose the resulting expression into five distinct components (① through ⑤). This isolates the effects of the objective’s smoothness, the momentum weights, and the variance: T2 ≤ 2
K−1 X
ηk
K−1 X
ηk
k k X Y
(1 − αj )αt (L0 + L1 E [∥∇f (xt )∥∗ ]) exp (L1 Rt ηt−δt ) Rt ηt−δt
t=1 j=t+1
k=0
+ 2ρσ
(1 − αj ) (L0 + L1 E [∥∇f (xt )∥∗ ]) exp(L1 ηt−1 )ηt−1
t=1 j=t
k=0
+2
k Y k X
K−1 X
v u k k uX Y t ηk (1 − αj )2 α2 t
k=0
t=0 j=t+1
=①+②+③+④+⑤, where ①, ②, ③, ④ and ⑤ are defined as ① := 2L0
K−1 X
ηk
K−1 X
ηk
K−1 X
ηk
K−1 X k=0
⑤ := 2ρσ
K−1 X
(1 − αj )αt exp (L1 Rt ηt−δt ) Rt ηt−δt ;
k Y k X
(1 − αj ) exp(L1 ηt−1 )ηt−1 E [∥∇f (xt )∥∗ ] ;
t=1 j=t
k=0
④ := 2L1
k k Y X
t=1 j=t+1
k=0
③ := 2L1
(1 − αj ) exp(L1 ηt−1 )ηt−1 ;
t=1 j=t
k=0
② := 2L0
k Y k X
ηk
k k Y X
(1 − αj )αt exp (L1 Rt ηt−δt ) Rt ηt−δt E [∥∇f (xt )∥∗ ] ;
t=1 j=t+1
v u k k uX Y ηk t (1 − αj )2 α2 . t
k=0
t=0 j=t+1
We now proceed to bound these terms individually. Bounding ①. The first term, ①, represents the deterministic error discounted by our momentum parameter. By substituting our stepsize ηk and momentum parameter αj = j −1/2 , we can bound the nested products. We then control the resulting harmonic series by upper-bounding it with a logarithmic integral: 2
① ≤ 2L0 η̂ exp(L1 η̂) ≤ 2L0 η̂ 2 exp(L1 η̂) 2
K−1 k Y k XX
(1 − j −1/2 )t−3/4 (k + 1)−3/4
k=0 t=1 j=t K−1 X k=0
≤ 36L0 η̂ exp(L1 η̂)
18 k+1
Z K 1
1 dt + 1 t
2
≤ 36L0 η̂ exp(L1 η̂) (log K + 1) . 26
!
Bounding ②. The second term tracks the error impacted by the delay parameter δt ≤ Rt . Recalling that Rt = 1/αt for t ≥ 1, we have αt Rt = 1. We first bound the exponential term: 1 exp (L1 Rt ηt−δt ) = exp L1 ηt−δt αt √ t ≤ exp L1 η̂ (t − δt + 1)3/4 √ t √ ≤ exp L1 η̂ . (t − t + 1)3/4 √
t Since (t−√t+1) 3/4 ≤ 1.03 for all t ≥ 1, we can upper-bound this as:
exp L1 η̂
√ t √ ≤ exp(1.03L1 η̂). (t − t + 1)3/4
To continue bounding ②, we need to bound the following term: k k Y X
(1 − αj )αt Rt ηt−δt = η̂
t=1 j=t+1
≤ η̂ ≤ η̂
k k Y X
1
3
t=1 j=t+1 k k X Y
(1 − j − /2 )(t − Rt + 1)− /4 1
3
t=1 j=t+1 k k Y X
(1 − j − /2 )(t − 1
√
t + 1)− /4 3
t=1 j=t+1
≤ 2η̂
k k X Y
(1 − j − /2 )t− /4 1
3
t=1 j=t+1
(∗)
≤ 4η̂ exp ≤
(1 − j − /2 )(t − δt + 1)− /4
1 1 − 1/2
(k + 1) /2− /4 1
3
36η̂ , (k + 1)1/4
where in (∗) we used inequality (5) from Lemma A.1. Using the previous inequality, we bound the resulting sum similar to ①: ② ≤ 2L0 exp(1.03L1 η̂) ≤ 2L0 exp(1.03L1 η̂)
K−1 X
ηk
k=0
(1 − αj )ηt−δt
t=1 j=t+1
k=0 K−1 X
k k X Y
ηk ·
36η̂ (k + 1)1/4
2
≤ 72L0 η̂ exp(1.03L1 η̂)(log K + 1). 3/4
Bounding ③. For ③, we swap the order of summation. Notice that ηt−1 = ηt (t+1) t3/4 Rearranging the sums yields: ③=2
K−1 X k=1
≤2
K−1 X
ηk
k Y k X
(1 − αj )L1 E [∥∇f (xt )∥∗ ] exp(L1 ηt−1 )ηt−1
t=1 j=t
ηt−1 exp(L1 ηt−1 )L1 E [∥∇f (xt )∥∗ ]
t=1
K−1 X k=t
27
ηk
k Y
(1 − αj ).
j=t+1
≤ 23/4 ηt .
For t ≥ 1, using the empty-product convention and inequality (4) from Lemma A.1, we have K−1 X
ηk
k Y
(1 − αj ) = ηt +
j=t+1
k=t
≤
K−1 X
ηk
k Y
(1 − αj )
j=t+1
k=t+1
η̂ 1 1 + 3η̂t− /4 ≤ 4η̂t− /4 . (t + 1)3/4
Therefore, ③ ≤ 8L1 ≤ 8L1 Bounding ④. K−1 X
k Y
K−1 X
η̂ ηt−1 exp(L1 ηt−1 ) 1/4 E [∥∇f (xt )∥∗ ] t t=1
K−1 X
η̂ ηk−1 exp(L1 ηk−1 ) 1/4 E [∥∇f (xk )∥∗ ] . k k=1
By the inequality (4) from Lemma A.1, (1 − αj )ηk = ηt +
k=t j=t+1
K−1 X
k Y
(1 − αj )ηk
k=t+1 j=t+1
= η̂(t + 1)− /4 + η̂ 3
k Y
K−1 X
(1 − j − /2 )(k + 1)− /4 ≤ 4η̂t− /4 . 1
3
1
k=t+1 j=t+1
We apply the same sum-swapping technique and the previous inequality: K−1 k X K−1 X Y ηt−δt ④ = 2L1 (1 − αj )ηk ηt−δt exp L1 E [∥∇f (xt )∥∗ ] αt t=1 k=t j=t+1 √ K−1 X η̂ t E [∥∇f (xt )∥∗ ] ≤ 8L1 ηt−δt exp L1 η̂ 3/4 1/4 (t − δt + 1) t t=1 √ K−1 X t η̂ √ ≤ 8L1 ηt−δt exp L1 η̂ E [∥∇f (xt )∥∗ ] 1/4 3/4 t (t − t + 1) t=1 K−1 X 1 η̂ ≤ 8L1 ηk−δk exp 2L1 η̂ 1/4 E [∥∇f (xk )∥∗ ] , 1/4 k k k=1 √
3/4
t 2 where in the last inequality we used (t−√t+1) 3/4 ≤ 1/4 . t
Bounding ⑤. The variance term ⑤ involves the standard summation of squared momentum weights. Expanding this geometric-like sum gives the classical logarithmic accumulation: v u k K−1 k k Y Y X uX (1 − αj )2 α2 + (1 − αj )2 α2 ⑤ = 2ρσ ηk t t
t=1 j=t+1
k=0
= 2ρσ
K−1 X k=1 K−1 X
0
j=1
v u k k k Y uX Y ηk t (1 − j −1/2 )t−1 + (1 − j −1/2 )2 + 2ρσ η̂ t=1 j=t+1
s
j=1
1 1 √ + 2ρσ η̂ 1 − 1/2 k+1 k=1 s K−1 X η̂ 18 ≤ 2ρσ + 2ρσ η̂ 3/4 (k + 1) (k + 1)1/2 k=1
≤ 2ρσ
ηk
2 exp
≤ 11ρσ η̂(log K + 1),
where in the first inequality we used the inequality (5) from Lemma A.1. 28
Combining the terms.
We now piece together T1 , T2 , and T3 . Notice that K−1
T3 =
1 X L1 ηk2 eL1 ηk E [∥∇f (xk )∥∗ ] . 2 k=0
Combining T3 with ③ and ④, we can group the coefficients in front of E [∥∇f (xk )∥∗ ]. 3/4
Bounding the ratios of sequential stepsizes (k+1) k3/4 the coefficients into a single envelope:
T3 + ③ + ④ ≤
K−1 X
L1 ηk
k=0
3/4
(k+1) ≤ 23/4 and (k−δ ≤ 23/4 , we consolidate +1)3/4 k
1 −1/4 exp 3L1 η̂(k + 1) E [∥∇f (xk )∥∗ ] . 33η̂ (k + 1)1/4
Substituting this back into the main inequality gives: K−1 X k=0
3 ηk E [∥∇f (xk )∥∗ ] ≤ ∆0 + L0 η̂ 2 eL1 η̂ + 36L0 η̂ 2 eL1 η̂ (log K + 1) 2 + 72L0 η̂ 2 e1.03L1 η̂ (log K + 1) + 11ρσ η̂(log K + 1) K−1 X η̂ −1/4 + 33L1 ηk exp 3L η̂(k + 1) E [∥∇f (xk )∥∗ ] . 1 (k + 1)1/4 k=0
Thresholding trick and absorbing the sum. To handle the recursive expected gradient term on the right-hand side, we split the summation into two parts: a transient phase (k < k0 ) and a stable phase (k ≥ k0 ). We define the threshold iteration k0 such that the multiplicative factor becomes less than or equal to 1/2, allowing us to absorb it into the left-hand side. We require: 33L1
1 η̂ −1/4 exp 3L η̂(k + 1) ≤ . 1 2 (k + 1)1/4
Because x exp(cx) is monotonically increasing, we can secure this bound for all k ≥ k0 if we conservatively choose: k0 := ⌈(70L1 η̂)4 ⌉ − 1. PK−1 For iterations k ≥ k0 , we subtract 12 k=k0 ηk E [∥∇f (xk )∥∗ ] from both sides. We are left with the task of explicitly bounding the transient sum from k = 0 to k0 − 1. First, we bound the local travel distance (radius) rk = ∥xk − x0 ∥. Summing the stepsizes and approximating via an integral: ! Z k k h ik X 1 1 1/4 1 1 ≤ η̂ 1 + dt = η̂ 1 + 4t = 4η̂k /4 − 3η̂ ≤ 4η̂k /4 . rk ≤ ητ = η̂ 3/4 3/4 1 τ t 1 τ =0 τ =1 k−1 X
Next, we bound E [∥∇f (xk )∥∗ ] using the generalized (L0 , L1 )-smoothness condition and the triangle inequality: ∥∇f (xk )∥∗ ≤ ∥∇f (xk ) − ∇f (x0 )∥∗ + ∥∇f (x0 )∥∗ ≤ (L0 + L1 ∥∇f (x0 )∥∗ ) exp (L1 ∥xk − x0 ∥) ∥xk − x0 ∥ + ∥∇f (x0 )∥∗ ≤ L0 rk eL1 rk + L1 rk eL1 rk ∥∇f (x0 )∥∗ + ∥∇f (x0 )∥∗ . 29
Substituting this gradient bound into our transient sum, we decompose it into three distinct sums to bound individually: kX 0 −1
η̂ −1/4 E [∥∇f (xk )∥∗ ] ≤ exp 3L η̂(k + 1) 1 (k + 1)1/4 k=0 kX 0 −1 η̂ 2 3L1 η̂ 33L1 L0 rk eL1 rk exp 1/4 (k + 1) (k + 1) k=1 | {z } 33L1 ηk
:=⑥
+ 33L1
kX 0 −1 k=1
η̂ 2 exp (k + 1)
|
3L1 η̂ (k + 1)1/4 {z
3L1 η̂ (k + 1)1/4 {z
L1 rk eL1 rk E [∥∇f (x0 )∥∗ ] }
:=⑦
+ 33L1
kX 0 −1 k=0
η̂ 2 exp (k + 1)
|
E [∥∇f (x0 )∥∗ ] . }
:=⑧
Bounding ⑥. We substitute the bound rk ≤ 4η̂k 1/4 − 3η̂ and factor out the constants. Noting that 1 k /4 1 k+1 ≤ (k+1)3/4 , we approximate the summation with an integral up to k0 : ⑥ ≤ 33L1
kX 0 −1 k=1
η̂ 2 exp (k + 1)
3L1 η̂ (k + 1)1/4
L0 (4η̂k /4 − 3η̂) exp(L1 (4η̂k /4 − 3η̂)) 1
1
kX 0 −1
1 1 exp(4L1 η̂k /4 ) 3/4 (k + 1) k=1 Z k0 1 1/4 ≤ 132L1 L0 η̂ 3 exp(4L1 η̂k0 ) dt 3/4 t 1 ≤ 132L1 L0 η̂ 3
1/4
1/4
≤ 528L1 L0 η̂ 3 exp(4L1 η̂k0 )k0 . 1/4
Since k0 ≤ 70L1 η̂, we obtain the coarse uniform bound: 2
2
⑥ ≤ 36960L21 η̂ 2 e280η̂ L1 L0 η̂ 2 . Bounding ⑦. The second term mimics the structure of the first but scales with L1 and the initial gradient norm. Following the exact same integral approximation: ⑦ ≤ 33L21 η̂ 2
kX 0 −1 k=1
1 3 1 1/4 − 3 E [∥∇f (x0 )∥∗ ] (4η̂k /4 − 3η̂) exp L1 η̂ + 4k (k + 1) (k + 1)1/4
kX 0 −1
1 1 exp(L1 η̂ · 4k /4 )E [∥∇f (x0 )∥∗ ] 3/4 (k + 1) k=1 Z k0 1 1/4 ≤ 132L21 η̂ 3 exp(4L1 η̂k0 ) dt E [∥∇f (x0 )∥∗ ] 3/4 t 1
≤ 132L21 η̂ 3
1/4
1/4
≤ 528L21 η̂ 3 exp(4L1 η̂k0 )k0 E [∥∇f (x0 )∥∗ ] . Again, evaluating at the threshold k0 , we obtain: ⑦ ≤ 36960L31 η̂ 4 exp(280L21 η̂ 2 )E [∥∇f (x0 )∥∗ ] . Bounding ⑧. The third term lacks the rk multiplier, rendering it a pure harmonic series weighted by an exponential. Bounding the exponential globally by its value at k = 0 and applying the standard 30
Pk0 −1
harmonic bound
k=0
1 k+1 ≤ log(k0 ) + 1:
⑧ ≤ 33L1 η̂ 2 exp(3L1 η̂)
kX 0 −1 k=0
1 E [∥∇f (x0 )∥∗ ] (k + 1)
2
≤ 33L1 η̂ exp(3L1 η̂)(log k0 + 1)E [∥∇f (x0 )∥∗ ]
≤ 132L1 η̂ 2 exp(3L1 η̂) max{1, log⌈70L1 η̂⌉}E [∥∇f (x0 )∥∗ ] .
PK−1 Final assembly of the bound. Combining T1 , T2 , T3 , moving the 12 k=k0 ηk E [∥∇f (xk )∥∗ ] to the left-hand side, and inserting our explicit bounds for the transient phase (⑥, ⑦, ⑧), we yield the global inequality: K−1
1 X 3 ηk E [∥∇f (xk )∥∗ ] ≤ ∆0 + L0 η̂ 2 eL1 η̂ + 36L0 η̂ 2 eL1 η̂ (log K + 1) 2 2 k=0
+ 72L0 η̂ 2 e1.03L1 η̂ (log K + 1) + 11ρσ η̂(log K + 1) 2
2
+ 36960L21 η̂ 2 e280η̂ L1 L0 η̂ 2 2 2
+ 36960L31 η̂ 4 e280L1 η̂ E [∥∇f (x0 )∥∗ ] + 132L1 η̂ 2 e3L1 η̂ max{1, log⌈70L1 η̂⌉}E [∥∇f (x0 )∥∗ ] . By Lemma A.2, we have 2
∥∇f (x0 )∥∗ ≥ 4∆0 ≥ L0 + L1 ∥∇f (x0 )∥∗
(
∥∇f (x0 )∥2∗ , 2L0 ∥∇f (x0 )∥∗ , 2L1
0 ∥∇f (x0 )∥∗ ≤ L L1 ;
otherwise.
Thus, we have
L0 ∥∇f (x0 )∥∗ ≤ max 8L1 ∆0 , L1
.
Therefore, we obtain K−1
3 1 X ηk E [∥∇f (xk )∥∗ ] ≤ ∆0 + L0 η̂ 2 eL1 η̂ + 36L0 η̂ 2 eL1 η̂ (log K + 1) 2 2 k=0
+ 72L0 η̂ 2 e1.03L1 η̂ (log K + 1) + 11ρσ η̂(log K + 1) 2
2
+ 36960L21 η̂ 2 e280η̂ L1 L0 η̂ 2 2 2
2 2
+ 295680L41 η̂ 4 e280L1 η̂ ∆0 + 36960L21 η̂ 2 e280L1 η̂ L0 η̂ 2 + 1056L21 η̂ 2 e3L1 η̂ max{1, log⌈70L1 η̂⌉}∆0 + 132e3L1 η̂ max{1, log⌈70L1 η̂⌉}L0 η̂ 2 . Standard smoothness.
In this case we have L1 = 0. Thus, taking η̂ = η we have
K−1
1 X ηk E [∥∇f (xk )∥∗ ] ≤ ∆0 + 11ρση(log K + 1) + 228L0 η 2 (log K + 1) . 2 k=0
PK−1 Dividing by k=0 ηk , collecting the terms and consolidating coefficients using Big-O notation, we arrive at the envelope Ψ(L0 , 0) = ∆0/η + ρσ + L0 η. Finally, bounding the minimum gradient by the weighted average over K steps: min k∈{0,...,K−1}
1
K−1 X
k=0
ηk k=0
E [∥∇f (xk )∥∗ ] ≤ PK−1
31
ηk E [∥∇f (xk )∥∗ ] ≤ O
Ψ(L0 , 0) log K K 1/4
.
Generalized smoothness: L1 is unknown. In this case we have L1 > 0. Thus, taking η̂ = η/17, we obtain K−1 X k=0
2 2 2 2 ηk E [∥∇f (xk )∥∗ ] ≤ O eL1 η ∆0 + ρσ η̂(log K + 1) + eL1 η L0 η̂ 2 (log K + 1) .
PK−1 Dividing by k=0 ηk and collecting the terms, we arrive at the envelope Ψ(L0 , L1 ). Finally, bounding the minimum gradient by the weighted average over K steps: K−1 X 1 Ψ(L0 , L1 ) log K min E [∥∇f (xk )∥∗ ] ≤ PK−1 , ηk E [∥∇f (xk )∥∗ ] ≤ O k∈{0,...,K−1} K 1/4 k=0 ηk k=0 which completes the proof. Generalized smoothness: L1 is known. Thus, taking η̂ = 1/17L1 , we obtain K−1 X σ L0 ηk E [∥∇f (xk )∥∗ ] ≤ O ∆0 + ρ (log K + 1) + 2 (log K + 1) . L1 L1 k=0
PK−1 Dividing by k=0 ηk and collecting the terms, we arrive at the envelope Ψ(L0 , L1 ) = L1 ∆0 + 0 ρσ + L L1 . Finally, bounding the minimum gradient by the weighted average over K steps: min k∈{0,...,K−1}
1 E [∥∇f (xk )∥∗ ] ≤ PK−1 k=0
K−1 X
ηk k=0
ηk E [∥∇f (xk )∥∗ ] ≤ O
which completes the proof.
32
Ψ(L0 , L1 ) log K K 1/4
,
D
Time complexity analysis under fixed computation times
We now prove the time complexity bounds under Assumption 2.4. We start with the case of a fixed delay threshold. Lemma 4.2 (Duration of R updates [Maranjyan et al., 2025d, Lemma 4.1]). Let the delay threshold of Algorithm 1 be fixed as Rk ≡ R. Under Assumption 2.4, the time needed to complete any R consecutive updates is at most R t(R) := 2 min Hm 1 + . m m∈[n] As an immediate corollary, if the threshold remains fixed at R for the first K updates, then partitioning these updates into ⌈K/R⌉ consecutive blocks and applying Lemma 4.2 yields K T (K) ≤ t(R) . R This blockwise estimate is the key ingredient in the proof of Theorem 4.3. D.1
Proof of Theorem 4.3
We restate the theorem for convenience. Theorem 4.3 (Fixed-parameter total time complexity). Under the assumptions and parameter choices of Theorem 4.1, under Assumption 2.4, the time complexity of reaching an ε-stationary point satisfies √ L0 ∆0 L1 ∆0 L0 ∆0 L0 ∆0 (ρσ)2 L21 ∆0 (ρσ)2 (ρσ)2 Tε = O min Hm + + + + + . ε2 ε ε mε4 mL0 ε2 mε2 m∈[n] Proof. Define L0 ∆0 L1 ∆0 , Bε := , ε2 ε L0 ∆0 (ρσ)2 L21 ∆0 (ρσ)2 (ρσ)2 := := Cε := , D , E . ε ε ε4 L0 ε2 ε2 and let K ε := Aε + Bε + Cε + Dε + Eε . By Theorem 4.1, it suffices to take Kε = cK ε for a universal constant c ≥ 1. Aε :=
Recall that α=
1,
σ = 0,
min 1,
1/2 1/2 ∆0 L 0 1/2 ρσKε
,
σ > 0,
R=
1 . α
We partition the Kε updates into ⌈Kε/R⌉ consecutive blocks and apply Lemma 4.2. This gives Kε 1 Tε ≤ t(R) = 2 ⌈Kε α⌉ min Hm 1 + . R mα m∈[n] Up to constants hidden in O(·), we can ignore the ceiling, so the last inequality gives Kε . Tε = O min Hm Kε α + m m∈[n] It remains to bound Kε α. Since Kε = cK ε , it suffices to bound K ε α. If σ = 0, then Cε = Dε = Eε = 0 and α = 1, so K ε α = Aε + Bε . If σ > 0, then α ≤ 1 and (Aε + Bε )α ≤ Aε + Bε . Moreover, since c ≥ 1 and Cε , Dε , Eε are terms in K ε , we have Cε , Dε , Eε ≤ Kε , and √ √ L0 ∆0 L0 ∆0 Cε α ≤ C ε √ ≤ Cε √ = Aε , ρσ Kε ρσ Cε √ √ L 0 ∆0 L0 ∆0 Dε α ≤ Dε √ ≤ Dε √ = Bε , ρσ Kε ρσ Dε √ √ √ L0 ∆0 L0 ∆0 L0 ∆0 Eε α ≤ E ε √ ≤ Eε √ = =: Gε . ε ρσ Kε ρσ Eε 33
Hence K ε α = O(Aε + Bε + Gε ), and therefore Kε α = O(Aε + Bε + Gε ). Also, Kε Aε + B ε + C ε + D ε + E ε Cε + D ε + E ε =O = O Aε + B ε + , m m m where the last equality uses m ≥ 1. Combining the last three bounds gives the claimed bound. The parameter-agnostic case is slightly more delicate, since the threshold does not remain fixed throughout the run. To handle this, we extend the same blockwise argument to a non-decreasing threshold sequence. Lemma D.1 (Time for non-decreasing adaptive threshold). Let K ≥ 1, and let Rk ≥ 1 be a non-decreasing delay threshold. For each integer r ∈ {⌊R0 ⌋, . . . , ⌊RK−1 ⌋}, let Kr := {k ∈ {0, . . . , K − 1} | ⌊Rk ⌋ = r}
denote the corresponding block of iterations, and define |Kr | Nr := . r Then the time needed by Algorithm 1 to complete its first K updates satisfies ⌊RK−1 ⌋
T (K) ≤
X
Nr t(r) .
r=⌊R0 ⌋
Proof. Split the iterations {0, 1, . . . , K − 1} into the blocks Kr . Since Rk is non-decreasing, all iterations in each Kr appear next to each other in time. For a fixed r, partition Kr into |Kr | Nr = r consecutive groups, each containing at most r iterations. Throughout each such group, ⌊Rk ⌋ = r, and hence Rk ≥ r. Replacing the actual threshold by the fixed threshold r can only make the rule more conservative. If a group contains exactly r iterations, then Lemma 4.2 shows that it is completed within time at most t(r). If the last group contains fewer than r iterations, its duration is still at most t(r), since completing fewer than r updates cannot take longer than completing r updates under the same threshold. Therefore, the total time spent on all iterations in Kr is at most Nr t(r). Summing this bound over all values of r gives ⌊RK−1 ⌋
T (K) ≤
D.2
X
Nr t(r) .
r=⌊R0 ⌋
Proof of Lemma 4.5
We now specialize Lemma D.1 to the square-root threshold schedule. Lemma√4.5 (Time complexity for square-root threshold). Let the delay threshold be R0 = 1 and Rk = k for k ≥ 1. Under Assumption 2.4, for any K ≥ 1, the time needed to complete K iterations satisfies √ K T (K) = O min Hm K+ . m m∈[n] Proof. Apply Lemma D.1 with R0 = 1 and Rk = k ≥ 1, the condition ⌊Rk ⌋ = r is equivalent to
√
k for k ≥ 1. For any integer r ≥ 1 and any
r2 ≤ k < (r + 1)2 .
The extra iteration k = 0 belongs to K1 . Therefore, for every r ≥ 1, |Kr | ≤ 2r + 2 . 34
Hence, Nr = Let rmax := ⌊RK−1 ⌋. Since rmax ≤ T (K) ≤ Since
P
rX max r=1
√
|Kr | 2r + 2 ≤ ≤4. r r
K, using Lemma D.1 and Lemma 4.2, we obtain rX max
Nr t(r) ≤
r minm ar,m ≤ minm
4 t(r) = 8
r=1
rX max
r min Hm 1 + . m m∈[n] r=1
P
r ar,m for nonnegative ar,m , it follows that
T (K) ≤ 8 min Hm m∈[n]
= 8 min m∈[n]
We conclude that
rX max
1+
r=1
rmax (rmax + 1) Hm rmax + Hm 2m
T (K) = O
r m
min Hm
m∈[n]
This proves the claim.
35
√
K K+ m
.
.
E
Time complexity under time-varying computation rates
Assumption 2.4 assumes that each worker has a constant gradient-computation time throughout the run. This leads to explicit time complexity bounds in terms of the harmonic means Hm , but it does not capture clusters whose worker speeds change over time because of interruptions or fluctuating load [Dutta et al., 2018, Li et al., 2020, Maranjyan et al., 2025c]. In this section we therefore move to a more general model with time-varying computation rates. E.1
Universal computation model
We use the universal computation model of Tyurin [2024]. Assumption E.1 (Universal computation model). For each worker i ∈ [n], there exists a nonnegative function pi : R+ → R+ that is continuous almost everywhere. For any 0 ≤ T1 ≤ T2 , the number of stochastic gradients completed by worker i during the interval [T1 , T2 ] is $Z % T2
Ni (T1 , T2 ) :=
pi (s) ds . T1
This assumption includes Assumption 2.4 as a special case: if pi (s) = 1/τi for all s ≥ 0, then T2 − T1 Ni (T1 , T2 ) = . τi E.2
Fixed-parameter schedule
We begin with the fixed-parameter regime, which is the counterpart of Theorem 4.3 under Assumption E.1. Unlike in the fixed computation model, the duration of a block of R updates now depends on the starting time of that block. For this reason, the resulting time bound is recursive rather than explicit. Lemma E.2 (Duration of R updates under the universal computation model [Maranjyan et al., 2025d, Lemma 5.1]). Assume Assumption E.1. Let the delay threshold of Algorithm 1 be fixed as Rk ≡ R. Assume that an iteration starts at time T0 . Then the next R updates of Algorithm 1 are completed within elapsed time % ) ( $ Z n X 1 T0 +t pi (s) ds ≥ R . t(R; T0 ) := min t ≥ 0 : 4 T0 i=1 Theorem E.3 (Fixed-parameter time complexity under the universal computation model). Under the assumptions and parameter choices of Theorem 4.1, let Kε be the number of iterations needed to reach an ε-stationary point, and let R = 1/α be the corresponding fixed delay threshold. Then, under Assumption E.1, the time complexity of reaching an ε-stationary point satisfies Tε ≤ T⌈Kε /R⌉ , where the sequence {Tk }k≥0 is defined recursively by T0 := 0, Tk := Tk−1 + t(R; Tk−1 ),
k ≥ 1.
Proof. Partition the first Kε updates into ⌈Kε /R⌉ consecutive blocks of length at most R. Applying Lemma E.2 blockwise yields the stated recursion and therefore the bound Tε ≤ T⌈Kε /R⌉ .
36
Here Kε denotes the iteration complexity from Theorem 4.1, namely, the number of updates needed to reach an ε-stationary point. Thus the theorem converts the iteration bound into a time bound by grouping these Kε updates into blocks of size R and tracking their completion times through the recursion above. In the classical smooth case L1 = 0, the iteration complexity from Theorem 4.1 has the same ε-dependence as Ringmaster ASGD, up to the norm-equivalence factor ρ. In the Euclidean setting, where ρ = 1, this recovers the Ringmaster ASGD scaling and hence the corresponding optimal time complexity. We next extend Lemma E.2 to non-decreasing threshold sequences. E.3
Adaptive threshold schedule
We now turn to the parameter-agnostic schedule, where the threshold increases with the iteration counter. We first derive a generic recursive bound for √ arbitrary non-decreasing thresholds and then specialize it to the square-root rule R0 = 1 and Rk = k for k ≥ 1. Lemma E.4 (Time for general adaptive threshold under the universal computation model). Let Rk be a non-decreasing delay threshold. For each integer r ∈ {⌊R0 ⌋, . . . , ⌊RK−1 ⌋}, let Kr := {k | ⌊Rk ⌋ = r}
denote the block of iterations where the threshold is constant. For every such r, partition Kr into |Kr | := Nr r consecutive sub-blocks, each containing at most r iterations. P⌊RK−1 ⌋ Let B := r=⌊R Nr , and enumerate all these sub-blocks in chronological order. Denote by rj the 0⌋ threshold value attached to the j-th sub-block, where j = 1, . . . , B. Define the sequence {Tj }B j=0 recursively by T0 := 0, Tj := Tj−1 + t(rj ; Tj−1 ),
j = 1, . . . , B.
Then the time needed to complete the first K iterations of Algorithm 1 satisfies T (K) ≤ TB . Proof. Partition the iterations {0, 1, . . . , K − 1} into disjoint blocks Kr according to the integer value of ⌊Rk ⌋. Since Rk is non-decreasing, all iterations in each Kr appear next to each other in time. Fix some r and consider one of its sub-blocks. By construction, this sub-block contains at most r consecutive iterations, and throughout this sub-block the threshold satisfies r ≤ Rk < r + 1.
Replacing the actual threshold by r can only make the rule more conservative on this sub-block. Hence Lemma E.2 implies that, if this sub-block starts at time T , then all its iterations are completed within time at most t(r; T ). Now enumerate all sub-blocks in chronological order and let rj be the threshold associated with the j-th one. Starting from time T0 = 0, the first sub-block is completed by time at most T1 = T0 + t(r1 ; T0 ). Repeating the same argument inductively, if the first j − 1 sub-blocks are completed by time at most Tj−1 , then the j-th sub-block is completed by time at most Tj = Tj−1 + t(rj ; Tj−1 ). After all
⌊RK−1 ⌋
B=
X
r=⌊R0 ⌋
|Kr | r
sub-blocks are processed, all first K iterations have been completed. Therefore, T (K) ≤ TB .
37
This is the universal-model analogue of Lemma D.1: it reduces the analysis of a non-decreasing threshold sequence to a recursion over threshold-constant sub-blocks. Lemma √ E.5 (Square-root threshold under the universal computation model). Let R0 = 1 and Rk = k for k ≥ 1. Under Assumption E.1, the time needed to complete the first K iterations of Algorithm 1 is bounded by T (K) ≤ S3⌊√K⌋+1 , √ 3⌊ K⌋+1
where the sequence {Sj }j=0
is defined recursively as
S0 := 0, j−1 Sj := Sj−1 + t max 1, ; Sj−1 , 3 Proof. We apply Lemma E.4 √ with R0 = 1 and Rk = every k ≥ 1, the condition ⌊ k⌋ = r is equivalent to
√
√ j = 1, . . . , 3⌊ K⌋ + 1.
k for k ≥ 1. For every integer r ≥ 1 and
r2 ≤ k < (r + 1)2 ,
The extra iteration k = 0 belongs to K1 . Hence
|Kr | ≤ 2r + 2
for all r ≥ 1. Therefore N1 ≤ 4, while Nr ≤ 3 for all r ≥ 2. √ Let rmax := ⌊ K⌋. Hence the total number of sub-blocks is at most B ≤ 4 + 3(rmax − 1) = 3rmax + 1.
Moreover, after the first threshold level, each new level contributes at most three additional sub-blocks, so the threshold attached to the j-th sub-block is at most j−1 max 1, . 3 Therefore, using the monotonicity of t(r; T ) in r, the recursion from Lemma E.4 is dominated by S0 := 0, j−1 ; Sj−1 , Sj := Sj−1 + t max 1, 3
j = 1, . . . , 3rmax + 1,
and the resulting completion time satisfies T (K) ≤ S3rmax +1 = S3⌊√K⌋+1 . Theorem E.6 (Parameter-agnostic time complexity under the universal computation model). Consider Algorithm 1 with √ 1 η α0 = 1, αk = √ ηk = (k ≥ 1), R0 = 1, Rk = k (k ≥ 1), (k + 1)3/4 k Let
Ψ4 ε4 be the iteration complexity from Theorem 4.4. Then, under Assumption E.1, the time complexity of reaching an ε-stationary point satisfies
e Kε = O
Tε ≤ S3⌊√Kε ⌋+1 = SO(Ψ e 2 /ε2 ) , e hides logarithmic factors in 1/ε. The sequence {Sj } is defined recursively by Here and below, O(·) S0 := 0,
j−1 Sj := Sj−1 + t max 1, ; Sj−1 . 3 38
Proof. By Theorem 4.4, the method reaches an ε-stationary point within 4 Ψ e Kε = O ε4 iterations. Applying Lemma E.5 with K = Kε yields Tε ≤ S3⌊√Kε ⌋+1 .
Since
√
E.4
Comparison with Ringmaster ASGD
e 2 /ε2 ), the claim follows. Kε = O(Ψ
We briefly compare the recursive time bounds obtained for Ringmaster ASGD and Ringmaster LMO under the universal computation model. Because the duration function t(r; T ) depends on both the block size r and the starting time T , this comparison is necessarily weaker than in the fixed computation model. For Ringmaster ASGD, the iteration complexity in Maranjyan et al. [2025d] is 48L∆σ 2 K̄ := . ε4 Take the fixed threshold R :=
lσm
. ε2 Applying the fixed-threshold bound above with this choice of R gives (R)
TεRingmaster ≤ TNRingmaster , where
K̄ NRingmaster := R
(R)
and the sequence {Tk } is defined by (R)
T0
:= 0,
(R) (R) (R) Tk := Tk−1 + t(R; Tk−1 ).
Since K̄ = O(σ 2 L∆/ε4 ) and R = Θ(σ/ε2 ), we obtain σL∆ . NRingmaster = O ε2 For our method, the iteration complexity from Theorem 4.4 is 4 e Ψ , Kε = O ε4 and Theorem E.6 yields the time bound Tε ≤ SNours , where and the sequence {Sj } is defined by
p Nours := 3⌊ Kε ⌋ + 1
S0 := 0,
j−1 Sj := Sj−1 + t max 1, ; Sj−1 . 3
Since
√
e 2 /ε2 ), this gives Kε = O(Ψ e Nours = O 39
Ψ2 ε2
.
Thus, from the point of view of the ε-dependence, the two methods are already quite close: both require on the order of 1/ε2 recursion steps, up to logarithmic factors for our method. Their problemdependent prefactors are different, so this should be interpreted only as an ε-scaling comparison. The difference lies in how these blocks are formed. Ringmaster ASGD uses the same threshold R in every block, while our method uses the increasing thresholds p p p 1, 1, 1, 1, 2, 2, 2, 3, 3, 3, . . . , ⌊ Kε ⌋, ⌊ Kε ⌋, ⌊ Kε ⌋. Because the quantity t(r; T ) depends on both the threshold r and the starting time T , the universal computation model does not by itself give a constant-factor comparison between (R)
TNRingmaster
and
SNours .
Indeed, the worker speeds may change arbitrarily with time, so two recursions with a similar number of steps need not be within a constant factor of each other. Therefore, in the fully general universal model, the clean conclusion is that both methods involve about 1/ε2 recursion steps, but we cannot claim that our bound matches the Ringmaster ASGD bound up to a universal constant without extra assumptions on the computation-rate functions {pi }ni=1 . If one wants a genuine constant-factor time comparison, a sufficient extra condition is that block durations of the relevant size remain stable under shifts of the starting time. Let p Rε⋆ := ⌊ Kε ⌋, and assume that there exists a constant C⋆ ≥ 1 such that t(Rε⋆ ; T ) ≤ C⋆ t(Rε⋆ ; 0)
Since
j ≤ Rε⋆ 3
for all T ≥ 0.
for j = 1, . . . , Nours ,
and t(r; T ) is nondecreasing in r, every increment in the recursion for {Sj } satisfies j ; Sj−1 ≤ t(Rε⋆ ; Sj−1 ) ≤ C⋆ t(Rε⋆ ; 0). t 3 Therefore,
Ψ2 p t( K ; 0) . ε ε2 In other words, under this localized start-time stability assumption, our adaptive schedule has the same time scaling as a fixed-threshold scheme at the natural block size Rε⋆ , up to constant and logarithmic factors. e Tε ≤ SNours ≤ C⋆ Nours t(Rε⋆ ; 0) = O
40