Slow is Fast: 1. Is the End of Shampoo Actually Muon?

By 苏剑林 | September 28, 2026

As is well known, the Muon optimizer has been popular for quite some time. While it hasn't completely replaced Adam, it has certainly become an unignorable and highly competitive optimizer. Meanwhile, Shampoo-style optimizers have not been "left behind" over the past two years. Beyond the original Shampoo, variants like KL-Shampoo, OKLS, PSGD, and SOAP have emerged, all claiming better results than Muon.

However, this article will pour some cold water on that: all Shampoo-like preconditioner optimizers based on Kronecker products likely have Muon as their ideal limit—the expensive preconditioning matrices we construct might mostly be "wasted effort."

Average Convergence

We start from the average convergence conclusion obtained in "Making Alchemy More Scientific (10): Breaking and Mending the Monotonicity Hypothesis": under the assumption of convexity, the update rule $\boldsymbol{\theta}_{t+1} = \boldsymbol{\theta}_t - \eta_t \boldsymbol{H}_t^{-1}\boldsymbol{g}_t$ satisfies the following inequality:

\begin{equation}\frac{\sum_{t=1}^T \eta_t \mathbb{E}[L(\boldsymbol{\theta}_t) - L(\boldsymbol{\theta}^*)]}{\sum_{t=1}^T \eta_t} \leq \mathbb{E}\left[\frac{R^2 \Lambda_T}{2\sum_{t=1}^T \eta_t} + \frac{\sum_{t=1}^T \eta_t^2 \Vert\boldsymbol{g}_t\Vert_{\boldsymbol{H}_t^{-1}}^2}{2\sum_{t=1}^T \eta_t}\right]\end{equation}

where

\begin{equation}\Lambda_T \triangleq \lambda_{\max}(\boldsymbol{H}_1) + \sum_{t=2}^T [\lambda_{\max}(\boldsymbol{H}_t - \boldsymbol{H}_{t-1})]_+\end{equation}

I won't elaborate on the meanings of the other symbols here; please refer to the original article. Intuitively, this convergence upper bound is divided into two terms: the first term represents how fast the descent is, and the second term represents how much the convergence process fluctuates. In actual training, we need to choose appropriate $\eta_t$ and $\boldsymbol{H}_t$ to balance these two terms and make them as small as possible.

That is, we assume: by minimizing the convergence upper bound, a better optimizer can be obtained. We have precedents for this approach: the article "Making Alchemy More Scientific (9): Classic Adaptive Gradient Algorithms" based on the assumptions $\boldsymbol{H}_t\succeq \boldsymbol{H}_{t-1}$ and $\eta_t$ being constant, replaced $\Lambda_T$ with its upper bound $\newcommand{tr}{\mathop{\text{tr}}}\tr(\boldsymbol{H}_T)$, and then minimized the upper bound to get $\boldsymbol{H}^* = \left(\sum_{t=1}^T \boldsymbol{g}_t\boldsymbol{g}_t^{\top}\right)^{1/2}$, leading to the AdaGrad optimizer.

And this assumption is the core starting point running through the derivation of this article.

Slowest Descent

Now let's switch perspectives: suppose there is an existing, fixed matrix $\boldsymbol{H}_t$, and we transform it into $\boldsymbol{O}_t^{\top}\boldsymbol{H}_t\boldsymbol{O}_t$, where $\boldsymbol{O}_t$ is an orthogonal matrix. We hope to choose an appropriate $\boldsymbol{O}_t$ such that

\begin{equation}\boldsymbol{\theta}_{t+1} = \boldsymbol{\theta}_t - \eta_t (\boldsymbol{O}_t^{\top}\boldsymbol{H}_t\boldsymbol{O}_t)^{-1}\boldsymbol{g}_t\end{equation}

becomes a better optimizer. How do we judge if it's good? According to the viewpoint in the previous section, we want to reduce the convergence upper bound. Note that orthogonal transformations do not change eigenvalues, so $\lambda_{\max}(\boldsymbol{O}_1^{\top}\boldsymbol{H}_1\boldsymbol{O}_1)=\lambda_{\max}(\boldsymbol{H}_1)$. Assuming $\boldsymbol{O}_t$ varies slowly with $t$, then

\begin{equation}\lambda_{\max}(\boldsymbol{O}_t^{\top}\boldsymbol{H}_t\boldsymbol{O}_t - \boldsymbol{O}_{t-1}^{\top}\boldsymbol{H}_{t-1}\boldsymbol{O}_{t-1}) \approx \lambda_{\max}(\boldsymbol{O}_t^{\top}\boldsymbol{H}_t\boldsymbol{O}_t - \boldsymbol{O}_t^{\top}\boldsymbol{H}_{t-1}\boldsymbol{O}_t) = \lambda_{\max}(\boldsymbol{H}_t - \boldsymbol{H}_{t-1})\end{equation}

In this way, $\Lambda_T$ remains approximately invariant under the substitution $\boldsymbol{H}_t\to \boldsymbol{O}_t^{\top}\boldsymbol{H}_t\boldsymbol{O}_t$. So if we only want to optimize $\boldsymbol{O}_t$, the objective function is reduced to the second term:

\begin{equation}\frac{\sum_{t=1}^T \eta_t^2 \Vert\boldsymbol{g}_t\Vert_{(\boldsymbol{O}_t^{\top}\boldsymbol{H}_t\boldsymbol{O}_t)^{-1}}^2}{2\sum_{t=1}^T \eta_t}\end{equation}

Since $\boldsymbol{O}_t$ for each $t$ can be different, we can minimize $\Vert\boldsymbol{g}_t\Vert_{(\boldsymbol{O}_t^{\top}\boldsymbol{H}_t\boldsymbol{O}_t)^{-1}}^2$ point-by-point, which can be written as:

\begin{equation}\min_{\boldsymbol{O}_t}\langle\boldsymbol{g}_t, (\boldsymbol{O}_t^{\top}\boldsymbol{H}_t\boldsymbol{O}_t)^{-1}\boldsymbol{g}_t\rangle\qquad\text{s.t.}\qquad \boldsymbol{O}_t^{\top}\boldsymbol{O}_t = \boldsymbol{I}\end{equation}

Where $\boldsymbol{\varphi}_t =(\boldsymbol{O}_t^{\top}\boldsymbol{H}_t\boldsymbol{O}_t)^{-1}\boldsymbol{g}_t$ is the update direction of the optimizer. We know $\mathcal{L}(\boldsymbol{\theta} - \eta_t \boldsymbol{\varphi}_t)\approx \mathcal{L}(\boldsymbol{\theta}) - \eta_t\langle\boldsymbol{g}_t,\boldsymbol{\varphi}_t\rangle$. Thus, another meaning of the above formula is to minimize the loss descent per step under a first-order approximation, i.e., the slowest descent.

The Debate Between Fast and Slow

Readers of articles like "Muon Sequel: Why We Chose to Try Muon?", "Steepest Descent on Manifolds: 1. SGD + Hypersphere", and "Above MuP: 2. Linear Layers and Steepest Descent" know that many optimizers are derived based on the principle of steepest descent (fastest descent). Does the "slowest descent" here conflict with that?

In fact, these represent two completely different ways of looking at optimizers. Steepest descent is a first-order approximate, local greedy solution, expecting that taking the optimal step at each moment leads to a globally better solution; it represents an optimistic assumption. Minimizing the convergence upper bound represents a global, pessimistic assumption, hoping that even in the worst-case scenario, we can still achieve certain convergence properties.

The reason why the latter results in the counter-intuitive form of slowest descent is that we have restricted the optimization space to orthogonal transformations, where the first term of the convergence upper bound is approximately negligible, leaving only the second term. This also reminds us that slowest descent shouldn't be considered for every optimization space; this is just a special case under orthogonal space, and the feasibility in other situations needs specific analysis.

The role of slowest descent here is to ensure that while descending, the optimization trajectory is smoother, thus allowing for a larger learning rate, or moving faster and more accurately at the same learning rate. If using the language of steepest descent, one would include assumptions about stability within the constraints of the update amount (as Muon does by choosing the spectral norm).

Furthermore, for stochastic optimization scenarios, this term is also related to noise—decomposing $\boldsymbol{g}_t$ into "true value + random noise," then $\mathbb{E}[\Vert\boldsymbol{g}_t\Vert_{\boldsymbol{H}_t^{-1}}^2]$ will yield a variance term. Therefore, slowest descent also means reducing the impact of noise on the optimization process, making optimization more accurate.

Matrix Optimization

The above sections discussed general principles. Now we formally combine them with matrix optimization. Consider matrix parameters $\boldsymbol{W}\in\mathbb{R}^{n\times m}$ with gradient $\boldsymbol{G}$. The flattened vector is $\boldsymbol{g}=\mathop{\text{vec}}(\boldsymbol{G})\in\mathbb{R}^{nm}$. Here $\boldsymbol{H}$ is an $nm\times nm$ matrix, which is difficult to think about.

For matrix scenarios, we consider preconditioner matrices that can be expressed as Kronecker products. However, rather than using the language of Kronecker products, I intend to understand it directly from a linear operator perspective—treating $\boldsymbol{H}^{-1}\boldsymbol{g}$ directly as:

\begin{equation}\boldsymbol{H}^{-1}\boldsymbol{g} = \mathop{\text{vec}}(\boldsymbol{P}\boldsymbol{G}\boldsymbol{Q})\end{equation}

That is, the optimizer's update rule is $\boldsymbol{W} - \eta\boldsymbol{P}\boldsymbol{G}\boldsymbol{Q}$ (written with Kronecker products as $\boldsymbol{H}^{-1} = \boldsymbol{Q}\otimes\boldsymbol{P}$), where $\boldsymbol{P}\in\mathbb{R}^{n\times n}$ and $\boldsymbol{Q}\in\mathbb{R}^{m\times m}$ are two semi-positive definite preconditioning matrices. $\boldsymbol{P}$ and $\boldsymbol{Q}$ can be related to the current gradient or even historical gradients, so this covers Muon and various Shampoo variants (non-momentum versions).

Now, fixing matrices $\boldsymbol{P}, \boldsymbol{G}, \boldsymbol{Q}$, we hope to choose appropriate orthogonal matrices $\boldsymbol{X}\in\mathbb{R}^{n\times n}$ and $\boldsymbol{Y}\in\mathbb{R}^{m\times m}$ such that the new update $\boldsymbol{X}^{\top}\boldsymbol{P}\boldsymbol{X}\boldsymbol{G}\boldsymbol{Y}^{\top}\boldsymbol{Q}\boldsymbol{Y}$ (written with Kronecker products as $\boldsymbol{O} = \boldsymbol{Y}\otimes\boldsymbol{X}$) becomes a better optimizer. The guiding principle is the "orthogonal slowest descent" discussed earlier:

\begin{equation}\min_{\boldsymbol{X},\boldsymbol{Y}} \tr(\boldsymbol{G}^{\top}\boldsymbol{X}^{\top}\boldsymbol{P}\boldsymbol{X}\boldsymbol{G}\boldsymbol{Y}^{\top}\boldsymbol{Q}\boldsymbol{Y})\qquad\text{s.t.}\qquad \boldsymbol{X}^{\top}\boldsymbol{X} = \boldsymbol{I}_n,\quad \boldsymbol{Y}^{\top}\boldsymbol{Y} = \boldsymbol{I}_m\end{equation}

Preliminary Solution

The complete problem is slightly complex. For simplicity, let's first solve the special case where $\boldsymbol{Q}=\boldsymbol{I}_m$ to get a feel for it. The problem simplifies to:

\begin{equation}\min_{\boldsymbol{X}} \tr(\boldsymbol{G}^{\top}\boldsymbol{X}^{\top}\boldsymbol{P}\boldsymbol{X}\boldsymbol{G})\qquad\text{s.t.}\qquad \boldsymbol{X}^{\top}\boldsymbol{X} = \boldsymbol{I}_n\end{equation}

Noting that $\tr(\boldsymbol{G}^{\top}\boldsymbol{X}^{\top}\boldsymbol{P}\boldsymbol{X}\boldsymbol{G})=\tr(\boldsymbol{G}\boldsymbol{G}^{\top}\boldsymbol{X}^{\top}\boldsymbol{P}\boldsymbol{X})$, this is exactly the von Neumann trace inequality problem derived in "The Rearrangement Inequality and its Generalizations". The minimum value is:

\begin{equation}\sum_{i=1}^n\lambda_i(\boldsymbol{G}\boldsymbol{G}^{\top}) \lambda_{n+1-i}(\boldsymbol{P})\end{equation}

$\lambda_i(\cdot)$ denotes the $i$-th eigenvalue of the corresponding matrix sorted from largest to smallest. Let the SVD of $\boldsymbol{G}$ be $\boldsymbol{U}\boldsymbol{\Sigma}\boldsymbol{V}^{\top}$ (full version, where $\boldsymbol{U}, \boldsymbol{V}$ are square orthogonal matrices and singular values are sorted descending). Then the eigenvalue decomposition of $\boldsymbol{G}\boldsymbol{G}^{\top}$ is $\boldsymbol{U}\boldsymbol{\Sigma}\boldsymbol{\Sigma}^{\top}\boldsymbol{U}^{\top}$. The $\boldsymbol{X}_*$ that achieves the minimum satisfies $\boldsymbol{X}_*^{\top}\boldsymbol{P}\boldsymbol{X}_* = \boldsymbol{U} \boldsymbol{\Lambda}^{\uparrow}_{\boldsymbol{P}}\boldsymbol{U}^{\top}$, where $\boldsymbol{\Lambda}^{\uparrow}_{\boldsymbol{P}}$ is the diagonal matrix formed by sorting the eigenvalues of $\boldsymbol{P}$ from smallest to largest. The optimal update is then:

\begin{equation}\boldsymbol{X}_*^{\top}\boldsymbol{P}\boldsymbol{X}_*\boldsymbol{G} = (\boldsymbol{U} \boldsymbol{\Lambda}^{\uparrow}_{\boldsymbol{P}}\boldsymbol{U}^{\top})( \boldsymbol{U}\boldsymbol{\Sigma}\boldsymbol{V}^{\top}) = \boldsymbol{U} (\boldsymbol{\Lambda}^{\uparrow}_{\boldsymbol{P}}\boldsymbol{\Sigma})\boldsymbol{V}^{\top}\end{equation}

When I first saw this result, I was quite shocked. It tells us that no matter what the original $\boldsymbol{P}$ was, after orthogonal slowest descent, it can at most retain its eigenvalues—the characteristic matrix completely disappears. The final update is $\boldsymbol{G}$'s left and right singular matrices sandwiching a diagonal matrix $\boldsymbol{\Lambda}^{\uparrow}_{\boldsymbol{P}}\boldsymbol{\Sigma}$. This diagonal matrix is the product of the eigenvalues of $\boldsymbol{P}$ and the singular values of $\boldsymbol{G}$ in reverse order. This method makes the overall singular values more uniform, consistent with Muon's philosophy.

This result suggests that no matter how expensive a preconditioner we make, the best result is just to "rush toward Muon." Next, we will see that this conclusion still applies to the general case where $\boldsymbol{Q}\neq\boldsymbol{I}_m$.

General Case

Let $\boldsymbol{A} = \boldsymbol{G}^{\top}\boldsymbol{X}^{\top}\boldsymbol{P}\boldsymbol{X}\boldsymbol{G}$ and $\boldsymbol{B} = \boldsymbol{G}\boldsymbol{Y}^{\top}\boldsymbol{Q}\boldsymbol{Y}\boldsymbol{G}^{\top}$. By the cyclic property of the trace:

\begin{equation}\tr(\boldsymbol{G}^{\top}\boldsymbol{X}^{\top}\boldsymbol{P}\boldsymbol{X}\boldsymbol{G}\boldsymbol{Y}^{\top}\boldsymbol{Q}\boldsymbol{Y}) = \tr(\boldsymbol{A}\boldsymbol{Y}^{\top}\boldsymbol{Q}\boldsymbol{Y}) = \tr(\boldsymbol{B}\boldsymbol{X}^{\top}\boldsymbol{P}\boldsymbol{X})\end{equation}

Let the optimal solutions be $\boldsymbol{X}_*, \boldsymbol{Y}_*$, and their corresponding $\boldsymbol{A}, \boldsymbol{B}$ be $\boldsymbol{A}_*, \boldsymbol{B}_*$. By fixing $\boldsymbol{X}_*$ and solving for $\boldsymbol{Y}_*$, and vice versa, according to the results of the previous section, we can write the stationary conditions:

\begin{equation}\boldsymbol{X}_*^{\top}\boldsymbol{P}\boldsymbol{X}_* = \boldsymbol{U}_{\boldsymbol{B}} \boldsymbol{\Lambda}_{\boldsymbol{P}}^{\uparrow} \boldsymbol{U}_{\boldsymbol{B}}^{\top},\qquad \boldsymbol{Y}_*^{\top}\boldsymbol{Q}\boldsymbol{Y}_* = \boldsymbol{U}_{\boldsymbol{A}} \boldsymbol{\Lambda}_{\boldsymbol{Q}}^{\uparrow} \boldsymbol{U}_{\boldsymbol{A}}^{\top}\end{equation}

where $\boldsymbol{U}_{\boldsymbol{A}}\boldsymbol{\Lambda}_{\boldsymbol{A}}^{\downarrow}\boldsymbol{U}_{\boldsymbol{A}}^{\top}$ and $\boldsymbol{U}_{\boldsymbol{B}}\boldsymbol{\Lambda}_{\boldsymbol{B}}^{\downarrow}\boldsymbol{U}_{\boldsymbol{B}}^{\top}$ are the eigenvalue decompositions of $\boldsymbol{A}_*$ and $\boldsymbol{B}_*$ respectively ("$\downarrow$" indicates descending order). Our task is to analyze the properties of $\boldsymbol{\Phi}=\boldsymbol{X}_*^{\top}\boldsymbol{P}\boldsymbol{X}_*\boldsymbol{G}\boldsymbol{Y}_*^{\top}\boldsymbol{Q}\boldsymbol{Y}_*$. The first conclusion is identical to the previous section: $\boldsymbol{\Phi}$ has the form $\boldsymbol{U}\boldsymbol{D}\boldsymbol{V}^{\top}$, where $\boldsymbol{U}\boldsymbol{\Sigma}\boldsymbol{V}^{\top}$ is the SVD of $\boldsymbol{G}$ and $\boldsymbol{D}$ is some diagonal matrix. The proof is not difficult; calculating directly from the stationary conditions:

\begin{equation}\boldsymbol{G}^{\top}\boldsymbol{\Phi} = \boldsymbol{U}_{\boldsymbol{A}}(\boldsymbol{\Lambda}_{\boldsymbol{A}}^{\downarrow}\boldsymbol{\Lambda}_{\boldsymbol{Q}}^{\uparrow})\boldsymbol{U}_{\boldsymbol{A}}^{\top},\qquad \boldsymbol{\Phi}\boldsymbol{G}^{\top} = \boldsymbol{U}_{\boldsymbol{B}}(\boldsymbol{\Lambda}_{\boldsymbol{B}}^{\downarrow}\boldsymbol{\Lambda}_{\boldsymbol{P}}^{\uparrow})\boldsymbol{U}_{\boldsymbol{B}}^{\top} \label{eq:G-Phi}\end{equation}

This shows that both $\boldsymbol{G}^{\top}\boldsymbol{\Phi}$ and $\boldsymbol{\Phi}\boldsymbol{G}^{\top}$ are symmetric matrices. Since $\boldsymbol{U}, \boldsymbol{V}$ are square orthogonal matrices, any $\boldsymbol{\Phi}$ can be uniquely written as $\boldsymbol{\Phi}=\boldsymbol{U}\boldsymbol{M}\boldsymbol{V}^{\top}$ (i.e., $\boldsymbol{M}=\boldsymbol{U}^{\top}\boldsymbol{\Phi}\boldsymbol{V}$). Substituting this gives $\boldsymbol{G}^{\top}\boldsymbol{\Phi}=\boldsymbol{V}\boldsymbol{\Sigma}^{\top}\boldsymbol{M}\boldsymbol{V}^{\top}$ and $\boldsymbol{\Phi}\boldsymbol{G}^{\top}=\boldsymbol{U}\boldsymbol{M}\boldsymbol{\Sigma}^{\top}\boldsymbol{U}^{\top}$. Their symmetry is equivalent to the symmetry of $\boldsymbol{\Sigma}^{\top}\boldsymbol{M}$ and $\boldsymbol{M}\boldsymbol{\Sigma}^{\top}$, which in component form is:

\begin{equation}\sigma_i M_{i,j} = \sigma_j M_{j,i},\qquad \sigma_j M_{i,j} = \sigma_i M_{j,i}\end{equation}

From this we get $(\sigma_i^2-\sigma_j^2)M_{i,j}=0$. Assuming the singular values of $\boldsymbol{G}$ are distinct (if not, treat it as the limit of a sequence of matrices), we conclude that $M_{i,j}=0$ for $i\neq j$, meaning $\boldsymbol{M}$ is a diagonal matrix $\boldsymbol{D}$. This result shows once again that even in the general case, the characteristic matrix of the preconditioning matrices is eliminated by the orthogonal slowest descent, leaving at most the eigenvalues.

Towards Equilibrium

Regrettably, we cannot easily write a closed-form solution for $\boldsymbol{D}$, but we can qualitatively deduce: the diagonal of $\boldsymbol{D}$ has a tendency toward equilibrium. Specifically, since $\boldsymbol{\Phi} = \boldsymbol{U}\boldsymbol{D}\boldsymbol{V}^{\top}$, we can write:

\begin{equation}\boldsymbol{G}^{\top}\boldsymbol{\Phi} = \boldsymbol{V}(\boldsymbol{\Sigma}^{\top}\boldsymbol{D})\boldsymbol{V}^{\top},\qquad \boldsymbol{\Phi}\boldsymbol{G}^{\top} = \boldsymbol{U}(\boldsymbol{D}\boldsymbol{\Sigma}^{\top})\boldsymbol{U}^{\top}\end{equation}

Comparing with Eq. $\eqref{eq:G-Phi}$, and based on the uniqueness of the eigenvalue decomposition (Note: Uniqueness depends on distinct eigenvalues; otherwise, counterexamples could be constructed, but as this is a heuristic derivation, we won't dwell on details), we can say $\boldsymbol{U}_{\boldsymbol{A}}, \boldsymbol{U}_{\boldsymbol{B}}$ and $\boldsymbol{V}, \boldsymbol{U}$ differ at most by a permutation matrix:

\begin{equation}\boldsymbol{U}_{\boldsymbol{A}} = \boldsymbol{V}\boldsymbol{\Pi}_{\boldsymbol{Q}},\qquad \boldsymbol{U}_{\boldsymbol{B}} = \boldsymbol{U}\boldsymbol{\Pi}_{\boldsymbol{P}}\end{equation}

where $\boldsymbol{\Pi}_{\boldsymbol{Q}}, \boldsymbol{\Pi}_{\boldsymbol{P}}$ are permutation matrices. Then:

\begin{equation}\begin{aligned} \boldsymbol{\Phi} = \boldsymbol{X}_*^{\top}\boldsymbol{P}\boldsymbol{X}_*\boldsymbol{G}\boldsymbol{Y}_*^{\top}\boldsymbol{Q}\boldsymbol{Y}_* = &\, (\boldsymbol{U}_{\boldsymbol{B}} \boldsymbol{\Lambda}_{\boldsymbol{P}}^{\uparrow} \boldsymbol{U}_{\boldsymbol{B}}^{\top}) \boldsymbol{G}(\boldsymbol{U}_{\boldsymbol{A}} \boldsymbol{\Lambda}_{\boldsymbol{Q}}^{\uparrow} \boldsymbol{U}_{\boldsymbol{A}}^{\top}) \\ = &\, (\boldsymbol{U}\boldsymbol{\Pi}_{\boldsymbol{P}} \boldsymbol{\Lambda}_{\boldsymbol{P}}^{\uparrow} \boldsymbol{\Pi}_{\boldsymbol{P}}^{\top}\boldsymbol{U}^{\top}) (\boldsymbol{U}\boldsymbol{\Sigma}\boldsymbol{V}^{\top}) (\boldsymbol{V}\boldsymbol{\Pi}_{\boldsymbol{Q}} \boldsymbol{\Lambda}_{\boldsymbol{Q}}^{\uparrow} \boldsymbol{\Pi}_{\boldsymbol{Q}}^{\top}\boldsymbol{V}^{\top}) \\ = &\, \boldsymbol{U}\Big((\boldsymbol{\Pi}_{\boldsymbol{P}} \boldsymbol{\Lambda}_{\boldsymbol{P}}^{\uparrow} \boldsymbol{\Pi}_{\boldsymbol{P}}^{\top})\boldsymbol{\Sigma}(\boldsymbol{\Pi}_{\boldsymbol{Q}} \boldsymbol{\Lambda}_{\boldsymbol{Q}}^{\uparrow} \boldsymbol{\Pi}_{\boldsymbol{Q}}^{\top})\Big)\boldsymbol{V}^{\top} \end{aligned}\end{equation}

Thus $\tr(\boldsymbol{G}^{\top}\boldsymbol{\Phi})=\tr\big(\boldsymbol{\Sigma}^{\top}(\boldsymbol{\Pi}_{\boldsymbol{P}} \boldsymbol{\Lambda}_{\boldsymbol{P}}^{\uparrow} \boldsymbol{\Pi}_{\boldsymbol{P}}^{\top})\boldsymbol{\Sigma}(\boldsymbol{\Pi}_{\boldsymbol{Q}} \boldsymbol{\Lambda}_{\boldsymbol{Q}}^{\uparrow} \boldsymbol{\Pi}_{\boldsymbol{Q}}^{\top})\big)$. Note that every term inside the trace is a diagonal matrix. To find its minimum, we are essentially looking for a matching relationship between the eigenvalues of $\boldsymbol{P}, \boldsymbol{Q}, \boldsymbol{G}^{\top}\boldsymbol{G}$ such that the sum of their products is minimized. This is essentially a three-sequence generalization of the rearrangement inequality. From "The Rearrangement Inequality and its Generalizations", we know the ideal result is that the component-wise products are equal after matching, which corresponds to $\boldsymbol{D}\propto\boldsymbol{\Sigma}^{-1}$ (taking the inverse of the diagonal elements, preserving shape).

Of course, since we are only seeking the minimum through matching, it's hard to achieve equality for all products. We generally expect it to be an intermediate state between $\boldsymbol{\Sigma}$ (worst case) and $\boldsymbol{\Sigma}^{-1}$ (best case). Their geometric mean is exactly the identity matrix! So this points once again toward the orthogonalization operation that makes all singular values uniform!

Extended Remarks

It must be stated that I have no intention of boosting Muon at the expense of Shampoo. On the contrary, I spent a lot of time testing various Shampoo variants recently, with the initial goal of finding an optimizer better than Muon. Unfortunately, after trying for a long time, I found that no Shampoo variant consistently outperformed Muon. Recently, while deriving convergence results for adaptive optimizers, I combined the two thoughts, resulting in the conclusion of this article.

Frankly, the conclusion of this article relies on the "orthogonal slowest descent" framework, which itself depends on a series of assumptions, each of which can be questioned. For instance, using the convergence upper bound to guide optimizer design is a pessimistic, worst-case criterion; a smaller upper bound does not guarantee faster actual progress. Furthermore, the assumption of slowly changing $\boldsymbol{O}_t$ and the first-order approximations used throughout are quite coarse and cannot be accurately quantified further.

More importantly, although minimizing the convergence upper bound appears more global, it is essentially also a greedy solution. Strictly speaking, the proposition itself might not hold! Because it assumes a given sequence of $\boldsymbol{g}_t$ to derive other parameters. However, if the learning rate or other parameters are modified, subsequent gradients will change; there is no such thing as a fixed, pre-given sequence of all gradients.

However, if we accept these informalities and the principle of orthogonal slowest descent, we obtain the preceding heuristic conclusion. it tells us: we only need to focus on adjusting the spectrum (singular values) of the gradient or momentum, without paying a high price for complex preconditioning matrices. This is at least engineering-friendly. Whether the Shampoo route is effective is one thing, but the engineering cost of its extra state variables and calculations is certainly high.

Summary

Starting from the average convergence upper bound of convex optimization, this article proposes the principle of "orthogonal slowest descent" and applies it to Kronecker product-based preconditioned optimizers (Muon, Shampoo family). We find that regardless of the original preconditioning matrices, orthogonal slowest descent "washes away" their characteristic matrices, leaving at most the eigenvalues to correct the singular value spectrum of the gradient (momentum) to balance it. All of this is highly homogeneous with Muon—as suggested by the title.