Fast Gauss-Newton Methods

• Fast Gauss-Newton (FGN) methods are optimization techniques that approximate Hessian curvature using reduced-cost Jacobian surrogates, preserving key convergence behavior.
• FGN leverages strategies such as local refinement, curvature compression, and structure exploitation to efficiently solve problems in computer vision, inverse problems, and deep learning.
• These methods balance computational speed with accuracy by selectively discarding costly components of the full Gauss-Newton update while maintaining fidelity to the underlying residual map.

Searching arXiv for papers on “Fast Gauss-Newton” and related variants to ground the article. Fast Gauss-Newton (FGN) denotes a family of computationally accelerated Gauss-Newton or generalized Gauss-Newton constructions that preserve Jacobian-based curvature information while reducing the dominant bottleneck of standard second-order updates. In recent literature, the name has been used for several non-equivalent mechanisms: a differentiable single-step depth-refinement layer in multi-view stereo, a gradient-only Gauss-Newton method for PDE-constrained inverse problems, and a curvature approximation for multiclass softmax cross-entropy, among others (Yu et al., 2020, Cherry et al., 6 Mar 2026, Korbit et al., 7 May 2026). This suggests that FGN is best viewed as a design pattern rather than a single canonical algorithm.

1. Terminological scope and major instantiations

The term is applied across computer vision, inverse problems, deep learning, signal processing, and control. What is shared is the attempt to retain the favorable convergence behavior of Gauss-Newton while replacing expensive Jacobian–Hessian machinery, dense linear algebra, or repeated hypothesis searches by lower-cost surrogates.

Setting	Fast mechanism	Representative paper
Multi-view stereo	Single differentiable Gauss-Newton depth-refinement layer after sparse-to-dense propagation	(Yu et al., 2020)
PDE-constrained inverse problems	No additional PDE solves beyond gradient evaluation	(Cherry et al., 6 Mar 2026)
Multiclass cross-entropy	Drop the PSD within-competitor covariance term of the multiclass GGN	(Korbit et al., 7 May 2026)
Motion optimization	Reformulate costs so the GN Hessian captures curvature efficiently	(Ratliff et al., 2016)
Low-rank signal estimation	Modified weighted GN through direct variable projection	(Zvonarev et al., 2021)
Nonlinear optimal control	GN-accelerated PANOC step cast as an LQR and solved by Riccati recursion	(Pas et al., 2022)

A plausible implication is that the adjective “fast” is relative to the ambient bottleneck of each field. In some papers it means fewer PDE solves, in others one-step local refinement, structured linear algebra, reduced class dependence, or lower memory footprint.

2. Shared computational pattern

Classical Gauss-Newton replaces the full Hessian of a nonlinear least-squares objective by a Jacobian Gram matrix. In the PDE-constrained formulation, for example, standard Gauss-Newton uses

$$H_k^{\rm GN}=\sum_{i=1}^N J_i(m_k)^\top J_i(m_k)+\nabla^2 R(m_k),$$

with step

$p_k^{\rm GN}=-(H_k^{\rm GN})^{-1}\nabla F(m_k),$

followed by globalization through line search or trust region (Cherry et al., 6 Mar 2026). Fast variants typically leave this outer logic intact but alter how the curvature is represented, assembled, or solved.

One recurring strategy is local refinement around a good initialization. In Fast-MVSNet, the Gauss-Newton layer is applied only after sparse high-resolution depth prediction and learned propagation, so the update is a local correction rather than a global depth search (Yu et al., 2020). Another strategy is curvature compression: the multiclass cross-entropy FGN keeps only the true-vs-rest curvature term and discards a PSD covariance term over competitors (Korbit et al., 7 May 2026). A third strategy is reuse of already-paid-for information: the gradient-only PDE method reconstructs a Gauss-Newton Jacobian from per-source losses and gradients already computed during ordinary gradient evaluation (Cherry et al., 6 Mar 2026). A fourth is structure exploitation: block-diagonal Kronecker factorizations, Sherman-Morrison-Woodbury reductions, hierarchical-matrix compression, direct variable projection, and Riccati recursions all replace generic dense solves by domain-specific linear algebra (Botev et al., 2017, Ren et al., 2019, Chen et al., 2019, Zvonarev et al., 2021, Pas et al., 2022).

These constructions are fast for different reasons, but they all preserve the central Gauss-Newton principle: curvature is approximated from first-order sensitivity information rather than from full second derivatives.

3. Representative mathematical constructions

In multiclass softmax cross-entropy, the full generalized Gauss-Newton matrix admits the exact decomposition

$$H^{\mathrm{GGN}} = H^{\mathrm{FGN}} + p_\dagger\, J_{-\star}^\top \bigl(\operatorname{diag}(\rho)-\rho\rho^\top\bigr) J_{-\star},$$

with

$$H^{\mathrm{FGN}} = p_\star p_\dagger\, J_s^\top J_s.$$

Here the loss and gradient are unchanged, and the approximation enters only at the curvature level. The discarded term is PSD, so $H^{\mathrm{FGN}}$ is a PSD under-approximation of the multiclass GGN, and it is exact for binary classification (Korbit et al., 7 May 2026).

In large-scale PDE-constrained inverse problems, the fast construction starts from

$\Phi(m)=\frac12\sum_{i=1}^N \|r_i(m)\|^2$

and introduces scalar residual norms $\rho_i(m)=\|r_i(m)\|=\sqrt{2\phi_i(m)}$. Because

$$\nabla \rho_i(m)=\frac{\nabla \phi_i(m)}{\sqrt{2\phi_i(m)}},$$

the Jacobian

$$J^{\rm GO}(m)= \begin{bmatrix} \nabla \rho_1(m)^\top\ \vdots\ \nabla \rho_N(m)^\top \end{bmatrix}$$

can be built from quantities already available from gradient computation. The resulting approximation

$$H_k^{\rm GO}=J^{\rm GO}(m_k)^\top J^{\rm GO}(m_k)+\nabla^2R(m_k)$$

eliminates the extra PDE solves that standard GN-CG would require for Jacobian-vector products (Cherry et al., 6 Mar 2026).

In Fast-MVSNet, the final refinement stage minimizes a multi-view feature discrepancy

$p_k^{\rm GN}=-(H_k^{\rm GN})^{-1}\nabla F(m_k),$0

with residuals $p_k^{\rm GN}=-(H_k^{\rm GN})^{-1}\nabla F(m_k),$1, Jacobians

$p_k^{\rm GN}=-(H_k^{\rm GN})^{-1}\nabla F(m_k),$2

and single-step Gauss-Newton increment

$p_k^{\rm GN}=-(H_k^{\rm GN})^{-1}\nabla F(m_k),$3

The layer is differentiable, avoids hypothesis sampling, and in practice converges quickly with only a single step of update (Yu et al., 2020).

These examples show that FGN is often loss-specific or model-specific. The acceleration is not usually a generic optimizer wrapper; it is encoded in the algebra of the residual map.

4. Structure-exploiting accelerations across domains

In motion optimization, the key result is asymptotic rather than purely heuristic: for a broad class of discrete-time trajectory objectives,

$p_k^{\rm GN}=-(H_k^{\rm GN})^{-1}\nabla F(m_k),$4

with rate $p_k^{\rm GN}=-(H_k^{\rm GN})^{-1}\nabla F(m_k),$5, where $p_k^{\rm GN}=-(H_k^{\rm GN})^{-1}\nabla F(m_k),$6 (Ratliff et al., 2016). This motivates reformulating rigid-body kinetic energy as a squared velocity norm in a 12D Euclidean task space,

$p_k^{\rm GN}=-(H_k^{\rm GN})^{-1}\nabla F(m_k),$7

so that $p_k^{\rm GN}=-(H_k^{\rm GN})^{-1}\nabla F(m_k),$8. The stated benefit is fast and accurate curvature calculation requiring only as much computation as a single evaluation of the traditional representation (Ratliff et al., 2016).

In deep learning, several fast GN lines coexist. One practical approach uses a block-diagonal approximation to the Gauss-Newton matrix for feedforward networks, with Kronecker-factored blocks that make per-layer inversion cheap and competitive with first-order methods (Botev et al., 2017). A different line gives entry-wise evaluation of the Gauss-Newton Hessian for multilayer perceptrons with a precomputation step and an importance-sampling step, reducing the cost of a sampled entry from $p_k^{\rm GN}=-(H_k^{\rm GN})^{-1}\nabla F(m_k),$9 to $$H^{\mathrm{GGN}} = H^{\mathrm{FGN}} + p_\dagger\, J_{-\star}^\top \bigl(\operatorname{diag}(\rho)-\rho\rho^\top\bigr) J_{-\star},$$0, and enabling hierarchical-matrix approximations (Chen et al., 2019). A third line uses subsampled Gauss-Newton or Fisher matrices together with Levenberg-Marquardt damping and Sherman-Morrison-Woodbury reductions, with convergence to a stationary point proved for the semi-stochastic full-gradient variant (Ren et al., 2019).

In low-rank signal estimation, the modified weighted Gauss-Newton method is implemented through direct variable projection onto a space of low-rank signals, with asymptotic iteration cost

$$H^{\mathrm{GGN}} = H^{\mathrm{FGN}} + p_\dagger\, J_{-\star}^\top \bigl(\operatorname{diag}(\rho)-\rho\rho^\top\bigr) J_{-\star},$$1

and it can handle missing values without increasing the computational cost (Zvonarev et al., 2021). In nonlinear optimal control, Gauss-Newton directions are embedded into PANOC; for optimal control problems the GN step becomes an equality-constrained finite-horizon LQR solved by Riccati recursion, and the reported benchmark showed the method was more than twice as fast as the traditional L-BFGS variant of PANOC (Pas et al., 2022).

5. Exactness, convergence, and accuracy regimes

FGN methods are not uniformly approximate in the same sense. Some are exact in special regimes. The multiclass cross-entropy construction is exact for binary classification because the conditional competitor covariance vanishes when there is only one competitor (Korbit et al., 7 May 2026). In motion optimization, the claim is asymptotic exactness: as time discretization becomes small, the true Hessian limits to the Gauss-Newton Hessian quickly for the covered class of objective terms (Ratliff et al., 2016).

Other variants are explicitly local. The q-Gauss-Newton method for unconstrained nonlinear optimization is presented as much faster than q-steepest descent algorithms because it approximates the matrix of q-second order derivatives with the first-order q-Jacobian matrix, but its convergence is assured only when the initial guess is close enough to the solution (Protic et al., 2021). In phase retrieval, the Gauss-Newton algorithm combines a modified spectral initialization with re-sampled local updates and proves quadratic convergence in the real case under Gaussian measurements; the number of iterations required to reach accuracy $$H^{\mathrm{GGN}} = H^{\mathrm{FGN}} + p_\dagger\, J_{-\star}^\top \bigl(\operatorname{diag}(\rho)-\rho\rho^\top\bigr) J_{-\star},$$2 is $$H^{\mathrm{GGN}} = H^{\mathrm{FGN}} + p_\dagger\, J_{-\star}^\top \bigl(\operatorname{diag}(\rho)-\rho\rho^\top\bigr) J_{-\star},$$3 under the stated sampling assumptions (Gao et al., 2016).

Still other fast GN methods emphasize globalization. The derivative-free DFO-GN method replaces quadratic residual models by linear interpolation models inside a trust-region framework and retains $$H^{\mathrm{GGN}} = H^{\mathrm{FGN}} + p_\dagger\, J_{-\star}^\top \bigl(\operatorname{diag}(\rho)-\rho\rho^\top\bigr) J_{-\star},$$4 worst-case complexity to approximate first-order stationarity (Cartis et al., 2017). The flexible regularized Gauss-Newton family for solving nonlinear equations proves global convergence bounds and, under non-degeneracy assumptions, local quadratic convergence; its stochastic extension attains sub-linear or linear convergence under weak growth and Polyak-Lojasiewicz conditions (Yudin et al., 2021).

The literature therefore does not support a single universal convergence statement for FGN. Theoretical behavior ranges from exactness in special low-dimensional regimes, to asymptotic exactness, to local quadratic convergence, to globally convergent trust-region variants.

6. Trade-offs, limitations, and recurring misconceptions

A common misconception is that “fast” means a universally stronger optimizer. The papers show a more qualified picture. The multiclass softmax variant deliberately drops a PSD part of the curvature and is therefore an under-approximation of the full GGN; it stays closest to the full softmax GGN when competitor mass is concentrated or damping is large, and deviates as the dropped within-competitor covariance grows (Korbit et al., 7 May 2026). The PDE gradient-only method captures only the curvature in the data-space subspace spanned by the $$H^{\mathrm{GGN}} = H^{\mathrm{FGN}} + p_\dagger\, J_{-\star}^\top \bigl(\operatorname{diag}(\rho)-\rho\rho^\top\bigr) J_{-\star},$$5 scalar residual norms, and its effectiveness depends on the rank and conditioning of $$H^{\mathrm{GGN}} = H^{\mathrm{FGN}} + p_\dagger\, J_{-\star}^\top \bigl(\operatorname{diag}(\rho)-\rho\rho^\top\bigr) J_{-\star},$$6 (Cherry et al., 6 Mar 2026).

A second misconception is that FGN necessarily removes the need for a good initial estimate. In fact, several prominent variants remain strongly local. Fast-MVSNet applies Gauss-Newton only after sparse high-resolution depth prediction and learned propagation (Yu et al., 2020). The q-Gauss-Newton method requires an initial guess close enough to the optimal solution (Protic et al., 2021). The phase-retrieval algorithm proves quadratic convergence only after a modified spectral initializer places the iterate inside the local basin, and its proof further relies on resampling to restore independence between iterates and measurements (Gao et al., 2016).

A third misconception is that FGN always avoids expensive automatic-differentiation primitives. Some variants do, but others merely relocate the cost. The multiclass softmax method solves a row-space system matrix-free by conjugate gradient using JVPs and VJPs of the scalar margin map (Korbit et al., 7 May 2026). Subsampled neural-network GN methods still require backpropagation-based Jacobian products and damped linear solves (Ren et al., 2019). Optimal-control FGN remains efficient because the GN step is recast as an LQR and solved by Riccati recursion, not because second-order structure disappears (Pas et al., 2022).

The broad picture is therefore heterogeneous. Fast Gauss-Newton methods preserve the favorable geometry of Gauss-Newton only insofar as a particular reformulation, approximation, or structural reduction remains faithful to the underlying residual map. Where that fidelity is high, FGN can behave like a practical second-order method; where it is weak, the method can revert toward first-order behavior or depend heavily on damping, conditioning, and initialization.