Paper deep dive
Muon Is Theoretically Wrong For Convolutions, But Empirically Effective
Thibaut Boissin, Thomas Massena, Mathieu Serrurier, Franck Mamalet
Intelligence
Status: succeeded | Model: Gemma-4-26B-A4B | Prompt: intel-v1 | Confidence: 96%
Last extracted: 10/7/2026, 8:42:37 AM
Summary
This paper investigates the application of the Muon optimizer to convolutional neural networks. It highlights that standard Muon implementations reshape 4D convolutional kernels into 2D matrices, which breaks the theoretical geometric alignment required for proper orthogonalization. The authors propose Convolutional Newton-Schulz (Conv-NS), a method that approximates the polar factor directly in convolutional operator geometry while preserving kernel support. Experiments on CIFAR-10 and ImageNet show that Conv-NS achieves comparable accuracy to the reshape-based Muon, despite the latter being theoretically misaligned. The study suggests that exact convolutional orthogonalization may overconstrain updates, explaining why the theoretically correct method does not outperform the heuristic.
Entities (11)
Relation Signals (7)
Convolutional Newton-Schulz → evaluatedon → CIFAR-10
confidence 99% · achieve comparable accuracy on CIFAR-10 and ImageNet classification tasks
Convolutional Newton-Schulz → evaluatedon → ImageNet
confidence 99% · achieve comparable accuracy on CIFAR-10 and ImageNet classification tasks
Convolutional Newton-Schulz → comparedwith → Muon
confidence 98% · Conv-NS and reshape-based Muon are both computationally efficient and achieve comparable accuracy
ResNet-50 → usedin → ImageNet
confidence 98% · a standard ResNet-50... is trained on ImageNet-1k
Convolutional Newton-Schulz → proposedby → Thibaut Boissin
confidence 95% · To investigate this, we formalize... and introduce Convolutional Newton-Schulz (Conv-NS)... Thibaut Boissin
Muon → uses → Reshape-Kernel Orthogonalization
confidence 95% · Standard implementations reshape these tensors into matrices... RKO(G K ) = unflat PF(flat(G K ))
Reshape-Kernel Orthogonalization → critique → Muon
confidence 90% · a shortcut which breaks the theoretical understanding behind Muon
Cypher Suggestions (0)
No Cypher suggestions yet.
Abstract
Abstract:Muon, an optimizer known for its efficiency, has a clear interpretation for matrix-valued updates, but convolutional kernels are stored as four-dimensional tensors. Standard implementations reshape these tensors into matrices, a shortcut which breaks the theoretical understanding behind Muon. To investigate this, we formalize the corresponding optimization objective directly in convolutional operator geometry and introduce Convolutional Newton-Schulz (Conv-NS), which approximates the polar factor in this geometry while preserving kernel support. When applied in fast training experiments, Conv-NS and reshape-based Muon are both computationally efficient and achieve comparable accuracy on CIFAR-10 and ImageNet classification tasks. However, as one could expect a theoretically aligned Conv-NS to outperform reshape-based Muon, we investigate this mismatch between practice and theoretical understanding, with the hypothesis that exact convolutional orthogonalization may overconstrain updates. These findings highlight Muon's strong practical performance while opening directions for its further development on convolutions. Our code is publicly available at \href{this https URL}{github conv-muon}.
Tags
Links
- Source: https://arxiv.org/abs/2610.07103v1
- Canonical: https://arxiv.org/abs/2610.07103v1
Trouble viewing inline? Open PDF directly →
Full Text
79,759 characters extracted from source content.
Expand or collapse full text
MUON IS THEORETICALLY WRONG FOR CONVOLU- TIONS, BUT EMPIRICALLY EFFECTIVE Thibaut Boissin 1,2 Thomas Massena 2,3 Mathieu Serrurier 2 Franck Mamalet 1 1 Institut de Recherche Technologique Saint-Exupery, France 2 IRIT, France 3 Innovation & Research Division, SNCF, France ABSTRACT Muon, an optimizer known for its efficiency, has a clear interpretation for matrix- valued updates, but convolutional kernels are stored as four-dimensional tensors. Standard implementations reshape these tensors into matrices, a shortcut which breaks the theoretical understanding behind Muon. To investigate this, we formalize the corresponding optimization objective directly in convolutional operator geome- try and introduce Convolutional Newton-Schulz (Conv-NS), which approximates the polar factor in this geometry while preserving kernel support. When applied in fast training experiments, Conv-NS and reshape-based Muon are both computa- tionally efficient and achieve comparable accuracy on CIFAR-10 and ImageNet classification tasks. However, as one could expect a theoretically aligned Conv-NS to outperform reshape-based Muon, we investigate this mismatch between prac- tice and theoretical understanding, with the hypothesis that exact convolutional orthogonalization may overconstrain updates. These findings highlight Muon’s strong practical performance while opening directions for its further development on convolutions. Our code is publicly available at github conv-muon. 1INTRODUCTION A recent line of work in deep learning optimization departs from element-wise update rules such as SGD and AdamW and instead exploits the matrix structure of neural network parameters. This perspective has led to influential methods ranging from preconditioned second-order optimizers such as Shampoo (Gupta et al., 2018; Anil et al., 2020), which accumulate Kronecker-factored curvature estimates, to first-order methods such as Muon (Jordan et al., 2024b; Bernstein, 2025). Muon replaces the gradient of a linear layer by an approximation of its semi-orthogonal polar factor, effectively promoting well-conditioned, near-orthogonal updates at the operator level. This optimizer has demonstrated competitive performance across tasks such as transformer pretraining (Jordan et al., 2024a; Karpathy, 2025), image diffusion model training (Schaipp, 2025), and large language model training (AI et al., 2025; Jordan et al., 2024a; Liu et al., 2025). We provide a more detailed discussion of related work in Section A. While this construction is natural for linear layers, convolutional layers do not fit trivially into the same framework. In this paper, we study the use of Muon on 2D convolutions. A convolutional layer is parameterized by a compact kernelK ∈ R c out ×c in ×k h ×k w , but acts on activations through a much larger structured linear operatorT (K), typically Toeplitz or doubly circulant depending on the boundary conditions, that need not be materialized during the forward pass to avoid redundant computations. In practice, current implementations apply Muon by reshaping the 4D kernel into a matrix of size R c out ×(c in k h k w ) and orthogonalizing this reshaped matrix. However, orthogonalizing the reshaped kernel does not generally correspond to orthogonalizing the equivalent linear transformation matrix (see fig. 1, upper part). The resulting update is therefore misaligned with the operator geometry of the layer being optimized. We formalize this mismatch in section 2 and introduce the first convolution-aware construction that directly targets the polar factor of the induced operator while preserving the finite support of the kernel. In this paper, we use orthogonalization to refer to procedures that drive the singular values of the convolutional operator toward one (as seen in fig. 1, right). 1 arXiv:2610.07103v1 [cs.LG] 5 Oct 2026 Figure 1: Problem framing for Muon on convolutions. Muon’s current implementation (option 1) reshape the update and applyPF, however this result in a non orthogonal convolution as observed on the singular value distibution. However there is no trivial solution since the direct application ofPF to the equivalent linear operator leads to a an orthogonal matrix, but the structure of the convolution is lost in the process. The3× 3convolution withc i = 2,c o = 1was sampled from a gaussian distribution. Independently from the less generally understood question of generalization (Wilson et al., 2017; Balles & Hennig, 2018; Luo et al., 2019), we run benchmarks to prioritize speed of convergence (section 5). Surprisingly, reshape-based Muon, despite its structural mismatch, and our convolution- aware alternative Conv-NS achieve comparable performance on these benchmarks, while both substantially outperforming SGD with momentum in compute-constrained settings. These results raise the following question: why does more faithful convolutional orthogonalization not improve optimization performance? The insights we provide are the following: •Standard Muon overlooks convolutional geometry. Applying Muon to convolutions typically involves reshaping each kernel into a matrix before orthogonalization. This does not generally orthogonalize the convolution operator. •Convolutional orthogonalization is computationally practical. We introduce Conv-NS, which preserves kernel support and empirically approaches convolutional orthogonality with near-minimal gradient modification. •Better orthogonalization does not consistently improve training. Conv-NS and reshape- based Muon achieve competitive performance on CIFAR-10 and ImageNet, while increasing Conv-NS iterations can degrade training performance. • Extending Muon to convolutions requires more than polar factorization. We revisit the optimization principle underlying Muon: choosing an update that maximizes predicted loss reduction while limiting its effect on layer outputs. Reshaping does not enforce this limit, while convolutional orthogonality imposes a stronger constraint. •These limitations open directions for convolutional Muon. Both implementations already achieve strong accuracy–runtime trade-offs in the tested settings, motivating updates that better balance predicted loss reduction and control of layer outputs. 2MUON FOUNDATIONS AND THE CONVOLUTION MISMATCH Notations For any matricesX ∈ R m×n , we denote by∥X∥ F its Frobenius norm, by⟨X,Y⟩ F the associated inner product, and by∥X∥ op its operator (or spectral) norm. The Stiefel manifold is denoted whenm ≤ n, bySt(m,n) := X ∈ R m×n : X ⊤ = I 1 , whereIdenotes them× m identity matrix. When the dimensions are clear from the context, we simply writeSt. We denote the squared Frobenius distance between two matrices X,Y ∈ R m×n , byL frob (X,Y ) =∥X − Y∥ 2 F . 2 2.1ABOUT THE MUON UPDATE RULE A dense hidden layer acts on activations via a weight matrixW ∈ R d out ×d in . Given a momentum or gradient matrixG, the Muon optimizer (Jordan et al., 2024b) computes the update direction by approximately orthogonalize the gradient matrix. IfG = U ΣV ⊤ is the singular value decomposition, then the Muon update is: ∆ ⋆ W =−η PF(G) =−ηUV ⊤ (1) From an optimization perspective, the polar factorPF(G)can equivalently be viewed as the projection ofGonto the Stiefel manifoldSt, and thus defined by a constrained minimization over the objective L frob (G,Q): PF(G) = argmin Q∈St L frob (G,Q)(2) Since computing an exact SVD at every step is computationally prohibitive, Muon practically approximatesPF(G)using a small number of Newton-Schulz iterations (Bj ̈ orck & Bowie, 1971; Amsel et al., 2025; Grishina et al., 2024) 2 . For an appropriately scaled matrix X , the iteration: X ← (1 + β)X − βX ⊤ X(3) drives the orthogonality penalty∥X ⊤ − I∥ 2 F iteratively to zero, yielding the polar factor. 2.2HOW STANDARD MUON HANDLES CONVOLUTIONS Convolutional layers do not fit cleanly into this dense picture. The trainable parameters and their gradient updates are finite-support 4D kernels:G K ∈M kernel , whereM kernel = R c out ×c in ×k h ×k w Unfortunately, materializing the equivalent linear transformation to the convolution is infeasible at scale. Therefore, current implementations of Muon sidestep this by using a structural shortcut: they treat the 4D kernel as a 2D matrix by flattening its spatial dimensions, apply the exact same dense Newton-Schulz orthogonalization routine, and unflatten the result: Denoting the flattening and unflattening operators asflat : R c out ×c in ×k h ×k w → R c out ×(c in k h k w ) andunflat : R c out ×(c in k h k w ) → R c out ×c in ×k h ×k w , the Reshape-Kernel Orthogonalization (RKO) update is: RKO(G K ) = unflat PF(flat(G K )) (4) This formulation effectively enforces an orthogonality constraint on the flattened kernel matrix. We can formalize the surrogate objective implicitly minimized by RKO asL flat : L flat (∆ K ) =∥ flat(∆ K ) flat(∆ K ) ⊤ − I∥ 2 F (5) By drivingL flat to zero, RKO makes the flattened matrix perfectly semi-orthogonal. While compu- tationally convenient,flat(G K )ignores the true spatial and channel interaction of the convolution, optimizing a mismatched geometric surrogate. For example, it cannot distinguish a3× 3convolution withc in = 8andc out = 8from a pointwise1× 1convolution withc in = 72andc out = 8, even though those objects act fundamentally differently on activations. Notably, Wang et al. (2020); Li et al. (2019) show that RKO does not lead to an orthogonal linear operator, as illustrated in the upper part of fig. 1. 3OUR CORRECTION: THE CONVOLUTIONAL POLAR FACTOR To correctly adapt Muon to convolutions, we must work in the space in which the layer actually acts. For a given input activation map, convolution by a kernelKis formally equivalent to a structured linear operatorT (K)∈ R (c out h out w out )×(c in h in w in ) , typically a doubly circulant or Toeplitz matrix depending on boundary conditions (Achour et al., 2022). The equivalent Muon update for a convolution layer would therefore be the polar factor (Equation (2)) of the induced operator:PF(T (G K )) ∈ R (c out h out w out )×(c in h in w in ) . As settled for dense matrices, the polar factor of this polar factorPF(T (G K ))could be approximated by iterative algorithms 1 the condition for column orthogonal matrices n≤ m is X ⊤ X = I , with I of size n× n 2 Variants with fifth order matrix multiplication is often used 3 (Equation (3)). But, first applying Newton-Schultz on large matrices likeT (G K )would be costly, and the resulting polar factor generally loses its structural sparsity:PF(T (G K ))cannot be converted back into a validM kernel update (lower part of fig. 1). Therefore, a principled definition of the convolutional polar factor requires an extra constraint: the solution operator must remain in the subspace of fixed-support convolutionsM kernel : ConvPF(G K ) = argmin ∆ K ∈M kernel T (∆ K )∈St L frob (T (G K ),T (∆ K )).(6) Closeness to the Stiefel manifold is assessed through the semi-orthogonality residual ∥T (∆ K )T (∆ K ) ⊤ − I∥ 2 . Wang et al. (2020); Qi et al. (2020); Achour et al. (2022) have shown that this analogously translates to a spatial self-convolution identity on the underlying kernel. We denote this kernel-level orthogonality distance asL orth : L orth (∆ K ) =∥ conv ∆ K , ∆ K , padding = P, stride = s − I r0 ∥ 2 F (7) = 1 h out w out ∥T (∆ K )T (∆ K ) ⊤ − I∥ 2 F ,(8) whereI r0 is the identity impulse kernel (zero everywhere except at the central position where it is the identity matrix over channels) andh out w out is the number of spatial locations. The eq. (8), formally demonstrated in Achour et al. (2022), establishes the exact equivalence 3 . TheConvPF(G K )can thus be expressed using theL orth : ConvPF(G K ) =argmin ∆ K ∈M kernel L orth (∆ K )=0 L frob (T (G K ),T (∆ K )).(9) Thus the convolutional polar factor is the equivalent for convolution operators of the target update of Muon for dense matrices. However, neither the SVD decomposition of the operatorT (G K ), nor the RKO heuristic of Muon is able to compute or even approximate ConvPF(G K ). 4CONVOLUTIONAL NEWTON-SCHULZ: A CONVOLUTION-AWARE ALTERNATIVE We introduce Convolutional Newton-Schulz (Conv-NS), an efficient convolution-aware iterative algorithm that converges toConvPF(G K ). Conv-NS serves a dual purpose: it provides a convenient method to orthogonalize the convolution gradient matrix that can be used as a gradient update, evaluated in section 5, and a principled reference for measuring how RKO differs from the true convolutional polar factor. We exploit this second role in section 6 to characterize RKO on real training gradients. The proposed algorithm enables to compute the ground-truth comparison that the existing literature lacks. We propose an iterative projected gradient descent scheme to jointly drive to orthogonality and constrain to theM kernel manifold. After each iterationtof the order-3 Newton-Schulz (also known as Bj ̈ orck) algorithm (Bj ̈ orck & Bowie, 1971), a structural projectionΠ supp is applied to force the kernel back ontoM kernel : T ( e ∆ (t) K )← (1 + β)T (∆ (t) K )− βT (∆ (t) K )T (∆ (t) K ) ⊤ T (∆ (t) K )(10) ∆ (t+1) K ← αΠ supp e ∆ (t) K ,(11) whereβ t is the orthogonalization step size, andαa normalization factor to avoid expansion (justifica- tion of this parameter is given in section B). Without the projection (Equation (11)), the kernel size will grow at each iteration, and would converge to PF(T (G K )) as described in Section 3. In practice, these iterations would require huge matrix multiplications with high computational redundancy due to the doubly circulant or Toeplitz structure ofT (∆ (t) K ). We instead propose to 3 For simplicity, we assume the row-orthogonal (wide) casec out ≤ c in s 2 , analogous to our dense formulation in section 2.1. In the column-orthogonal case,L orth (∆ K )cannot reach0(its global minimum is bounded by c out − c in s 2 ), but its minimizer nonetheless yields a semi-orthogonal convolutional operator. 4 Algorithm 1 Conv-NS Require:Gradient kernelG K , iterationsT, step size β 1: ∆ (0) K ← Normalize(G K ) 2: for t = 0,...,T − 1 do 3: C (t) ← ∆ (t) K ⊛ (∆ (t) K ) ⊤ ⊛ ∆ (t) K 4: e ∆ (t) K ← (1 + β)∆ (t) K − βC (t) 5:∆ (t+1) K ← αΠ supp e ∆ (t) K 6: end for 7: return ∆ (T ) K Figure 2: Conv-NS orthogonalizes convo- lutional updates while preserving kernel support. Each iteration applies a Newton– Schulz step through block convolutions, then crops the result to the original kernel support. 0.000.250.500.751.001.251.501.752.00 Singular value 0 1 2 3 4 5 Density Random 3×3 kernel, 32 channels Random kernel Conv-NS (5 iters) — 1.96 ms Conv-NS (10 iters) — 3.63 ms Conv-NS (25 iters) — 8.79 ms RKO (5 iters) — 0.06 ms Target = 1 Median over 25 runs · CPU, 1 thread · float32 Figure 3: Conv-NS converges to convolutional or- thogonality. Increasing Conv-NS iterations concen- trates the spectrum of singular values around one, whereas RKO, used in standard Muon implementa- tions, leaves it dispersed. All methods start from the same random 3× 3 kernel. use block-convolution to perform the NS-3 scheme directly on the kernels∆ (t) K . Following Li et al. (2019); Boissin et al. (2025), the block convolution operator⊛combines two kernels into a single kernel representing their composition (without materializing their full convolution matrices). For A∈ R c int ×c in ×k A h ×k A w andB ∈ R c out ×c int ×k B h ×k B w , their block convolution is defined with respect to the induced linear operators as: T (B⊛ A) =T (B)T (A) with B⊛ A∈ R c out ×c in ×(k A h +k B h −1)×(k A w +k B w −1) .(12) We also use the definition of the transposed convolution kernelT (∆ ⊤ K ) =T (∆ K ) ⊤ , which can be computed efficiently in kernel space. The complete Convolutional Newton-Schultz algorithm is provided as a kernel formulation in Algorithm 1. Lines 3 and 4 are equivalent to the order-3 Newton-Schulz iteration. Importantly, the block convolution expands the spatial support of the operands, soC (t) kernel size is generally larger than∆ (t) K one. The projectionΠ supp crops the expanded kernel by dropping the values outside the original k h × k w support, to converge to convPF(T (G K )). In Section C, we provide the full proofs, together with an explicit gradient formula, showing that the Conv-NS algorithm converges to a minimum ofL orth . section B provide more details about the exact implementation and how strided and grouped convolutions are covered. In Section D, we also propose convolutional adaptations of the main grafting factors commonly applied to Muon updates. Computing a ground truth for ConvPF. SinceConvPFminimizes gradient modification under convolutional orthogonality (eq. (9)), we construct a numerical reference by relaxing this constraint into a penalty: min ∆ K ∈M kernel λL orth (∆ K ) + (1− λ)L frob (∆ K ,G K ),0≤ λ < 1.(13) Here,G K is the optimizer momentum rescaled to unit operator norm to separate orthogonalization- induced modification from initial normalization. Varyingλtrades fidelity to the original gradient against convolutional orthogonality. Since this objective is nonconvex, we estimate the orthogonality– fidelity frontier using Adam with simulated annealing noise injection for each weighting, retaining the best result across multiple trials (section K.1). This provides a numerical reference for evaluating Conv-NS, rather than a certificate of global optimality. Empirical validation:Figure 3 compares the singular-value spectra of the RKO Muon update and Conv-NS after different numbers of iterations. It shows that Conv-NS progressively approaches an orthogonal convolution operator. We compare Conv-NS against this numerical reference on synthetic 5 non feasible zone suboptimal zone 10 1 10 0 10 1 10 2 L orth 10 1 10 0 10 1 10 2 L frob 2 3 4 5 10 25 100 1 3 5 100 Pareto Conv-NS Muon (RKO) ConvPF( K) K Figure 4: Conv-NS approaches convolutional orthogonality with near-minimal gradient mod- ification. Conv-NS follows the numerically estimated orthogonality–fidelity frontier: stronger orthogonality requires progressively larger departures from the original gradient along this frontier. RKO remains near the low-modification end, with higher orthogonality error. Metrics are averaged across convolutional layers and three stages of CIFAR-10 training. kernels (Figure 8) and real training gradients (Figure 4). For the latter, we collect snapshots at 1%, 50%, and 99% of a CIFAR-10 training run using the architecture of section 5. We evaluate Conv-NS at iteration counts2, 3, 4, 5, 10, 25, 100against the per-layer frontiers and average metrics across layers and stages. Conv-NS remains close to the estimated frontier, approaching orthogonality with near-minimal gradient modification among the solutions found numerically. Increasing the iteration count moves along this trade-off toward stronger orthogonality and larger gradient modification. Per-stage results are reported in section F. 5EXPERIMENTS: CIFAR-10 AND IMAGENET CLASSIFICATION In the following, we compare Conv-NS Muon against the standard reshape-based Muon on several image-classification benchmarks with convolutional architectures. Experimental setup.We evaluate Conv-NS in two classification settings that stress different aspects of convolutional optimization. First, CIFAR-10 AirBench94 (Jordan, 2024) is a speed-optimized benchmark tailored to reach 94% accuracy in less than 4 seconds on a single RTX4090. Second, a standard ResNet-50 with atimm-style recipe shortened to 30 epochs is trained on ImageNet-1k (Deng et al., 2009). In addition to these standard short-training budgets, we evaluate an extreme 1-epoch setting on both datasets. This 1-epoch regime is specifically designed to isolate early optimization efficiency and highlight differences before network capacity becomes the limiting factor. On CIFAR- 10, reshape-based Muon and SGD use the original recipe’s heavily optimized settings, while Conv-NS is tuned with a 1000-configuration budget. On ImageNet, reshape-based Muon and Conv-NS receive equal 48-configuration grid-search budgets. Full experimental details are reported in section E; results are summarized in table 1. Both Muon-Reshape and Muon-ConvNS are efficient.The results highlight that Muon variants drastically outperform SGD in short training settings. This is particularly prevalent in the extreme 1-epoch ImageNet setting, where both Muon (RKO) and the proposed Muon-Conv-NS methods achieve an accuracy above52%, showcasing a massive gap compared to SGD’s30.47%. This is, to our knowledge, the best result that can be achieved without techniques such as replay buffers (Hu et al., 2021). A deeper comparison with these techniques can be found in section J. However, as training gets longer, performance becomes increasingly capped by the network’s capacity rather than the optimizer’s speed. This dynamic explains why the performance gap between Muon and SGD narrows in the 8-epoch and 30-epoch settings, and it also contextualizes why SGD remains highly competitive for very long training schedules, such as the 600-epoch runs used in modern ResNet baselines (Wightman et al., 2021), which prioritize generalization. 6 Dataset (Budget)OptimizerVal. Acc. (%)Train Acc. (%)Runtime CIFAR-10 (8 epochs)Muon - reshape94.02± 0.1497.66± 0.313.074± 0.005s 8 layers CNN, bs = 2000Muon - Conv-NS93.98± 0.1497.17± 0.383.608± 0.018s 1xRTX4090SGD91.25± 0.0993.00± 0.352.975± 0.009s CIFAR-10 (1 epoch)Muon - reshape87.18± 0.00286.46± 0.010.509± 0.003s 8 layers CNN, bs = 250Muon - Conv-NS87.75± 0.00285.44± 0.020.559± 0.003s 1xRTX4090SGD84.99± 0.00283.52± 0.020.496± 0.004s ImageNet (30 epochs)Muon - reshape74.8260.511 h 54m 06s±12m 51s ResNet-50, bs = 2048Muon - Conv-NS74.5358.701 h 49m 20s±1m 42s 4xH100SGD70.2548.031 h 33m 56s±8m 07s ImageNet (1 epoch)Muon - reshape54.0528.1014m 57s±4m 57s ResNet-50, bs = 512Muon - Conv-NS52.9526.7413m 13s±44s 1xH100SGD30.4710.1013m 00s±3m 46s Table 1: CIFAR-10 and ImageNet classification results. Conv-NS essentially matches reshape- based Muon on CIFAR-10 AirBench94 and remains competitive on ImageNet-1k, while reshape- based Muon keeps a slight overall edge. Train and validation accuracies are shown alongside runtimes for both the full budget (8/30 epochs) and the 1-epoch setting. Conv-NS closes the theory gap with competitive efficiency. Conv-NS successfully fixes the theoretical operator mismatch while achieving performance that is on par with the reshape-based heuristic. On CIFAR-10, Conv-NS matches the accuracy of reshape-based Muon (93.98%vs. 94.02%). On ImageNet, it achieves a competitive74.53%accuracy, only slightly below reshape- based Muon’s74.82%. In terms of runtime, the principled approach is highly practical: it introduces only a marginal overhead on the small CIFAR-10 network and actually runs slightly faster than the reshape baseline on the larger ResNet-50 architecture. This demonstrates that closing the theory- practice gap does not require sacrificing computational efficiency. Additional comparisons with other unpublished baselines are reported in section K.2. 6WHAT SHOULD CONVOLUTIONAL ORTHOGONALIZATION OPTIMIZE? The near-parity with Conv-NS observed on both CIFAR-10 and ImageNet raises a central question: Why does improving operator orthogonalization in Muon for convolution updates fail to outperform the existing heuristic? More accurate orthogonalization can hurt training.We sweep the number of Conv-NS iterations in the one-epoch CIFAR-10 setting of section 5. Under training-loss selection, validation accuracy peaks at 10 iterations (87.56%), then decreases at 25 and 100 iterations (86.42%and85.13%); training loss also worsens (Table 2). The same qualitative trend holds under validation- and training- accuracy-based selection (Table 6). In the eight-epoch setting, Figure 5 shows that increasing the iteration count substantially reducesL orth throughout training, yet the 100-iteration setting generally learns more slowly than intermediate iteration counts. Thus, more accurate convolutional orthogonalization does not consistently improve training performance. Learning-rate sweeps for each iteration count are provided in section H. Iterations (Conv-NS)151025100 Val. Acc. (%)84.9386.3587.5686.4285.13 Train Acc. (%)86.3286.7287.4486.5286.32 Train Loss1.16931.15191.13521.15711.1821 Muon LR0.5790.5130.4870.2000.094 Muon Momentum0.7990.7240.2500.4060.706 Table 2: One-epoch CIFAR-10 sweep results. We report the single top configuration per sweep under training-loss selection; each configuration line is already averaged over repeated runs with the same hyperparameters. The best one-epoch operating point occurs at 10 Conv-NS iterations. 7 Lower orthogonality errors do not guarantee better training. We train the same CIFAR-10 network for eight epochs with RKO and Conv-NS usingn ∈ 1, 5, 10, 25, 100, with five runs per setting. Conv-NS settings share a learning rate of0.52, momentum0.6, and damping0.99; RKO uses learning rate0.24and the same momentum. We track training and validation accuracy alongsideL orth ,L flat , andL frob , averaged across convolutional layers. RKO’s performance raises the question of whether its flattened-kernel objectiveL flat is better suited to training thanL orth . However, Figure 5 shows that high-iteration Conv-NS reduces both objectives and can achieve lower L flat than RKO without improving validation accuracy. Thus, lower values of either orthogonality objective alone do not explain the performance ordering. Full residual comparisons are reported in section G. Since Conv-NS empirically approaches convolutional orthogonality with near-minimal gradient modification (Section 4), these results motivate examining whether theConvPFtarget itself imposes unnecessary constraints. Convolutional structure limits the LMO justification for polar factorization. Several works in the literature have proposed explanations for the empirical performance of Muon and gradient orthogonalization. One of the most prominent interpretations relates Muon to approximately solving a Linear Minimization Oracle (LMO) over a spectral-norm ball of radius η 4 (Pethick et al., 2025): ∆ ⋆ W = arg min ∥∆ W ∥ op ≤η ⟨G, ∆ W ⟩ F .(14) This formulation measures the size of the update through the operator norm (spectral norm), thereby constraining the worst-case change in the layer output (Bernstein & Newhouse, 2025). By duality be- tween the spectral and nuclear norms, the closed-form solution is proportional to the semi-orthogonal polar factor of the gradient (Bernstein & Newhouse, 2024; Newhouse, 2025). For completeness, we provide the proof of the optimality of PF(G) in section L. For convolutional layers, the LMO must instead be formulated in the operator space in which the layer acts. The minimization is therefore restricted to structured linear operators induced by fixed-support kernels, with G K , ∆ K ∈M kernel : ∆ ⋆ K =arg min ∆ K ∈M kernel ∥T (∆ K )∥ op ≤1 ⟨T (G K ),T (∆ K )⟩ F .(15) The connection between the convolutional LMO in eq. (15) and the convolutional polar factor ConvPF(G K )follows the same proof as for dense matrices (Section L). As a first indication, ConvPF(G K )satisfies the unit-ball constraint∥T (∆)∥ op ≤ 1by construction, unlike the RKO solution (section 2.2). Moreover, decomposingL frob (T (G K ),T (∆ K ))shows thatConvPF(G K ) solves the LMO over the Stiefel-constrained feasible set (Section M): − ConvPF(G K )∈ argmin ∆ K ∈M kernel T (∆ K )∈St ⟨T (G K ),T (∆ K )⟩ F (16) The main difference between Equation (15) and Equation (16) lies in the minimization set which are respectively∥T (∆)∥ op ≤ 1 andT (∆ K )∈ St. Finally, in the dense case, polar-factor optimality follows fromCH(St)(the Convex Hull ofSt) being the operator-norm unit ball. Does this remain true under the kernel-size constraint in eq. (16)? Indeed, the optimumT (ConvPF(G K ))lies both inT (M kernel ) =T (∆), ∆∈M kernel , and on the Stiefel manifoldSt. Therefore, the convex hull of their intersection is contained in the intersection of their convex hulls: CH(T (M kernel )∩ St)⊆T (M kernel )∩ CH(St) =∆∈M kernel such that∥T (∆)∥ op ≤ 1 (17) The right-hand side is exactly the feasible set of the convolutional LMO in eq. (15). Since the inclusion is generally strict,ConvPF(G K )is optimal over a more constrained set than the operator-norm unit ball restricted toT (M kernel ). This distinction offers a possible explanation for the previous results: approachingConvPFmore accurately enforces orthogonality constraints beyond what the convolutional LMO requires. The 4 In the following, η may be omitted and replaced by the unit ball for simplicity 8 012345678 epoch 0.6 0.8 Training accuracy 012345678 epoch 0.4 0.6 0.8 Validation accuracy 012345678 epoch 10 −1 10 0 10 1 10 2 orth (mean) 012345678 epoch 10 2 3 × 10 1 4 × 10 1 6 × 10 1 2 × 10 2 flat (mean) 012345678 epoch 10 5 Frob (mean) flatten-muon conv-muon (n=1) conv-muon (n=5) conv-muon (n=10) conv-muon (n=25) conv-muon (n=100) Figure 5: Better convolutional and flattened orthogonality does not guarantee better training. Increasing Conv-NS iterations lowers bothL orth andL flat , yet high iteration counts slow early learning and increase gradient modification. Curves show means over five eight-epoch CIFAR-10 runs; shaded bands indicate one standard deviation. accompanying increase in gradient modification may therefore hinder early optimization, even as bothL orth andL flat decrease (Table 2 and fig. 5). Both approaches therefore have limitations: reshape-based Muon does not enforce the convolutional operator-norm constraint (under constraint), while Conv-NS targets the stronger orthogonality con- straint (over constraint). Neither is guaranteed to solve the convolutional LMO; a principled extension of Muon to convolutions requires more than computing the convolutional polar factor. 7CONCLUSION In this work, we revisited Muon for convolutional layers from the perspective of the induced convolution operator. We showed that the standard reshape-based implementation does not, in general, orthogonalize this operator. We introduced Convolutional Newton-Schulz (Conv-NS), a convolution-aware procedure that pre- serves kernel support and converges toward the convolutional polar factor. Despite enforcing better operator-level orthogonality, Conv-NS only achieves performance comparable to reshape-based Muon on CIFAR-10 and ImageNet. Moreover, increasing the number of Conv-NS iterations could degrade performance, suggesting that exact convolutional orthogonality may overconstrain the update. We finally relate this observation to the LMO interpretation of Muon: the convolutional polar factor is optimal under the Stiefel constraint, but not necessarily the minimizer of the full operator-norm LMO. This suggests that adapting Muon to convolutions requires balancing operator-level control with sufficient freedom in the update. Finally, future work could explore more various training tasks like generative modeling or self- supervised representation learning with Conv-NS; along with its extensions for non-2D convolutional layers. Indeed, real training gradient analysis on classification tasks seems to indicate that reshape- based Muon empirically finds a compromise between gradient fidelity and orthogonality, this is not the case on synthetic gradients of different distributions (see fig. 8). This distribution-dependent optimality indicates that reshape-based Muon could be suboptimal on different training tasks with different gradient distributions. 9 ACKNOWLEDGMENTS This work was carried out within the DEEL project, 5 which is part of IRT Saint Exup ́ ery and the ANITI AI cluster. The authors acknowledge the financial support from DEEL’s Industrial and Academic Members and the France 2030 program - Grant agreements n°ANR-10-AIRT-01 and n°ANR-23-IACL-0002. This work was granted access to the HPC resources of IDRIS under the allocation 2025-AD011016381 made by GENCI. REFERENCES El Mehdi Achour, Franc ̧ois Malgouyres, and Franck Mamalet. Existence, stability and scalability of orthogonal convolutional neural networks. Journal of Machine Learning Research, 23(347), 2022. Essential AI, Ishaan Shah, Anthony M. Polloreno, Karl Stratos, Philip Monk, Adarsh Chaluvaraju, Andrew Hojel, Andrew Ma, Anil Thomas, Ashish Tanwer, Darsh J. Shah, Khoi Nguyen, Kurt Smith, Michael Callahan, Michael Pust, Mohit Parmar, Peter Rushton, Platon Mazarakis, Ritvik Kapila, Saurabh Srivastava, Somanshu Singla, Tim Romanski, Yash Vanjani, and Ashish Vaswani. Practical Efficiency of Muon for Pretraining, May 2025. Noah Amsel, David Persson, Christopher Musco, and Robert M. Gower. The Polar Express: Optimal Matrix Sign Methods and their Application to the Muon Algorithm. October 2025. Rohan Anil, Vineet Gupta, Tomer Koren, Kevin Regan, and Yoram Singer. Scalable second order optimization for deep learning. arXiv preprint arXiv:2002.09018, 2020. Lukas Balles and Philipp Hennig. Dissecting Adam: The Sign, Magnitude and Variance of Stochastic Gradients. In Proceedings of the 35th International Conference on Machine Learning, p. 404–413. PMLR, July 2018. ISSN: 2640-3498. Jeremy Bernstein. Deriving muon, 2025. URLhttps://jeremybernste.in/writing/ deriving-muon. Jeremy Bernstein and Laker Newhouse. Old Optimizer, New Norm: An Anthology. December 2024. URL https://openreview.net/forum?id=ux18f5nOpD. Jeremy Bernstein and Laker Newhouse. Modular Duality in Deep Learning. In Proceedings of the 42nd International Conference on Machine Learning, p. 3920–3930. PMLR, October 2025. URLhttps://proceedings.mlr.press/v267/bernstein25a.html. ISSN: 2640- 3498. ̊ Ake Bj ̈ orck and Clazett Bowie. An iterative algorithm for computing the best estimate of an orthogonal matrix. SIAM Journal on Numerical Analysis, 8(2):358–364, 1971. Thibaut Boissin, Franck Mamalet, Thomas Fel, Agustin Martin Picard, Thomas Massena, and Mathieu Serrurier. An adaptive orthogonal convolution scheme for efficient and flexible cnn architectures. In International Conference on Machine Learning, 2025. J ́ er ˆ ome Bolte, Aris Daniilidis, Olivier Ley, and Laurent Mazet. Characterizations of lojasiewicz inequalities and applications. Transactions of the American Mathematical Society, 362(6):3319– 3363, 2010. Ali Ebrahimpour Boroojeny, Matus Telgarsky, and Hari Sundaram. Spectrum extraction and clipping for implicitly linear layers. In International Conference on Artificial Intelligence and Statistics, p. 2971–2979. PMLR, 2024. Franz Louis Cesista. Frequency domain muon for convolutional neural networks: Simplified, 2026. URL https://leloykun.github.io/ponder/freqmuon/. 5 https://w.deel.ai/ 10 Blaise Delattre, Quentin Barth ́ elemy, Alexandre Araujo, and Alexandre Allauzen. Efficient bound of lipschitz constant for convolutional layers by gram iteration. In International Conference on Machine Learning, p. 7513–7532. PMLR, 2023. Blaise Delattre, Quentin Barth ́ elemy, and Alexandre Allauzen. Spectral norm of convolutional layers with circular and zero paddings. arXiv preprint arXiv:2402.00240, 2024. Jia Deng, Wei Dong, Richard Socher, Li-Jia Li, Kai Li, and Li Fei-Fei. ImageNet: A large-scale hierarchical image database. In IEEE Conference on Computer Vision and Pattern Recognition (CVPR), p. 248–255, 2009. Sara Dragutinovi ́ c and Rajesh Ranganath. To Use or not to Use Muon: How Simplicity Bias in Optimizers Matters. March 2026. Kevin Frans, Pieter Abbeel, and Sergey Levine. What really matters in matrix-whitening optimizers? In OPT 2025: Optimization for Machine Learning, 2025. Ekaterina Grishina, Mikhail Gorbunov, and Maxim Rakhuba. Tight and efficient upper bound on spectral norm of convolutional layers. In European Conference on Computer Vision, p. 19–34. Springer, 2024. Vineet Gupta, Tomer Koren, and Yoram Singer. Shampoo: Preconditioned stochastic tensor optimiza- tion. In International Conference on Machine Learning, p. 1842–1850. PMLR, 2018. Huiyi Hu, Ang Li, Daniele Calandriello, and Dilan Gorur. One pass imagenet. In NeurIPS 2021 Workshop on ImageNet: Past, Present, and Future, 2021. Keller Jordan. 94% on CIFAR-10 in 3.29 seconds on a single GPU, 2024. URLhttp://arxiv. org/abs/2404.00498. 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. URLhttps://github.com/KellerJordan/ modded-nanogpt. 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. URLhttps: //kellerjordan.github.io/posts/muon/. Michel Journ ́ e, Yurii Nesterov, Peter Richt ́ arik, and Rodolphe Sepulchre. Generalized power method for sparse principal component analysis. Journal of Machine Learning Research, 11(2), 2010. Andrej Karpathy. nanochat: The best chatgpt that $100 can buy, 2025. URLhttps://github. com/karpathy/nanochat. Ji-HaKim.Frequency-domainmuonforconvfilters-orthogonalizing theoperator,2026.URLhttps://jiha-kim.github.io/posts/ frequency-domain-muon-for-conv-filters/. Qiyang Li, Saminul Haque, Cem Anil, James Lucas, Roger B Grosse, and J ̈ orn-Henrik Jacobsen. Preventing gradient attenuation in lipschitz constrained convolutional networks. Advances in neural information processing systems, 32, 2019. Jingyuan Liu, Jianlin Su, Xingcheng Yao, Zhejun Jiang, Guokun Lai, Yulun Du, Yidao Qin, Weixin Xu, Enzhe Lu, Junjie Yan, Yanru Chen, Huabin Zheng, Yibo Liu, Shaowei Liu, Bohong Yin, Weiran He, Han Zhu, Yuzhi Wang, Jianzhou Wang, Mengnan Dong, Zheng Zhang, Yongsheng Kang, Hao Zhang, Xinran Xu, Yutao Zhang, Yuxin Wu, Xinyu Zhou, and Zhilin Yang. Muon is Scalable for LLM Training, February 2025. URLhttp://arxiv.org/abs/2502.16982. arXiv:2502.16982 [cs]. Liangchen Luo, Yuanhao Xiong, and Yan Liu. Adaptive gradient methods with dynamic bound of learning rate. In International Conference on Learning Representations, 2019. URLhttps: //openreview.net/forum?id=Bkg3g2R9FX. 11 Laker Newhouse. Understanding muon, 2025. URL https://lakernewhouse.com/muon. Thomas Pethick, Wanyun Xie, Kimon Antonakopoulos, Zhenyu Zhu, Antonio Silveti-Falls, and Volkan Cevher. Training Deep Learning Models with Norm-Constrained LMOs. In Proceedings of the 42nd International Conference on Machine Learning, p. 49069–49104. PMLR, October 2025. ISSN: 2640-3498. Omead Pooladzandi. Orthogonalizing all layers for muon optimizer via psgd’s matrixizer, 2026. URL https://github.com/opooladz/matrixizer-cifar10/tree/master. Haozhi Qi, Chong You, Xiaolong Wang, Yi Ma, and Jitendra Malik. Deep isometric learning for visual recognition. In International conference on machine learning, p. 7824–7835. PMLR, 2020. Artem Riabinin, Egor Shulgin, Kaja Gruntkowska, and Peter Richt ́ arik. Gluon: Making Muon & Scion Great Again! (Bridging Theory and Practice of LMO-based Optimizers for LLMs). June 2025. Fabian Schaipp. Optimization benchmark for diffusion models on dynamical systems. In EurIPS 2025 Workshop on Principles of Generative Modeling (PriGM), 2025. URLhttps://openreview. net/forum?id=demP9pzBDX. Hanie Sedghi, Vineet Gupta, and Philip M. Long. The singular values of convolutional layers. In International Conference on Learning Representations, 2019. URLhttps://openreview. net/forum?id=rJevYoA9Fm. Alexandra Senderovich, Ekaterina Bulatova, Anton Obukhov, and Maxim Rakhuba. Towards Practical Control of Singular Values of Convolutional Layers. Advances in Neural Information Processing Systems, 35:10918–10930, December 2022. Egor Shulgin, Sultan AlRashed, Francesco Orabona, and Peter Richt ́ arik. Beyond the ideal: Analyzing the inexact muon update. In The Third Conference on Parsimony and Learning (Recent Spotlight Track), 2026. URL https://openreview.net/forum?id=veVNVUSieC. Bhavya Vasudeva, Puneesh Deora, Yize Zhao, Vatsal Sharan, and Christos Thrampoulidis. How muon’s spectral design benefits generalization: A study on imbalanced data.In C. Vondrick, B. Hariharan, C. Raffel, L. Pinto, D. Yang, and A. Faust (eds.), In- ternational Conference on Learning Representations, volume 2026, p. 142276–142311, 2026. URLhttps://proceedings.iclr.c/paper_files/paper/2026/file/ e6231c5f46598cfd09f1970524e0436-Paper-Conference.pdf. Nikhil Vyas, Depen Morwani, Rosie Zhao, Itai Shapira, David Brandfonbrener, Lucas Janson, and Sham Kakade. Soap: Improving and stabilizing shampoo using adam for language modeling. In International Conference on Learning Representations, volume 2025, p. 93423–93444, 2025. Jiayun Wang, Yubei Chen, Rudrasis Chakraborty, and Stella X Yu. Orthogonal convolutional neural networks. In Proceedings of the IEEE/CVF conference on computer vision and pattern recognition, p. 11505–11515, 2020. Kaiyue Wen, David Hall, Tengyu Ma, and Percy Liang. Fantastic Pretraining Optimizers and Where to Find Them. In International Conference on Learning Representations, 2026. Ross Wightman, Hugo Touvron, and Herv ́ e J ́ egou. Resnet strikes back: An improved training procedure in timm. In NeurIPS 2021 Workshop on ImageNet: Past, Present, and Future, 2021. Ashia C Wilson, Rebecca Roelofs, Mitchell Stern, Nati Srebro, and Benjamin Recht. The Marginal Value of Adaptive Gradient Methods in Machine Learning. In Advances in Neural Information Processing Systems, volume 30. Curran Associates, Inc., 2017. Tian Xie, Haoming Luo, Haoyu Tang, Yiwen Hu, Jason Klein Liu, Qingnan Ren, Yang Wang, Wayne Xin Zhao, Rui Yan, Bing Su, Chong Luo, and Baining Guo. Controlled LLM Train- ing on Spectral Sphere, January 2026. URLhttp://arxiv.org/abs/2601.08393. arXiv:2601.08393 [cs]. 12 ARELATED WORK A.1MUON AND MATRIX-WHITENING OPTIMIZATION The development of Linear Minimization Oracle (LMO) and matrix-whitening optimizers, such as Muon (Jordan et al., 2024b; Bernstein, 2025) and SOAP (Vyas et al., 2025; Wen et al., 2026), aims to improve upon standard element-wise algorithms like AdamW. Frans et al. (2025) argue that the gains of matrix-whitening optimizers are not explained solely by spectral normalization: in their Transformer language-modeling experiments, variance-adapted variants such as Adam, SOAP, and AdaMuon outperform sign-descent counterparts such as Signum, SPlus, and Muon. Applying such optimizers alters the implicit bias of the network, fundamentally shifting convergence toward max-margin solutions with respect to the spectral norm rather than prioritizing simple dominant components (Dragutinovi ́ c & Ranganath, 2026; Vasudeva et al., 2026). While this structural bias often benefits generalization on complex distributions, controlling the resulting weight drift remains challenging, leading to variants like Spectral Sphere Optimizer (SSO) (Xie et al., 2026) which enforces strict module-wise spectral properties. Despite reports of substantial compute-time frontier breakthroughs and improved large-batch efficiency (AI et al., 2025), recent rigorous scaling laws show that the speedup relative to AdamW diminishes at larger parameter scales depending on careful hyperparameter tuning (Wen et al., 2026). Furthermore, theoretical stepsizes are often bridged to empirical reality by accounting for inexact matrix updates (Shulgin et al., 2026) and layer-wise geometry (Riabinin et al., 2025). A.2COMMUNITY EXPLORATIONS FOR CONVOLUTIONAL ORTHOGONALIZATION The structural gap between 4D kernels and 2D spectral optimizers has spurred recent community efforts to adapt Muon for CNNs, notably Freq-Muon (Kim, 2026; Cesista, 2026), which orthogonal- izes per-frequency channel-mixing matrices via the DFT, and Matrixizer (Pooladzandi, 2026), which searches for an optimal 2D factorization of the gradient tensor before applying Newton–Schulz. These implementations highlight an active demand for a CNN-compatible Muon, yet no peer-reviewed work has formalized these approaches or provided rigorous evaluations. We describe both methods in detail and report appendix-only comparisons in section K. BIMPLEMENTATION OF CONV-NS Algorithm 2 Convolutional Newton–Schulz (Conv-NS) Require: Kernel K, stride s = (s h ,s w ), groups g Require: N = 12, β = 0.5, α = 0.91, ε = 10 −2 1: (c o ,c i /g,k h ,k w )← shape(K) 2: K ← SPECTRALRESCALE(K)▷ Three Gram iterations 3: p← s h ⌊(k h − 1)/s h ⌋, s w ⌊(k w − 1)/s w ⌋ 4: for t = 1,...,N do 5: B ← PACKGROUPS(K,g)▷ Only for grouped convolutions 6: A← CONV2D(B,K; s,p,g) 7: C ← CONVTRANSPOSE2D(A,K; s,p,g, 0)▷ Performs Convolution and cropping 8: C ← UNPACKGROUPS(C,g)▷ Only for grouped convolutions 9: K ← α (1 + β)K− βC 10: end for 11: return K/α We initialize the kernel by spectral rescaling using three Gram iterations (Delattre et al., 2023; 2024). Group packing rearranges the kernel so that grouped convolutions can process its output filters as a batch; unpacking restores the original layout. The convolution–transposed-convolution pair evaluates the cubic Newton–Schulz correction on the original kernel support, with cropping implemented implicitly by the padding choices. We use a constantβ = 0.5: increasingβcan destabilize the iteration, whereas decreasing it slows convergence. Although the unprojected polynomial preserves the operator-norm unit ball for0 < β ≤ 0.5, support cropping need not do so. Damping byα = 0.91 13 provides a stability margin intended to keep the iterates inside this ball. The final division byα removes the last damping factor. The block-convolution derivation concerns ungrouped, unit-stride convolutions; the grouped and strided extensions were checked by unit tests. Efficiency. Cropping preserves the original kernel support, giving each iteration a fixed compu- tational and memory cost for fixed kernel dimensions, stride, and group count. The Conv-NS loop operates directly on real-valued kernels and requires no FFT. It therefore avoids complex-valued frequency-domain arithmetic and transforms on the activation-sized grid, whose spatial dimensions are typically much larger than the kernel support. CCONVERGENCE ANALYSIS The convolutional Newton–Schulz update should be interpreted in two stages. Without the support projection, it is exactly the classical Bj ̈ orck / Newton–Schulz iteration written in kernel coordinates. Once the projection onto fixed spatial support is inserted, the method is no longer exact Newton– Schulz on the full convolution operator; it becomes a projected gradient method for a support- constrained orthogonality objective. This objective coincides with theL orth surrogate of Achour et al. (2022) in the stride-aware setting. Exact block-convolution form of Newton–Schulz whenS = 1. Assume circular padding and unit stride. Let b K (t+1) = (1 + β) b K (t) − β b K (t) ⊛ ( b K (t) ) ⊤ ⊛ b K (t) , where the support of b K (t) is allowed to grow with t. Define X (t) :=T circ ( b K (t) ). Proposition C.1 (Unprojected block-convolution Newton–Schulz). For every t, X (t+1) = (1 + β)X (t) − βX (t) (X (t) ) ⊤ X (t) . Hence the unprojected kernel iteration is exactly the classical Bj ̈ orck / Newton–Schulz iteration on the induced convolution operator. In particular, if∥X (0) ∥ 2 ≤ 1and0 < β ≤ 1 2 , thenX (t) converges to the polar factor PF(X (0) ). Proof. Using the block-convolution identities, T circ (A⊛ B) =T circ (A)T circ (B), T circ (K ⊤ ) =T circ (K) ⊤ , we obtain T circ b K (t) ⊛ ( b K (t) ) ⊤ ⊛ b K (t) = X (t) (X (t) ) ⊤ X (t) . ApplyingT circ to the kernel recursion yields the claimed matrix iteration. The convergence statement is then the standard Bj ̈ orck / Newton–Schulz result. Implicit Newton–Schulz as gradient descent on the orthogonality objective. LetM kernel denote the linear space of kernels with fixed spatial supportk h ×k w , and letΠ supp be the orthogonal projection onto that space (center crop). In the row-orthogonal case, define Φ bc (K) := 1 4 ∥K⊛ K ⊤ − I∥ 2 F , K ∈M kernel . The column-orthogonal case is identical after replacing K⊛ K ⊤ − I by K ⊤ ⊛ K− I . Proposition C.2 (Gradient formula). For every K ∈M kernel , ∇Φ bc (K) = Π supp (K⊛ K ⊤ − I)⊛ K = Π supp K⊛ K ⊤ ⊛ K− K . Consequently, Π supp (1 + β)K− βK⊛ K ⊤ ⊛ K = K− β∇Φ bc (K). 14 Proof. Let H ∈M kernel . Since D(K⊛ K ⊤ )[H] = H⊛ K ⊤ + K⊛ H ⊤ , we have DΦ bc (K)[H] = 1 2 K⊛ K ⊤ − I, H⊛ K ⊤ + K⊛ H ⊤ . Now use that the adjoint of convolution byK ⊤ is convolution byKfollowed by the center crop onto the support of H : ⟨H⊛ K ⊤ ,Z⟩ =⟨H, Π supp (Z⊛ K)⟩. Since K⊛ K ⊤ − I is self-transpose, the two terms coincide, hence DΦ bc (K)[H] = Π supp (K⊛ K ⊤ − I)⊛ K ,H . This identifies the gradient. The final identity follows because Π supp (K) = K. Remark C.3 (WhyL orth is the right strided objective). For strideS > 1, exact Newton–Schulz on the actual layer operator must be written for A S (K) := S S T (K), and its cubic term is A S (K)A S (K) ⊤ A S (K) = S S T (K)T (K) ⊤ S ⊤ S S S T (K), which is not the same asS S T (K⊛K ⊤ ⊛K). Therefore, in the strided case, the correct fixed-support objective is not the plain unstrided block-convolution residual but the stride-aware regularizerL orth of (Achour et al., 2022), whose gradient is implemented in practice by a standard convolution followed by a transposed convolution. The formula above is recovered exactly when S = 1. Projected descent on the support-constrainedL orth objective.In the row-orthogonal case define Φ S (K) := 1 4 L orth (K) = 1 4 ∥conv(K,K, padding = P, stride = S)− I r0 ∥ 2 F , with P =⌊(k− 1)/S⌋S. The practical update can then be written abstractly as K (t+1) = K (t) − β t ∇Φ S (K (t) ), or equivalently, in implementation form, K (t+1) = Π supp K (t) − β t G(K (t) ) , whereG(K (t) )is obtained by one direct convolution and one transposed convolution. WhenS = 1, this reduces to the block-convolution formula above. Theorem C.4 (Convergence to a stationary point). Assume that the iterates remain in a bounded setD ⊂M kernel . SinceΦ S is a polynomial on the finite-dimensional spaceM kernel , its gradient is Lipschitz onD; let L be such a Lipschitz constant. If the step sizes satisfy 0 < β ≤ β t ≤β < 2 L , then Φ S (K (t+1) )≤ Φ S (K (t) )− c∥∇Φ S (K (t) )∥ 2 F , c := β 1− L β 2 > 0. Consequently, ∞ X t=0 ∥∇Φ S (K (t) )∥ 2 F <∞, ∇Φ S (K (t) )→ 0, and every accumulation point of (K (t) ) t≥0 is a critical point of the support-constrained problem min K∈M kernel Φ S (K). Moreover, becauseΦ S is a polynomial (hence a Kurdyka–Łojasiewicz function), the whole sequence (K (t) ) t≥0 converges to a single critical point K ⋆ of Φ S . 15 Proof. By the descent lemma for L-smooth functions, Φ S (K− β∇Φ S (K))≤ Φ S (K)− β∥∇Φ S (K)∥ 2 F + Lβ 2 2 ∥∇Φ S (K)∥ 2 F . Applying this with K = K (t) and β = β t gives Φ S (K (t+1) )≤ Φ S (K (t) )− β t 1− Lβ t 2 ∥∇Φ S (K (t) )∥ 2 F , hence the announced inequality. Summing overtproves that P t ∥∇Φ S (K (t) )∥ 2 F <∞ , and therefore ∇Φ S (K (t) )→ 0. Criticality of accumulation points follows by continuity of∇Φ S . Finally, the KL property of polynomial functions upgrades subsequential convergence to convergence of the whole bounded sequence. Corollary C.5 (Operator-level consequence). Under circular boundary conditions andSN ≥ 2k− 1, the results summarized in the LORTH analysis imply thatL orth controls both the Frobenius and spectral orthogonality errors of the induced convolution operator. Therefore, ifΦ S (K (t) )→ 0, then the operator-level orthogonality residuals vanish. More generally, if the terminal valueΦ S (K (T ) )is small, then the spectral residual is small, and the induced convolution is approximately isometric. Thus the projected convolutional Newton–Schulz iteration should be viewed as a support-constrained orthogonality solver: when an orthogonal kernel exists, it can converge to such a solution; otherwise it converges to the best stationary approximation inside the fixed-support kernel class. Initial conditions. To guarantee the convergence of the Newton–Schulz iteration, the initial op- erator norm of the convolution must be strictly bounded, satisfying∥T (K)∥ 2 ≤ 1.0. Although a simple Frobenius normalization of the kernel tensor is sufficient to achieve this, employing a tighter spectral estimate accelerates convergence. The accurate evaluation of the Lipschitz constant for convolutional layers is an active area of research in its own right: initial approaches leveraged the Fast Fourier Transform (FFT) to compute a lower bound for layers with circular padding (Sedghi et al., 2019). This foundation was subsequently refined by methods offering improved computational efficiency (Senderovich et al., 2022), support for diverse padding strategies (Grishina et al., 2024), and the ability to extract broader regions of the operator spectrum (Boroojeny et al., 2024). In our implementation, we rely on the techniques established by Delattre et al. (2023) and Delattre et al. (2024), which allow us to efficiently compute a certifiable upper bound on the spectral norm across various padding modes. Finally, while the target implies that∥T (K)∥ 2 ≤ 1.0, the projection step, coupled with numerical imprecision, can yield an iterate with∥T (K (t) )∥ 2 ≥ 1.0. To avoid this, we add a damping factor; when β = 0.5, a damping factor of 0.9 yields good results. Convergence guarantees: Without the support projection, the block-convolution Newton-Schulz iteration is exactly the classical Bj ̈ orck / Newton–Schulz iteration on the induced operatorT (∆ K ), and therefore converges to its polar factor whenever∥T (∆ K )∥ 2 ≤ 1andβ ≤ 0.5. Once the projectionΠ supp is inserted, the method is no longer exact Newton–Schulz on the full operator; instead, it becomes projected gradient descent on a support-constrained orthogonality objective Φ S = ∥T (∆)T (∆) T − I∥ 2 that coincides with theL orth surrogate of Achour et al. (2022) in the strided setting (proposition C.2). BecauseΦ S is a polynomial on the finite-dimensional kernel space, its gradient is Lipschitz and it satisfies the Kurdyka–Łojasiewicz property: the iterates are guaranteed to converge to a single critical point ofΦ S , with P t ∥∇Φ S (∆ (t) K )∥ 2 F < ∞(Bolte et al., 2010). In particular, when an orthogonal kernel exists within the fixed-support class, the method can converge to it; otherwise it converges to the best stationary approximation. The initial condition∥T (∆ K )∥ 2 ≤ 1 is enforced via a spectral-norm upper bound computed efficiently following Delattre et al. (2023; 2024), and a damping factor handles residual numerical drift across iterations. DCORRECTION OF THE GRAFTING FACTOR: Typical implementations of Muon multiply the orthogonal update by a scalar factor based on the matrix dimensions. While there is no strict consensus, these grafting factors allow balancing learning between layers and help with learning rate transferability. Following the logic of this paper, these 16 Jordan et al. (2024b)Liu et al. (2025)Bernstein & Newhouse (2025) Dense matrix q max 1, n m 0.2 p max(n,m) p n m (n× m) Reshaped r max 1, c out c in k w k h 0.2 p max (c out ,c in k w k h ) 1 k w k h q c out c in kernel Adaptation r max 1, c out c in s 2 0.2 q max c in i w i h ,c out i w s i h s q c out c in s 2 forT (∆ K ) Table 3: Grafting factors for Muon updates. The adaptation forT (∆ K )applies the dense-matrix formula using the dimensions of the induced convolution operator,n = c out i h s i w s andm = c in i h i w . factors should be computed according to the dimensions of the induced operatorT (∆ K )rather than the dimensions of the reshaped kernelK. Table 3 presents different scaling conventions and their adaptation. Note that Bernstein & Newhouse (2025) uses a different logic than the other two by proposing distinct formulas for dense and convolutional layers, where the dense formula does not rely on the kernel sizes k w and k h . ECIFAR-10 AND IMAGENET CLASSIFICATION DETAILS CIFAR-10 AirBench94. The CIFAR-10 experiment uses the AirBench94 benchmark (Jordan, 2024), a speed-optimized 8-epoch training recipe designed to reach 94% accuracy in only a few seconds on a single GPU. The network is the compact AirBench CNN with three convolutional groups (64→256→256 channels) built from3×3stride-1 convolutions. Runs use batch size 2000, 8 training epochs, 3 whitening epochs, horizontal-flip augmentation, and translation augmentation with padding 2. Conv-NS and reshape-based Muon use the same architecture, training schedule, augmentation pipeline, batch size, and non-convolutional optimizer groups; only the convolutional orthogonalization routine inside Muon differs. The SGD baseline replaces Muon entirely. For reshape-based Muon and SGD, we use the original AirBench author’s benchmark-specific optimized settings. Because no reference Conv-NS setting existed for this benchmark, we tune Conv-NS with a 1000-configuration search budget. The selected Conv-NS configuration uses Muon learning rate0.52, Muon momentum0.6, zero Muon weight decay, bias learning rate0.053, head learning rate0.67, Adam weight-decay scale6×10 −6 , and SGD momentum0.85for the remaining parameter groups. For the convolutional orthogonalization itself, we useβ init = 0.25,β end = 0.5, 25 Conv-NS iterations, damping 0.99, γ = 0.01, and ε = 0.01. Conv-NS, reshape-based Muon, and SGD are each reported over 25 independent runs. We report both ordinary validation accuracy and test-time-augmented validation accuracy. On a single RTX4090, the non-optimized PyTorch implementation of Conv-NS runs in4.45seconds per run, compared with 3.68seconds for reshape-based. Runtime with torch compilation enabled are reported in main paper. ImageNet-1k ResNet-50.The ImageNet experiment uses ImageNet-1k (Deng et al., 2009) with a ResNet-50 backbone and atimm-style training recipe shortened to 30 epochs. Training uses batch size 512 across 4 H100 GPUs (so global batch size is 2048), a linear learning-rate schedule with no warmup, TA-Wide augmentation, mixup (α=0.1), and CutMix (α=1.0). We compare theflatten strategy (reshape-based Muon) against the bjorck strategy (Conv-NS). Both strategies receive the same grid-search budget over learning rate∈ 0.01, 0.02, 0.05, 0.1, momentum∈ 0.8, 0.9, 0.95, weight decay∈ 2×10 −5 , 10 −4 , and cautious weight decay∈ true, false. This gives 48 configurations per strategy and 96 runs in total. Both strategies reach their best validation top-1 accuracy at lr= 0.02, momentum= 0.95, weight decay= 2×10 −5 , with cautious weight decay enabled. 17 FPER-STAGE PARETO ANALYSIS ON REAL TRAINING GRADIENTS Section 6 reports the orthogonality–fidelity Pareto front averaged over three training stages. Figure 6 shows the same analysis separately for gradient snapshots captured at the beginning, middle, and end of CIFAR-10 training. The key observation—that reshape-based Muon (RKO) lies on or very near the Pareto front for real gradients—holds consistently across all three stages, confirming that this is not an artifact of a particular training phase. 10 2 4 × 10 1 6 × 10 1 2 × 10 2 L orth 10 4 10 3 10 2 10 1 10 0 10 1 10 2 K ′ K 2F Begin 2 3 4 5 10 25 100 10 2 4 × 10 1 6 × 10 1 2 × 10 2 L orth Middle 2 3 4 5 10 25 100 10 2 4 × 10 1 6 × 10 1 2 × 10 2 L orth End 2 3 4 5 10 25 100 ParetoMuonNS3-conv Figure 6: Per-stage Pareto trade-off on real training gradients. Each panel shows the Pareto front computed on gradient snapshots from a single training stage (begin, middle, end). Reshape-based Muon (RKO) lies on the Pareto front at all stages, and Conv-NS traces the front as the iteration count increases. GHIERARCHICAL OPTIMIZATION OF ORTHOGONALITY OBJECTIVES This section provides the full numerical results comparing the orthogonality residuals of different operators, as discussed in the main text. Table 4 presents the results for synthetic, randomly initialized kernels, while Table 5 shows the results for real gradient snapshots captured during training. On synthetic kernels, table 4 illustrates an asymmetric relationship between the two orthogonality objectives. RKO nearly eliminates the flattened-kernel residual (L flat = 0.02), but leaves substantial convolutional error (L orth = 14.40). Rescaling RKO to unit operator norm reduces this error to 9.13while sacrificing flattened orthogonality, showing that spectral normalization alone does not resolve the mismatch. In contrast, ConvPF achieves low residuals on both objectives (0.10and 0.02, respectively), while ten Conv-NS iterations yield intermediate values (3.38and1.11). These results show that, in this setting, convolutional orthogonalization can also satisfy the flattened-kernel objective, whereas optimizing the latter alone is insufficient to ensure convolutional orthogonality. NodeTrialsLipschitzL orth L flat K102.31± 0.0736.57± 1.537.11± 0.69 RKO(∆K)101.83± 0.0314.40± 0.400.02± 0.00 RKO(∆K)/σ(RKO(∆K))101.00± 0.009.13± 0.147.85± 0.20 ConvPF (∆K)101.07± 0.020.10± 0.030.02± 0.00 ConvNS 10 (∆K)101.14± 0.013.38± 0.091.11± 0.04 Table 4: Orthogonality hierarchy on synthetic data. Highlighting the difference between the convolutional polar factor and the reshape-based orthogonalization (RKO). By measuring both the operator-level orthogonality residualL orth (Polar factor objective) and the flattened-kernel orthogonality residualL flat (RKO objective) we can observe a hierarchy between the two: the convolutional polar factor achieves both operator-level and flattened orthogonality, while RKO only satisfies the latter. On real training gradients, table 5 confirms the mismatch between flattened and convolutional orthogonality across all three training stages: RKO maintains a low flattened-kernel residual, yet its operator norm remains close to 3, compared with 1.83 in the synthetic setting. ConvPF consistently achieves near-zero convolutional error, while ten Conv-NS iterations provide an intermediate level of orthogonalization. However, the relationship between the two objectives is less clear-cut than on synthetic kernels: ConvPF yieldsL flat ≈ 29, substantially above RKO’s0.22–0.51. Thus, although 18 StageNodeTrialsLipschitzL orth L flat Begin∆K6137.34± 47.945.2e+071.1e+07 Begin∆K/σ(∆K)61.00± 0.00152.45± 119.32191.04± 104.36 Begin RKO(∆K)62.99± 0.04461.04± 304.040.51± 0.78 Begin RKO(∆K)/σ(RKO(∆K))61.00± 0.00118.74± 103.28151.18± 81.79 Begin ConvPF (∆K)61.08± 0.070.09± 0.0829.26± 62.44 Begin ConvNS 10 (∆K)61.28± 0.09108.20± 101.75137.14± 87.38 Middle∆K6119.15± 47.634.7e+071.2e+07 Middle∆K/σ(∆K)61.00± 0.00150.68± 117.47189.18± 102.69 Middle RKO(∆K)62.94± 0.12471.96± 287.200.22± 0.11 Middle RKO(∆K)/σ(RKO(∆K))61.00± 0.00117.73± 101.64149.60± 80.70 Middle ConvPF (∆K)61.05± 0.070.04± 0.0729.25± 62.47 Middle ConvNS 10 (∆K)61.26± 0.1071.22± 61.3087.57± 58.03 End∆K6126.52± 50.666.0e+071.5e+07 End∆K/σ(∆K)61.00± 0.00150.98± 117.89189.51± 103.03 End RKO(∆K)62.95± 0.11457.84± 281.820.22± 0.12 End RKO(∆K)/σ(RKO(∆K))61.00± 0.00117.51± 101.34149.55± 80.27 End ConvPF (∆K)61.04± 0.050.08± 0.0729.25± 62.45 End ConvNS 10 (∆K)61.26± 0.1072.59± 64.1089.26± 59.25 Table 5: Orthogonality hierarchy on real gradients. The same comparison as in table 4, evaluated on real gradient snapshots from a CIFAR-10 training run across different stages (begin, middle, end). The hierarchy between the objectives remains consistent with the synthetic setting. the convolutional-orthogonality ordering persists throughout training, these measurements do not establish that ConvPF also minimizes the flattened-kernel residual on real gradients. HEFFECT OF ITERATION COUNT IN THE 8-EPOCH SETTING In Figure 7 the 8-epoch CIFAR-10 AirBench setting is reused to evaluate the effect of different iteration counts. For a given iteration count, the learning rate is swept along the momentum. Five runs are conducted per learning rate, and the momentum with the best average (test time augmentation) validation accuracy across learning rates is selected for each iteration count; shaded bands show ±1standard deviation over these five runs. While low iteration budget exhibits a similar behavior observed by (Shulgin et al., 2026): a larger number of iterations improves performance with a shift in optimal lr and momentum, which can be interpreted as the result of a faster training. 10 1 10 0 Muon learning rate 86.0% 88.0% 90.0% 92.0% 94.0% TTA validation accuracy Selected momentum 1 iterations (m=0.95) 5 iterations (m=0.95) 10 iterations (m=0.9) 25 iterations (m=0.9) 100 iterations (m=0.9) Figure 7: Effect of Conv-NS iteration count on CIFAR-10 AirBench. TTA validation accuracy over the Muon learning-rate sweep for different numbers of Newton–Schulz iterations. For each iteration count, momentum is selected by the best average TTA validation accuracy across learning rates; shaded bands show±1 standard deviation over five runs. 19 IONE-EPOCH SELECTION CRITERIA Robustness to the selection criterion.To verify that the non-monotonic effect of Conv-NS iteration count is not an artifact of the hyperparameter-selection metric, we reselect the best configurations using either TTA (test time augmentation) validation accuracy or training accuracy (Table 6). In both cases, performance improves up to 10 iterations and degrades at 25 and 100 iterations, despite substantial shifts in the selected learning rate. Moreover, the same trend under training-accuracy selection shows that the degradation is already present at the optimization level rather than arising solely from generalization. These results indicate that the Conv-NS iteration count is not merely an approximation-accuracy parameter: pushing the update closer to exact convolutional orthogonality can itself hinder training. Top by TTA validation accuracy Iterations (Conv-NS)151025100 TTA Val. Acc. (%)85.3686.6687.7586.6185.90 Val. Acc. (%)84.6785.9186.8485.5884.97 Train Acc. (%)84.4486.2085.4485.4485.76 Muon LR0.4330.4320.3810.1810.143 Top by training accuracy Iterations (Conv-NS)151025100 TTA Val. Acc. (%)84.9386.4687.5586.0585.13 Val. Acc. (%)84.2285.7186.6985.2983.94 Train Acc. (%)86.3287.3288.1287.0086.32 Muon LR0.5790.1940.4410.1580.094 Table 6: Alternate one-epoch CIFAR-10 selection criteria. These blocks report the single top configuration per sweep under validation-oriented and train-accuracy-oriented selection criteria; each configuration line is already averaged over repeated runs with the same hyperparameters. JCOMPARISON WITH ONE PASS IMAGENET One Pass ImageNet (OPIN; Hu et al. (2021)) studies ResNet-50 training from a randomly ordered ImageNet stream in which each training example arrives once. Its one-pass constraint permits retaining and revisiting examples in a limited replay buffer. The proposed error-prioritized replay method preferentially samples examples with larger predictive errors and applies importance weights to correct the sampling bias. Withkreplay minibatches per incoming minibatch, the effective training budget isE eff = k + 1epochs. OPIN Table 2 reports results at 2, 4, 6, and 9 effective epochs, using buffers containing 1%, 5%, or 10% of the dataset, alongside full-storage multi-epoch baselines. We compare these results with our one-epoch ImageNet experiments, which use no replay and therefore haveE eff = 1. To distinguish data-processing cost from optimization steps, we report both budgets relative to a 90-epoch run with batch size 128. For dataset sizeNand global batch sizeB, our batch size of 512 gives approximately 1.11% of the example-processing budget and 0.278% of the update budget. These quantities describe a single training run and exclude hyperparameter search. The computation proxy does not account for optimizer overhead and should not be interpreted as measured FLOPs or wall-clock time. Both Muon variants exceed all OPIN results at two effective epochs, including the full-storage baseline (46.1%), while processing half as many training examples and taking approximately eight times fewer optimizer updates, without a replay buffer. They remain below the four-effective-epoch replay results (55.5–57.2%). Thus, the comparison highlights strong early learning under a smaller data-processing and update budget. Because the training recipes differ, it is a cross-paper comparison rather than a controlled optimizer ablation; the resource ratios alone do not establish an equivalent wall-clock speedup. 20 10 4 10 2 10 0 10 2 10 4 L orth 10 4 10 2 10 0 10 2 K K 2F 2 3 4 5 10 25 50 100 Pareto Freq-Muon Matrixizer Muon NS3-conv ConvPF( K) K Figure 8: Pareto trade-off on synthetic Gaussian kernels. The grey curve is the empirical Pareto front computed via simulated annealing. Conv-NS iterates (blue dots, labelled by iteration count) trace the front from the rescaled SGD corner toward the convolutional polar factor as the number of iterations increases. Reshape-based Muon (RKO) lies far inside the suboptimal zone, reflecting a structural misalignment with the operator-level objective on random kernels. Matrixizer and Freq- Muon are far off-scale. KCOMMUNITY BASELINE COMPARISONS We compare Conv-NS against two recent community-driven baselines -the only publicly available attempts at convolution-aware Muon to date -neither of which has been peer-reviewed or formally evaluated. These comparisons are included for context, but remain appendix material because the main empirical claims center on Conv-NS, reshape-based Muon, and SGD. Freq-Muon (Kim, 2026; Cesista, 2026). Under circular boundary conditions, a convolution operator is block-circulant and is therefore diagonalized by the spatial DFT intoM 2 independent channel-mixing matrices ˆ K[ω]∈ C c out ×c in , one per frequency binω. Freq-Muon applies the standard Muon NS5 polar projection independently to each ˆ K[ω]and recovers the spatial kernel via inverse FFT. When the FFT grid matches the full input spatial dimensions (M = i h = i w ), this procedure yields an exact orthogonalization of the circular convolution operator. In practice, however,Mis set to a small fixed value (e.g.M =8) to keep costs tractable, breaking the exact operator correspondence and turning the method into a surrogate. Additionally, the inverse FFT result is cropped to the original k h × k w support, discarding spatial taps beyond the kernel size. Matrixizer (Pooladzandi, 2026). Instead of the standardc out × (c in k h k w )flattening used in reshape-based Muon, Matrixizer searches over all permutations of the gradient tensor’s dimensions Validation top-1 accuracy (%) Method Effective epochs Computation proxy Relative updates Batch size Updates (approx.) No replay Replay 1% Replay 5% Replay 10% Multi-epoch 100% OPIN22/902/9012820,018–44.745.145.746.1 OPIN44/904/9012840,036–55.557.157.259.0 OPIN66/906/9012860,055–58.961.362.264.1 OPIN99/909/9012890,082–59.363.265.068.2 Muon (reshape), ours11/90 1/3605122,50254.05– Muon (Conv-NS), ours11/90 1/3605122,50252.95– SGD, our recipe11/90 1/3605122,50230.47– Table 7: Comparison of OPIN Table 2 with the one-epoch ImageNet results in our Table 1. All experiments use ResNet-50. Replay column percentages indicate buffer capacity relative to the dataset; the multi-epoch baseline stores the full dataset. Our runs require no replay storage. Both budget ratios use a 90-epoch, batch-128 reference. Update counts are estimates usingN = 1,281,167, rounded to the nearest integer; exact counts depend on final-batch handling. The OPIN counts interpret each reported replay step as a separate minibatch update. 21 and all possible split points to find the 2D factorization that minimizesm 2 + n 2 , i.e. the total preconditioner cost (following the PSGD logic of Pooladzandi (2026)). The standard Muon NS5 iteration is then applied to this optimally shaped matrix. While the resulting reshape is typically better balanced than the default flattening, the method remains fundamentally reshape-based: it does not account for stride, groups, dilation, or the structured embedding∆K 7→T (∆K)that governs how the convolution acts on data. K.1SYNTHETIC ANALYSIS OF THE ORTHOGONALITY–FIDELITY TRADE-OFF Computing the empirical Pareto front. To evaluate how different orthogonalization strategies navigate the inherent trade-off between operator-level orthogonality (L orth ) and fidelity to the original kernel (L frob = ∥∆K ′ − ∆K∥ 2 F ), we first construct a ground-truth Pareto front on synthetic data. We initialize a random Gaussian convolution kernel withc in = 16,c out = 16, andk h = k w = 3. Because the scalarized objectiveλL orth + (1− λ)L frob offers no convexity guarantees, we rely on a simulated annealing approach to explore the landscape. We sample 50 scalar weightingsλ and optimize each configuration using Adam with simulated annealing noise injection over 1500 iterations. The best result is kept across multiple trials to avoid local minima, yielding an empirical Pareto front that serves as a rigorous baseline. Evaluating baselines on the synthetic landscape. Against this optimal trade-off curve, we overlay the operating points of Convolutional Newton–Schulz (Conv-NS) at iteration counts 2, 3, 4, 5, 10, 25, 50, 100, reshape-based Muon (RKO), Matrixizer (Pooladzandi, 2026), and Freq- Muon (Cesista, 2026). The results are reported in fig. 8. Conv-NS traces the Pareto front.As the iteration count grows from 2 to 100, the Conv-NS operat- ing points slide along the simulated-annealing Pareto front, progressing from the low-modification region (highL orth , lowL frob ) toward the exact convolutional polar factor (lowL orth , highL frob ). This confirms that the projected iteration efficiently navigates the support-constrained landscape and that the iteration count serves as a user-controllable dial for the orthogonality–fidelity trade-off. Muon (RKO) is misaligned with the operator-level objective.Unlike its behavior on real training gradients (fig. 4), reshape-based Muon sits far above the Pareto front on synthetic kernels. It modifies the kernel (L frob ≈ 2.3) without improving operator-level orthogonality (L orth ≈ 15.6, compared to∼10.2for the unmodified kernel). Rather than a systematic degradation, this reflects a structural misalignment: on unstructured random Gaussian data, the reshape surrogate pushes the kernel in a direction that does not coincide with the true operator-level objective. A contributing factor is operator-norm inflation—orthogonalizing the flattened matrix can increase the spectral norm of the induced convolution—though the magnitude of this effect varies across initializations. Importantly, as shown in section 6, the structure of real training gradients corrects this misalignment, allowing RKO to fall surprisingly close to the Pareto front in practice. Matrixizer and Freq-Muon are far from optimal. Matrixizer (L orth ≈ 157) and Freq-Muon (L orth ≈ 3.7× 10 4 ) lie orders of magnitude above the Pareto front, indicating that the approximations underlying these methods—Kronecker factorization and partial-spectrum FFT orthogonalization, respectively—introduce substantial artifacts on this objective. Practical take-away. Even a modest iteration budget (e.g. 5 iterations) already places Conv-NS near the Pareto front at a useful orthogonality level (L orth ≈ 4.9), offering a practical compromise between compute cost and operator-level fidelity. In contrast, no single-shot baseline achieves a comparable position on the front. K.2CIFAR-10 COMMUNITY BASELINE RESULTS All Muon-family methods: Conv-NS, reshape-based Muon, Matrixizer, and Freq-Muon, share the same architecture, training schedule, and augmentation pipeline; only the orthogonalization routine inside the optimizer differs. Muon - Conv-NS and Muon - reshape are each reported over 25 independent runs (results slightly differs from table 1 which report statistics over 10 runs); Matrixizer, 22 MethodVal. acc. (%)TTA acc. (%)Runs Conv-NS (ours)93.17± 0.1293.98± 0.0925 Muon (reshape)93.10± 0.1493.99± 0.1425 Matrixizer92.9693.731 Freq-Muon92.5593.061 SGD90.4791.111 Table 8: CIFAR-10 AirBench94 community baseline comparison. Mean±standard deviation over 25 runs for Conv-NS and Muon (reshape); single-run baselines are preliminary and shown as context. Freq-Muon, and SGD are preliminary single-run results that will be extended before the final version, results are reported form their repositories in table 8. Conv-NS and reshape-based Muon achieve near-identical CIFAR-10 AirBench94 performance and both outperform the single-run SGD, Freq-Muon, and Matrixizer baselines. Because Matrixizer and Freq-Muon are preliminary community implementations without peer-reviewed standardized evaluations, we use these results as diagnostic context rather than as a central comparison in the main paper. LLMO DENSE MATRIX PROOF The objective of Muon is to find for a given layer of weightWthe optimal weight update that maximise the loss decrease under a constraint of the maximal variation of the output of the layer. Mathematically, the formulation of the loss decrease is based on the first-order Taylor approximation. L(W + ∆ W ) =L(W ) +⟨G, ∆ W ⟩ F + o(∥∆ W ∥).(18) While the constraint on the maximum variation of the layer output is given by the operator norm: sup ∥x∥ 2 =1 ∥Wx− (W + ∆ W )x∥ 2 =∥∆ W ∥ op (19) This lead to a LMO formulation in the space of the dense matrices: ∆ ⋆ W = arg min ∥∆ W ∥ op ≤η ⟨G, ∆ W ⟩ F The solution of this LMO is proportional to the the polar factor of the gradient defined by PF(G) = argmin Q L frob (G,Q) subject to ∥Q ⊤ − I∥ 2 F = 0 or withm≤ n, using the definition of the Stiefel ManifoldSt(m,n) = Q∈ R m×n : Q T = Id , noted St. PF(G) = argmin Q∈St L frob (G,Q) Proof. First,∀Q∈ St, the Frobenius norm is constant and equal to m: ∥Q∥ 2 F = tr(Q T ) = m. Therefore, ∥Q− G∥ 2 F =∥Q∥ 2 F +∥G∥ 2 F − 2⟨G,Q⟩ F = m +∥G∥ 2 F − 2⟨G,Q⟩ F Thus PF (G)∈ argmax Q∈St ⟨G,Q⟩ F , or equivalently−ηPF (G)∈ argmin Q∈ηSt ⟨G,Q⟩ F . But in the previous equation,the minimization is overηSt,and not the set ∆∈ R m×n :∥∆∥ op ≤ η. Since, Linear minimization does not distinguish a set from its convex hull, noted CH, so: −ηPF (G)∈argmin Q∈η CH(St(m,n)) ⟨G,Q⟩ F 23 . An important property of the Stiefel Manifold for dense matrices is that: CH(St(m,n)) = ∆∈ R m×n :∥∆∥ op ≤ 1 The proof is given in Journ ́ e et al. (2010) Thus −ηPF (G)∈ argmin ∥Q∥ op ≤η ⟨G,Q⟩ F . MCONVOLUTIONAL LMO PROOF This section presents the convolutional LMO construction. Given a convolutional layer with a kernel K, the first-order Taylor approximation is similar to the dense formulatin (Equation (18)): L(K + ∆ K ) =L(K) +⟨G K , ∆ K ⟩ F + o(∥∆ K ∥).(20) But the constraint on the minimum variation has to be evaluated in the space of the structured linear operatorT (K): sup ∥x∥ 2 =1 ∥T (K)x−T (K + ∆ K )x∥ 2 =∥T (∆ K )∥ op (21) The convolutional LMO formulation is thus: ∆ ⋆ K =arg min ∆ K ∈M kernel ∥T (∆ K )∥ op ≤η ⟨G K , ∆ K ⟩ F Furthermore, under circular padding and fixed spatial resolution, the Frobenius distance between two convolutional operators is proportional to the Frobenius distance between their underlying kernels: ⟨T (G K ),T (∆ K )⟩ F ∝⟨G K , ∆ K ⟩ F .(22) This lead to the convolutional LMO formulation provided in Equation (15). The link with the convolutional polar factorConvPF(G K )(Equation (6)) follows the same proof as for dense matrices: Proof. First,∀∆ K ∈M kernel such thatT (∆ K )∈ St, the Frobenius norm ofT (∆ K ) is constant. Therefore, ∥T (∆ K )−T (G K )∥ 2 F =∥T (∆ K )∥ 2 F +∥T (G K )∥ 2 F −2⟨T (G K ),T (∆ K )⟩ F = cste−2⟨T (G K ),T (∆ K )⟩ F Thus ConvPF(G K )∈ argmax ∆ K ∈M kernel T (∆ K )∈St ⟨T (G K ),T (∆ K )⟩ F , or equivalently−η ConvPF(G)∈ argmin ∆ K ∈M kernel T (∆ K )∈ηSt ⟨T (G K ),T (∆ K )⟩ F . Unfortunately the Convex Hull argument is no more an equality in the convolution setup since: CH(T (M kernel )∩ St)⊆T (M kernel )∩ CH(St) =∆∈M kernel such that∥T (∆)∥ op ≤ 1 The convolutional polar factorConvPFis thus the optimum on a subspace of the unit ball, inducing an overconstraint. 24